Dear Dmitry,
I successfully found the eigenvalues of your matrix using GSL 2.6 (see
attached code). The output of the program is:
---
$ ./eig
eval[1] = 0.000000000000e+00 + i 0.000000000000e+00
eval[2] = 0.000000000000e+00 + i 0.000000000000e+00
eval[3] = 0.000000000000e+00 + i 0.000000000000e+00
eval[4] = 4.300443739734e-16 + i 0.000000000000e+00
eval[5] = -6.000000000000e+00 + i 0.000000000000e+00
eval[6] = 6.000000000000e+00 + i 0.000000000000e+00
eval[7] = 6.000000000000e+00 + i 0.000000000000e+00
eval[8] = -6.000000000000e+00 + i 0.000000000000e+00
eval[9] = 5.377642775528e-17 + i 0.000000000000e+00
eval[10] = 3.000000000000e+00 + i 0.000000000000e+00
eval[11] = -3.000000000000e+00 + i 0.000000000000e+00
eval[12] = 6.000000000000e+00 + i 0.000000000000e+00
eval[13] = -6.000000000000e+00 + i 0.000000000000e+00
eval[14] = 3.000000000000e+00 + i 0.000000000000e+00
eval[15] = -2.912706738243e-16 + i 0.000000000000e+00
eval[16] = 6.000000000000e+00 + i 0.000000000000e+00
eval[17] = -6.000000000000e+00 + i 0.000000000000e+00
eval[18] = -3.000000000000e+00 + i 0.000000000000e+00
---
Can you clarify which version of GSL you are using? Also you might want
to print out your matrix m to make sure it is initialized correctly. I
did not try to debug your code.
Best,
Patrick
On 12/29/20 7:44 AM, Dmitry Cheshkov wrote:
> Dear Sirs!
>
> I have to calculate the eigenvalues of the following 18x18 real
> non-symmetric square matrix:
>
> 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0
> 0 0 0 0 0 0 0 0 0 0 0 0 0 2 0 0 0 0
> 0 0 0 0 0 0 0 0 0 0 0 0 2 0 1 0 0 0
> 0 0 0 0 0 0 0 0 0 0 0 3 0 0 0 0 0 3
> 0 0 0 0 0 0 0 0 0 0 6 0 0 0 0 0 6 0
> 0 0 0 0 0 0 0 0 0 6 0 3 0 0 0 6 0 3
> 0 0 0 0 0 0 0 0 0 0 0 0 0 0 2 0 0 0
> 0 0 0 0 0 0 0 0 0 0 0 0 0 4 0 0 0 0
> 0 0 0 0 0 0 0 0 0 0 0 0 4 0 2 0 0 0
> 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
> 0 0 0 0 2 0 0 0 0 0 0 0 0 0 0 0 0 0
> 0 0 0 2 0 1 0 0 0 0 0 0 0 0 0 0 0 0
> 0 0 3 0 0 0 0 0 3 0 0 0 0 0 0 0 0 0
> 0 6 0 0 0 0 0 6 0 0 0 0 0 0 0 0 0 0
> 6 0 3 0 0 0 6 0 3 0 0 0 0 0 0 0 0 0
> 0 0 0 0 0 2 0 0 0 0 0 0 0 0 0 0 0 0
> 0 0 0 0 4 0 0 0 0 0 0 0 0 0 0 0 0 0
> 0 0 0 4 0 2 0 0 0 0 0 0 0 0 0 0 0 0
>
> I wrote the folowing code:
>
> #include <iostream>
> #include <fstream>
> #include <gsl/gsl_math.h>
> #include <gsl/gsl_eigen.h>
>
> using namespace std;
>
> int main(void)
> {
>
> ifstream istr("matrix");
> int N = 18;
> gsl_matrix* m = gsl_matrix_alloc(N, N);
> gsl_vector_complex* eval = gsl_vector_complex_alloc(N);
> gsl_eigen_nonsymm_workspace* w = gsl_eigen_nonsymm_alloc(N);
> for (int i = 0; i < N; i++)
> for (int j = 0; j < N; j++)
> istr >> m->data[N * j + i];
> gsl_eigen_nonsymm(m, eval, w);
> for (int i = 0; i < N; i++)
> cout << eval->data[2 * i] << '\t' << eval->data[2 * i + 1] << endl;
> cout << endl;
> return 0;
>
> }
>
> If the file "matrix" contains the above mentioned matrix, I have the
> following error:
>
> gsl: francis.c:209: ERROR: maximum iterations reached without finding
> all eigenvalues
> Default GSL error handler invoked.
>
> While Wolfram Mathematica works well:
> In[3]:= Eigenvalues[matrix]
> Out[3]= {-6, -6, -6, -6, 6, 6, 6, 6, -3, -3, 3, 3, 0, 0, 0, 0, 0, 0}
>
> Is it possible to calculate this eigenvalues using GSL library and how?
>
>
> With best regards,
> Dmitry Cheshkov
>
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <gsl/gsl_math.h>
#include <gsl/gsl_eigen.h>
#include <gsl/gsl_vector.h>
#include <gsl/gsl_matrix.h>
int
main()
{
const double data[] = {
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 0, 1, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 3, 0, 0, 0, 0, 0, 3,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 6, 0, 0, 0, 0, 0, 6, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 6, 0, 3, 0, 0, 0, 6, 0, 3,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 4, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 4, 0, 2, 0, 0, 0,
0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 2, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 3, 0, 0, 0, 0, 0, 3, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 6, 0, 0, 0, 0, 0, 6, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
6, 0, 3, 0, 0, 0, 6, 0, 3, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 4, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 4, 0, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 };
const size_t N = 18;
gsl_matrix_const_view m = gsl_matrix_const_view_array(data, N, N);
gsl_matrix * A = gsl_matrix_alloc(N, N);
gsl_vector_complex * eval = gsl_vector_complex_alloc(N);
gsl_eigen_nonsymm_workspace * w = gsl_eigen_nonsymm_alloc(N);
size_t i;
gsl_matrix_memcpy(A, &m.matrix);
gsl_eigen_nonsymm(A, eval, w);
for (i = 0; i < N; ++i)
{
gsl_complex ei = gsl_vector_complex_get(eval, i);
printf("eval[%zu] = %.12e + i %.12e\n",
i + 1, GSL_REAL(ei), GSL_IMAG(ei));
}
gsl_matrix_free(A);
gsl_vector_complex_free(eval);
gsl_eigen_nonsymm_free(w);
return 0;
}
