scilab ヤコビ法について、errorがでます

質問

ヤコビ法で固有値、固有ベクトルを求めるプログラムを色々調べながら作っていたのですが、errorが出てしまい手詰まりになりました。あまり、scilab自体詳しくなく、多々間違いがあるかもしれませんが、学校でどうしても必要なため、どうか、解答の方をよろしくお願いします。

//初期化
A=[2 1 3;1 0 1;3 1 0]

Pinf=eye(3,3)
//ヤコビ法：Aの対角成分が固有値、 Pinfの行ベクトルが固有ベクトルになる
for i=1:21,
i
if abs(A(1,2))>abs(A(1,3)) then
theta=atan(2*A(1,2))/(A(1,2)-A(1,3))/2
P=[cos(theta),sin(theta);-sin(theta),cos(theta)];
elseif abs(A(1,3))>abs(A(2,3)) then
theta=atan(2*A(1,3))/(A(1,3)-A(2,3))/2
P=[cos(theta),sin(theta);-sin(theta),cos(theta)];
else theta=atan(2*A(2,3)/(A(2,3)-A(1,2)))/2
P=[cos(theta),sin(theta);-sin(theta),cos(theta)];
end
Pinf=P*Pinf
A=P*A*P'
,end

現在、下から4行目のendを入れなければ、end または 「else が抜けています...」とerrorが出てしまい、また、endを入れると、その下の　Pinf=P*Pinf　の所で、「一貫性がない掛け算です。」とerrorがでます。
以上が、現在困っている状況になります。
よろしくおねがいします。

回答No.3

f272 · ベストアンサー

あなたが書いたものからif文に関係ありそうなところを抜き出すと
if abs(A(1,2)) > abs(A(1,3)) then if abs(A(1,2)) > abs(A(2,3) then
elseif abs(A(1,2)) > abs(A(1,3)) then if abs(A(2,3)) > abs(A(1,2) then
elseif abs(A(1,3))>abs(A(2,3)) then if abs(A(1,3)) > abs(A(1,2)) then
else
end
となるけど，これで正しい文法に従っていると思うか？
多分，意図しているのはこれであろう。

if abs(A(1,2)) >= abs(A(1,3)) & abs(A(1,2)) >= abs(A(2,3)) then
elseif abs(A(2,3)) >= abs(A(1,3)) & abs(A(2,3)) >= abs(A(1,2)) then
elseif abs(A(1,3)) >= abs(A(2,3)) & abs(A(1,3)) >= abs(A(1,2)) then
else //ここには来ない
end

お礼

ありがとうございます。
おかげさまで、無事、プログラムを完成させることができました。

回答No.2

f272

少し修正。

> Pinfの行ベクトルが固有ベクトルになる

ということですから，「Pinf=P*Pinf」で問題なかった。
でも，Pinfの列ベクトルが固有ベクトルになるようにした方がいいんじゃないかなあ？

お礼

一応、考え直した結果、ある程度進展はしたと想うのですが・・・また、「end またはelseが抜けています」とerrorが出てしまいました・・・。何度も、申し訳ないのですが、もし、よろしければ、また間違えているところの指摘をお願いします。
↓がプログラムです。

//初期化
A=[2 1 3;1 0 1;3 1 0]

Pinf=eye(3,3)
//ヤコビ法：Aの対角成分が固有値、 Pinfの行ベクトルが固有ベクトルになる
for i=1:20,
i
if abs(A(1,2)) > abs(A(1,3)) then if abs(A(1,2)) > abs(A(2,3) then theta=atan(2*A(1,2))/(A(1,1)-A(2,2))/2,
P=[cos(theta),sin(theta),0;-sin(theta),cos(theta),0;0,0,1];
elseif abs(A(1,2)) > abs(A(1,3)) then if abs(A(2,3)) > abs(A(1,2) then theta=atan(2*A(2,3))/(A(2,2)-A(3,3))/2,
P=[1,0,0;0,cos(theta),sin(theta);0,-sin(theta),cos(theta)]
elseif abs(A(1,3))>abs(A(2,3)) then if abs(A(1,3)) > abs(A(1,2)) then theta=atan(2*A(1,3))/(A(1,1)-A(3,3))/2,
P=[cos(theta),0,sin(theta);0,1,0;-sin(theta),0,cos(theta)];
else theta=atan(2*A(1,2)/(A(1,1)-A(2,2)))/2,
P=[1,0,0;0,cos(theta),sin(theta);0,-sin(theta),cos(theta)];
end
Pinf=P*Pinf
A=P*A*P'
,end

回答No.1

f272

ヤコビ法なら非対角要素で絶対値最大のものを求めるのだが，そうなっていない。
3*3行列の計算をしたいはずだが，Pは2*2行列になっている。
そのうえthetaの計算がヤコビ法に従っていない。
Pinfの更新の計算もヤコビ法に従っていない。
もうどうしようもないくらいおかしなことばかりです。

お礼

ご指摘、ありがとうございます。
指摘を参考にもう一度考え直してみます。