さまざまな言語で数値計算
Only Do What Only You Can Do
反復法
反復法を使って, 最大固有値と, 対応する固有ベクトルを求める
$ n \times n $ の正方行列 $ \boldsymbol{A} $ と $ n $ 次元のベクトル $ \boldsymbol{x} $ について
$ \boldsymbol{Ax} = \lambda \boldsymbol{x} $ (ただし $ \boldsymbol{x} \neq 0 $ ) が成り立つとき
$ \lambda $ を固有値, $ \boldsymbol{x} $ を固有ベクトルという.
最初に適当なベクトル $ \boldsymbol{x}_0 $ から始めて $ \boldsymbol{x}_{k+1} = \boldsymbol{Ax}_{k} $ を反復すると
$ \boldsymbol{x}_{k} $ は行列 $ \boldsymbol{A} $ の最大固有値に対応する固有ベクトルに収束する.
固有値はレイリー(Rayleigh)商
- $ \lambda = \displaystyle \frac{\boldsymbol{x}_k^T \boldsymbol{x}_{k+1}}{\boldsymbol{x}_k^T \boldsymbol{x}_k} $
べき乗法、累乗法とも言う.
VBScript
Option Explicit Private Const N = 3 'ベキ乗法で最大固有値を求める Call Main Private Sub Main() Dim a: a = Array(_ Array(5.0, 4.0, 1.0, 1.0), _ Array(4.0, 5.0, 1.0, 1.0), _ Array(1.0, 1.0, 4.0, 2.0), _ Array(1.0, 1.0, 2.0, 4.0)) Dim x: x = Array(1.0 ,0.0 ,0.0 ,0.0) 'ベキ乗法 Dim lambda: lambda = power(a, x) WScript.StdOut.WriteLine WScript.StdOut.WriteLine "eigenvalue" WScript.StdOut.WriteLine Right(Space(14) & FormatNumber(lambda, 10, -1, 0, 0), 14) WScript.StdOut.WriteLine "eigenvector" Call disp_vector(x) End Sub 'ベキ乗法 Private Function power(ByVal a, ByRef x0) Dim i, j, k Dim lambda: lambda = 0.0 '正規化 (ベクトル x0 の長さを１にする) Call normarize(x0) Dim e0: e0 = 0.0 For i = 0 To N e0 = e0 + x0(i) Next For k = 1 To 200 '１次元配列を表示 WScript.StdOut.Write Right(Space(3) & k, 3) & vbTab Call disp_vector(x0) '行列の積 x1 = A × x0 Dim x1: x1 = Array(0.0, 0.0, 0.0, 0.0) For i = 0 To N For j = 0 To N x1(i) = x1(i) + a(i)(j) * x0(j) Next Next '内積 Dim p0: p0 = 0.0 Dim p1: p1 = 0.0 For i = 0 To N p0 = p0 + x1(i) * x1(i) p1 = p1 + x1(i) * x0(i) Next '固有値 lambda = p0 / p1 '正規化 (ベクトル x1 の長さを１にする) Call normarize(x1) '収束判定 Dim e1: e1 = 0.0 For i = 0 To N e1 = e1 + x1(i) Next If Abs(e1 - e0) < 0.00000000001 Then Exit For x0 = x1 e0 = e1 Next power = lambda End Function '１次元配列を表示 Private Sub disp_vector(ByVal row) Dim col For Each col In row WScript.StdOut.Write Right(Space(14) & FormatNumber(col, 10, -1, 0, 0), 14) & vbTab Next WScript.StdOut.WriteLine End Sub '正規化 (ベクトルの長さを１にする) Private Sub normarize(ByRef x()) Dim s: s = 0.0 Dim i For i = 0 To N s = s + x(i) * x(i) Next s = Sqr(s) For i = 0 To N x(i) = x(i) / s Next End Sub
Z:\>cscript //nologo Z:\1101.vbs 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
JScript
var N = 4 // ベキ乗法で最大固有値を求める main() function main() { var a = [[5.0, 4.0, 1.0, 1.0], [4.0, 5.0, 1.0, 1.0], [1.0, 1.0, 4.0, 2.0], [1.0, 1.0, 2.0, 4.0]] var x = [1.0 ,0.0 ,0.0 ,0.0] // ベキ乗法 var lambda = power(a, x) WScript.StdOut.WriteLine() WScript.StdOut.WriteLine("eigenvalue") WScript.StdOut.WriteLine((" " + lambda.toFixed(10)).slice(-14)) WScript.StdOut.WriteLine("eigenvector") disp_vector(x) } // ベキ乗法 function power(a, x0) { var lambda = 0.0 // 正規化 (ベクトル x0 の長さを１にする) normarize(x0) var e0 = 0.0 for (var i = 0; i < N; i++) e0 += x0[i] for (var k = 1; k <= 200; k++) { // １次元配列を表示 WScript.StdOut.Write((" " + k).slice( -3) + "\t") disp_vector(x0) // 行列の積 x1 = A × x0 var x1 = [0.0 ,0.0 ,0.0 ,0.0] for (var i = 0; i < N; i++) for (var j = 0; j < N; j++) x1[i] += a[i][j] * x0[j] // 内積 var p0 = 0.0 var p1 = 0.0 for (var i = 0; i < N; i++) { p0 += x1[i] * x1[i] p1 += x1[i] * x0[i] } // 固有値 lambda = p0 / p1 // 正規化 (ベクトル x1 の長さを１にする) normarize(x1) // 収束判定 var e1 = 0.0 for (var i = 0; i < N; i++) e1 += x1[i] if (Math.abs(e0 - e1) < 0.00000000001) break for (var i = 0; i < N; i++) x0[i] = x1[i] e0 = e1 } return lambda } // １次元配列を表示 function disp_vector(row) { for (var i = 0; i < N; i++) WScript.StdOut.Write((" " + row[i].toFixed(10)).slice(-14) + "\t") WScript.StdOut.WriteLine() } // 正規化 (ベクトルの長さを１にする) function normarize(x) { var s = 0.0 for (var i = 0; i < N; i++) s += x[i] * x[i] s = Math.sqrt(s) for (var i = 0; i < N; i++) x[i] /= s }
Z:\>cscript //nologo Z:\1101.js 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
PowerShell
set-variable -name N -value 3 -option constant # １次元配列を表示 function disp_vector($row) { foreach ($col in $row) { Write-Host ([String]::Format("{0,14:F10}", $col)) -nonewline } Write-Host } # 正規化 (ベクトルの長さを１にする) function normarize($x) { $s = 0.0 foreach ($i in 0..$N) { $s += $x[$i] * $x[$i] } $s = [Math]::Sqrt($s) foreach ($i in 0..$N) { $x[$i] /= $s } } # ベキ乗法 function power($a, $x0) { $lambda = 0.0 # 正規化 (ベクトル x0 の長さを１にする) normarize $x0 $e0 = 0.0 foreach ($i in 0..$N) { $e0 += $x0[$i] } foreach ($k in 1..200) { # １次元配列を表示 Write-Host ([String]::Format("{0,3:D}`t", $k)) -nonewline disp_vector $x0 # 行列の積 x1 = A × x0 $x1 = (0.0, 0.0, 0.0, 0.0) foreach ($i in 0..$N) { foreach ($j in 0..$N) { $x1[$i] += $a[$i][$j] * $x0[$j] } } # 内積 $p0 = 0.0 $p1 = 0.0 foreach ($i in 0..$N) { $p0 += $x1[$i] * $x1[$i] $p1 += $x1[$i] * $x0[$i] } # 固有値 $lambda = $p0 / $p1 # 正規化 (ベクトル x1 の長さを１にする) normarize $x1 # 収束判定 $e1 = 0.0 foreach ($i in 0..$N) { $e1 += $x1[$i] } if ([Math]::Abs($e0 - $e1) -lt 0.00000000001) { break } foreach ($i in 0..$N) { $x0[$i] = $x1[$i] } $e0 = $e1 } $lambda } # ベキ乗法で最大固有値を求める $a = (5.0, 4.0, 1.0, 1.0), (4.0, 5.0, 1.0, 1.0), (1.0, 1.0, 4.0, 2.0), (1.0, 1.0, 2.0, 4.0) $x = (1.0 ,0.0 ,0.0 ,0.0) # ベキ乗法 $lambda = (power $a $x) Write-Host Write-Host "eigenvalue" Write-Host ([String]::Format("{0,14:F10}", $lambda)) Write-Host "eigenvector" disp_vector $x
Z:\>powershell -file Z:\1101.ps1 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
Perl
use constant N => 3; # ベキ乗法で最大固有値を求める main(); sub main { my @a = ([5.0, 4.0, 1.0, 1.0], [4.0, 5.0, 1.0, 1.0], [1.0, 1.0, 4.0, 2.0], [1.0, 1.0, 2.0, 4.0]); my @x = (1.0, 0.0, 0.0, 0.0); # ベキ乗法 Power iteration $lambda = power(\@a, \@x); print "\neigenvalue\n"; printf("%14.10f", $lambda); print "\neigenvector\n"; disp_vector(\@x); } # ベキ乗法 Power iteration sub power { my ($a, $x0) = @_; my $lambda = 0.0; # 正規化 (ベクトル x0 の長さを１にする) normarize(\@$x0); my $e0 = 0.0; for $i (0..N) { $e0 += $$x0[$i]; } for $k (1..200) { # １次元配列を表示 printf("%3d\t", $k); disp_vector(\@$x0); # 行列の積 x1 = A × x0 my @x1; for $i (0..N) { for $j (0..N) { $x1[$i] += $$a[$i][$j] * $$x0[$j]; } } # 内積 inner product my $p0 = 0.0; my $p1 = 0.0; for $i (0..N) { $p0 += $x1[$i] * $x1[$i]; $p1 += $x1[$i] * $$x0[$i]; } # 固有値 eigenvalue $lambda = $p0 / $p1; # 正規化 (ベクトル x1 の長さを１にする) normarize(\@x1); # 収束判定 my $e1 = 0.0; for $i (0..N) { $e1 += $x1[$i]; } last if (abs($e0 - $e1) < 0.00000000001); for $i (0..N) { $$x0[$i] = $x1[$i]; } $e0 = $e1; } $lambda; } # １次元配列を表示 sub disp_vector { my ($x) = @_; for $i (0..N) { printf("%14.10f\t", $$x[$i]); } print "\n"; } # 正規化 (ベクトルの長さを１にする) sub normarize { my ($x) = @_; my $s = 0.0; for $i (0..N) { $s += $$x[$i] * $$x[$i]; } $s = sqrt($s); for $i (0..N) { $$x[$i] /= $s; } }
Z:\>perl Z:\1101.pl 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
PHP
<?php define("N", 3); # ベキ乗法で最大固有値を求める main(); function main() { $a = [[5.0, 4.0, 1.0, 1.0], [4.0, 5.0, 1.0, 1.0], [1.0, 1.0, 4.0, 2.0], [1.0, 1.0, 2.0, 4.0]]; $x = [1.0, 0.0, 0.0, 0.0]; # ベキ乗法 $lambda = power($a, $x); print "\neigenvalue\n"; printf("%14.10f", $lambda); print "\neigenvector\n"; disp_vector($x); } # ベキ乗法 function power(&$a, &$x0) { $lambda = 0.0; # 正規化 (ベクトル x0 の長さを１にする) normarize($x0); $e0 = 0.0; foreach (range(0, N) as $i) { $e0 += $x0[$i]; } foreach (range(1, 200) as $k) { # １次元配列を表示 printf("%3d\t", $k); disp_vector($x0); # 行列の積 x1 = A × x0 $x1 = [0.0, 0.0, 0.0, 0.0]; foreach (range(0, N) as $i) { foreach (range(0, N) as $j) { $x1[$i] += $a[$i][$j] * $x0[$j]; } } # 内積 $p0 = 0.0; $p1 = 0.0; foreach (range(0, N) as $i) { $p0 += $x1[$i] * $x1[$i]; $p1 += $x1[$i] * $x0[$i]; } # 固有値 $lambda = $p0 / $p1; # 正規化 (ベクトル x1 の長さを１にする) normarize($x1); # 収束判定 $e1 = 0.0; foreach (range(0, N) as $i) { $e1 += $x1[$i]; } if (abs($e0 - $e1) < 0.00000000001) break; foreach (range(0, N) as $i) { $x0[$i] = $x1[$i]; } $e0 = $e1; } return $lambda; } # １次元配列を表示 function disp_vector($row) { foreach ($row as $col) printf("%14.10f\t", $col); print "\n"; } # 正規化 (ベクトルの長さを１にする) function normarize(&$x) { $s = 0.0; foreach (range(0, N) as $i) { $s += $x[$i] * $x[$i]; } $s = sqrt($s); foreach (range(0, N) as $i) { $x[$i] /= $s; } } ?>
Z:\>php Z:\1101.php 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
Python
# coding: Shift_JIS import math N = 4 # １次元配列を表示 def disp_vector(row): for col in row: print "%14.10f\t" % col, print "" # 正規化 (ベクトルの長さを１にする) def normarize(x): s = 0.0 for i in range(0, N, 1): s += x[i] * x[i] s = math.sqrt(s) for i in range(0, N, 1): x[i] /= s # ベキ乗法 def power(a, x0): lambda0 = 0.0 # 正規化 (ベクトル x0 の長さを１にする) normarize(x0) e0 = 0.0 for i in range(0, N, 1): e0 += x0[i] for k in range(1, 200, 1): # １次元配列を表示 print "%3d\t" % k, disp_vector(x0) # 行列の積 x1 = A × x0 x1 = [ 0.0, 0.0, 0.0, 0.0] for i in range(0, N, 1): for j in range(0, N, 1): x1[i] += a[i][j] * x0[j] # 内積 p0 = 0.0 p1 = 0.0 for i in range(0, N, 1): p0 += x1[i] * x1[i] p1 += x1[i] * x0[i] # 固有値 lambda0 = p0 / p1 # 正規化 (ベクトル x1 の長さを１にする) normarize(x1) # 収束判定 e1 = 0.0 for i in range(0, N, 1): e1 += x1[i] if (abs(e0 - e1) < 0.00000000001): break for i in range(0, N, 1): x0[i] = x1[i] e0 = e1 return lambda0 # ベキ乗法で最大固有値を求める a = [[5.0, 4.0, 1.0, 1.0], [4.0, 5.0, 1.0, 1.0], [1.0, 1.0, 4.0, 2.0], [1.0, 1.0, 2.0, 4.0]] x = [1.0 ,0.0 ,0.0 ,0.0] # ベキ乗法 lambda0 = power(a, x) print "" print "eigenvalue" print "%14.10f" % lambda0 print "eigenvector" disp_vector(x)
Z:\>python Z:\1101.py 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
Ruby
# coding: Shift_JIS N = 3 # １次元配列を表示 def disp_vector(row) row.each do |col| printf("%14.10f\t", col) end puts "" end # 正規化 (ベクトルの長さを１にする) def normarize(x) s = 0.0 (0..N).each do |i| s += x[i] * x[i] end s = Math.sqrt(s) (0..N).each do |i| x[i] /= s end end # ベキ乗法 def power(a, x0) lambda = 0.0 # 正規化 (ベクトル x0 の長さを１にする) normarize(x0) e0 = 0.0 (0..N).each do |i| e0 += x0[i] end (1..200).each do |k| # １次元配列を表示 printf("%3d\t", k) disp_vector(x0) # 行列の積 x1 = A × x0 x1 = [ 0.0, 0.0, 0.0, 0.0] (0..N).each do |i| (0..N).each do |j| x1[i] += a[i][j] * x0[j] end end # 内積 p0 = 0.0 p1 = 0.0 (0..N).each do |i| p0 += x1[i] * x1[i] p1 += x1[i] * x0[i] end # 固有値 lambda = p0 / p1 # 正規化 (ベクトル x1 の長さを１にする) normarize(x1) # 収束判定 e1 = 0.0 (0..N).each do |i| e1 += x1[i] end break if ((e0 - e1).abs < 0.00000000001) (0..N).each do |i| x0[i] = x1[i] end e0 = e1 end lambda end # ベキ乗法で最大固有値を求める a = [[5.0, 4.0, 1.0, 1.0], [4.0, 5.0, 1.0, 1.0], [1.0, 1.0, 4.0, 2.0], [1.0, 1.0, 2.0, 4.0]] x = [1.0 ,0.0 ,0.0 ,0.0] # ベキ乗法 lambda0 = power(a, x) puts "" puts "eigenvalue" printf("%14.10f\n", lambda0) puts "eigenvector" disp_vector(x)
Z:\>ruby Z:\1101.rb 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
Groovy
N = 3 // ベキ乗法で最大固有値を求める main() def main() { def a = [[5.0, 4.0, 1.0, 1.0], [4.0, 5.0, 1.0, 1.0], [1.0, 1.0, 4.0, 2.0], [1.0, 1.0, 2.0, 4.0]] as double[][] def x = [1.0 ,0.0 ,0.0 ,0.0] as double[] // ベキ乗法 lambda = power(a, x) println() println("eigenvalue") printf("%14.10f\n" , lambda) println("eigenvector") disp_vector(x) } // ベキ乗法 def power(a, x0) { lambda = 0.0 // 正規化 (ベクトル x0 の長さを１にする) normarize(x0) e0 = 0.0 for (i in 0..N) e0 += x0[i] for (k in 1..200) { // １次元配列を表示 printf("%3d\t" , k) disp_vector(x0) // 行列の積 x1 = A × x0 x1 = [0.0 ,0.0 ,0.0 ,0.0] for (i in 0..N) for (j in 0..N) x1[i] += a[i][j] * x0[j] // 内積 p0 = 0.0 p1 = 0.0 for (i in 0..N) { p0 += x1[i] * x1[i] p1 += x1[i] * x0[i] } // 固有値 lambda = p0 / p1 // 正規化 (ベクトル x1 の長さを１にする) normarize(x1) // 収束判定 e1 = 0.0 for (i in 0..N) e1 += x1[i] if (Math.abs(e0 - e1) < 0.00000000001) break for (i in 0..N) x0[i] = x1[i] e0 = e1 } lambda } // １次元配列を表示 def disp_vector(row) { for (col in row) { printf("%14.10f\t" , col) } println() } // 正規化 (ベクトルの長さを１にする) def normarize(x) { s = 0.0 for (i in 0..N) s += x[i] * x[i] s = Math.sqrt(s) for (i in 0..N) x[i] /= s }
Z:\>gvy Gvy1101.gvy 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
Pascal
program Pas1101(arg); {$MODE delphi} uses SysUtils, Math; const N = 3; type TwoDimArray = array [0..N,0..N] of Double; // １次元配列を表示 procedure disp_vector(row:array of Double); var i:Integer; begin for i := Low(row) to High(row) do write(format('%14.10f'#9, [row[i]])); writeln(); end; // 正規化 (ベクトルの長さを１にする) procedure normarize(var x:array of Double); var s:Double = 0.0; i:Integer; begin for i := Low(x) to High(x) do s := s + x[i] * x[i]; s := Sqrt(s); for i := Low(x) to High(x) do x[i] := x[i] / s; end; // ベキ乗法 function power(a:TwoDimArray; var x0:array of Double):Double; var lambda:Double = 0.0; p0, p1, e0, e1:Double; i, j, k:Integer; x1:array [0..N] of Double = (0.0, 0.0, 0.0, 0.0); begin // 正規化 (ベクトル x0 の長さを１にする) normarize(x0); e0 := 0.0; for i := Low(x0) to High(x0) do e0 := e0 + x0[i]; for k := 1 to 100 do begin // １次元配列を表示 write(format('%3d'#9, [k])); disp_vector(x0); // 行列の積 x1 = A × x0 for i := Low(x1) to High(x1) do x1[i] := 0.0; for i := Low(x1) to High(x1) do begin for j := Low(x0) to High(x0) do x1[i] := x1[i] + a[i][j] * x0[j]; end; // 内積 p0 := 0.0; p1 := 0.0; for i := Low(x0) to High(x0) do begin p0 := p0 + x1[i] * x1[i]; p1 := p1 + x1[i] * x0[i]; end; // 固有値 lambda := p0 / p1; // 正規化 (ベクトル x1 の長さを１にする) normarize(x1); // 収束判定 e1 := 0.0; for i := Low(x0) to High(x0) do e1 := e1 + x1[i]; if Abs(e1 - e0) < 0.00000000001 then break; for i := Low(x0) to High(x0) do x0[i] := x1[i]; e0 := e1; end; power := lambda; end; // ベキ乗法で最大固有値を求める var a :TwoDimArray = ((5.0, 4.0, 1.0, 1.0) ,(4.0, 5.0, 1.0, 1.0) ,(1.0, 1.0, 4.0, 2.0) ,(1.0, 1.0, 2.0, 4.0)); x :array [0..N] of Double = (1.0, 0.0, 0.0, 0.0); lambda:Double; begin // ベキ乗法 lambda := power(a, x); writeln(); writeln('eigenvalue'); writeln(format('%14.10f', [lambda])); writeln('eigenvector'); disp_vector(x); end.
Z:\>fpc -v0 -l- Pas1101.pp Z:\>Pas1101 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
Ada
with TEXT_IO, Ada.Integer_Text_IO, Ada.Long_Float_Text_IO, Ada.Numerics.Long_Elementary_Functions; use TEXT_IO, Ada.Integer_Text_IO, Ada.Long_Float_Text_IO, Ada.Numerics.Long_Elementary_Functions; procedure Ada1101 is N : Constant Integer := 3; type Long_Float_Array is array (0..N) of Long_Float; type Long_Float_TwoDimArray is array (0..N, 0..N) of Long_Float; a:Long_Float_TwoDimArray := ((5.0, 4.0, 1.0, 1.0) ,(4.0, 5.0, 1.0, 1.0) ,(1.0, 1.0, 4.0, 2.0) ,(1.0, 1.0, 2.0, 4.0)); x:Long_Float_Array := (1.0, 0.0, 0.0, 0.0); -- １次元配列を表示 procedure disp_vector(row:Long_Float_Array) is begin for i in row'Range loop Put(row(i), Fore=>3, Aft=>10, Exp=>0); Put(Ascii.HT); end loop; New_Line; end disp_vector; -- 正規化 (ベクトルの長さを１にする) procedure normarize(x:in out Long_Float_Array) is s:Long_Float; begin s := 0.0; for i in x'Range loop s := s + x(i) * x(i); end loop; s := Sqrt(s); for i in x'Range loop x(i) := x(i) / s; end loop; end normarize; -- ベキ乗法 procedure power(a:Long_Float_TwoDimArray; x0:in out Long_Float_Array; lambda:out Long_Float) is p0, p1, e0, e1:Long_Float; x1:Long_Float_Array := (0.0, 0.0, 0.0, 0.0); begin -- 正規化 (ベクトル x0 の長さを１にする) normarize(x0); e0 := 0.0; for i in x0'Range loop e0 := e0 + x0(i); end loop; for k in 1..100 loop -- １次元配列を表示 Put(k, Width=> 3); Put(Ascii.HT); disp_vector(x0); -- 行列の積 x1 = A × x0 for i in x1'Range loop x1(i) := 0.0; end loop; for i in x1'Range loop for j in x1'Range loop x1(i) := x1(i) + a(i, j) * x0(j); end loop; end loop; -- 内積 p0 := 0.0; p1 := 0.0; for i in x0'Range loop p0 := p0 + x1(i) * x1(i); p1 := p1 + x1(i) * x0(i); end loop; -- 固有値 lambda := p0 / p1; -- 正規化 (ベクトル x1 の長さを１にする) normarize(x1); -- 収束判定 e1 := 0.0; for i in x0'Range loop e1 := e1 + x1(i); end loop; if Abs(e1 - e0) < 0.00000000001 then exit; end if; for i in x0'Range loop x0(i) := x1(i); end loop; e0 := e1; end loop; end power; -- ベキ乗法で最大固有値を求める lambda:Long_Float := 0.0; begin -- ベキ乗法 power(a, x, lambda); Put_Line(""); Put_Line("eigenvalue"); Put(lambda, Fore=>3, Aft=>10, Exp=>0); New_Line; Put_Line("eigenvector"); disp_vector(x); end Ada1101;
xxxxxx@yyyyyy /Z $ gnatmake Ada1101.adb xxxxxx@yyyyyy /Z $ Ada1101 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
VB.NET
Option Explicit Module VB1101 Private Const N As Integer = 3 'ベキ乗法で最大固有値を求める Public Sub Main() Dim a(,) As Double = {{5.0, 4.0, 1.0, 1.0}, _ {4.0, 5.0, 1.0, 1.0}, _ {1.0, 1.0, 4.0, 2.0}, _ {1.0, 1.0, 2.0, 4.0}} Dim x() As Double = {1.0 ,0.0 ,0.0 ,0.0} 'ベキ乗法 Dim lambda As Double = power(a, x) Console.WriteLine() Console.WriteLine("eigenvalue") Console.WriteLine(String.Format("{0,14:F10}", lambda)) Console.WriteLine("eigenvector") disp_vector(x) End Sub 'ベキ乗法 Private Function power(ByVal a(,) As Double, ByRef x0() As Double) As Double Dim lambda As Double = 0.0 '正規化 (ベクトル x0 の長さを１にする) normarize(x0) Dim e0 As Double = 0.0 For i As Integer = 0 To N e0 += x0(i) Next For k As Integer = 1 To 100 '１次元配列を表示 Console.Write(String.Format("{0,3:D}{1}", k, vbTab)) disp_vector(x0) '行列の積 x1 = A × x0 Dim x1() As Double = {0.0, 0.0, 0.0, 0.0} For i As Integer = 0 To N For j As Integer = 0 To N x1(i) += a(i, j) * x0(j) Next Next '内積 Dim p0 As Double = 0.0 Dim p1 As Double = 0.0 For i As Integer = 0 To N p0 += x1(i) * x1(i) p1 += x1(i) * x0(i) Next '固有値 lambda = p0 / p1 '正規化 (ベクトル x1 の長さを１にする) normarize(x1) '収束判定 Dim e1 As Double = 0.0 For i As Integer = 0 To N e1 += x1(i) Next If Math.Abs(e1 - e0) < 0.00000000001 Then Exit For For i As Integer = 0 To N x0(i) = x1(i) Next e0 = e1 Next Return lambda End Function '１次元配列を表示 Private Sub disp_vector(ByVal row() As Double) For Each col As Double In row Console.Write(String.Format("{0,14:F10}{1}", col, vbTab)) Next Console.WriteLine() End Sub '正規化 (ベクトルの長さを１にする) Private Sub normarize(ByRef x() As Double) Dim s As Double = 0.0 For i As Integer = 0 To N s += x(i) * x(i) Next s = Math.Sqrt(s) For i As Integer = 0 To N x(i) /= s Next End Sub End Module
Z:\>vbc -nologo VB1101.vb Z:\>VB1101 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
C#
using System; public class CS1101 { private const int N = 4; // ベキ乗法で最大固有値を求める public static void Main() { double[,] a = {{5.0, 4.0, 1.0, 1.0}, {4.0, 5.0, 1.0, 1.0}, {1.0, 1.0, 4.0, 2.0}, {1.0, 1.0, 2.0, 4.0}}; double[] x = {1.0 ,0.0 ,0.0 ,0.0}; // ベキ乗法 double lambda = power(a, x); Console.WriteLine(); Console.WriteLine("eigenvalue"); Console.WriteLine(string.Format("{0,14:F10}", lambda)); Console.WriteLine("eigenvector"); disp_vector(x); } // ベキ乗法 private static double power(double[,] a, double[] x0) { double lambda = 0.0; // 正規化 (ベクトル x0 の長さを１にする) normarize(x0); double e0 = 0.0; for (int i = 0; i < N; i++) e0 += x0[i]; for (int k = 1; k <= 200; k++) { // １次元配列を表示 Console.Write(string.Format("{0,3:D}\t", k)); disp_vector(x0); // 行列の積 x1 = A × x0 double[] x1 = new double[N]; for (int i = 0; i < N; i++) for (int j = 0; j < N; j++) x1[i] += a[i, j] * x0[j]; // 内積 double p0 = 0.0; double p1 = 0.0; for (int i = 0; i < N; i++) { p0 += x1[i] * x1[i]; p1 += x1[i] * x0[i]; } // 固有値 lambda = p0 / p1; // 正規化 (ベクトル x1 の長さを１にする) normarize(x1); // 収束判定 double e1 = 0.0; for (int i = 0; i < N; i++) e1 += x1[i]; if (Math.Abs(e0 - e1) < 0.00000000001) break; for (int i = 0; i < N; i++) x0[i] = x1[i]; e0 = e1; } return lambda; } // １次元配列を表示 private static void disp_vector(double[] row) { foreach (double col in row) Console.Write(string.Format("{0,14:F10}\t", col)); Console.WriteLine(); } // 正規化 (ベクトルの長さを１にする) private static void normarize(double[] x) { double s = 0.0; for (int i = 0; i < N; i++) s += x[i] * x[i]; s = Math.Sqrt(s); for (int i = 0; i < N; i++) x[i] /= s; } }
Z:\>csc -nologo CS1101.cs Z:\>CS1101 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
Java
import java.lang.*; public class Java1101 { private static final int N = 4; // ベキ乗法で最大固有値を求める public static void main(String []args) { double[][] a = {{5.0, 4.0, 1.0, 1.0}, {4.0, 5.0, 1.0, 1.0}, {1.0, 1.0, 4.0, 2.0}, {1.0, 1.0, 2.0, 4.0}}; double[] x = {1.0 ,0.0 ,0.0 ,0.0}; // ベキ乗法 double lambda = power(a, x); System.out.println(); System.out.println("eigenvalue"); System.out.println(String.format("%14.10f", lambda)); System.out.println("eigenvector"); disp_vector(x); } // ベキ乗法 private static double power(double[][] a, double[] x0) { double lambda = 0.0; // 正規化 (ベクトル x0 の長さを１にする) normarize(x0); double e0 = 0.0; for (int i = 0; i < N; i++) e0 += x0[i]; for (int k = 1; k <= 100; k++) { // １次元配列を表示 System.out.print(String.format("%3d\t", k)); disp_vector(x0); // 行列の積 x1 = A × x0 double[] x1 = new double[N]; for (int i = 0; i < N; i++) for (int j = 0; j < N; j++) x1[i] += a[i][j] * x0[j]; // 内積 double p0 = 0.0; double p1 = 0.0; for (int i = 0; i < N; i++) { p0 += x1[i] * x1[i]; p1 += x1[i] * x0[i]; } // 固有値 lambda = p0 / p1; // 正規化 (ベクトル x1 の長さを１にする) normarize(x1); // 収束判定 double e1 = 0.0; for (int i = 0; i < N; i++) e1 += x1[i]; if (Math.abs(e0 - e1) < 0.00000000001) break; for (int i = 0; i < N; i++) x0[i] = x1[i]; e0 = e1; } return lambda; } // １次元配列を表示 private static void disp_vector(double[] row) { for (double col: row) System.out.print(String.format("%14.10f\t", col)); System.out.println(); } // 正規化 (ベクトルの長さを１にする) private static void normarize(double[] x) { double s = 0.0; for (int i = 0; i < N; i++) s += x[i] * x[i]; s = Math.sqrt(s); for (int i = 0; i < N; i++) x[i] /= s; } }
Z:\>javac Java1101.java Z:\>java Java1101 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
C++
#include <iostream> #include <iomanip> #include <math.h> #include <vector> using namespace std; const int N = 4; // ベキ乗法 double power(double a[N][N], double x0[N]); // １次元配列を表示 void disp_vector(double row[N]); // 正規化 (ベクトルの長さを１にする) void normarize(double x[N]); // ベキ乗法で最大固有値を求める int main() { double a[N][N] = {{5.0, 4.0, 1.0, 1.0}, {4.0, 5.0, 1.0, 1.0}, {1.0, 1.0, 4.0, 2.0}, {1.0, 1.0, 2.0, 4.0}}; double x[N] = {1.0 ,0.0 ,0.0 ,0.0}; // ベキ乗法 double lambda = power(a, x); cout << endl << "eigenvalue" << endl; cout << setw(14) << fixed << setprecision(10) << lambda; cout << endl << "eigenvector" << endl; disp_vector(x); return 0; } // ベキ乗法 double power(double a[N][N], double x0[N]) { double lambda = 0.0; // 正規化 (ベクトル x0 の長さを１にする) normarize(x0); double e0 = 0.0; for (int i = 0; i < N; i++) e0 += x0[i]; for (int k = 1; k <= 200; k++) { // １次元配列を表示 cout << setw(3) << k << "\t"; disp_vector(x0); // 行列の積 x1 = A × x0 double x1[N] = {0.0 ,0.0 ,0.0 ,0.0}; for (int i = 0; i < N; i++) for (int j = 0; j < N; j++) x1[i] += a[i][j] * x0[j]; // 内積 double p0 = 0.0; double p1 = 0.0; for (int i = 0; i < N; i++) { p0 += x1[i] * x1[i]; p1 += x1[i] * x0[i]; } // 固有値 lambda = p0 / p1; // 正規化 (ベクトル x1 の長さを１にする) normarize(x1); // 収束判定 double e1 = 0.0; for (int i = 0; i < N; i++) e1 += x1[i]; if (fabs(e0 - e1) < 0.00000000001) break; for (int i = 0; i < N; i++) x0[i] = x1[i]; e0 = e1; } return lambda; } // １次元配列を表示 void disp_vector(double row[N]) { for (int i = 0; i < N; i++) cout << setw(14) << fixed << setprecision(10) << row[i] << "\t"; cout << endl; } // 正規化 (ベクトルの長さを１にする) void normarize(double x[N]) { double s = 0.0; for (int i = 0; i < N; i++) s += x[i] * x[i]; s = sqrt(s); for (int i = 0; i < N; i++) x[i] /= s; }
Z:\>bcc32 -q CP1101.cpp cp1101.cpp: Z:\>CP1101 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
Objective-C
#import <Foundation/Foundation.h> #import <math.h> #define N 4 // ベキ乗法 double power(double a[N][N], double x0[N]); // １次元配列を表示 void disp_vector(double row[N]); // 正規化 (ベクトルの長さを１にする) void normarize(double x[N]); // ベキ乗法で最大固有値を求める int main() { double a[N][N] = {{5.0, 4.0, 1.0, 1.0}, {4.0, 5.0, 1.0, 1.0}, {1.0, 1.0, 4.0, 2.0}, {1.0, 1.0, 2.0, 4.0}}; double x[N] = {1.0 ,0.0 ,0.0 ,0.0}; // ベキ乗法 double lambda = power(a, x); printf("\neigenvalue\n"); printf("%14.10f\n", lambda); printf("eigenvector\n"); disp_vector(x); return 0; } // ベキ乗法 double power(double a[N][N], double x0[N]) { int i, j, k; double lambda = 0.0; // 正規化 (ベクトル x0 の長さを１にする) normarize(x0); double e0 = 0.0; for (i = 0; i < N; i++) e0 += x0[i]; for (k = 1; k <= 200; k++) { // １次元配列を表示 printf("%3d\t", k); disp_vector(x0); // 行列の積 x1 = A × x0 double x1[N] = {0.0 ,0.0 ,0.0 ,0.0}; for (i = 0; i < N; i++) for (j = 0; j < N; j++) x1[i] += a[i][j] * x0[j]; // 内積 double p0 = 0.0; double p1 = 0.0; for (i = 0; i < N; i++) { p0 += x1[i] * x1[i]; p1 += x1[i] * x0[i]; } // 固有値 lambda = p0 / p1; // 正規化 (ベクトル x1 の長さを１にする) normarize(x1); double e1 = 0.0; // 収束判定 for (i = 0; i < N; i++) e1 += x1[i]; if (fabs(e0 - e1) < 0.00000000001) break; for (i = 0; i < N; i++) x0[i] = x1[i]; e0 = e1; } return lambda; } // １次元配列を表示 void disp_vector(double row[N]) { int i; for (i = 0; i < N; i++) printf("%14.10f\t", row[i]); printf("\n"); } // 正規化 (ベクトルの長さを１にする) void normarize(double x[N]) { int i; double s = 0.0; for (i = 0; i < N; i++) s += x[i] * x[i]; s = sqrt(s); for (i = 0; i < N; i++) x[i] /= s; }
xxxxxx@yyyyyy /Z $ gcc -o OC1101 OC1101.m -lobjc -lgnustep-base -I $INCLUDE -L $LIB $CFLAGS xxxxxx@yyyyyy /Z $ OC1101 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
D
import std.stdio; import std.math; const int N = 4; // ベキ乗法で最大固有値を求める void main(string[] args) { double a[N][N] = [[5.0, 4.0, 1.0, 1.0], [4.0, 5.0, 1.0, 1.0], [1.0, 1.0, 4.0, 2.0], [1.0, 1.0, 2.0, 4.0]]; double x[N] = [1.0 ,0.0 ,0.0 ,0.0]; // ベキ乗法 double lambda = power(a, x); writefln("\neigenvalue"); writefln("%14.10f", lambda); writefln("eigenvector"); disp_vector(x); } // ベキ乗法 double power(double[N][N] a, ref double[N] x0) { double lambda = 0.0; // 正規化 (ベクトル x0 の長さを１にする) normarize(x0); double e0 = 0.0; foreach (i; 0..N) e0 += x0[i]; foreach (k; 1..200) { // １次元配列を表示 writef("%3d\t", k); disp_vector(x0); // 行列の積 x1 = A × x0 double[N] x1 = [0.0 ,0.0 ,0.0 ,0.0]; foreach (i; 0..N) foreach (j; 0..N) x1[i] += a[i][j] * x0[j]; // 内積 double p0 = 0.0; double p1 = 0.0; foreach (i; 0..N) { p0 += x1[i] * x1[i]; p1 += x1[i] * x0[i]; } // 固有値 lambda = p0 / p1; // 正規化 (ベクトル x1 の長さを１にする) normarize(x1); double e1 = 0.0; // 収束判定 foreach (i; 0..N) e1 += x1[i]; if (fabs(e0 - e1) < 0.00000000001) break; foreach (i; 0..N) x0[i] = x1[i]; e0 = e1; } return lambda; } // １次元配列を表示 void disp_vector(double[N] row) { foreach(col; row) writef("%14.10f\t", col); writefln(""); } // 正規化 (ベクトルの長さを１にする) void normarize(ref double[N] x) { double s = 0.0; foreach (i; 0..N) s += x[i] * x[i]; s = sqrt(s); foreach (i; 0..N) x[i] /= s; }
Z:\>dmd D1101.d Z:\>D1101 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
Go
package main import "fmt" import "math" const N = 4 // ベキ乗法で最大固有値を求める func main() { var a [N][N]float64 = [N][N]float64{{5.0, 4.0, 1.0, 1.0}, {4.0, 5.0, 1.0, 1.0}, {1.0, 1.0, 4.0, 2.0}, {1.0, 1.0, 2.0, 4.0}} var x = []float64 {1.0 ,0.0 ,0.0 ,0.0} // ベキ乗法 var lambda float64 = power(a, x) fmt.Println("\neigenvalue") fmt.Printf("%14.10f\n", lambda) fmt.Println("eigenvector") disp_vector(x) } // ベキ乗法 func power(a[N][N]float64, x0[]float64) float64 { var lambda float64 = 0.0 // 正規化 (ベクトル x0 の長さを１にする) normarize(x0) var e0 float64 = 0.0 for i := 0; i < N; i++ { e0 += x0[i] } for k := 1; k < 200; k++ { // １次元配列を表示 fmt.Printf("%3d\t", k) disp_vector(x0) // 行列の積 x1 = A × x0 var x1 = []float64 {0.0 ,0.0 ,0.0 ,0.0} for i := 0; i < N; i++ { for j := 0; j < N; j++ { x1[i] += a[i][j] * x0[j] } } // 内積 var p0 float64 = 0.0 var p1 float64 = 0.0 for i := 0; i < N; i++ { p0 += x1[i] * x1[i] p1 += x1[i] * x0[i] } // 固有値 lambda = p0 / p1 // 正規化 (ベクトル x1 の長さを１にする) normarize(x1) var e1 float64 = 0.0 // 収束判定 for i := 0; i < N; i++ { e1 += x1[i] } if (math.Abs(e0 - e1) < 0.00000000001) { break } for i := 0; i < N; i++ { x0[i] = x1[i] } e0 = e1 } return lambda } // １次元配列を表示 func disp_vector(row[]float64) { for _, col := range row { fmt.Printf("%14.10f\t", col) } fmt.Println("") } // 正規化 (ベクトルの長さを１にする) func normarize(x[]float64) { var s float64 = 0.0 for i := 0; i < N; i++ { s += x[i] * x[i] } s = math.Sqrt(s) for i := 0; i < N; i++ { x[i] /= s } }
Z:\>go run GO1101.go 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
Scala
object Scala1101 { val N = 3 // ベキ乗法で最大固有値を求める def main(args: Array[String]) { var a:Array[Array[Double]] = Array(Array(5.0, 4.0, 1.0, 1.0) ,Array(4.0, 5.0, 1.0, 1.0) ,Array(1.0, 1.0, 4.0, 2.0) ,Array(1.0, 1.0, 2.0, 4.0)) var x:Array[Double] = Array(1.0 ,0.0 ,0.0 ,0.0) // ベキ乗法 val (lambda, x1) = power(a, x) println() println("eigenvalue") println("%14.10f".format(lambda)) println("eigenvector") disp_vector(x1) } // ベキ乗法 def power(a:Array[Array[Double]], x:Array[Double]): (Double, Array[Double]) = { var x0:Array[Double] = normarize(x) var x1:Array[Double] = Array(0.0 ,0.0 ,0.0 ,0.0) var lambda:Double = 0.0 var e:Double = 0.0 var k = 1 do { e = x0.sum // １次元配列を表示 print("%3d\t".format(k)) disp_vector(x0) k += 1 // 行列の積 x1 = A × x0 x1 = Array(0.0 ,0.0 ,0.0 ,0.0) for (i <- 0 to N) for (j <- 0 to N) x1(i) += a(i)(j) * x0(j) // 内積 var (t1, t2) = (x0 zip x1).map(t => t match {case (f, s) => (s * s, s * f)}).unzip // 固有値 lambda = t1.sum / t2.sum // 正規化 (ベクトル x1 の長さを１にする) x1 = normarize(x1) x0 = x1.clone } while (k <= 100 && Math.abs(e - x1.sum) >= 0.00000000001) return (lambda, x0) } // １次元配列を表示 def disp_vector(row:Array[Double]) = { row.foreach { col => print("%14.10f\t".format(col)) } println() } // 正規化 (ベクトルの長さを１にする) def normarize(x:Array[Double]):Array[Double] = { var s:Double = x.map(n => n * n).sum s = Math.sqrt(s) x.map(n => n / s) } }
Z:\>scala Scala1101.scala 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 10.0000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660
F#
module Fs1101 open System let N = 3 // １次元配列を表示 let disp_vector (row:float[]) = row |> Array.iter (fun x -> printf "%14.10f" x) printfn "" // 正規化 (ベクトルの長さを１にする) let normarize(x:float[]):float[] = let s = x |> Array.map (fun n -> n * n) |> Array.sum x |> Array.map(fun n -> n / Math.Sqrt(s)) // ベキ乗法 let power (a:float[][]) (x:float[]) : float * float[] = let mutable x0:float[] = normarize x let mutable x1:float[] = [|0.0 ;0.0 ;0.0 ;0.0|] let mutable lambda = 0.0 let mutable e = 0.0 let mutable k = 1 let mutable sw = true while sw || (k <= 100 && Math.Abs(e - (x1 |> Array.sum)) >= 0.00000000001) do if sw then sw <- false e <- x0 |> Array.sum // １次元配列を表示 printf "%3d\t" k disp_vector x0 k <- k + 1 // 行列の積 x1 = A × x0 x1 <- [|0.0 ;0.0 ;0.0 ;0.0|] for i in [0..N] do for j in [0..N] do x1.[i] <- x1.[i] + a.[i].[j] * x0.[j] // 内積 let t1, t2 = Array.zip x0 x1 |> Array.map(fun t -> match t with | (f, s) -> (s * s, s * f)) |> Array.unzip // 固有値 lambda <- (t2 |> Array.sum) / (t1 |> Array.sum) // 正規化 (ベクトル x1 の長さを１にする) x1 <- normarize x1 x0 <- Array.copy x1 (lambda, x0) // ベキ乗法で最大固有値を求める let mutable a:float[][] = [| [|5.0; 4.0; 1.0; 1.0|]; [|4.0; 5.0; 1.0; 1.0|]; [|1.0; 1.0; 4.0; 2.0|]; [|1.0; 1.0; 2.0; 4.0|] |] let mutable x:float[] = [|1.0 ;0.0 ;0.0 ;0.0|] // ベキ乗法 let (lambda, x1) = power a x printfn "" printfn "eigenvalue" printfn "%14.10f" lambda printfn "eigenvector" disp_vector x1 exit 0
Z:\>fsi --nologo --quiet Fs1101.fs 1 1.0000000000 0.0000000000 0.0000000000 0.0000000000 2 0.7624928517 0.6099942813 0.1524985703 0.1524985703 3 0.6745979924 0.6589096670 0.2353248811 0.2353248811 4 0.6517382447 0.6501601860 0.2761602732 0.2761602732 5 0.6421032522 0.6419452154 0.2963188771 0.2963188771 6 0.6373267026 0.6373108932 0.3063082594 0.3063082594 7 0.6349074765 0.6349058954 0.3112772079 0.3112772079 8 0.6336860412 0.6336858831 0.3137548428 0.3137548428 9 0.6330719648 0.6330719490 0.3149919005 0.3149919005 10 0.6327640473 0.6327640457 0.3156099832 0.3156099832 省略 26 0.6324555367 0.6324555367 0.3162277566 0.3162277566 27 0.6324555344 0.6324555344 0.3162277613 0.3162277613 28 0.6324555332 0.6324555332 0.3162277637 0.3162277637 29 0.6324555326 0.6324555326 0.3162277648 0.3162277648 30 0.6324555323 0.6324555323 0.3162277654 0.3162277654 31 0.6324555322 0.6324555322 0.3162277657 0.3162277657 32 0.6324555321 0.6324555321 0.3162277659 0.3162277659 33 0.6324555321 0.6324555321 0.3162277659 0.3162277659 34 0.6324555321 0.6324555321 0.3162277660 0.3162277660 35 0.6324555320 0.6324555320 0.3162277660 0.3162277660 eigenvalue 0.1000000000 eigenvector 0.6324555320 0.6324555320 0.3162277660 0.3162277660