GSL eigenvalues order

I am using the functions gsl_eigen_nonsymm and/or gsl_eigen_symm from the GSL library to find the eigenvalues of an L x L matrix M[i][j] which is also a function of time t = 1,....,N so i have M[i][j][t] to get the eigenvalues for every t i allocate an L x L matrix E[i][j] = M[i][j][t] and diagonalize it for every t.

The problem is that the program gives the eigenvalues at different order after some iteration. For example (L = 3) if at t = 0 i get eigen[t = 0] = {l1,l2,l3}(0) at t = 1 i may get eigen[t = 1] = {l3,l2,l1}(1) while i need to always have {l1,l2,l3}(t) To be more concrete: consider the matrix M (t) ) = {{0,t,t},{t,0,2t},{t,t,0}} the eigenvalues will always be (approximatevly) l1 = -1.3 t , l2 = -t , l3 = 2.3 t When i tried to diagonalize it (with the code below) i got several times a swap in the result of the eigenvalues. Is there a way to prevent it? I can't just sort them by magnitude i need them to be always in the same order (whatever it is) a priori. (the code below is just an example to elucidate my problem)

EDIT: I can't just sort them because a priori I don't know their value nor if they reliably have a structure like l1<l2<l3 at every time due to statistical fluctuations, that's why i wanted to know if there is a way to make the algorithm behave always in the same way so that the order of the eigenvalues is always the same or if there is some trick to make it happen.

Just to be clearer I'll try to re-describe the toy problem I presented here. We have a matrix that depends on time, I, maybe naively, expected to just get lambda_1(t).....lambda_N(t), instead what I see is that the algorithm often swaps the eigenvalues at different times, so if at t = 1 I've got ( lambda_1,lambda_2,lambda_3 )(1) at time t = 2 (lambda_2,lambda_1,lambda_3)(2) so if for instance I wanted to see how lambda_1 evolves in time I can't because the algorithm mixes the eigenvalues at different times. The program below is just an analytical-toy example of my problem: The eigenvalues of the matrix below are l1 = -1.3 t , l2 = -t , l3 = 2.3 t but the program may give me as an output(-1.3,-1,2.3)(1), (-2,-2.6,4.6)(2), etc As previously stated, I was wondering then if there is a way to make the program order the eigenvalues always in the same way despite of their actual numerical value, so that i always get the (l1,l2,l3) combination. I hope it is more clear now, please tell me if it is not.

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <gsl/gsl_linalg.h>
#include <gsl/gsl_eigen.h>
#include <gsl/gsl_sort_vector.h>

main() {
    int L = 3, i, j, t;
    int N = 10;
    double M[L][L][N];
    gsl_matrix *E = gsl_matrix_alloc(L, L);
    gsl_vector_complex *eigen = gsl_vector_complex_alloc(L);
    gsl_eigen_nonsymm_workspace *  w = gsl_eigen_nonsymm_alloc(L);

    for(t = 1; t <= N; t++) {
        M[0][0][t-1] = 0;
        M[0][1][t-1] = t;
        M[0][2][t-1] = t;
        M[1][0][t-1] = t;
        M[1][1][t-1] = 0;
        M[1][2][t-1] = 2.0 * t;
        M[2][1][t-1] = t;
        M[2][0][t-1] = t;
        M[2][2][t-1] = 0;

        for(i = 0; i < L; i++) {
            for(j = 0; j < L; j++) {
                gsl_matrix_set(E, i, j, M[i][j][t - 1]);
            }
        }

        gsl_eigen_nonsymm(E, eigen, w); /*diagonalize E which is M at t fixed*/
        printf("#%d\n\n", t);

        for(i = 0; i < L; i++) {
            printf("%d\t%lf\n", i, GSL_REAL(gsl_vector_complex_get(eigen, i)))
        }
        printf("\n");
   }
}

Solution

I don't think you can do what you want. As t changes the output changes.

My original answer mentioned ordering on the pointers but looking at the data structure it won't help. When the eigenvalues have been computed the values are stored in E. You can see them as follows.

gsl_eigen_nonsymm(E, eigen, w);
double *mdata = (double*)E->data;
printf("mdata[%i]    \t%lf\n", 0, mdata[0]);
printf("mdata[%i]    \t%lf\n", 4, mdata[4]);
printf("mdata[%i]    \t%lf\n", 8, mdata[8]);

The following code is how the the data in the eigenvector is layed out...

double *data  = (double*)eigen->data;
for(i = 0; i < L; i++) {
  printf("%d n  \t%zu\n", i, eigen->size);
  printf("%d    \t%lf\n", i, GSL_REAL(gsl_vector_complex_get(eigen, i)));
  printf("%d r  \t%lf\n", i, data[0]);
  printf("%d i  \t%lf\n", i, data[1]);
  printf("%d r  \t%lf\n", i, data[2]);
  printf("%d i  \t%lf\n", i, data[3]);
  printf("%d r  \t%lf\n", i, data[4]);
  printf("%d i  \t%lf\n", i, data[5]);
}

If, and you can check this when you see the order change, the order of the data in mdata changes AND the order in data changes then the algorithm does not have a fixed order ie you cannot do what you're asking it to do. If the order does not change in mdata and it changes in data then you have a solution but I really doubt that will be the case.