[vtkusers] Jacobi eigenvectors incorrect?

[Bill Lorensen]

I know zip about matlab. Bryn seems to have success. Maybe he can help.

Bill

On Fri, Oct 30, 2009 at 8:02 AM, David Doria <daviddoria+vtk at gmail.com> wrote:
> On Wed, Oct 28, 2009 at 5:57 AM, Fred Fred <stan1313 at hotmail.fr> wrote:
>> BTW, which is the easiest way to compute the real eigenvalues of a
>> non-symetric matrix in VTK? Of course such eigenvalues may not exist but in
>> case they exist, I wonder which is the better way to use.
>>
>>> Date: Tue, 27 Oct 2009 08:39:34 -0400
>>> From: bill.lorensen at gmail.com
>>> To: daviddoria+vtk at gmail.com
>>> CC: vtkusers at vtk.org
>>> Subject: Re: [vtkusers] Jacobi eigenvectors incorrect?
>>>
>>> According to the documentation in vtkMath.cxx, the code
>>> vtkJacobiN finds the solution of an nxn real SYMMETRIC matrix.
>>> Can you try your tests with a symmetric matrix and compare the results
>>> between vtk and matlab.
>>>
>>> Bill
>
> I am still not getting the same evecs even with a symmetric matrix:
>
> #include "vtkMath.h"
>
> int main (int argc, char *argv[])
> {
>  double *a[3];
>  double a0[3];
>  double a1[3];
>  double a2[3];
>  a[0] = a0; a[1] = a1; a[2] = a2;
>
>  //set the matrix to all zeros
>  for(unsigned int i = 0; i < 3; i++)
>  {
>    a0[i] = a1[i] = a2[i] = 0.0;
>  }
>
>  /*
>  1 2 3
>  2 5 9
>  3 9 8
>  */
>  a0[0] = 1;
>  a0[1] = 2;
>  a0[2] = 3;
>  a1[0] = 2;
>  a1[1] = 5;
>  a1[2] = 9;
>  a2[0] = 3;
>  a2[1] = 9;
>  a2[2] = 8;
>
>  // Extract eigenvectors from covariance matrix
>  double *v[3];
>  double v0[3];
>  double v1[3];
>  double v2[3];
>
>  v[0] = v0; v[1] = v1; v[2] = v2;
>  double eigval[3];
>  vtkMath::Jacobi(a,eigval,v);
>  vtkstd::cout << "eigvals: " << eigval[0] << " " << eigval[1] << " "
> << eigval[2] << vtkstd::endl;
>  vtkstd::cout << "v0: " << v0[0] << " " << v0[1] << " " << v0[2] <<
> vtkstd::endl;
>  vtkstd::cout << "v1: " << v1[0] << " " << v1[1] << " " << v1[2] <<
> vtkstd::endl;
>  vtkstd::cout << "v2: " << v2[0] << " " << v2[1] << " " << v2[2] <<
> vtkstd::endl;
>
>  vtkstd::cout << "Evec corresponding to smallest eval: " << v2[0] <<
> " " << v2[1] << " " << v2[2] << vtkstd::endl;
>  return 0;
> }
>
> Here is the matlab:
>
>
> a=[1 2 3; 2 5 9; 3 9 8];
> [v d] = eig(a);
>
> The columns of V are the eigenvectors of a:
>
> v =
>
>    0.1372    0.9644    0.2259
>    0.7382   -0.2516    0.6259
>   -0.6605   -0.0809    0.7464
>
> The output of my vtk c++:
> Evec corresponding to smallest eval: 0.746449 0.0808535 -0.660513
>
> I see this vector reversed in the bottom row of the matlab eigenvector
> matrix - is this a coincidence or am I misinterpreting the output from
> the Jacobi function and just getting these things in the wrong order?
> This is exactly what I was talking about when I suggested that we
> either reimplement or make a nice wrapper so that users can be sure
> what is going on - setting (r,c) type entries in a "matrix" and then
> getting "vectors" back instead of ** style arrays.
>
> Thanks,
>
> David
> _______________________________________________
> Powered by www.kitware.com
>
> Visit other Kitware open-source projects at http://www.kitware.com/opensource/opensource.html
>
> Please keep messages on-topic and check the VTK FAQ at: http://www.vtk.org/Wiki/VTK_FAQ
>
> Follow this link to subscribe/unsubscribe:
> http://www.vtk.org/mailman/listinfo/vtkusers
>