0

我正在为我的一个班级做作业,我应该使用我选择的程序(我选择 Matlab)编写代码,以使用四阶龙格-库塔方法求解贝塞尔函数微分方程。作为参考,贝塞尔函数 DE 是:

x^2*(J_n)''+x*(J_n)'+(x^2-n^2)*J_n=0。

我可以通过以下方式将其分成两个耦合的一阶 DE:

(J_n)'=Z_n 和

(Z_n)'+(1/x)*Z_n+[(x^2-n^2)/x^2]*J_n=0。

在此任务之前,我没有使用 Matlab 或任何其他编程语言的经验。我知道 Matlab 有“ode45”命令,但我应该自己编写代码,而不是依赖于 Matlab 的命令。到目前为止,我一直在研究 Bessel 函数的 n=0 情况,但是当我尝试绘制函数时,我一直收到错误消息。我说的当前错误是:“'double' 类型的输入参数的未定义函数或方法 'J'。” 但我不知道如何解决这个错误,也不知道我的代码是否正确。有人可以告诉我哪里出了问题或者编写此代码的正确方法是什么?

h=0.01;             %step size

J_0(1)=1;           %initial condition for J_0

Z_0(1)=1;           %initial condition for Z_0-This value should be zero
                    %but Matlab gives me an error. To fix this, I input
                    %Z_0(1)-1 to use the correct value for Z_0(1).

x(1)=0.001;         %first value of x

dZ(Z_0,J_0)=(-1/x)*(Z_0-1)-J_0;

for i=[1:1:10]

    dZ1=(-1/x)*(Z_0-1)-J_0;
    dJ1=(Z_0(1)-1)*h;

    dZ2=(-1/x)*(Z_0-1+0.5*h)-(J_0+0.5*h*dJ1);
    dJ2=((Z_0(1)-1)+dZ1)*h;

    dZ3=(-1/x)*(Z_0-1+0.5*h)-(J_0+0.5*h*dJ2);
    dJ3=((Z_0(1)-1)+dZ1+dZ2)*h;

    dZ4=(-1/x)*(Z_0-1+h)-(J_0+h*dJ3);
    dJ4=((Z_0(1)-1)+dZ1+dZ2+dZ3)*h;

    J(i+1)=J(i)+(h/6)*(dJ1+2*dJ2+2*dJ3+dJ4);
end

plot(J_0);

提前感谢您的帮助

4

1 回答 1

0

你的问题是在线:

J(i+1)=J(i)+(h/6)*(dJ1+2*dJ2+2*dJ3+dJ4);

在赋值运算符的右侧,您使用J以前从未设置过的变量i取值 1。对我来说看起来像是一个错字(应该J_0改为吗?)

另外,在循环中i计算你的dJdZ东西时不要忘记你的索引。for

于 2013-04-13T22:38:29.880 回答