回覆列表
  • 1 # 我是阿嘛

    function [Y] = RK45(t,X,f,h)

    K1=f(t,X);

    K2=f(t+h/2,X+h/2*K1);

    K3=f(t+h/2,X+h/2*K2);

    K4=f(t+h,X+h*K3);

    Y=X+h/6*(K1+2*K2+2*K3+K4);

    end

    以上是4階龍格庫塔法的程式碼:

    自己寫函式,存為f.m

    function dxdt = f (t,x)

    dxdt(1)=exp(x(1)*sin(t))+x(2);

    dxdt(2)=exp(x(2)*cos(t))+x(1); % x(1)是你的f,x(2)是你的g

    dxdt=dxdt(:);

    end

    自己給出t0,x0,h的值(初始時間,初值,步長)

    如果求t0到t1的軌跡的話:給個例子如下

    t0=0;t1=5;h=0.02;x0=[-1;-1];

    T=t0:h:t1;X=zeros(length(x0),length(T));X(:,1)=x0;

    for j=1:length(T)-1

    X(:,j+1)=RK45(T(j),X(:,j),@(t,x) f(t,x),h);

    end

    plot(T,X(1,:));

    hold on;

    plot(T,X(2,:),"r");

    具體引數自己設定

  • 中秋節和大豐收的關聯?
  • 叫你早睡的人和陪你熬夜的人,你會選擇哪一個?