GSL Eigenwerte bestellen

Question

ich die Funktionen bin mit gsl_eigen_nonsymm und/oder gsl_eigen_symm aus der GSL-Bibliothek die Eigenwerte einer L x L Matrix M[i][j] die auch t = 1,....,N eine Funktion der Zeit zu finden, so dass ich M[i][j][t] haben die Eigenwerte für jeden ti zu bekommen ordnen Sie eine L x L-Matrix E[i][j] = M[i][j][t] zu und diagonalisieren Sie sie für jedes t.

Das Problem ist, dass das Programm die Eigenwerte nach einiger Iteration in unterschiedlicher Reihenfolge gibt. Zum Beispiel (L = 3), wenn bei t = 0 i eigen[t = 0] = {l1,l2,l3}(0) bei t = 1 bekomme ich eigen[t = 1] = {l3,l2,l1}(1) bekommen kann, während ich immer haben müssen {l1,l2,l3}(t) konkreter werden: Betrachten Sie die Matrix M (t)) = {{0,t,t},{t,0,2t},{t,t,0}} die Eigenwerte immer (approximatevly) sein l1 = -1.3 t , l2 = -t , l3 = 2.3 t Als ich versuchte, diagonalisieren es (mit dem Code unten) habe ich mehrere Male einen Swap im Ergebnis der Eigenwerte bekommen. Gibt es einen Weg, das zu verhindern? Ich kann sie nicht einfach nach Größe sortieren, ich brauche sie immer in der gleichen Reihenfolge (was auch immer es ist) a priori. (der folgende Code ist nur ein Beispiel, um mein Problem zu erhellen)

EDIT: Ich kann sie nicht einfach sortieren, weil ich a priori ihren Wert nicht kenne oder wenn sie zu jeder Zeit zuverlässig eine Struktur wie l1<l2<l3 haben Statistische Fluktuationen, deshalb wollte ich wissen, ob es einen Weg gibt, den Algorithmus immer gleich zu verhalten, so dass die Reihenfolge der Eigenwerte immer gleich ist oder wenn es einen Trick gibt, um das zu erreichen.

Nur um klarer zu sein, ich werde versuchen, das Spielzeugproblem, das ich hier vorgestellt habe, neu zu beschreiben. Wir haben eine Matrix, die von der Zeit abhängt, ich, vielleicht naiv, erwarte nur lambda_1(t).....lambda_N(t) zu bekommen, stattdessen sehe ich, dass der Algorithmus oft die Eigenwerte zu verschiedenen Zeiten tauscht, also wenn ich zum Beispiel bei t = 1 I've got (lambda_1,lambda_2,lambda_3)(1) at time t = 2 (lambda_2,lambda_1,lambda_3)(2) sehen wollte, wie Lambda_1 entwickelt sich in der Zeit kann ich nicht, weil der Algorithmus die Eigenwerte zu verschiedenen Zeiten mischt. Das Programm unten ist nur ein analytisches Spielzeug Beispiel für mein Problem: Die Eigenwerte der Matrix unten sind l1 = -1.3 t , l2 = -t , l3 = 2.3 t, aber das Programm kann mir als Ausgabe geben (-1.3,-1,2.3)(1), (-2,-2.6,4.6)(2), etc Wie zuvor gesagt, ich frage mich dann, ob es eine Möglichkeit gibt, das Programm zu bestellen die Eigenwerte sind immer gleich, trotz ihres tatsächlichen numerischen Wertes, so dass ich immer die Kombination (l1, l2, l3) bekomme. Ich hoffe es ist jetzt klarer, bitte sag mir ob es nicht so ist.

#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"); 
    } 
} 

Answer 1

Ich glaube nicht, können Sie tun, was Sie wollen. Als t ändert sich der Ausgang.

Meine ursprüngliche Antwort erwähnt Reihenfolge auf den Zeigern, aber mit Blick auf die Datenstruktur wird es nicht helfen. Wenn die Eigenwerte berechnet wurden, werden die Werte in E gespeichert. Sie können sie wie folgt sehen.

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]); 

Der folgende Code ist, wie die Daten in dem Eigenvektor heraus gelegt ...

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]); 
} 

