function [X,Y]=odeRK4(str,a,b,ya,n)
% 四阶龙格-库塔方法
%f 函数
%a和b 求解区域的端点
%ya 初始条件y(a)
%n 求解步数
%f 如果是文件，调用[x,y]=odeRK4(@f,a,b,ya,n).
%f 如果是匿名函数，调用[x,y]=odeRK4(f,a,b,ya,n).
h=(b-a)/n;
Y=zeros(1,n+1);
X=a:h:b;
Y(1)=ya;
f_xy= @(x,y)str;
for j=1:n
k1=f_xy(X(j),Y(j));
k2=f_xy(X(j)+h/2,Y(j)+h/2*k1);
k3=f_xy(X(j)+h/2,Y(j)+h/2*k2);
k4=f_xy(X(j)+h,Y(j)+h*k3);
Y(j+1)=Y(j)+h/6*(k1+k4)+h/3*(k2+k3);
end
for k=1:n+1
fprintf('x[%d]=%f\ty[%d]=%f\n',k-1,X(k),k-1,Y(k));
end
end
clear all
figure('name','四阶龙格-库塔', 'NumberTitle','off') %建立名字是四阶龙格-库塔的图形窗口
uicontrol(gcf,'style','text', 'position',[50 350 350 50], 'string','常微分方程组初值问题数值解', 'fontsize',16);
% 建立名为常微分方程组初值问题数值解的求解的静态编辑框
text_fun1= uicontrol(gcf, 'Style','Text','String','常微分方程表达式: ','FontSize',14,'Position',[55 278 160 25],'HorizontalAlignment','Left');
edit_fun1= uicontrol(gcf, 'Style','Edit','String','-x*x*y*y*y','Position',[220 278 140 25],'HorizontalAlignment','left', 'FontSize',14);
text_fun2= uicontrol(gcf, 'Style','Text','String','区间左端点: ','FontSize',14,'Position',[55 248 160 25],'HorizontalAlignment','Left');
edit_fun2= uicontrol(gcf,'style','Edit','String','0', 'position',[220 248 140 25],'HorizontalAlignment','left', 'FontSize',14);
text_fun3= uicontrol(gcf, 'Style','Text','String','区间右端点: ','FontSize',14,'Position',[55 218 160 25],'HorizontalAlignment','Left');
edit_fun3= uicontrol(gcf,'style','Edit','String','5', 'position',[220 218 140 25],'HorizontalAlignment','left', 'FontSize',14);
text_fun4= uicontrol(gcf, 'Style','Text','String','初值: ','FontSize',14,'Position',[55 188 160 25],'HorizontalAlignment','Left');
edit_fun4= uicontrol(gcf,'style','Edit','String','1', 'position',[220 188 140 25],'HorizontalAlignment','left', 'FontSize',14);
text_fun5= uicontrol(gcf, 'Style','Text','String','求解步数: ','FontSize',14,'Position',[55 158 160 25],'HorizontalAlignment','Left');
edit_fun5= uicontrol(gcf,'style','Edit','String','20', 'position',[220 158 140 25],'HorizontalAlignment','left', 'FontSize',14);
uicontrol(gcf,'style','push', 'position',[55 108 150 30], 'string','开始运行', 'fontsize',10, 'call',['odefun =num2str(get(edit_fun1,''string''));a = str2num(get(edit_fun2,''string'')),b = str2num(get(edit_fun3,''string'')),ya = str2num(get(edit_fun4,''string'')),n = str2num(get(edit_fun5,''string''));','[x,y]=odeRK4(odefun,a,b,ya,n)'])
% 建立名为开始运行的静态编辑框
push_Quit = uicontrol(gcf, 'Style','Push','String','Quit','value',0,'FontSize',10, 'Pos',[450 50 100 30], 'Call', 'fun_quit');
运行报错
在赋值 A(:) = B 中，A 和 B 中的元素数目必须相同。
出错 odeRK4 (line 19)
Y(j+1)=Y(j)+h/6*(k1+k4)+h/3*(k2+k3);
计算 UIControl Callback 时出错
