MATLABEuler法解常微分方程
欧拉法(euler)求解常微分方程的matlab程序及案例

欧拉法(euler)求解常微分方程的matlab程序及案例欧拉方法是最初用于求解常微分方程的数值方法之一,它是一种显式的一步法,具有易于实施的优点,特别适合初学者使用。
本文将介绍欧拉法的原理和使用MATLAB求解常微分方程的具体方法,同时给出一个简单的实例进行说明。
一、欧拉法原理考虑一个一阶常微分方程y'=f(t,y),欧拉法的基本思想是将时间步长Δt均分成n个小步长,从y(t0)开始依次计算每个时刻的值,得到一列估计值y1, y2, …, yn。
欧拉法的计算公式为:(1)y1=y(t0+Δt)=y(t0)+Δtf(t0, y0)(2)y2=y(t0+2Δt)=y(t0+Δt)+Δtf(t0+Δt, y1)(3)yn=y(t0+nΔt)=y(t0+(n-1)Δt)+Δtf(t0+(n-1)Δt, yn-1)可以看出,欧拉法的核心在于利用已知的t和y计算f(t,y),从而获得y的逼近值。
但是需要注意的是,步长Δt越小,计算所需的时间和内存就越多,而精度却并不一定提高。
因此在实际应用中需要结合具体问题选择合适的步长。
二、MATLAB求解常微分方程的具体方法(1)定义常微分方程我们以一个简单的例子开始,考虑求解y'=1-y,y(0)=0.5在[0,1]区间内的积分。
首先定义匿名函数dydt,将其传到ode45中求解:dydt=@(t,y)1-y;[t,y]=ode45(dydt,[0 1],0.5);plot(t,y,'-o')运行以上代码可以得到结果,其中plot函数用于绘制图像。
但是,由于求解过程中计算机执行到ode45函数时可能需要很长时间,因此需要更快捷的方法。
(2)利用欧拉法求解方程欧拉法求解方程首先需要定义步长Δt,这里设Δt为0.1。
定义起始值y=[0.5]、时间向量t=0:Δt:1,然后计算列向量y的估计值:t=0:0.1:1;y=zeros(size(t));y(1)=0.5;for n=1:length(t)-1y(n+1)=y(n)+0.1*(1-y(n));endplot(t,y,'-o')以上代码的执行结果与前面的ode45方法相同,但是速度更快。
matlab中的向后euler方法

文章标题:深入理解Matlab中的向后Euler方法在数值计算中,求解常微分方程(ODE)是一个十分常见也是重要的任务。
其中,向后Euler方法是一种常用的数值求解方法之一。
在Matlab中,我们可以通过调用内置函数ode15s来使用向后Euler方法来求解ODE。
在本文中,我们将深入探讨Matlab中使用向后Euler方法求解ODE的原理和应用,并通过具体例子来展示其在实际工程中的价值。
1. 向后Euler方法的原理和特点向后Euler方法是一种隐式数值求解方法,其基本思想是将微分方程转化为差分方程,然后通过迭代的方式逼近实际解。
与前向Euler方法相比,向后Euler方法具有更好的稳定性和收敛性,特别在处理刚性ODE的时候表现更为突出。
在Matlab中,我们可以使用ode15s函数来调用向后Euler方法进行数值求解,其使用形式为[y, t] =ode15s(@(t,y)odefun,tspan,y0,options)。
2. 向后Euler方法在实际工程中的应用在实际工程中,ODE求解是一个非常重要的任务。
在控制系统设计中,经常需要求解系统的状态方程以完成系统设计和性能评估;在仿真和建模中,也需要对系统进行数值求解以获得系统的动态响应。
在这些应用中,向后Euler方法常常被用来求解刚性ODE,以获得更为准确和稳定的结果。
3. 使用Matlab进行向后Euler方法的数值求解在Matlab中,使用向后Euler方法进行数值求解非常方便。
通过调用ode15s函数,我们可以通过简单的参数设置来实现对ODE的求解。
下面通过一个简单的例子来展示如何使用Matlab中的向后Euler方法来求解ODE。
考虑一个一阶常微分方程dy/dt = -2y,初始条件为y(0) = 1。
我们可以使用Matlab中的ode15s函数来对其进行数值求解。
具体代码如下:```matlabfunction yprime = odefun(t,y)yprime = -2*y;end[t, y] = ode15s(@odefun, [0 10], 1);plot(t, y, '-')xlabel('t');ylabel('y(t)');title('Solution of dy/dt = -2y using Backward Euler method');```通过运行以上代码,我们可以得到dy/dt = -2y的数值解。
用MATLAB解常微分方程