Wenn und Sie können dies überprüfen, wenn Sie die Bestellung ändern zu sehen, die Reihenfolge der Daten in mdata ändert UND die Reihenfolge in data ändert dann der Algorithmus hat keine feste Reihenfolge, dh Sie können nicht tun, was Sie es tun wollen. Wenn sich der Auftrag nicht in mdata ändert und es sich in data ändert, dann haben Sie eine Lösung, aber ich bezweifle wirklich, dass das der Fall sein wird.

Answer 2

Ihre Frage macht keinen Sinn. Eigenwerte haben für sie keine inhärente Ordnung. Es klingt für mich so, als würden Sie Eigenwerte von M_t definieren, die etwas mit L_1 (M_t), ..., L_n (M_t) zu tun haben, und dann verfolgen, wie sie sich in der Zeit ändern. Angenommen, Ihr Prozess, der M_t antreibt, ist kontinuierlich, so werden auch Ihre Eigenwerte sein.Mit anderen Worten, sie ändern sich nicht signifikant, wenn Sie kleine Änderungen an M_t vornehmen. Wenn Sie also eine Reihenfolge definieren, indem Sie L_1 < L_2 ... < L_n erzwingen, ändert sich diese Reihenfolge für kleine Änderungen in t nicht. Wenn zwei Eigenwerte gekreuzt sind, müssen Sie eine Entscheidung darüber treffen, wie Sie Änderungen zuweisen. Wenn Sie "zufällige Fluktuationen" haben, die größer sind als der typische Abstand zwischen Ihren Eigenwerten, wird dies im Wesentlichen unmöglich.

Hier ist eine andere Möglichkeit, Eigenvektoren zu verfolgen, die sich als besser erweisen könnten. Angenommen, Ihre Eigenvektoren sind v_i, mit den Komponenten v_ij. Was Sie tun, ist zuerst Ihre Eigenvektoren zu "normalisieren", so dass v_i1 nicht negativ ist, d.h. einfach das Vorzeichen jedes Eigenvektors entsprechend umkehren. Dies definiert eine Ordnung auf Ihren Eigenwerten durch eine Ordnung auf v_i1, der ersten Komponente jedes Eigenvektors. Auf diese Weise können Sie immer noch Eigenwerte verfolgen, die sich überschneiden. Wenn sich Ihre Eigenvektoren jedoch auf die erste Komponente kreuzen, sind Sie in Schwierigkeiten.

Answer 3

die Dokumentation nach, kehren diese Funktionen ungeordnete:

https://www.gnu.org/software/gsl/manual/html_node/Real-Symmetric-Matrices.html

diese Funktion kann die Eigenwerte des realen symmetrischen Matrix A. zusätzlichen Arbeitsbereiches der entsprechenden Größe berechnet in w bereitgestellt werden müssen. Der diagonale und der untere dreieckige Teil von A werden während der Berechnung zerstört, aber der strenge obere dreieckige Teil wird nicht referenziert. Die Eigenwerte werden im Vektor eval gespeichert und sind ungeordnet.

Auch die Funktionen, die Ergebnisse geordnete Rückkehr tun dies durch einfache Aufstiegs-/Abstiegs Größe:

https://www.gnu.org/software/gsl/manual/html_node/Sorting-Eigenvalues-and-Eigenvectors.html

Diese Funktion sortiert gleichzeitig die gespeicherten Eigenwerte im Vektor eval und die entsprechenden realen Eigenvektoren gespeichert in den Spalten der Matrix evec in aufsteigender oder absteigender Reihenfolge entsprechend dem Wert des Parameters sort_type, wie oben gezeigt.

Wenn Sie sich für die zeitliche Entwicklung des Eigenwerts suchen, wie genau das tun Sie wurden für die zeitabhängige Darstellungen, zum Beispiel zu tun und lösen:

lambda_1(t).....lambda_N(t) 

Für Ihre einfache Zeit-as -scalar Beispiel

l1 = -1.3 t , l2 = -t , l3 = 2.3 t 

Sie haben buchstäblich eine Parametrisierung aller möglichen Lösungen und weil man sie Identifikatoren zugewiesen haben ln Sie nicht in der Frage der Entartung führen. Selbst wenn M[i][j] nichtlineare Funktionen von t sind, sollte es keine Rolle spielen, da das System selbst linear ist und Lösungen rein durch die charakteristische Gleichung berechnet werden (die t konstant halten wird, während für Lambda gelöst wird).