当前位置:文档之家› MATlAB 数值求解ODE方程

MATlAB 数值求解ODE方程

MATlAB 数值求解ODE方程
MATlAB 数值求解ODE方程

线性代数方程组数值解法及MATLAB实现综述

线性代数方程组数值解法及MATLAB实现综述廖淑芳20122090 数计学院12计算机科学与技术1班(职教本科)一、分析课题 随着科学技术的发展,提出了大量复杂的数值计算问题,在建立电子计算机成为数值计算的主要工具以后,它以数字计算机求解数学问题的理论和方法为研究对象。其数值计算中线性代数方程的求解问题就广泛应用于各种工程技术方面。因此在各种数据处理中,线性代数方程组的求解是最常见的问题之一。关于线性代数方程组的数值解法一般分为两大类:直接法和迭代法。 直接法就是经过有限步算术运算,可求的线性方程组精确解的方法(若计算过程没有舍入误差),但实际犹如舍入误差的存在和影响,这种方法也只能求得近似解,这类方法是解低阶稠密矩阵方程组级某些大型稀疏矩阵方程组的有效方法。直接法包括高斯消元法,矩阵三角分解法、追赶法、平方根法。 迭代法就是利用某种极限过程去逐步逼近线性方程组精确解的方法。迭代法具有需要计算机的存储单元少,程序设计简单,原始系数矩阵在计算过程始终不变等优点,但存在收敛性级收敛速度问题。迭代法是解大型稀疏矩阵方程组(尤其是微分方程离散后得到的大型方程组)的重要方法。迭代法包括Jacobi法SOR法、SSOR法等多种方法。 二、研究课题-线性代数方程组数值解法 一、直接法 1、Gauss消元法 通过一系列的加减消元运算,也就是代数中的加减消去法,以使A对角线以下的元素化为零,将方程组化为上三角矩阵;然后,再逐一回代求解出x向量。

1.1消元过程 1. 高斯消元法(加减消元):首先将A 化为上三角阵,再回代求解。 11121121222212n n n n nn n a a a b a a a b a a a b ?? ? ? ? ??? (1)(1)(1)(1)(1)11 121311(2)(2)(2)(2)222322 (3)(3)(3)3333()()000 00 n n n n n nn n a a a a b a a a b a a b a b ?? ? ? ? ? ? ?? ? 步骤如下: 第一步:1 11 1,2,,i a i i n a -? +=第行第行 11121121222212 n n n n nn n a a a b a a a b a a a b ?? ? ? ? ??? 111211(2)(2)(2)2222 (2)(2)(2)2 00n n n nn n a a a b a a b a a b ?? ? ? ? ??? 第二步:(2)2 (2)222,3, ,i a i i n a -?+=第行第行 111211(2)(2)(2)2222 (2)(2)(2)2 00n n n nn n a a a b a a b a a b ?? ? ? ? ??? 11121311(2)(2)(2)(2)222322 (3)(3)(3)33 33(3)(3)(3)3 00000 n n n n nn n a a a a b a a a b a a b a a b ?? ? ? ? ? ? ?? ? 类似的做下去,我们有: 第k 步:() ()k ,1, ,k ik k kk a i i k n a -?+=+第行第行。 n -1步以后,我们可以得到变换后的矩阵为: 11121311(2)(2)(2)(2)222322 (3)(3)(3)3333()()00000 n n n n n nn n a a a a b a a a b a a b a b ?? ? ? ? ? ? ?? ?

Matlab求解微分方程组及偏微分方程组

第四讲 Matlab 求解微分方程(组) 理论介绍:Matlab 求解微分方程(组)命令 求解实例:Matlab 求解微分方程(组)实例 实际应用问题通过数学建模所归纳得到的方程,绝大多数都是微分方程,真正能得到代数方程的机会很少.另一方面,能够求解的微分方程也是十分有限的,特别是高阶方程和偏微分方程(组).这就要求我们必须研究微分方程(组)的解法:解析解法和数值解法. 一.相关函数、命令及简介 1.在Matlab 中,用大写字母D 表示导数,Dy 表示y 关于自变量的一阶导数,D2y 表示y 关于自变量的二阶导数,依此类推.函数dsolve 用来解决常微分方程(组)的求解问题,调用格式为: X=dsolve(‘eqn1’,’eqn2’,…) 函数dsolve 用来解符号常微分方程、方程组,如果没有初始条件,则求出通解,如果有初始条件,则求出特解. 注意,系统缺省的自变量为t 2.函数dsolve 求解的是常微分方程的精确解法,也称为常微分方程的符号解.但是,有大量的常微分方程虽然从理论上讲,其解是存在的,但我们却无法求出其解析解,此时,我们需要寻求方程的数值解,在求常微分方程数值解方面,MATLAB 具有丰富的函数,我们将其统称为solver ,其一般格式为: [T,Y]=solver(odefun,tspan,y0) 说明:(1)solver 为命令ode45、ode23、ode113、ode15s 、ode23s 、ode23t 、ode23tb 、ode15i 之一. (2)odefun 是显示微分方程'(,)y f t y =在积分区间tspan 0[,]f t t =上从0t 到f t 用初始条件0y 求解. (3)如果要获得微分方程问题在其他指定时间点012,,,,f t t t t 上的解,则令 tspan 012[,,, ]f t t t t =(要求是单调的). (4)因为没有一种算法可以有效的解决所有的ODE 问题,为此,Matlab 提供了多种求解器solver ,对于不同的ODE 问题,采用不同的solver.

