- 締切済み
Octaveでのオイラー法とルンゲクッタ法について
講義ノートを見ながら数値解析ソフトOctaveの操作を勉強しています。 読んでいる講義ノートにオイラー法とルンゲ・クッタ法でf''(x)+f'(x)+x=0の微分方程式を解く箇所があるのですが、書いている通りに打ち込むとx=0の直線になってしまいます。lsodeで解くとうまくいくのですが。学び始めたばかりの初心者でまったく分からないので、、どうすればいいか教えてください。以下が講義ノートに書いてあったプログラムです。xdot2,euler,rk4は定義された関数で、rk4eulerがメインプログラムです。 % xdot2.m function xdot2 = xdot2(x,t) xdot2=[x(2) -x(2)-x(1)]; % euler.m function x=euler(F,x,t) [m n]=size(t); %t のサイズ(=ステップ数)を抽出 for k=1:m-1 %tのサイズ-1回繰り返す。 tau=t(k+1)-t(k); x(k+1,:)=x(k,:)+tau.*feval(F,x(k,:),t(k)); %上の式がオイラー法の公式 end % rk4.m function x=rk4(F,x,t) [m n]=size(t); for k=1:m-1 %t のサイズ-1 だけステップを繰り返す。 tau = t(k+1)-t(k); %刻み幅 F1=feval(F,x(k,:),t(k)); F2=feval(F,x(k,:)+1/2*tau.*F1,t(k)+1/2*tau); F3=feval(F,x(k,:)+1/2*tau.*F2,t(k)+1/2*tau); F4=feval(F,x(k,:)+tau*F3,t(k)+tau); x(k+1,:)= x(k,:)+tau/6*(F1 +2*F2 +2*F3 +F4); %以上が4次のRK法の公式 end %rk4euler.m %オイラー、RK、lsode の比較 t=[0:0.2:10]’; x0=[0 1]; %定義域、初期値の定義 R1=euler('xdot2',x0,t); x1=R1(:,1); %オイラー法 R2= rk4('xdot2',x0,t); x2=R2(:,1); %RK 法 R3=lsode('xdot2',x0,t); x3=R3(:,1); %lsode plot(t,x1,'rx-.',t,x2,'bx',t,x3,'ko'); xlabel('t'); ylabel('x'); legend('Euler','RK4','lsode'); title('rk4euler1.m');
- みんなの回答 (1)
- 専門家の回答
みんなの回答
- f272
- ベストアンサー率46% (8469/18132)
おかしなところは rk4euler.mの2行目で t=[0:0.2:10]’; x0=[0 1]; %定義域、初期値の定義 に変な文字「’」がある。 これを「'」に変えれば,ちゃんと動く。