|
|
> No.3419[元記事へ]
島村1243さんへのお返事です。
> Amとφを入力値とし、
> Am=0.1~1.0 step 0.1
> φ=PI/18（パラメータで、後でいろいろと変える）
> 下記の連立方程式を満足する出力値Edとθを求めるプログラムをご教示ください。
これで解を求められると思います。
!2元非線形連立方程式の解
!ニュートン・ラフソン（Newton-Raphson）法
! 1変数の場合、Xi+1=Xi-f(Xi)/f'(Xi) 変形して、Xi+1=Xi+(-f(Xi)/f'(Xi))
! ここで、漸化式の増分部分⊿x=-f(Xi)/f'(Xi)は、f'(Xi)*⊿x=-f(Xi)となる。
! 行列表記で解釈すると、Jx=bとなる。
!
!これを多変数に拡張する。
FUNCTION f(Ed,θ)
LET IS=SQR((Am/SQR(2)*Ed*SIN(φ))^2+(Vs-Am/SQR(2)*Ed*COS(φ))^2)/Xs
LET f=Ed^2/R -Vs*IS*COS(θ)
END FUNCTION
FUNCTION g(Ed,θ)
LET IS=SQR((Am/SQR(2)*Ed*SIN(φ))^2+(Vs-Am/SQR(2)*Ed*COS(φ))^2)/Xs
LET g=(Am/SQR(2)*Ed*SIN(θ+φ)-Xs*IS) -Vs*SIN(θ)
END FUNCTION
LET Am=0.5 !0＜Am≦1
LET φ=PI/18 !0＜φ＜π/2
LET R=100
LET Xs=PI
LET Vs=100
LET d=0.01 !∂
LET xi=0 !初期値
LET yi=0
FOR i=1 TO 50 !漸化式を収束させる
DIM x(2) !⊿x、⊿y
!連立１次方程式Jx=bを得る
DIM J(2,2)
LET fxy=f(xi,yi) !前進差分で近似して求める（平均変化率）
LET J(1,1)=(f(xi+d,yi)-fxy)/d !J=(∂f/∂x ∂f/∂y) ヤコビ(Jacobi)行列
LET J(1,2)=(f(xi,yi+d)-fxy)/d ! (∂g/∂x ∂g/∂y)
LET gxy=g(xi,yi)
LET J(2,1)=(g(xi+d,yi)-gxy)/d
LET J(2,2)=(g(xi,yi+d)-gxy)/d
DIM b(2)
LET b(1)=-fxy !b=(-f)
LET b(2)=-gxy ! (-g)
!---------- ※ガウスの消去法などで２元連立１次方程式を解く
DIM Ji(2,2) !逆行列を求める Full BASIC版
MAT Ji=INV(J) !正則なら
MAT x=Ji*b
!----------
LET xi=xi+x(1) !Xi+1=Xi+⊿x
LET yi=yi+x(2) !Yi+1=Yi+⊿y
PRINT xi,yi !確認 ∥Xi+1 - Xi∥≦ε(1+∥Xi+1∥)
NEXT i
PRINT f(xi,yi),g(xi,yi) !確認 ∥f(xi)∥<δ
END
|
|