?一階問題舉例：

高階問題舉例 ：

?

?

常微分方程數(shù)值解：向前歐拉方法之一階問題

clc,clear,close all;
a=0;%初始時刻
b=2*pi;%結(jié)束時刻
n=100;%離散點數(shù)量
x0=0;%初值
h=(b-a)/n;%步長
x=x0 + [0:n]*h;%離散點數(shù)組
fun=inline('sin(x)+y','x','y');
y0=1;

%計算
y(1)=y0;
for i=1:n?????? %迭代
??? f=feval(fun,x(i),y(i));
??? y(i+1)=y(i)+h*f;
end
plot(x,y,'o');

hold on

n=1000;
h=(b-a)/n;%步長
x=x0 + [0:n]*h;%離散點數(shù)組
for i=1:n
??? f=feval(fun,x(i),y(i));
??? y(i+1)=y(i)+h*f;
end
plot(x,y,'-');

常微分方程數(shù)值解：向前歐拉方法之高階問題

clc,clear,close all;
a=1;%初始時刻
b=10;%結(jié)束時刻
n=100;%離散點數(shù)量
x0=1;%初值
h=(b-a)/n;%步長
x=x0 + [0:n]*h;%離散點數(shù)組
y0=[1,10,30];

%計算
y(1,:)=y0';
for i=1:n?????? %迭代
??? f=feval('fun_3rd',x(i),y(i,:));%x是標(biāo)量，y是向量
%feval()函數(shù)執(zhí)行指定的函數(shù)。
%將想要執(zhí)行的函數(shù)以及相應(yīng)的參數(shù)一起作為feval()的參數(shù)，
%feval()的輸出等于想要執(zhí)行的函數(shù)的輸出。
??? f=f';
??? y(i+1,:)=y(i,:)+h*f;
end
plot(x,y(:,1),'o');%所有的行的第一列

龍格庫塔方法（精確度高）：一階問題

clc,clear,close all; a=0;%初始時刻 b=2*pi;%結(jié)束時刻 n=1000;%離散點數(shù)量 x0=0;%初值 h=(b-a)/n;%步長 x=x0 + [0:n]*h;%離散點數(shù)組 fun=inline('sin(x)+y','x','y'); y0=1;%計算 y(1)=y0; for i=1:n %迭代k1=feval(fun,x(i),y(i));k2=feval(fun,x(i)+h/2,y(i)+h*k1/2);k3=feval(fun,x(i)+h/2,y(i)+h*k2/2);k4=feval(fun,x(i)+h,y(i)+h*k3);y(i+1)=y(i)+h*(k1+2*k2+2*k3+k4)/6; end plot(x,y,'o');

龍格庫塔方法（精確度高）：高階問題

%龍格庫塔方法 %高階問題,內(nèi)置函數(shù) %[a,b]=ode45('fun_3rd',[1,10],[1 10 30])%x的范圍，y的初值 %plot(a,b(:,1)) clc,clear,close all; a=1;%初始時刻 b=10;%結(jié)束時刻 n=100;%離散點數(shù)量 x0=1;%初值 h=(b-a)/n;%步長 x=x0 + [0:n]*h;%離散點數(shù)組 y0=[1,10,30];%計算 y(1,:)=y0'; for i=1:n %迭代k1=feval('fun_3rd',x(i),y(i,:));k1=k1';k2=feval('fun_3rd',x(i)+h/2,y(i,:)+h*k1/2);k2=k2';k3=feval('fun_3rd',x(i)+h/2,y(i,:)+h*k2/2);k3=k3';k4=feval('fun_3rd',x(i)+h,y(i,:)+h*k3);k4=k4';y(i+1,:)=y(i,:)+h*(k1+2*k2+2*k3+k4)/6; endplot(x,y(:,1),'o');

改進(jìn)歐拉方法：一階問題

%向前歐拉法的改進(jìn)方法 %改進(jìn)歐拉方法 %一階問題 clc,clear,close all; a=0;%初始時刻 b=2*pi;%結(jié)束時刻 n=100;%離散點數(shù)量 x0=0;%初值 h=(b-a)/n;%步長 x=x0 + [0:n]*h;%離散點數(shù)組 fun=inline('sin(x)+y','x','y'); y0=1;%計算 y(1)=y0; for i=1:n %迭代k1=feval(fun,x(i),y(i));k2=feval(fun,x(i)+h,y(i)+h);y(i+1)=y(i)+h/2*(k1+k2); end plot(x,y,'o');

改進(jìn)歐拉方法：高階問題

%高階問題 clc,clear,close all; a=1;%初始時刻 b=10;%結(jié)束時刻 n=1000;%離散點數(shù)量 x0=1;%初值 h=(b-a)/n;%步長 x=x0 + [0:n]*h;%離散點數(shù)組 y0=[1,10,30];%計算 y(1,:)=y0'; for i=1:n %迭代f1=feval('fun_3rd',x(i),y(i,:));%x是標(biāo)量，y是向量f1=f1';f2=feval('fun_3rd',x(i)+h,y(i,:)+h);f2=f2';y(i+1,:)=y(i,:)+h/2*(f1+f2); end plot(x,y(:,1),'o');%所有的行的第一列

偏微分方程數(shù)值解：一維熱傳導(dǎo)方程

%偏微分方程數(shù)值解 %1.空間有多大 %2.時間有多長 %邊界條件和初值 %顯示向前歐拉算法 clc,clear,close all; A=0.5; xf=2;%邊界點 T=0.1;%最后到達(dá)的時刻 M=50;%空間上分段，空間離散個數(shù) N=100;%時間上分段，時間離散個數(shù) fun0=inline('sin(pi*x)','x');%初值 bx0=inline('0');%左邊界條件 bxf=inline('0');%右邊界條件 dx=(xf-0)/M;%空間步長 dt=(T-0)/N;%時間步長 x=[0:M]'*dx;%空間離散點 t=[0:N]'*dt;%時間離散點 r=A*dt/dx/dx; r1=1-2*r;for i=1:M+1 %離散初始時刻的值u(i,1)=feval(fun0,x(i)); endfor n=1:N+1 %邊界條件u([1,M+1],n)=[bx0(t(n));bxf(t(n))]; endfor k=1:Nfor i=2:Mu(i,k+1)=r*(u(i+1,k)+u(i-1,k))+r1*u(i,k);end endmesh(t,x,u)

總結(jié)

以上是生活随笔為你收集整理的微分方程数值解的全部內(nèi)容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網(wǎng)站內(nèi)容還不錯，歡迎將生活随笔推薦給好友。