2016-04-19 9 views
7

を注文するので、私はすべてのTIのための固有値を取得するM[i][j][t]を持っても、時間t = 1,....,Nの関数である長さ×L行列M[i][j]の固有値を見つけるために、GSLライブラリからgsl_eigen_nonsymmおよび/またはgsl_eigen_symm L x L行列E[i][j] = M[i][j][t]を割り当て、すべてのtに対して対角化する。GSLの固有値は、私は関数を使用しています

問題は、プログラムがいくつかの繰り返しの後で異なる順序で固有値を与えることです。私は対角化しようとしたときの固有値は常にl1 = -1.3 t , l2 = -t , l3 = 2.3 t(approximatevly)になりますマトリックスM (t)) = {{0,t,t},{t,0,2t},{t,t,0}}を考えてみます。t = 0で私はt = 1eigen[t = 0] = {l1,l2,l3}(0)を取得する場合、私は常に、より具体的に{l1,l2,l3}(t) を持っている必要がありながら、例(L = 3)iについてeigen[t = 1] = {l3,l2,l1}(1)を得ることができますそれは(以下のコードで)私は固有値の結果に数倍のスワップを得ました。それを防ぐ方法はありますか?私は大きさでそれらを並べ替えることはできません。私はいつも同じ順序(何でも)先験的にする必要があります。 (下のコードは私の問題を解明するための一例に過ぎません)

編集:先験的な私はその価値を知らないし、信頼性の高い構造のために毎回l1<l2<l3のように統計的なゆらぎのため、アルゴリズムが常に同じように動作するようにする方法があるかどうかを知りたかったのは、固有値の順序が常に同じであるように、またはそれを実現するためのトリックがある場合です。

ここで私が提示したおもちゃの問題を再説明しようとします。私たちは時間に依存する行列を持っています。私はおそらく純粋にlambda_1(t).....lambda_N(t)を得ると予想されていますが、アルゴリズムはしばしば異なる時刻に固有値を入れ替えることが多いので、例えばt = 1 I've got (lambda_1,lambda_2,lambda_3)(1) at time t = 2 (lambda_2,lambda_1,lambda_3)(2)にあるとすれば、アルゴリズムが異なる時間に固有値を混ぜ合わせるので、時間がたつにつれて進化する。下のプログラムは、私の問題の分析的なおもちゃの例です。以下の行列の固有値はl1 = -1.3 t , l2 = -t , l3 = 2.3 tですが、プログラムは私に出力として与えるかもしれません。(-1.3,-1,2.3)(1), (-2,-2.6,4.6)(2), etc先に述べたように、私はプログラムのオーダー固有値は実際の数値にもかかわらず常に同じ方法であるので、私は常に(l1、l2、l3)の組み合わせを得ます。私はそれがより明確であることを今願っています、そうでないか教えてください。

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

私はあなたがそれらを並べ替えることは望まないと分かりますが、gslにはソート機能もあります。私はあなたがすでにそれを知っていると思いますが、ただの場合... gsl_eigen_nonsymmv_sort –

+0

固有値のための「先験的な順序付け」はありません。あなたのソートを妨げるものを詳しく教えてください。 – Phillip

+0

私は質問を編集しましたが、基本的には、固有値が時々、その行動を「ランダムに」変化させるような統計的な変動のため、通常の方法で並べ替えることはできません。 – Fra

答えて

1

私はあなたが望むことはできないと思います。 tは出力の変更を変更します。

私の元の答えは、ポインタの順序を述べているが、データ構造を見ればそれは役に立ちません。固有値が計算されると、値はEに格納されます。次のように値が表示されます。

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

次のコードは、場合

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

...固有ベクトルのデータがレイアウトされている方法です、そしてあなたが順序の変更を参照するときには、データの順序を、これを確認することができますmdataに変更し、dataの順番を変更すると、アルゴリズムには一定の順序がありません。つまり、実行するように求めることはできません。注文がmdataで変更されず、dataに変更された場合は解決策がありますが、その場合は本当に疑問です。

+0

ごめんなさいもう少し詳しく教えていただけますか?ありがとうございます – Fra

+0

あなたの答えにいくつかの努力を置くか、単に答えないでください。これをコメントとして投稿できました。または、答えに詳細と詳細を記入してください。あなたの答えを編集するか、単に削除してください。この状態では役に立たない。 –

+0

@AshishAhujaツもっと努力してみてください。 – Harry

1

あなたの質問には意味がありません。固有値に固有の順序はありません。 M_tの固有値をL_1(M_t)、...、L_n(M_t)に類似したものに定義し、それらの時間変化を追跡したいと思うように聞こえます。あなたのプロセスがM_tを連続的に動かすと仮定すると、あなたの固有値もそうです。言い換えれば、M_tを少し変更しても、それらは大きく変化しません。したがって、L_1 < L_2 ... < L_nを適用して順序付けを定義すると、tの小さな変更に対してこの順序は変更されません。 2つの固有値が交差する場合は、変更を割り当てる方法を決定する必要があります。あなたの固有値間の典型的な距離よりも大きな「ランダムな変動」がある場合、これは本質的に不可能になります。

ここでは、固有ベクトルをトラッキングする別の方法があります。これを行うには、固有ベクトルがv_iで、コンポーネントv_ijであるとします。まず、v_i1が非負であるように固有ベクトルを「正規化」します。つまり、各固有ベクトルの符号を適切に反転させるだけです。これは、各固有ベクトルの最初の要素であるv_i1の順序付けによって、固有値の順序を定義します。このようにして、相互に交差する固有値を追跡することができます。しかし、あなたの固有ベクトルが最初のコンポーネントと交差すると、困ったことになります。

0

は、ドキュメントによれば、これらの機能は、順不同返す:

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

この機能は、適切なサイズの実対称行列Aの追加ワークスペースの固有値をWに提供されなければならない計算します。 Aの対角線と下三角部分は計算中に破棄されますが、厳密な上三角部分は参照されません。固有値はベクトルevalに格納され、順序付けられていません。

結果を命じ返すとしても機能し、シンプル昇順/降順の大きさによってそう:同時に

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

この関数は、ベクトルevalの中に保存された固有値と対応する実固有ベクトルをソート行列evecの列に、上に示したパラメータsort_typeの値に従って昇順または降順に格納されます。

あなたは、あなたがやっているように行い、時間依存の表現のために解決し、固有値の時間発展を探しているなら例えば:あなたの簡単なために

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

時間など-scalarたとえば、

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

あなたは、文字通りすべての可能な解決策のパラメータを持っていて、それらをあなたが縮退の問題に実行しないln識別子が割り当てられたので。 M[i][j]がtの非線形関数であっても、システム自体は線形であり、解は純粋に(ラムダを解く際にt定数を保持する)特性方程式によって計算されるので問題ではありません。

関連する問題