Re: [eigen] Re: efficient trace of product of symmetric matrices |

[ Thread Index | Date Index | More lists.tuxfamily.org/eigen Archives ]

*To*: eigen@xxxxxxxxxxxxxxxxxxx*Subject*: Re: [eigen] Re: efficient trace of product of symmetric matrices*From*: Gael Guennebaud <gael.guennebaud@xxxxxxxxx>*Date*: Tue, 21 Feb 2012 23:20:16 +0100*Authentication-results*: mr.google.com; spf=pass (google.com: domain of gael.guennebaud@xxxxxxxxx designates 10.42.162.194 as permitted sender) smtp.mail=gael.guennebaud@xxxxxxxxx; dkim=pass header.i=gael.guennebaud@xxxxxxxxx*Dkim-signature*: v=1; a=rsa-sha256; c=relaxed/relaxed; d=gmail.com; s=gamma; h=mime-version:in-reply-to:references:from:date:message-id:subject:to :content-type:content-transfer-encoding; bh=wXxI7nzY+UC6mAKH7Vmu96fsppUkoEZHo1pNnG16hv8=; b=D19+onOun4VIJHY5VU8uXoDP5HOQl1CBf2Y7t5FIKJ87Hmjt99eIUfbGJB0/ebDvnI HDOSZdUJw8sEHMrNVGp82+rRc5I98pHV+uVMG9nethVXWoMlVAkvuuXEqC9zQ5rL5fb4 hnfoH5XLKeJZmrfxyRZ6Nix4uozUQ2QTylBhw=

Hi, On Mon, Feb 20, 2012 at 11:28 PM, Patrick <pat05.mlst@xxxxxxxxx> wrote: > Hi everyone, > > I'm wondering whether it is currently possible to compute > > tr(A*B), where A \in Sym(n,n), B \in Sym(n,n) > > efficiently using the development branch of Eigen3. > Since A and B are symmetric and > > tr(A*B) = sum_ij a_ij * b_ij, > > the best way I can think of doing this is > > double tr=0; > for(size_t i=0; i<A.rows(); ++i) > { > tr += A(i,i) * B(i,i); > for(size_t j=0; j<i; ++j) > tr += 2 * A(i,j)*B(i,j); > } > > I noticed that the triangularView methods do not expose the array > interface, so am I correct in assuming I cannot vectorise this using > triangular views? what you need here is a lazy triangular * triangular product, and indeed there is no such thing yet. Nevertheless you can still write your inner loop using an Eigen expression and get vectorization if A is row major and B column major. > Furthermore, when trying to make a matrix symmetric by calling > > A.triangularView<StrictlyUpper>() = A.triangularView<StrictlyLower>(), > > I noticed that no assignment was taking place, i.e. the matrix A remained as > is. > Is this the expected behaviour? Yes that's the intended behavior. You have to transpose the right hand side: A.triangularView<StrictlyUpper>() = A.triangularView<StrictlyLower>().transpose(), Actually, you can even do: A.triangularView<StrictlyUpper>() = A.transpose(); The triangularView on the lhs behave like a writing mask. cheers, gael. > > Thanks. > > Kind regards, > > Patrick

**References**:

**Messages sorted by:**[ date | thread ]- Prev by Date:
**[eigen] Re: efficient trace of product of symmetric matrices** - Next by Date:
**Re: [eigen] Creating a EulerAngle object** - Previous by thread:
**[eigen] Re: efficient trace of product of symmetric matrices** - Next by thread:
**Re: [eigen] Status of AVX support**

Mail converted by MHonArc 2.6.19+ | http://listengine.tuxfamily.org/ |