用Matlab解微分方程

用Matlab 软件求解微分方程 1.解析解 (1)一阶微分方程 求21y dx dy +=的通解:dsolve('Dy=1+y^2','x') 求y x dx dy -+=21的通解:dsolve('Dy=1+x^2-y','x') 求?????=+=1 )0(12y y dx dy 的特解:dsolve('Dy=1+y^2',’y(0)=1’,'x') (2)高阶微分方程 求解???-='==-+'+''. 2)2(,2)2(,0)(222πππy y y n x y x y x 其中,21=n ,命令为: dsolve('x^2*D2y+x*Dy+(x^2-0.5^2)*y=0','y(pi/2)=2,Dy(pi/2)=-2/pi','x') 求042=-+'-'''x y y y 的通解,命令为: dsolve('D3y-2*Dy+y-4*x=0','x') 输出为: ans=8+4*x+C1*exp(x)+C2*exp(-1/2*(5^(1/2)+1)*x)+C3*exp(1/2*(5^(1/2)-1)*x) (3)一阶微分方程组 求???+-='+='). (3)(4)(),(4)(3)(x g x f x g x g x f x f 的通解:[f,g]=dsolve('Df=3*f+4*g','Dg=-4*f+3*g','x') 输出为: f =exp(3*x)*(cos(4*x)*C1+sin(4*x)*C2) g =-exp(3*x)*(sin(4*x)*C1-cos(4*x)*C2) 若再加上初始条件1)0(,0)0(==g f ,则求特解: [f,g]=dsolve('Df=3*f+4*g','Dg=-4*f+3*g','f(0)=0,g(0)=1','x') 输出为: f =exp(3*x)*sin(4*x) g =exp(3*x)*cos(4*x) 2.数值解 (1)一阶微分方程

Matlab求解微分方程(组)及偏微分方程(组)

第四讲Matlab求解微分方程(组) 理论介绍:Matlab求解微分方程(组)命令 求解实例:Matlab求解微分方程(组)实例 实际应用问题通过数学建模所归纳得到得方程,绝大多数都就是微分方程,真正能得到代数方程得机会很少、另一方面,能够求解得微分方程也就是十分有限得,特别就是高阶方程与偏微分方程(组)、这就要求我们必须研究微分方程(组)得解法:解析解法与数值解法、 一.相关函数、命令及简介 1、在Matlab中,用大写字母D表示导数,Dy表示y关于自变量得一阶导数,D2y 表示y关于自变量得二阶导数,依此类推、函数dsolve用来解决常微分方程(组)得求解问题,调用格式为: X=dsolve(‘eqn1’,’eqn2’,…) 函数dsolve用来解符号常微分方程、方程组,如果没有初始条件,则求出通解,如果有初始条件,则求出特解、 注意,系统缺省得自变量为t 2、函数dsolve求解得就是常微分方程得精确解法,也称为常微分方程得符号解、但就是,有大量得常微分方程虽然从理论上讲,其解就是存在得,但我们却无法求出其解析解,此时,我们需要寻求方程得数值解,在求常微分方程数值解方 面,MATLAB具有丰富得函数,我们将其统称为solver,其一般格式为: [T,Y]=solver(odefun,tspan,y0) 说明:(1)solver为命令ode45、ode23、ode113、ode15s、ode23s、ode23t、ode23tb、ode15i之一、 (2)odefun就是显示微分方程在积分区间tspan上从到用初始条件求解、 (3)如果要获得微分方程问题在其她指定时间点上得解,则令tspan(要求就是单调得)、 (4)因为没有一种算法可以有效得解决所有得ODE问题,为此,Matlab提供了多种求解器solver,对于不同得ODE问题,采用不同得solver、 表1 Matlab中文本文件读写函数

用MATLAB解常微分方程

实验四 求微分方程的解 一、问题背景与实验目的 实际应用问题通过数学建模所归纳而得到的方程,绝大多数都是微分方程,真正能得到代数方程的机会很少.另一方面,能够求解的微分方程也是十分有限的,特别是高阶方程和偏微分方程(组).这就要求我们必须研究微分方程(组)的解法,既要研究微分方程(组)的解析解法(精确解),更要研究微分方程(组)的数值解法(近似解). 对微分方程(组)的解析解法(精确解),Matlab 有专门的函数可以用,本实验将作一定的介绍. 本实验将主要研究微分方程(组)的数值解法(近似解),重点介绍 Euler 折线法. 二、相关函数(命令)及简介 1.dsolve ('equ1','equ2',…):Matlab 求微分方程的解析解.equ1、equ2、…为方程(或条件).写方程(或条件)时用 Dy 表示y 关于自变量的一阶导数,用用 D2y 表示 y 关于自变量的二阶导数,依此类推. 2.simplify(s ):对表达式 s 使用 maple 的化简规则进行化简. 例如: syms x simplify(sin(x)^2 + cos(x)^2) ans=1 3.[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 = combine 4.[T,Y] = solver(odefun,tspan,y 0) 求微分方程的数值解. 说明: (1) 其中的 solver 为命令 ode45、ode23、ode113、ode15s 、ode23s 、ode23t 、ode23tb 之一. (2) odefun 是显式常微分方程:?????==0 0)() ,(y t y y t f dt dy (3) 在积分区间 tspan =],[0f t t 上,从0t 到f t ,用初始条件0y 求解.

Matlab求解微分方程(组)及偏微分方程(组)

第四讲 Matlab 求解微分方程(组) 理论介绍:Matlab 求解微分方程(组)命令 求解实例:Matlab 求解微分方程(组)实例 实际应用问题通过数学建模所归纳得到的方程,绝大多数都是微分方程,真正能得到代数方程的机会很少.另一方面,能够求解的微分方程也是十分有限的,特别是高阶方程和偏微分方程(组).这就要求我们必须研究微分方程(组)的解法:解析解法和数值解法. 一.相关函数、命令及简介 1.在Matlab 中,用大写字母D 表示导数,Dy 表示y 关于自变量的一阶导数,D2y 表示y 关于自变量的二阶导数,依此类推.函数dsolve 用来解决常微分方程(组)的求解问题,调用格式为: X=dsolve(‘eqn1’,’eqn2’,…) 函数dsolve 用来解符号常微分方程、方程组,如果没有初始条件,则求出通解,如果有初始条件,则求出特解. 注意,系统缺省的自变量为t 2.函数dsolve 求解的是常微分方程的精确解法,也称为常微分方程的符号解.但是,有大量的常微分方程虽然从理论上讲,其解是存在的,但我们却无法求出其解析解,此时,我们需要寻求方程的数值解,在求常微分方程数值解方面,MATLAB 具有丰富的函数,我们将其统称为solver ,其一般格式为: [T,Y]=solver(odefun,tspan,y0) 说明:(1)solver 为命令ode45、ode23、ode113、ode15s 、ode23s 、ode23t 、ode23tb 、ode15i 之一. (2)odefun 是显示微分方程'(,)y f t y =在积分区间tspan 0[,]f t t =上从0t 到f t 用初始条件0y 求解. (3)如果要获得微分方程问题在其他指定时间点012,,, ,f t t t t 上的解,则令 tspan 012[,,,]f t t t t =(要求是单调的). (4)因为没有一种算法可以有效的解决所有的ODE 问题,为此,Matlab 提供

数理方程基于matlab的数值解法

数理方程数值解法与其在matlab软件上的实现张体强1026222 廖荣发1026226 [摘要] 数学物理方程的数值解在实际生活中越来越使用,首先基于偏微分数值解的思想上,通过matlab软件的功能,研究其数学物理方程的数值解,并通过对精确解和数值解进行对比,追究其数值解的可行性,在此,给出相关例子和程序代码,利于以后的再次研究和直接使用。 [关键字] 偏微分方程数值解matlab 数学物理方程的可视化 一:研究意义 在我们解数学物理方程,理论上求数学物理方程的定解有着多种解法,但是有许多定解问题却不能严格求解,只能用数值方法求出满足实际需要的近似解。而且实际问题往往很复杂,这时即便要解出精确解就很困难,有时甚至不可能,另一方面,在建立数学模型时,我们已作了很多近似,所以求出的精确解也知识推导出的数学问题的精确解,并非真正实际问题的精确解。因此,我们有必要研究近似解法,只要使所求得的近似解与精确解之间的误差在规定的范围内,则仍能满足实际的需要,有限差分法和有限元法是两种最常用的

求解数学物理方程的数值解法,而MATLAB 在这一方面具有超强的数学功能,可以用来求其解。 二:数值解法思想和步骤 2.1:网格剖分 为了用差分方法求解上述问题,将求解区域 {}(,)|01,01x t x t Ω=≤≤≤≤作剖分。将空间区间[0,1]作m 等分,将时 间[0,1]区间作n 等分,并记 1/,1/,,0,,0j k h m n x jh j m t k k n ττ===≤≤=≤≤。分别称h 和τ 为空间和 时间步长。用两簇平行直线,0,,0j k x x j m t t k n =≤≤=≤≤将Ω分割成矩形网格。 2.2:差分格式的建立 0u u t x ??-=??………………………………(1) 设G 是,x t 平面任一有界域,据Green 公式(参考数学物理方程第五章): ( )()G u u dxdt udt udx t x Γ??-=--??? ? 其中G Γ=?。于是可将(1)式写成积分守恒形式: ()0udt udx Γ --=? (2) 我们先从(2)式出发构造熟知的Lax 格式设网格如下图所示

数学应用软件作业6-用Matlab求解微分方程(组)的解析解和数值解

数学应用软件作业6-用Matlab 求解微分方程(组)的解析解和数值解

注:上机作业文件夹以自己的班级姓名学号命名,文件夹包括如下上机报告和Matlab程序。 上机报告模板如下: 佛山科学技术学院 上机报告 课程名称数学应用软件 上机项目用Matlab求解微分方程(组)的解析解和数值解 专业班级姓名学号 一. 上机目的 1.了解求微分方程(组)的解的知识。 2.学习Matlab中求微分方程的各种解的函数,如 dsolve命令、ode45函数等等,其中注意把方程化为新的方程的形式。 3.掌握用matlab编写程序解决求解微分方程的问 题。 二. 上机内容 1、求高阶线性齐次方程:y’’’-y’’-3y’+2y=0。 2、求常微分方程组

2 210cos,2 24,0 t t t dx dy x t x dt dt dx dy y e y dt dt = - = ? +-== ?? ? ?++== ?? 3、求解 分别用函数ode45和ode15s计算求解,分别画出图形,图形分别标注标题。 4、求解微分方程 ,1 )0( ,1 '= + + - =y t y y 先求解析解,在[0,1]上作图; 再用ode45求数值解(作图的图形用“o”表示),在同一副图中作图进行比较,用不同的颜色表示。 三. 上机方法与步骤 给出相应的问题分析及求解方法,并写出Matlab 程序,并有上机程序显示截图。 题1:直接用命令dsolve求解出微分方程的通解。 Matlab程序:

dsolve('D3y-D2y-3*Dy+2*y','x') 题2:将微分方程组改写为 5cos2exp(2) 5cos2exp(2) (0)2,(0)0 dx t t x y xt dy t t x y dt x y ? =+--- ? ? ? =-+-+- ? ? == ? ? ? , 再用命令dsolve求解微分方程的通解。 Matlab程序: 建立timu2.m如下: [x,y]=dsolve('Dx=5*cos(t)+2*exp(-2*t)-x-y','Dy=-5*cos(t)+2*exp(-2*t)+x-y ','x(0)=2,y(0)=0','t') x=simple(x) y=simple(y)

利用matlab编写S函数求解微分方程

利用matlab编写S函数求解微分方程自动化专业综合设计报告 自动化专业综合设计报告

函数求解微S编写设计题目:利用 matlab 分方程 自动化系统仿真实验室所在 实验室: 郭卫平 指导教师: 律迪迪学生姓名 200990519114 班级文自0921 学号 成绩评定: 自动化专业综合设计报告

一、设计目的 了解使用simulink的扩展工具——S-函数,s函数可以利用matlab的丰富资源,而不仅仅局限于simulink提供的模块,而用c或c++等语言写的s函数还可以实现对硬件端口的操作,还可以操作windows API 等的,它的魅力在于完美结合了simulink 框图简洁明快的特点和编程灵活方便的优点,提供了增强和扩展sinulink能力的强大机制,同时也是使用RTW实现实时仿真的关键。 二、设计要求 求解解微分方程 y'=y-2x/y 自动化专业综合设计报告 y(0)=1 要求利用matlab编写S函数求解 三、设计内容(可加附页) 【步骤1】获取状态空间表达式。

在matlab中输入 dsolve(‘Dy=y-2*x/y','y(0)=1', 'x') 得到 y=(2*x+1).^(1/2); 【步骤2】建立s函数的m文件。 利用21·用S函数模板文件。 以下是修改之后的模板文件sfuntmpl.m 的内容。 function [sys,x0,str,ts] = sfuntmpl(t,x,u,flag) %SFUNTMPL S-function M-file General template define you can With % M-file S-functions, you own ordinary differential system equations (ODEs), discrete % equations, and/or just about any type of algorithm to be used within a %

一维抛物线型方程数值解法(1)(附图及matlab程序)

一维抛物线偏微分方程数值解法(1) 解一维抛物线型方程(理论书籍可以参看孙志忠:偏微分方程数值解法) Ut-Uxx=0, 00) U(x,0)=e^x, 0<=x<=1, U(0,t)=e^t,U(1,t)=e^(1+t), 00) %不用解线性方程组,由下一层(时间层)的值就直接得到上一层的值 %m,n为x,t方向的网格数,例如(2-0)/0.01=200; %e为误差,p为精确解 u=zeros(n+1,m+1); x=0+(0:m)*h1; t=0+(0:n)*h2; for(i=1:n+1) u(i,1)=exp(t(i)); u(i,m+1)=exp(1+t(i)); end for(i=1:m+1) u(1,i)=exp(x(i)); end for(i=1:n+1) for(j=1:m+1) f(i,j)=0; end end r=h2/(h1*h1); %此处r=a*h2/(h1*h1);a=1 要求r<=1/2差分格式才稳定for(i=1:n) for(j=2:m) u(i+1,j)=(1-2*r)*u(i,j)+r*(u(i,j-1)+u(i,j+1))+h2*f(i,j); end end for(i=1:n+1) for(j=1:m+1) p(i,j)=exp(x(j)+t(i)); e(i,j)=abs(u(i,j)-p(i,j)); end end

数学模型之微分方程及其MATLAB求解

数学模型之微分方程及其MATLAB求解 ---卫星轨迹等经典例题求解分析1. 考虑初值问题画图 y'''?3y ''?y 'y = 0 y(0) = 0 y '(0) =1 y ' '(0) = ?1 2、 3、 【实验步骤与程序】 1. M -文件建立m函数文件

function y=f(t,x) y=[x(2);x(3);9*x(3)^2+x(1)*x(2)]; 求解微分方程,命令如下: x0=[0;1;-1]; [t,y]=ode45(@mm,[0,2.5],x0); plot(y(:,1),y(:,2)); figure(2); plot3(y(:,1),y(:,2),y(:,3))

2、M -文件建立m函数文件 function dx=appollo(t,x) mu=1/82.45; mustar=1-mu; r1=sqrt((x(1)+mu)^2+x(3)^2); r2=sqrt((x(1)-mustar)^2+x(3)^2); dx=[x(2) 2*x(4)+x(1)-mustar*(x(1)+mu)/r1^3-mu*(x(1)-mustar)/r2^3 x(4) -2*x(2)+x(3)-mustar*x(3)/r1^3-mu*x(3)/r2^3];

求解微分方程,命令如下: x0=[1.2;0;0;-1.04935751]; options=odeset('reltol',1e-8); [t,y]=ode45(@appollo,[0,20],x0,options); plot(y(:,1),y(:,3)) title('Appollo卫星运动轨迹') xlabel('x') ylabel('y')

一维抛物线偏微分方程数值 解法(3)(附图及matlab程序)

一维抛物线偏微分方程数值解法(3) 上一篇参看一维抛物线偏微分方程数值解法(2)(附图及matlab程序) 解一维抛物线型方程(理论书籍可以参看孙志忠:偏微分方程数值解法) Ut-Uxx=0, 00) U(x,0)=e^x, 0<=x<=1, U(0,t)=e^t,U(1,t)=e^(1+t), 0

r=h2/(h1*h1); for(i=1:n) %外循环,先固定每一时间层,每一时间层上解一线性方程组% a(1)=0;b(1)=1+r;c(1)=-r/2;d(1)=r/2* (u(i+1,1)+u(i,1))+h2*f(i,j)... +(1-r)*u(i,2)+r/2*u(i,3); for(k=2:m-2) a(k)=-r/2;b(k)=1+r;c(k)=- r/2;d(k)=h2*f(i,j)+r/2*u(i,k)+(1-r)... *u(i,k+1)+r/2*u(i,k+2); %输入部分系数矩阵,为0的矩阵元素不输入% end a(m-1)=-r/2;b(m-1)=1+r;d(m-1)=h2*f(i,j)+r/2* (u(i,m+1)+u(i+1,m+1)... )+r/2*u(i,m-1)+(1-r)*u(i,m); for(k=1:m-2) %开始解线性方程组消元过程 a(k+1)=-a(k+1)/b(k); b(k+1)=b(k+1)+a(k+1)*c(k); d(k+1)=d(k+1)+a(k+1)*d(k); end u(i+1,m)=d(m-1)/b(m-1); %回代过程% for(k=m-2:-1:1) u(i+1,k+1)=(d(k)-c(k)*u(i+1,k+2))/b(k); end end for(i=1:n+1) for(j=1:m+1) p(i,j)=exp(x(j)+t(i)); %p为精确解 e(i,j)=abs(u(i,j)-p(i,j));%e为误差 end end [u p e x t]=CN(0.1,0.005,10,200);surf(x,t,e); shading interp; >> xlabel('x');ylabel('t');zlabel('e'); >> title('误差曲面')

双曲方程基于matlab的数值解法

双曲型方程基于MATLAB 的数值解法 (数学1201,陈晓云,41262022) 一:一阶双曲型微分方程的初边值问题 0,01,0 1.(,0)cos(),0 1. (0,)(1,)cos(),0 1. u u x t t x u x x x u t u t t t ππ??-=≤≤≤≤??=≤≤=-=≤≤ 精确解为 ()t x cos +π 二:数值解法思想和步骤 2.1:网格剖分 为了用差分方法求解上述问题,将求解区域{}(,)|01,01x t x t Ω=≤≤≤≤作剖分。将空间区间[0,1]作m 等分,将时间[0,1]区间作n 等分,并记 1/,1/,,0,,0j k h m n x jh j m t k k n ττ===≤≤=≤≤。分别称h 和τ为空间和时 间步长。用两簇平行直线,0,,0j k x x j m t t k n =≤≤=≤≤将Ω分割成矩形网格。 2.2:差分格式的建立 0u u t x ??-=?? 2.2.1:Lax-Friedrichs 方法 对时间、空间采用中心差分使得 2h 1 1111)(2 1u u x u u u u u t u k j k j k j k j k j k j -+-++-= +=-= ????τ τ 则由上式得到Lax-Friedrichs 格式 1 11111()202k k k k k j j j j j u u u u u h τ+-+-+-+-+=

截断误差为 ()[]k k k j h j j R u L u Lu =- 1 11111()22k k k k k k k j j j j j j j u u u u u u u h t x τ+-+-+-+-??=+-+?? 23222 3 (),(0,0)26k k j j u u h O h j m k n t x ττ??= -=+≤≤≤≤?? 所以Lax-Friedrichs 格式的截断误差的阶式2()O h τ+ 令/s h τ=:则可得差分格式为 1111 11(),(0,0)222 k k k k k j j j j j s s u u u u u j m k n +--++=-+++≤≤≤≤ 0cos()(0)j j u x j m π=≤≤ 0cos(),cos(),(0)k k k m k u t u t k n ππ==-≤≤ 其传播因子为: ()()()e e G h i h i s h i h i σσσστσ---=-+e e 221, 化简可得: ()()()()()h s G h is h G στσσστ σsin 11,sin cos ,2 2 2--=-= 所以当1s ≤时,()1,≤τσG ,格式稳定。 * 2.2.2:LaxWendroff 方法 用牛顿二次插值公式可以得到LaxWendroff 的差分格式,在此不详细分析,它的截断误差为() h 2 2 +O τ ,是二阶精度;当2s ≤时,()1,≤τσG , 格式稳定。在这里主要用它与上面一阶精度的Lax-Friedrichs 方法进行简单对比。 2.3差分格式的求解

用matlab求解常微分方程

实验六 用matlab 求解常微分方程 1.微分方程的概念 未知的函数以及它的某些阶的导数连同自变量都由一已知方程联系在一起的方程称为微分方程。如果未知函数是一元函数,称为常微分方程。常微分方程的一般形式为 0),,",',,()(=n y y y y t F 如果未知函数是多元函数,成为偏微分方程。联系一些未知函数的一组微分方程组称为微分方程组。微分方程中出现的未知函数的导数的最高阶解数称为微分方程的阶。若方程中未知函数及其各阶导数都是一次的,称为线性常微分方程,一般表示为 )()(')()(1)1(1)(t b y t a y t a y t a y n n n n =++++-- 若上式中的系数n i t a i ,,2,1),( =均与t 无关,称之为常系数。 2.常微分方程的解析解 有些微分方程可直接通过积分求解.例如,一解常系数常微分方程1+=y dt dy 可化为 dt y dy =+1,两边积分可得通解为 1-=t ce y .其中c 为任意常数.有些常微分方程可用一些技巧,如分离变量法,积分因子法,常数变异法,降阶法等可化为可积分的方程而求得解析解. 线性常微分方程的解满足叠加原理,从而他们的求解可归结为求一个特解和相应齐次微分方程的通解.一阶变系数线性微分方程总可用这一思路求得显式解。高阶线性常系数微分方程可用特征根法求得相应齐次微分方程的基本解,再用常数变异法求特解。 一阶常微分方程与高阶微分方程可以互化,已给一个n 阶方程 ),,",',()1()(-=n n y y y t f y 设)1(21,,',-===n n y y y y y y ,可将上式化为一阶方程组 ?????????====-),,,,(''''2113221n n n n y y y t f y y y y y y y 反过来,在许多情况下,一阶微分方程组也可化为高阶方程。所以一阶微分方程组与高阶常微分方程的理论与方法在许多方面是相通的,一阶常系数线性微分方程组也可用特征根法求解。 3.微分方程的数值解法 除常系数线性微分方程可用特征根法求解,少数特殊方程可用初等积分法求解外,大部分微分方程无限世界,应用中主要依靠数值解法。考虑一阶常微分方程初值问题 ???=<<=000)()),(,()('y t y t t t t y t f t y f

MATLAB求解常微分方程数值解

利用MATLAB求解常微分方程数值解

目录 1. 内容简介 (1) 2. Euler Method(欧拉法)求解 (1) . 显式Euler法和隐式Euler法 (2) . 梯形公式和改进Euler法 (3) . Euler法实用性 (5) 3. Runge-Kutta Method(龙格库塔法)求解 (6) . Runge-Kutta基本原理 (6) . MATLAB中使用Runge-Kutta法的函数 (8) 4. 使用MATLAB求解常微分方程 (8) . 使用ode45函数求解非刚性常微分方程 (8) . 刚性常微分方程 (9) 5. 总结 (10) 参考文献 (11) 附录 (12) 1. 显式Euler法数值求解 (12) 2. 改进Euler法数值求解 (12) 3. 四阶四级Runge-Kutta法数值求解 (13) 4. 使用ode45求解 (14)

1.内容简介 把《高等工程数学》看了一遍,增加对数学内容的了解,对其中数值解法比较感兴趣,这大概是因为在其它各方面的学习和研究中经常会遇到数值解法的问题。理解模型然后列出微分方程,却对着方程无从下手,无法得出精确结果实在是让人难受的一件事情。 实际问题中更多遇到的是利用数值法求解偏微分方程问题,但考虑到先从常微分方程下手更为简单有效率,所以本文只研究常微分方程的数值解法。把一个工程实际问题弄出精确结果远比弄清楚各种细枝末节更有意思,因此文章中不追求非常严格地证明,而是偏向如何利用工具实际求解出常微分方程的数值解,力求将课程上所学的知识真正地运用到实际方程的求解中去,在以后遇到微分方程的时候能够熟练运用MATLAB得到能够在工程上运用的结果。 文中求解过程中用到MATLAB进行数值求解,主要目的是弄清楚各个函数本质上是如何对常微分方程进行求解的,对各种方法进行MATLAB编程求解,并将求得的数值解与精确解对比,其中源程序在附录中。最后考察MATLAB中各个函数的适用范围,当遇到实际工程问题时能够正确地得到问题的数值解。 2.Euler Method(欧拉法)求解 Euler法求解常微分方程主要包括3种形式,即显式Euler法、隐式Euler法、梯形公式法,本节内容分别介绍这3种方法的具体内容,并在最后对3种方法精度进行对比,讨论Euler法的实用性。 本节考虑实际初值问题 使用解析法,对方程两边同乘以得到下式

第二章非线性方程(组)的数值解法的matlab程序

本章主要介绍方程根的有关概念,求方程根的步骤,确定根的初始近似值的方法(作图法,逐步搜索法等),求根的方法(二分法,迭代法,牛顿法,割线法,米勒(M üller )法和迭代法的加速等)及其MATLAB 程序,求解非线性方程组的方法及其MATLAB 程序. 2.1 方程(组)的根及其MATLAB 命令 2.1.2 求解方程(组)的solve 命令 求方程f (x )=q (x )的根可以用MATLAB 命令: >> x=solve('方程f(x)=q(x)',’待求符号变量x ’) 求方程组f i (x 1,…,x n )=q i (x 1,…,x n ) (i =1,2,…,n )的根可以用MATLAB 命令: >>E1=sym('方程f1(x1,…,xn)=q1(x1,…,xn)'); ……………………………………………………. En=sym('方程fn(x1,…,xn)=qn(x1,…,xn)'); [x1,x2,…,xn]=solve(E1,E2,…,En, x1,…,xn) 2.1.3 求解多项式方程(组)的roots 命令 如果)(x f 为多项式,则可分别用如下命令求方程0)(=x f 的根,或求导数)('x f (见表 2-1). 2.1.4 求解方程(组)的fsolve 命令 如果非线性方程(组)是多项式形式,求这样方程(组)的数值解可以直接调用上面已经介绍过的roots 命令.如果非线性方程(组)是含有超越函数,则无法使用roots 命令,需要调用MATLAB 系统中提供的另一个程序fsolve 来求解.当然,程序fsolve 也可以用于多项式方程(组),但是它的计算量明显比roots 命令的大. fsolve 命令使用最小二乘法(least squares method )解非线性方程(组) (F X =)0 的数值解,其中X 和F (X )可以是向量或矩阵.此种方法需要尝试着输入解X 的初始值(向量或矩阵)X 0,即使程序中的迭代序列收敛,也不一定收敛到(F X =)0的根(见例2.1.8). fsolve 的调用格式: X=fsolve(F,X0) 输入函数)(x F 的M 文件名和解X 的初始值(向量或矩阵)X 0,尝试着解方程(组)

Matlab解微分方程(ODE+PDE)

常微分方程: 1 ODE解算器简介(ode**) 2 微分方程转换 3 刚性/非刚性问题(Stiff/Nonstiff) 4 隐式微分方程(IDE) 5 微分代数方程(DAE) 6 延迟微分方程(DDE) 7 边值问题(BVP) 偏微分方程(PDEs)Matlab解法 偏微分方程: 1 一般偏微分方程组(PDEs)的命令行求解 2 特殊偏微分方程(PDEs)的PDEtool求解 3 陆君安《偏微分方程的MATLAB解法 先来认识下常微分方程(ODE)初值问题解算器(solver) [T,Y,TE,YE,IE] = odesolver(odefun,tspan,y0,options) sxint = deval(sol,xint) Matlab中提供了以下解算器: 输入参数: odefun:微分方程的Matlab语言描述函数,必须是函数句柄或者字符串,必须写成Matlab

规范格式(也就是一阶显示微分方程组),这个具体在后面讲解 tspan=[t0 tf]或者[t0,t1,…tf]:微分变量的范围,两者都是根据t0和tf的值自动选择步长,只是前者返回所有计算点的微分值,而后者只返回指定的点的微分值,一定要注意对于后者tspan必须严格单调,还有就是两者数据存储时使用的内存不同(明显前者多),其它没有任何本质的区别 y0=[y(0),y’(0),y’’(0)…]:微分方程初值,依次输入所有状态变量的初值,什么是状态变量在后面有介绍 options:微分优化参数,是一个结构体,使用odeset可以设置具体参数,详细内容查看帮助 输出参数: T:时间列向量,也就是ode**计算微分方程的值的点 Y:二维数组,第i列表示第i个状态变量的值,行数与T一致 在求解ODE时,我们还会用到deval()函数,deval的作用就是通过结构体solution计算t 对应x值,和polyval之类的很相似! 参数格式如下: sol:就是上次调用ode**函数得道的结构体解 xint:需要计算的点,可以是标量或者向量,但是必须在tspan范围内 该函数的好处就是如果我想知道t=t0时的y值,不需要重新使用ode计算,而直接使用上次计算的得道solution就可以 [教程] 微分方程转换为一阶显示微分方程组方法 好,上面我们把Matlab中的常微分方程(ODE)的解算器讲解的差不多了,下面我们就具体开始介绍如何使用上面的知识吧! 现实总是残酷的,要得到就必须先付出,不可能所有的ODE一拿来就可以直接使用,因此,在使用ODE解算器之前,我们需要做的第一步,也是最重要的一步,借助状态变量将微分

线性代数方程组数值解法及MATLAB实现综述

线性代数方程组数值解法及MATLAB 实现综述 廖淑芳 20122090 数计学院 12计算机科学与技术1班(职教本科) 一、分析课题 随着科学技术的发展,提出了大量复杂的数值计算问题,在建立电子计算机成为数值计算的主要工具以后,它以数字计算机求解数学问题的理论和方法为研究对象。其数值计算中线性代数方程的求解问题就广泛应用于各种工程技术方面。因此在各种数据处理中,线性代数方程组的求解是最常见的问题之一。关于线性代数方程组的数值解法一般分为两大类:直接法和迭代法。 直接法就是经过有限步算术运算,可求的线性方程组精确解的方法(若计算过程没有舍入误差),但实际犹如舍入误差的存在和影响,这种方法也只能求得近似解,这类方法是解低阶稠密矩阵方程组级某些大型稀疏矩阵方程组的有效方法。直接法包括高斯消元法,矩阵三角分解法、追赶法、平方根法。 迭代法就是利用某种极限过程去逐步逼近线性方程组精确解的方法。迭代法具有需要计算机的存储单元少,程序设计简单,原始系数矩阵在计算过程始终不变等优点,但存在收敛性级收敛速度问题。迭代法是解大型稀疏矩阵方程组(尤其是微分方程离散后得到的大型方程组)的重要方法。迭代法包括Jacobi 法SOR 法、SSOR 法等多种方法。 二、研究课题-线性代数方程组数值解法 一、 直接法 1、 Gauss 消元法 通过一系列的加减消元运算,也就是代数中的加减消去法,以使A 对角线以下的元素化为零,将方程组化为上三角矩阵;然后,再逐一回代求解出x 向量。 1.1消元过程 1. 高斯消元法(加减消元):首先将A 化为上三角阵,再回代求解。 11121121222212n n n n nn n a a a b a a a b a a a b ?? ? ? ? ???L L M M O M M L (1)(1)(1)(1)(1)11121311(2)(2)(2)(2)222322(3)(3)(3)3333()()000000n n n n n nn n a a a a b a a a b a a b a b ?? ? ? ? ? ? ???L L L M M M O M M L 步骤如下:

相关主题
文本预览
相关文档 最新文档