实验四 求微分方程的解一、问题背景与实验目的实际应用问题通过数学建模所归纳而得到的方程,绝大多数都是微分方程,真正能得到代数方程的机会很少.另一方面,能够求解的微分方程也是十分有限的,特别是高阶方程和偏微分方程(组).这就要求我们必须研究微分方程(组)的解法,既要研究微分方程(组)的解析解法(精确解),更要研究微分方程(组)的数值解法(近似解).对微分方程(组)的解析解法(精确解),Matlab 有专门的函数可以用,本实验将作一定的介绍.本实验将主要研究微分方程(组)的数值解法(近似解),重点介绍 Euler 折线法.二、相关函数(命令)及简介1.dsolve ('equ1','equ2',…):Matlab 求微分方程的解析解.equ1、equ2、…为方程(或条件).写方程(或条件)时用 Dy 表示y 关于自变量的一阶导数,用用 D2y 表示 y 关于自变量的二阶导数,依此类推.2.simplify(s ):对表达式 s 使用 maple 的化简规则进行化简. 例如: syms xsimplify(sin(x)^2 + cos(x)^2) ans=13.[r,how]=simple(s):由于 Matlab 提供了多种化简规则,simple 命令就是对表达式 s 用各种规则进行化简,然后用 r 返回最简形式,how 返回形成这种形式所用的规则.例如: syms x[r,how]=simple(cos(x)^2-sin(x)^2) r = cos(2*x) how = combine4.[T,Y] = solver(odefun,tspan,y 0) 求微分方程的数值解. 说明:(1) 其中的 solver 为命令 ode45、ode23、ode113、ode15s 、ode23s 、ode23t 、ode23tb 之一.(2) odefun 是显式常微分方程:⎪⎩⎪⎨⎧==00)(),(y t y y t f dt dy(3) 在积分区间 tspan =],[0f t t 上,从0t 到f t ,用初始条件0y 求解.(4) 要获得问题在其他指定时间点 ,210,,t t t 上的解,则令 tspan = ],,,[,210f t t t t (要求是单调的). (5) 因为没有一种算法可以有效地解决所有的 ODE 问题,为此,Matlab 提供了多种求解器 Solver ,对于不同的ODE 问题,采用不同的Solver .(6) 要特别的是:ode23、ode45 是极其常用的用来求解非刚性的标准形式的一阶常微分方程(组)的初值问题的解的 Matlab 的常用程序,其中:ode23 采用龙格-库塔2 阶算法,用3 阶公式作误差估计来调节步长,具有低等的精度.ode45 则采用龙格-库塔4 阶算法,用5 阶公式作误差估计来调节步长,具有中等的精度.5.ezplot (x,y,[tmin,tmax]):符号函数的作图命令.x,y 为关于参数t 的符号函数,[tmin,tmax] 为 t 的取值范围.6.inline():建立一个内联函数.格式:inline('expr', 'var1', 'var2',…) ,注意括号里的表达式要加引号.例:Q = dblquad(inline('y*sin(x)'), pi, 2*pi, 0, pi)三、实验内容1. 几个可以直接用 Matlab 求微分方程精确解的例子:例1:求解微分方程22x xe xy dxdy-=+,并加以验证.求解本问题的Matlab 程序为:syms x y %line1 y=dsolve('Dy+2*x*y=x*exp(-x^2)','x') %line2 diff(y,x)+2*x*y-x*exp(-x^2) %line3 simplify(diff(y,x)+2*x*y-x*exp(-x^2)) %line4 说明:(1) 行line1是用命令定义x,y 为符号变量.这里可以不写,但为确保正确性,建议写上;(2) 行line2是用命令求出的微分方程的解:1/2*exp(-x^2)*x^2+exp(-x^2)*C1(3) 行line3使用所求得的解.这里是将解代入原微分方程,结果应该为0,但这里给出:-x^3*exp(-x^2)-2*x*exp(-x^2)*C1+2*x*(1/2*exp(-x^2)*x^2+exp(-x^2)*C1)(4) 行line4 用 simplify() 函数对上式进行化简,结果为 0, 表明)(x y y =的确是微分方程的解.例2:求微分方程0'=-+x e y xy 在初始条件e y 2)1(=下的特解,并画出解函数的图形.求解本问题的 Matlab 程序为: syms x yy=dsolve('x*Dy+y-exp(x)=0','y(1)=2*exp(1)','x')ezplot(y)微分方程的特解为:y=1/x*exp(x)+1/x* exp (1) (Matlab 格式),即xe e y x+=,解函数的图形如图 1:图1例3:求微分方程组⎪⎪⎩⎪⎪⎨⎧=--=++035y x dtdy e y x dtdx t在初始条件0|,1|00====t t y x 下的特解,并画出解函数的图形.求解本问题的 Matlab 程序为: syms x y t[x,y]=dsolve('Dx+5*x+y=exp(t)','Dy-x-3*y=0','x(0)=1','y(0)=0','t') simple(x); simple(y);ezplot(x,y,[0,1.3]);axis auto微分方程的特解(式子特别长)以及解函数的图形均略. 2. 用ode23、ode45等求解非刚性的标准形式的一阶常微分方程(组)的初值问题的数值解(近似解).例4:求解微分方程初值问题⎪⎩⎪⎨⎧=++-=1)0(2222y xx y dx dy 的数值解,求解范围为区间[0, 0.5].fun=inline('-2*y+2*x^2+2*x','x','y'); [x,y]=ode23(fun,[0,0.5],1); x'; y';plot(x,y,'o-') >> x' ans =0.0000 0.0400 0.0900 0.1400 0.1900 0.2400 0.2900 0.3400 0.3900 0.4400 0.4900 0.5000 >> y'ans =1.0000 0.9247 0.8434 0.7754 0.7199 0.6764 0.6440 0.6222 0.6105 0.6084 0.6154 0.6179 图形结果为图 2.图2例 5:求解描述振荡器的经典的 Ver der Pol 微分方程.7,0)0(',1)0(,0)1(222====+--μμy y y dt dy y dty d 分析:令,,121dt dx x y x ==则.)1(,1221221x x x dtdx x dt dx --==μ 先编写函数文件verderpol.m :function xprime = verderpol(t,x) global mu;xprime = [x(2);mu*(1-x(1)^2)*x(2)-x(1)]; 再编写命令文件vdp1.m : global mu; mu = 7; y0=[1;0][t,x] = ode45('verderpol',[0,40],y0); x1=x(:,1);x2=x(:,2); plot(t,x1)图形结果为图3.图33. 用 Euler 折线法求解前面讲到过,能够求解的微分方程也是十分有限的.下面介绍用 Euler 折线法求微分方程的数值解(近似解)的方法.Euler 折线法求解的基本思想是将微分方程初值问题⎪⎩⎪⎨⎧==00)(),,(y x y y x f dxdy化成一个代数方程,即差分方程,主要步骤是用差商hx y h x y )()(-+替代微商dx dy ,于是:⎪⎩⎪⎨⎧==-+)()),(,()()(00x y y x y x f h x y h x y k k k k 记)(,1k k k k x y y h x x =+=+,从而)(1h x y y k k +=+,则有1,,2,1,0).,(,),(1100-=⎪⎩⎪⎨⎧+=+==++n k y x hf y y h x x x y y k k k k k k 例 6:用 Euler 折线法求解微分方程初值问题⎪⎩⎪⎨⎧=+=1)0(,22y y x y dxdy 的数值解(步长h 取0.4),求解范围为区间[0,2].解:本问题的差分方程为1,,2,1,0).2),( ),(,,4.0,1,021100-=⎪⎪⎪⎩⎪⎪⎪⎨⎧+=+=+====++n k y x y y x f y x hf y y h x x h y x k k k k k k (其中: 相应的Matlab 程序见附录 1. 数据结果为:0 1.0000 0.4000 1.4000 0.8000 2.1233 1.2000 3.1145 1.6000 4.4593 2.0000 6.3074图形结果见图4:图4特别说明:本问题可进一步利用四阶 Runge-Kutta 法求解,读者可将两个结果在一个图中显示,并和精确值比较,看看哪个更“精确”?(相应的 Matlab 程序参见附录 2).四、自己动手1. 求微分方程0sin 2')1(2=-+-x xy y x 的通解.2. 求微分方程x e y y y x sin 5'2''=+-的通解.3. 求微分方程组⎪⎪⎩⎪⎪⎨⎧=-+=++00y x dtdy y x dtdx在初始条件0|,1|00====t t y x 下的特解,并画出解函数()y f x =的图形. 4. 分别用 ode23、ode45 求上述第 3 题中的微分方程初值问题的数值解(近似解),求解区间为[0,2]t ∈.利用画图来比较两种求解器之间的差异.5. 用 Euler 折线法求解微分方程初值问题⎪⎩⎪⎨⎧=-=1)0(,12'32y y x y y 的数值解(步长h 取0.1),求解范围为区间[0,2].6. 用四阶 Runge-Kutta 法求解微分方程初值问题⎩⎨⎧=-=1)0(,cos 'y x e y y x 的数值解(步长h 取0.1),求解范围为区间[0,3].四阶 Runge-Kutta 法的迭代公式为(Euler 折线法实为一阶 Runge-Kutta 法):1,,2,1,0),()2,2()2,2(),()22(6,),(342312143211100-=⎪⎪⎪⎪⎪⎪⎪⎩⎪⎪⎪⎪⎪⎪⎪⎨⎧++=++=++==++++=+==++n k hL y h x f L L h y h x f L L h y h x f L y x f L L L L L hy y h x x x y y k k k k k k k k k k k k 相应的 Matlab 程序参见附录 2.试用该方法求解第5题中的初值问题. 7. 用 ode45 方法求上述第 6 题的常微分方程初值问题的数值解(近似解),从而利用画图来比较两者间的差异.五、附录附录 1:(fulu1.m)clearf=sym('y+2*x/y^2'); a=0; b=2; h=0.4;n=(b-a)/h+1; x=0; y=1;szj=[x,y];for i=1:n-1y=y+h*subs(f,{'x','y'},{x,y});x=x+h;szj=[szj;x,y];endszjplot(szj(:,1),szj(:,2))附录2:(fulu2.m)clearf=sym('y-exp(x)*cos(x)');a=0;b=3;h=0.1;n=(b-a)/h+1;x=0;y=1;szj=[x,y];for i=1:n-1l1=subs(f,{'x','y'},{x,y});l2=subs(f,{'x','y'},{x+h/2,y+l1*h/2});l3=subs(f,{'x','y'},{x+h/2,y+l2*h/2});l4=subs(f,{'x','y'},{x+h,y+l3*h});y=y+h*(l1+2*l2+2*l3+l4)/6;x=x+h;szj=[szj;x,y];endszjplot(szj(:,1),szj(:,2))。
matlab欧拉法解常微分方程

matlab欧拉法解常微分方程matlab欧拉法解常微分方程欧拉法是解决微分方程数值计算的一种基本方法,是通过估算函数图像的变化来得到函数的近似值。
而matlab是一种强大的数值计算软件,也能轻易地实现欧拉法解常微分方程的计算。
步骤一:选择解题模型选择合适的数学模型很重要。
对于已经给定的微分方程,需要将它化为标准的形式。
例如,我们有如下的微分方程:y’ = 2y - 3,y(0) = 1将其化为标准的形式:dy/dx = 2y -3 将初始值y(0) = 1带入。
步骤二:确定计算步长确定计算步长h。
步长的大小与计算精度有直接关系,步长太小,计算量将很大,而精度较高;步长太大,精度较低。
步长的计算公式为:h = (b-a)/n其中,a和b是区间限制,n是初始步数。
步骤三:使用欧拉法求解微分方程根据欧拉法的公式,假设在t时刻函数y的值是y(t),求在下一个时刻t+h如何估算y值,公式为:y(t+h) = y(t)+ h * y'(t)将y'(t)=2y-3代入上式,得:y(t+h) = y(t)+ h* (2y(t)-3)接下来,根据初始值y(0) = 1,带入计算步长可得出一系列的近似值。
步骤四:绘制函数图像对于计算结果,应绘制出函数的近似图像。
通过matlab绘制y(x)的图像,也可以通过计算的数据进行近似曲线的绘制。
步骤五:测试计算结果通过计算结果与初始值进行比较,看算法是否正确和有效。
也可以将步长不同的计算结果进行比较,判断精度和计算效率的高低。
欧拉法解常微分方程在matlab中的使用,相较于手工计算,更具有高效、准确、方便的优势。
正因如此,在各类数学、物理、工程等领域中都有着广泛的应用。
matlab欧拉法求解微分方程例题

matlab欧拉法求解微分方程例题当使用Matlab的欧拉法(Euler's method)求解微分方程时,需要将微分方程转化为离散的差分方程。
下面以一个简单的一阶常微分方程为例来说明。
假设我们要求解以下的微分方程:dy/dx = x + y并给定初始条件 y(0) = 1。
我们可以使用欧拉法来逼近该微分方程的解。
在Matlab中,可以按照以下步骤进行求解:1. 定义步长和计算步数:h = 0.1; % 步长 N = 10; % 计算步数2. 初始化变量并给定初始条件:x = zeros(N+1, 1); % 初始化 x 数组 y = zeros(N+1, 1); % 初始化 y 数组 x(1) = 0; % 初始条件 x(0) = 0 y(1) = 1; % 初始条件 y(0) = 13. 使用欧拉法进行迭代计算:for i = 1:N x(i+1) = x(i) + h; % 计算下一个 x 值 y(i+1) = y(i) + h * (x(i) + y(i)); % 根据差分方程计算下一个 y 值end4. 绘制结果:plot(x, y, 'o-'); % 绘制连线图 xlabel('x'); ylabel('y'); title('Approximate Solution using Euler''s Method');完整的代码如下:h = 0.1; % 步长 N = 10; % 计算步数 x = zeros(N+1, 1); % 初始化 x 数组 y = zeros(N+1, 1); % 初始化 y 数组 x(1) = 0; % 初始条件 x(0) = 0 y(1) = 1; % 初始条件 y(0) = 1 for i = 1:N x(i+1) = x(i) + h; % 计算下一个 x 值y(i+1) = y(i) + h * (x(i) + y(i)); % 根据差分方程计算下一个 y 值 end plot(x, y, 'o-'); % 绘制连线图xlabel('x'); ylabel('y'); title('Approximate Solution using Euler''s Method');运行该代码,将得到近似的解曲线图,表示微分方程的数值解。
MATLAB常微分方程数值解——欧拉法、改进的欧拉法与四阶龙格库塔方法

MATLAB常微分⽅程数值解——欧拉法、改进的欧拉法与四阶龙格库塔⽅法MATLAB常微分⽅程数值解作者:凯鲁嘎吉 - 博客园1.⼀阶常微分⽅程初值问题2.欧拉法3.改进的欧拉法4.四阶龙格库塔⽅法5.例题⽤欧拉法,改进的欧拉法及4阶经典Runge-Kutta⽅法在不同步长下计算初值问题。
步长分别为0.2,0.4,1.0.matlab程序:function z=f(x,y)z=-y*(1+x*y);function R_K(h)%欧拉法y=1;fprintf('欧拉法:x=%f, y=%f\n',0,1);for i=1:1/hx=(i-1)*h;K=f(x,y);y=y+h*K;fprintf('欧拉法:x=%f, y=%f\n',x+h,y);endfprintf('\n');%改进的欧拉法y=1;fprintf('改进的欧拉法:x=%f, y=%f\n',0,1);for i=1:1/hx=(i-1)*h;K1=f(x,y);K2=f(x+h,y+h*K1);y=y+(h/2)*(K1+K2);fprintf('改进的欧拉法:x=%f, y=%f\n',x+h,y);endfprintf('\n');%龙格库塔⽅法y=1;fprintf('龙格库塔法:x=%f, y=%f\n',0,1);for i=1:1/hx=(i-1)*h;K1=f(x,y);K2=f(x+h/2,y+(h/2)*K1);K3=f(x+h/2,y+(h/2)*K2);K4=f(x+h,y+h*K3);y=y+(h/6)*(K1+2*K2+2*K3+K4);fprintf('龙格库塔法:x=%f, y=%f\n',x+h,y);end结果:>> R_K(0.2)欧拉法:x=0.000000, y=1.000000欧拉法:x=0.200000, y=0.800000欧拉法:x=0.400000, y=0.614400欧拉法:x=0.600000, y=0.461321欧拉法:x=0.800000, y=0.343519欧拉法:x=1.000000, y=0.255934改进的欧拉法:x=0.000000, y=1.000000改进的欧拉法:x=0.200000, y=0.807200改进的欧拉法:x=0.400000, y=0.636118改进的欧拉法:x=0.600000, y=0.495044改进的欧拉法:x=0.800000, y=0.383419改进的欧拉法:x=1.000000, y=0.296974龙格库塔法:x=0.000000, y=1.000000龙格库塔法:x=0.200000, y=0.804636龙格库塔法:x=0.400000, y=0.631465龙格库塔法:x=0.600000, y=0.489198龙格库塔法:x=0.800000, y=0.377225龙格库塔法:x=1.000000, y=0.291009>> R_K(0.4)欧拉法:x=0.000000, y=1.000000欧拉法:x=0.400000, y=0.600000欧拉法:x=0.800000, y=0.302400改进的欧拉法:x=0.000000, y=1.000000改进的欧拉法:x=0.400000, y=0.651200改进的欧拉法:x=0.800000, y=0.405782龙格库塔法:x=0.000000, y=1.000000龙格库塔法:x=0.400000, y=0.631625龙格库塔法:x=0.800000, y=0.377556>> R_K(1)欧拉法:x=0.000000, y=1.000000欧拉法:x=1.000000, y=0.000000改进的欧拉法:x=0.000000, y=1.000000改进的欧拉法:x=1.000000, y=0.500000龙格库塔法:x=0.000000, y=1.000000龙格库塔法:x=1.000000, y=0.303395注意:在步长h为0.4时,要将for i=1:1/h改为for i=1:0.8/h。
matlab用欧拉法求解常微分方程

matlab用欧拉法求解常微分方程在数学和科学领域中,常微分方程是一种非常有用的工具,用于描述许多自然和物理现象。
MATLAB是一种强大的数学软件,可以用来解决许多数学问题。
本文将介绍如何使用欧拉法在MATLAB中求解常微分方程。
欧拉法是一种基本的数值方法,用于近似解决微积分方程问题。
该方法使用离散时间步长,将微积分方程转换成差分方程,并不断迭代求解。
欧拉法的实现非常简单,因此它很适合初学者。
下面是使用欧拉法在MATLAB中求解常微分方程的步骤:1. 定义常微分方程以 y' = -0.5y + 3sin(t) 为例,我们先定义常微分方程。
在MATLAB中,可以使用 anonymous functions 实现:dydt = @(t,y) -0.5*y + 3*sin(t);2. 定义时间范围和时间步长我们需要定义时间范围和时间步长,以便在一定时间范围内求解差分方程。
在这个例子中,我们定义时间范围为 0 到 10,并定义时间步长为 0.1:tspan = [0 10];h = 0.1;3. 定义初始条件我们需要定义初始条件,即 y(0) 的值。
在这个例子中,我们假设 y(0) = 1:y0 = 1;4. 求解差分方程现在我们可以使用欧拉法求解差分方程了。
在MATLAB中,可以使用 odeEuler 函数(需要自己编写):[t,y] = odeEuler(dydt,tspan,y0,h);5. 可视化结果最后,我们可以将结果可视化,以便更好地理解求解过程。
在这个例子中,我们可以用 plot 函数将求解结果绘制出来:plot(t,y)xlabel('Time')ylabel('y(t)')title('Solution of y'' = -0.5y + 3sin(t) using Euler''s method')以上就是使用欧拉法在MATLAB中求解常微分方程的基本步骤。
euler方法求解常微分方程matlab

euler方法求解常微分方程matlab以euler方法求解常微分方程matlab常微分方程是数学中的重要分支之一,它描述了自然界和工程中的许多现象和过程。
求解常微分方程的方法有很多,其中一种常用的方法是欧拉方法。
欧拉方法是一种数值解常微分方程的方法,它通过将微分方程转化为差分方程,从而得到近似解。
在matlab中,我们可以使用欧拉方法来求解常微分方程。
下面我们将以一个具体的例子来说明如何使用matlab来求解常微分方程。
假设我们要求解的常微分方程是一阶线性常微分方程:dy/dx = f(x, y)其中f(x, y)是已知的函数。
我们需要给定一个初始条件y(x0) = y0,其中x0和y0是已知的常数。
我们需要定义函数f(x, y)。
在matlab中,我们可以使用匿名函数来定义函数。
例如,如果我们要求解的常微分方程是dy/dx = x^2 + y,那么我们可以定义函数f(x, y)如下:f = @(x, y) x^2 + y接下来,我们需要定义初始条件x0和y0。
假设x0 = 0,y0 = 1,我们可以定义初始条件如下:x0 = 0;y0 = 1;然后,我们需要定义步长h,即每一步的增量。
步长h越小,求解的结果越精确,但计算量也会增加。
在matlab中,我们可以使用input函数来让用户输入步长h。
例如,我们可以这样定义步长h:h = input('请输入步长h:');接下来,我们需要定义求解的区间。
假设我们要求解的区间是0到1,我们可以定义区间如下:a = 0;b = 1;然后,我们需要计算步数n。
步数n可以通过区间长度除以步长h 来得到。
在matlab中,我们可以使用ceil函数来向上取整。
例如,我们可以这样计算步数n:n = ceil((b - a) / h);接下来,我们需要定义一个数组x和一个数组y,用来存储每一步的计算结果。
我们可以使用zeros函数来创建这两个数组,并将初始条件存储在数组中。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
Euler法解常微分方程
Euler法解常微分方程算法:
Step 1 分别取积分上限、积分下限、步长
Step 2计算判断是否成立,成立转到Step 3,否则继续进行Step 4
Step 3 计算
Step 4
Euler法解常微分方程算程序:
function euler2(fun,y0,A,h)
%fun--y'
%y0---初值
%A----x取值范围
%a----x左区间端点值
%b----x右区间端点值
%h----给定步长
x=min(A);
b=max(A);
y=y0;
while x<b-h
b=y;
y=y+h*feval(fun,x,b)
x=x+h;
end
例:用Euler法计算下列初值问题(取步长h=0.2)
输入:fun=inline('-y-x*y^2')
euler2(fun,1,[0 0.6],0.2)
得到:
y =
0.8000
y =
0.6144
y =
0.4613
指导教师:年月日改进Euelr法解常微分方程
改进Euler法解常微分方程算法:
Step 1 分别取积分上限、积分下限、步长
Step 2 取一个以h为步长,a,b分别为左右端点的矩阵
Step 3 (1)做显性Euler预测
(2)将带入
Step 4计算判断是否成立,成立返回Step 3,否则继续进行Step 5 Step 5
改进Euler法解常微分方程算程序:
function gaijineuler2(fun,y0,A,h)
%fun--y'
%y0---初值
%A----x取值范围
%a----x左区间端点值
%b----x右区间端点值
%h----给定步长
a=min(A);
b=max(A);
x=a:h:b;
y(1)=y0;
for i=1:length(x)-1
w1=feval(fun,x(i),y(i));
y(i+1)=y(i)+h*w1;
w2=feval(fun,x(i+1),y(i+1));
y(i+1)=y(i)+h*(w1+w2)/2;
end
x=x'
y=y'
例:用改进Euler法计算下列初值问题(取步长h=0.25)
输入:fun=inline('-x*y^2')
gaijineuler2(fun,2,[0 5],0.25)
得到:
x =
0.2500
0.5000
0.7500
1.0000
1.2500
1.5000
1.7500
2.0000
2.2500
2.5000
2.7500
3.0000
3.2500
3.5000
3.7500
4.0000
4.2500
4.5000
4.7500
5.0000
y =
2.0000
1.8750
1.5939
1.2824
1.0096
0.7932
0.6282
0.5037
0.4097
0.3379
0.2824
0.2389
0.2043
0.1765
0.1538
0.1352
0.1196
0.1066
0.0955
0.0861
0.0779
指导教师:年月日禁止复制北京石油化工学院Onlyunited。