Re: [eigen] discrepancy in (triangular view + transpose) and self-adjointness
• To: eigen@xxxxxxxxxxxxxxxxxxx
• Subject: Re: [eigen] discrepancy in (triangular view + transpose) and self-adjointness
• From: Manoj Rajagopalan <rmanoj@xxxxxxxxx>
• Date: Wed, 14 Jul 2010 19:41:58 -0400
• Organization: EECS Dept., University of Michigan, Ann Arbor, MI, USA

```Attached compilable C++ program.

thanks,
Manoj

On Wednesday 14 July 2010 07:14:02 pm Benoit Jacob wrote:
> Hi,
>
> It would help if you could make a compilable test program, so we can
> quickly run stuff.
>
> Benoit
>
> 2010/7/14 Manoj Rajagopalan <rmanoj@xxxxxxxxx>:
> > Hi all,
> >
> > discrepancy?
> >
> >   Given:
> >        int const N = 5;
> >        MatrixXi A(N,N);
> >        A.setRandom();
> >
> >        TriangularView<MatrixXi,Lower> const A_l(A);
> >        TriangularView<MatrixXi,Lower|ZeroDiag> const A_zl(A);
> >
> > the following must be true:
> >
> >   A.selfadjointView<Lower>() == A_l + A_zl.transpose()
> >
> > Basically, I am adding the lower triangle of A to its own transpose,
> > without adding the diagonal a second time. So the following is calculated
> > correctly by Eigen:
> >
> >    MatrixXi C1(N,N), C2(N,N);
> >    C2 = A_l.toDenseMatrix() + A_zl.transpose().toDenseMatrix();
> >    assert(C1 == C2);
> >
> > But, the following code produces a discrepancy:
> >
> >    MatrixXi B(N,N);
> >    B.setRandom();
> >    C1 = A.selfadjointView<Lower>() * B;
> >    C2 = (A_l * B)   +   (A_zl.transpose() * B);
> >    assert(C1 == C2); // <--------  FALSE!
> >
> >     After a long day my brain is unable to figure out if I am doing
> > something extremely silly - it seems too trivial to be a bug. Can someone
> >
> > thanks,
> > Manoj

```
```#include <Eigen/Core>
#include <iostream>

using namespace Eigen;
using namespace std;

int main(void)
{
int const N = 5;
MatrixXi A(N,N), B(N,N), C1(N,N), C2(N,N), C_diff(N,N);
A.setRandom();
B.setRandom();
// C1 = A.selfadjointView<Lower>() * B;

TriangularView<MatrixXi,Lower> const A_l(A);
TriangularView<MatrixXi,Lower|ZeroDiag> const A_zl(A);
// C2 = A_l * B;
// C2 += A_zl.transpose() * B;
C2 = A_l.toDenseMatrix() + A_zl.transpose().toDenseMatrix();
C_diff = C1 - C2;

cout << "A =\n" << A << endl;
cout << "\nA.triangularView<Lower>() =\n" << A_l.toDenseMatrix() << endl;
cout << "\nA.triangularView<Lower|ZeroDiag>() =\n" << A_zl.toDenseMatrix() << endl;
cout << "C1 =\n" << C1 << endl;
cout << "C2 =\n" << C2 << endl;
cout << "C_diff =\n" << C_diff << endl;

assert(C1 == C2);
return 0;
}

```

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