文档库 最新最全的文档下载
当前位置:文档库 › 偏微分方程实验1

偏微分方程实验1

偏微分方程实验1
偏微分方程实验1

《偏微分方程数值解》

课程实验报告(一)

1 实验题目

给定初值问题(注:共两次上机时间4学时)

1

02

)0(42'≤≤??

?=--=x y x y y

其精确解为

12)(2+-=-x e x y x

。取h=0.1, 分别用显式Euler 法、隐式梯形法、二级二阶Runge-Kutta 、三级三阶Runge-Kutta 、四级四阶Runge-Kutta 计算数值解,并与精确解比较。

2 求解方法

1、每种方法的迭代公式: (1)显式Euler 法

p[i]=p[i-1]+h*(*f)(x[i-1],p[i-1]); (2)隐式梯形法

p[i]=p[i-1]+h*(*f)(x[i-1],p[i-1]);; (3)二级二阶Runge-Kutta

K1=(*f)(x[i-1],p[i-1]);

K2=(*f)(x[i-1]+h/2,p[i-1]+h*K1/2); p[i]=p[i-1]+h*K2; (4)三级三阶Runge-Kutta

K1=(*f)(x[i-1],p[i-1]);

K2=(*f)(x[i-1]+h/2,p[i-1]+h*K1/2); K3=(*f)(x[i-1]+h/2,p[i-1]+h*K2/2); K4=(*f)(x[i-1]+h,p[i-1]+h*K3); p[i]=p[i-1]+h*(K1+2*K2+2*K3+K4)/6; (5)四级四阶Runge-Kutta

K1=(*f)(x[i-1],p[i-1]);

K2=(*f)(x[i-1]+h/2,p[i-1]+h*K1/2); K3=(*f)(x[i-1]+h/2,p[i-1]+h*K2/2); K4=(*f)(x[i-1]+h,p[i-1]+h*K3);

p[i]=p[i-1]+h*(K1+2*K2+2*K3+K4)/6;

3 程序源代码

#include

#include

#define N 11

double fun(double x,double y) //微分方程{

return -2*y-4*x;

}

double p_fun(double x) //原函数

{

return exp(-2*x)-2*x+1;

}

double* Exact(double x0, double h) //精确解{

int i;

double *p=new double[N];

for(i=0;i

return p;

}

double* Euler(double h,double y0) //显式欧拉法{

int i;

double *p=new double[N];

p[0]=y0;

for(i=1;i

}

double* Y_Euler(double h,double y0) //隐式梯形法{

int i;

double *p=new double[N];

p[0]=y0;

double t;

for(i=1;i

{

t=p[i-1]+h*fun((i-1)*h,p[i-1]);

p[i]=p[i-1]+0.5*h*(fun((i-1)*h,p[i-1])+fun(i*h,t)); }

return p;

}

double* Two_R_K(double h,double y0) //二阶二级Runge-Kutta {

int i;

double *p=new double[N];

p[0]=y0;

double k1,k2;

for(i=1;i

{

k1=fun((i-1)*h,p[i-1]);

k2=fun((i-0.5)*h,p[i-1]+0.5*h*k1);

p[i]=p[i-1]+h*k2;

}

return p;

}

double* Three_R_K(double h,double y0) //三阶三级Runge-Kutta {

int i;

double *p=new double[N];

p[0]=y0;

double k1,k2,k3;

for(i=1;i

{

k1=fun((i-1)*h,p[i-1]);

k2=fun((i-0.5)*h,p[i-1]+0.5*h*k1);

k3=fun(i*h,p[i-1]+h*(-k1+2*k2));

p[i]=p[i-1]+h/6*(k1+4*k2+k3);

}

return p;

}

double* Four_R_K(double h,double y0) //四阶四级Runge-Kutta {

int i;

double *p=new double[N];

p[0]=y0;

double k1,k2,k3,k4;

for(i=1;i

{

k1=fun((i-1)*h,p[i-1]);

k2=fun((i-0.5)*h,p[i-1]+0.5*h*k1);

k3=fun((i-0.5)*h,p[i-1]+0.5*h*k2);

k4=fun(i*h,p[i-1]+h*k3);

p[i]=p[i-1]+h/6*(k1+2*k2+2*k3+k4);

}

return p;

}

//主函数

int main(void)

{

double x,*p0,*p1,*p2,*p3,*p4,*p5;

double h,y0=2;

int i;

//输入步长

printf("请输入步长:");

scanf("%lf",&h);

//调用函数求解

p0=Exact(0,h);

p1=Euler(h,y0);

p2=Y_Euler(h,y0);

p3=Two_R_K(h,y0);

p4=Three_R_K(h,y0);

p5=Four_R_K(h,y0);

//输出

printf("\n x 精确解欧拉方法隐式梯形法二阶R_K 三阶R_K 四阶R_K\n") ;

for(i=0;i

printf("%6f %10lf %10lf %10lf %10lf %10lf %10lf\n",h*i,p0[i],p1[i],p2[i],p3[i],p4[i],p 5[i]);

}

4 实验结果及分析

编译运行程序,得到下图

可以看出,欧拉方法的精确度为10?1,隐式梯形法和二级二阶Runge-Kutta 的精确度为10?2,三级三阶Runge-Kutta的精确度为10?3,四级四阶Runge-Kutta的精确度为10?5。

偏微分方程数值解实验报告

偏微分方程数值解实验报告

1、用有限元方法求下列边值问题的数值解:''()112x -y +y =2s i n ,0∈∈??∈(0,)?, 其中取1ν= 要求画出解曲面。迭代格式如下: 1221212111111111122142212n n n n n n j j j j j j n n n n n n j j j j j j V V V V V V h h V V V V V V h h τ++++++++++-+-??-()-()()-()??++?????? ??-+-+??=+??????

1、 %Ritz Galerkin方法求解方程 function u1=Ritz(x) %定义步长 h=1/100; x=0:h:1; n=1/h; a=zeros(n-1,1); b=zeros(n,1); c=zeros(n-1,1); d=zeros(n,1); %求解Ritz方法中内点系数矩阵 for i=1:1:n-1 b(i)=(1/h+h*pi*pi/12)*2; d(i)=h*pi*pi/2*sin(pi/2*(x(i)+h))/2+h*pi*pi/2*sin(pi/2*x(i+1))/2; end %右侧导数条件边界点的计算 b(n)=(1/h+h*pi*pi/12); d(n)=h*pi*pi/2*sin(pi/2*(x(i)+h))/2; for i=1:1:n-1 a(i)=-1/h+h*pi*pi/24; c(i)=-1/h+h*pi*pi/24; end %调用追赶法 u=yy(a,b,c,d) %得到数值解向量 u1=[0,u] %对分段区间做图 plot(x,u1) %得到解析解 y1=sin(pi/2*x); hold on plot(x,y1,'o') legend('数值解','解析解') function x=yy(a,b,c,d) n=length(b); q=zeros(n,1); p=zeros(n,1); q(1)=b(1); p(1)=d(1); for i=2:1:n

差分法求解偏微分方程MAAB

南京理工大学 课程考核论文 课程名称:高等数值分析 论文题目:有限差分法求解偏微分方程 姓名:罗晨 学号: 成绩: 有限差分法求解偏微分方程 一、主要内容 1.有限差分法求解偏微分方程,偏微分方程如一般形式的一维抛物线型方程:具体求解的偏微分方程如下: 2.推导五种差分格式、截断误差并分析其稳定性; 3.编写MATLAB程序实现五种差分格式对偏微分方程的求解及误差分析;

4.结论及完成本次实验报告的感想。 二、推导几种差分格式的过程: 有限差分法(finite-differencemethods )是一种数值方法通过有限个微分方程近似求导从而寻求微分方程的近似解。有限差分法的基本思想是把连续的定解区域用有限个离散点构成的网格来代替;把连续定解区域上的连续变量的函数用在网格上定义的离散变量函数来近似;把原方程和定解条件中的微商用差商来近似,积分用积分和来近似,于是原微分方程和定解条件就近似地代之以代数方程组,即有限差分方程组,解此方程组就可以得到原问题在离散点上的近似解。 推导差分方程的过程中需要用到的泰勒展开公式如下: ()2100000000()()()()()()()......()(()) 1!2!! n n n f x f x f x f x f x x x x x x x o x x n +'''=+-+-++-+-(2-1) 求解区域的网格划分步长参数如下: 11k k k k t t x x h τ ++-=?? -=?(2-2) 2.1古典显格式 2.1.1古典显格式的推导 由泰勒展开公式将(,)u x t 对时间展开得 2,(,)(,)( )()(())i i k i k k k u u x t u x t t t o t t t ?=+-+-?(2-3) 当1k t t +=时有 21,112,(,)(,)( )()(())(,)()() i k i k i k k k k k i k i k u u x t u x t t t o t t t u u x t o t ττ+++?=+-+-??=+?+?(2-4) 得到对时间的一阶偏导数 1,(,)(,)()=()i k i k i k u x t u x t u o t ττ+-?+?(2-5) 由泰勒展开公式将(,)u x t 对位置展开得 223,,21(,)(,)()()()()(())2!k i k i k i i k i i u u u x t u x t x x x x o x x x x ??=+-+-+-??(2-6) 当11i i x x x x +-==和时,代入式(2-6)得

一阶线性偏微分方程

第七章 一阶线性偏微分方程 研究对象 一阶线性齐次偏微分方程 0),,,(),,,() ,,,(2122121211=??++??+??n n n n n x u x x x X x u x x x X x u x x x X 1基本概念 1) 一阶线性齐次偏微分方程 形如 0),,,(),,,(),,,(2122121211=??++??+??n n n n n x u x x x X x u x x x X x u x x x X (7.1) 的方程,称为一阶线性齐次偏微分方程,其中n x x x ,,,21 是自变量,u 是n x x x ,,,21 的未知函数,n X X X ,,,21 是域n R D ?内的已知函数,并设n X X X ,,,21 在域D 内不同时为零。 2) 一阶拟线性偏微分方程 形如 );,,,();,,,();,,,(21211211z x x x Z x z z x x x Y x z z x x x Y n n n n n =??++?? (7.2) 的方程,称为一阶拟线性偏微分方程,其中Z Y Y Y n ;,,,21 是1+n 个变元z x x x n ;,,,21 的已知函数。n Y Y Y ,,,21 在其定义域1+?'n R D 内不同时为零。 所谓“拟线性”是指方程仅对未知函数的各个一阶偏导数是线性的,以下总设n Y Y Y ,,,21 和Z 在域D '内连续可微。 3) 特征方程组 常微分方程组 n n X dx X dx X dx === 2211 (7.3) 称为一阶线性齐次偏微分方程(7.1)的特征方程组。 常微分方程组

偏微分方程数值解实验报告

精品文档 偏微分方程数值解 上 机 实 验 报 告 (一)实验一 一、上机题目: 用线性元求解下列边值问题的数值解:

精品文档 ′′22?? ?? ??,0

精品文档 (二)实验二 四、上机题目: 求解 Helmholtz 方程的边值问题: u k 2u 1 ,于(0,1)*(0,1) u0,于1{ x0,0y1} U{0x1, y 1} 1{ x0,0y1} U{0x1, y1} u 0,于2{0x1, y 0} U { x1,0y1} n 其中 k=1,5,10,15,20 五、实验程序:

偏微分方程数值解课程设计

课程设计报告 课程:偏微分方程数值解学号: 姓名: 班级: 教师:

《偏微分方程数值解》 课程设计指导书 一.课程设计的目的 1.帮助掌握偏微分方程数值解相关知识。 2.理解偏微分方程数值解差分隐格式解决自由振动方程问题的方法。 3.锻炼编写程序代码的能力。 二.设计名称 差分法求自由振动问题的周期解。 三.设计要求 1.要求写出差分隐格式的理论方法。 2.要求编写matlab 程序,画出函数图形。 3.要求写出实验总结及心得体会。 四.设计题目 用差分法求自由振动问题的周期解: 2222000,,0|0,|sin (0,)(2,)t t u u x t t x u u x t u t u t π==???-=-∞<<∞>???? ??==??? =??? 要求用差分隐格式求解,其中14 θ= 。 五.设计细则 1.区域剖分: 构造上式的差分逼近,取空间步长h 和时间步长τ,用两族平行直线 ?? ?===±±=== ,2,1,0,, ,2,1,0,n n t t j jh x x n j τ 作矩形网格。 2.离散格式: 显格式: 于网点),(n j t x 用Taylor 展式,并整理方程得: ??? ?? ??--++=+-++==-+-++-121121102 10102100 )1(2)(),()()1()]()([2),(n j n j n j n j n j j j j j j j j u u r u u r u x x r x x r u x u τ?????

隐格式: 上述显格式并不是绝对稳定的差分格式,为了得到绝对稳定的差分格式,用第1-n 层、 n 层、1+n 层的中心差商的权平均去逼近xx u ,得到下列差分格式: ? ??? ?? ???+-++--++-=+-+-++==----+-++-+++-++-]22)21(2[2), ()()1()]()([2),(2111112112111112 211102 10102100h u u u h u u u h u u u a u u u x x r x x r u x u n j n j n j n j n j n j n j n j n j n j n j n j j j j j j j j θθθττ?????其中10≤≤θ是参数。当0=θ时就是显格式,而当4 1 =θ时可以证明该格式绝对稳定。 隐格式的矩阵形式是: ??? ???????? ???????????=??????????????????????????????????????????????? ?-+-+-+-+--+-+-+++122111121121 12222 222 2222221212121J J j n J n J n j n n z z z z z u u u u u r r r r r r r r r r r r θθθθ θθθθθ θ θθ 其中: 1 111111122]2()2)(21[(-----+-+-++-++--=n j n j n j n j n j n j n j n j j u u u u u u u u r z θθ 3.格式稳定性: 1)显格式: 显格式稳定的充分必要条件是:网格比1

偏微分方程数值及matlab实验报告.docx

偏微分方程数值实验报告八 实验题目:利用有限差分法求解 u ( x) u(x) f (x), u( 1) 0, u(1) 0. 真解为 u( x) e x 2 (1 x 2 ) 实现算法:对于两点边值问题 d 2u f , x l , dx 2 (1) u(a),u(b) , 其中 l ( a, b) (a b), f 为 l [ a,b] 上的连续函数, , 为给定常数 . 其相应的有限差分法的算法如下: 1.对求解区域做网格剖分,得到计算网格 .在这里我们对区间 l 均匀剖分 n 段,每个剖分单元 b a 的剖分步长记为 h . n 2.对微分方程中的各阶导数进行差分离散,得到差分方程 .运用的离散方法有: 方法一 :用待定系数和泰勒展开进行离散 d 2u( x i ) i 1 u( x i 1) i u( x i ) i 1 u( x i 1) d( x i )2 方法二:利用差商逼近导数 d 2u( x i ) u( x i 1 ) 2u( x i ) u( x i 1 ) d( x i )2 h 2 将(2) 带入 (1)可以得到 u(x i 1) 2u(x i ) u(x i 1 ) ) R i (u) , h 2 f ( x i 其中 R i (u) 为无穷小量,这时我们丢弃 R i (u) ,则有在 x i 处满足的计算公式: u(x i 1) 2u( x i ) u( x i 1 ) 1,..., n 1 h 2 f ( x i ), i 3.根据边界条件,进行边界处理 .由 (1)可得 u 0 , u n (2) (3) (4) 称(3)(4)为逼近 (1) 的差分方程,并称相应的数值解向量 U n 1 为差分解, u i 为 u( x i ) 的近似值 . 4.最后求解线性代数方程组,得到数值解向量U n 1 .

(完整版)偏微分方程的MATLAB解法

引言 偏微分方程定解问题有着广泛的应用背景。人们用偏微分方程来描述、解释或者预见各种自然现象,并用于科学和工程技术的各个领域fll。然而,对于广大应用工作者来说,从偏微分方程模型出发,使用有限元法或有限差分法求解都要耗费很大的工作量,才能得到数值解。现在,MATLAB PDEToolbox已实现对于空间二维问题高速、准确的求解过程。 偏微分方程 如果一个微分方程中出现的未知函数只含一个自变量,这个方程叫做常微分方程,也简称微分方程;如果一个微分方程中出现多元函数的偏导数,或者说如果未知函数和几个变量有关,而且方程中出现未知函数对几个变量的导数,那么这种微分方程就是偏微分方程。 常用的方法有变分法和有限差分法。变分法是把定解问题转化成变分问题,再求变分问题的近似解;有限差分法是把定解问题转化成代数方程,然后用计算机进行计算;还有一种更有意义的模拟法,它用另一个物理的问题实验研究来代替所研究某个物理问题的定解。虽然物理现象本质不同,但是抽象地表示在数学上是同一个定解问题,如研究某个不规则形状的物体里的稳定温度分布问题,由于求解比较困难,可作相应的静电场或稳恒电流场实验研究,测定场中各处的电势,从而也解决了所研究的稳定温度场中的温度分布问题。 随着物理科学所研究的现象在广度和深度两方面的扩展,偏微分方程的应用范围更广泛。从数学自身的角度看,偏微分方程的求解促使数学在函数论、变分法、级数展开、常微分方程、代数、微分几何等各方面进行发展。从这个角度说,偏微分方程变成了数学的中心。

一、MATLAB方法简介及应用 1.1 MATLAB简介 MATLAB是美国MathWorks公司出品的商业数学软件,用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言和交互式环境,主要包括MATLAB和Simulink两大部分。 1.2 Matlab主要功能 数值分析 数值和符号计算 工程与科学绘图 控制系统的设计与仿真 数字图像处理 数字信号处理 通讯系统设计与仿真 财务与金融工程 1.3 优势特点 1) 高效的数值计算及符号计算功能,能使用户从繁杂的数学运算分析中解脱出来; 2) 具有完备的图形处理功能,实现计算结果和编程的可视化; 3) 友好的用户界面及接近数学表达式的自然化语言,使学者易于学习和掌握; 4) 功能丰富的应用工具箱(如信号处理工具箱、通信工具箱等) ,

一阶偏微分方程基本知识

一阶偏微分方程基本知识 这一章我们来讨论一阶线性偏微分方程和一阶拟线性偏微分方程的解法,因为它们都可以化为常微分方程的首次积分问题,所以我们先来介绍常微分方程的首次积分。 1一阶常微分方程组的首次积分 1.1首次积分的定义 从第三章我们知道,n 阶常微分方程 ()()() 1,,'',',-=n n y y y x f y Λ, ( 1.1) 在变换 ()1'12,,,,n n y y y y y y -===L ( 1.2) 之下,等价于下面的一阶微分方程组 ()()()1 112221212,,,,,,,,,,,,,,. n n n n n dy f x y y y dx dy f x y y y dx dy f x y y y dx ?=?? ?=???? ?=??L L M M M M L ( 1.3) 在第三章中,已经介绍过方程组( 1.3)通解的概念和求法。但是除了常系数线性方程组外,求一般的( 1.3)的解是极其困难的。然而在某些情况下,可以使用所谓“可积组合”法求通积分,下面先通过例子说明“可积组合”法,然后介绍一阶常微分方程组“首次积分”的概念和性质,以及用首次积分方法来求解方程组( 1.3)的问题。先看几个例子。

例1 求解微分方程组 ()()22221,1.dx dy y x x y x y x y dt dt =-+-=--+- ( 1.4) 解:将第一式的两端同乘x ,第二式的两端同乘y ,然后相加,得到 ()() 12222-++-=+y x y x dt dy y dt dx x , ()()()2222221 12 d x y x y x y dt +=-++-。 这个微分方程关于变量t 和()22x y +是可以分离,因此不难求得其解为 122 2221C e y x y x t =+-+, ( 1.5) 1C 为积分常数。( 1.5)叫做( 1.4)的首次积分。 注意首次积分( 1.5)的左端(),,V x y t 作为x ,y ,和t 的函数并不等于常数;从上面的推导可见,当(),()x x t y y t ==时微分方程组( 1.4)的解时,(),,V x y t 才等于常数1C ,这里的常数1C 应随解而异。因为式( 1.4)是一个二阶方程组,一个首次积分( 1.5)不足以确定它的解。为了确定( 1.4)的解,还需要找到另外一个首次积分。 将第一式两端同乘y ,第二式两端同乘x ,然后用第一式减去第二式,得到 22y x dt dy x dt dx y +=-, 即 () 22y x dt dx y dt dy x +-=-, 亦即 1arctan -=?? ? ?? dt x y d 。 积分得

第一章 偏微分方程和一阶线性偏微分方程解

第一章 偏微分方程和一阶线性偏微分方程解 本章介绍典型的几个偏微分方程。给出了最简单的偏微分方程(一阶线性偏微分方程)解的特征线方法。 典型的偏微分方程:扩散方程t xx u ku =,t u k u =?;波动方程2tt xx u c u =,2tt u c u =?。这是本课程讨论的主要两类方程。 偏微分方程的各类边值条件也是本章讨论的一个重点。 §1.1 一维空间中的偏微分方程 例1 (刚性污染流的方程) 假设均匀直线管道中的水流含污染物质的线密度是(,)u x t (即x 处在时刻t 的污染物的密度) 。如果流速是c ,问题:(,)u x t 满足什么样的方程? 解 如图,在[,]x x x +?内的流体,经过时间t ?,一定处于[,]x c t x x c t +?+?+?。所含污染物应相同,即 (,)(,)x x x x c t x x c t u t d u t t d ξξξξ+?+?+?+?= +?? ? , 由此 (,)(,)u x t u x c t t t =+?+?, 从而, 0t x u cu +=。 【End 】 可见偏微分方程是一个至少为两元的函数及其偏导数所满足的方程。 例2 (扩散方程) 假设水流静止,在t ?时间内,流经x 处的污染物质(不计高阶无穷小)与该处浓度的方向导数(浓度变化)成正比,比例系数为k : ()x u dm t k dt ku dt x ?==?, 所以,在时间段12[,]t t 内,通过12[,]x x 的污染物为 2 1 2 1 [(,)(,)]t x x t k u x t u x t dt -?。 在时刻1t 和2t ,在12[,]x x 内的污染物分别为2 1 1(,)x x u x t dx ?和2 1 2(,)x x u x t dx ? ,由物质守恒定律 2 2 2 1 1 1 2 1 2 1 (,)(,)[(,)(,)]x x t x x x x t u x t dx u x t dx k u x t u x t dt -=-??? 由1t ,2t 的任意性,

偏微分方程数值及matlab实验报告

偏微分方程数值实验报告八 实验题目:利用有限差分法求解 . 0)1(,0)1(),()()(==-=+''-u u x f x u x u 真解为 ) 1()(22 x e x u x -=-实现算法:对于两点边值问题 , )(,)(,,d 22βα==∈=-b u a u l x f dx u (1) 其中),(b a l =f b a ),(<为],[b a l =上的连续函数,βα,为给定常数. 其相应的有限差分法的算法如下: 1.对求解区域做网格剖分,得到计算网格.在这里我们对区间l 均匀剖分n 段,每个剖分单元的剖分步长记为n a b h -= .2.对微分方程中的各阶导数进行差分离散,得到差分方程.运用的离散方法有:方法一:用待定系数和泰勒展开进行离散 )()()()(d ) (d 11112 2++--++≈i i i i i i i i x u x u x u x x u ααα方法二:利用差商逼近导数 2 112 2) ()(2)()(d )(d h x u x u x u x x u i i i i i -++-≈(2) 将(2)带入(1)可以得到 )()() ()(2)(2 11u R x f h x u x u x u i i i i i +=+-- -+, 其中)(u R i 为无穷小量,这时我们丢弃)(u R i ,则有在i x 处满足的计算公式: 1,...,1)() ()(2)(2 11-==+-- -+n i x f h x u x u x u i i i i ,(3) 3.根据边界条件,进行边界处理.由(1)可得 β α==n u u ,0(4) 称(3)(4)为逼近(1)的差分方程,并称相应的数值解向量1-n U 为差分解,i u 为)(i x u 的近似值.4.最后求解线性代数方程组,得到数值解向量1 -n U .

偏微分方程实验1

《偏微分方程数值解》 课程实验报告(一)

1 实验题目 给定初值问题(注:共两次上机时间4学时) 1 02 )0(42'≤≤?? ?=--=x y x y y 其精确解为 12)(2+-=-x e x y x 。取h=0.1, 分别用显式Euler 法、隐式梯形法、二级二阶Runge-Kutta 、三级三阶Runge-Kutta 、四级四阶Runge-Kutta 计算数值解,并与精确解比较。 2 求解方法 1、每种方法的迭代公式: (1)显式Euler 法 p[i]=p[i-1]+h*(*f)(x[i-1],p[i-1]); (2)隐式梯形法 p[i]=p[i-1]+h*(*f)(x[i-1],p[i-1]);; (3)二级二阶Runge-Kutta K1=(*f)(x[i-1],p[i-1]); K2=(*f)(x[i-1]+h/2,p[i-1]+h*K1/2); p[i]=p[i-1]+h*K2; (4)三级三阶Runge-Kutta K1=(*f)(x[i-1],p[i-1]); K2=(*f)(x[i-1]+h/2,p[i-1]+h*K1/2); K3=(*f)(x[i-1]+h/2,p[i-1]+h*K2/2); K4=(*f)(x[i-1]+h,p[i-1]+h*K3); p[i]=p[i-1]+h*(K1+2*K2+2*K3+K4)/6; (5)四级四阶Runge-Kutta K1=(*f)(x[i-1],p[i-1]);

K2=(*f)(x[i-1]+h/2,p[i-1]+h*K1/2); K3=(*f)(x[i-1]+h/2,p[i-1]+h*K2/2); K4=(*f)(x[i-1]+h,p[i-1]+h*K3); p[i]=p[i-1]+h*(K1+2*K2+2*K3+K4)/6; 3 程序源代码 #include #include #define N 11 double fun(double x,double y) //微分方程{ return -2*y-4*x; } double p_fun(double x) //原函数 { return exp(-2*x)-2*x+1; } double* Exact(double x0, double h) //精确解{ int i; double *p=new double[N]; for(i=0;i

(整理)一阶线性偏微分方程.

第七章 一阶线性偏微分方程 例7-1 求方程组 ()()()yz B A Cdz xz A C Bdy yz C B Adx -=-=- 通积分,其中C B A ,,为互不 相等的常数。 解 由第一个等式可得 xyz ydy A C B xyz xdx C B A -=-, 即有 0=---ydy A C B xdx C B A , 两边积分得方程组的一个首次积分 122,C y A C B x C B A z y x Φ=---= ),(。 由第二个等式可得 xyz zdz B A C xyz ydy A C B -=-, 即有 0=---zdz B A C ydy A C B , 两边积分得方程组的另一个首次积分 222,C z B A C y A C B z y x Ψ=---= ),(。 由于,雅可比矩阵 ? ???? ?????------=????? ???? ????ψ??ψ??ψ ??Φ??Φ ??Φ ?=?ψΦ?z B A C y A C B y A C B x C B A y y x z y x z y x 002),,(),( 的秩为2,这两个首次积分相互独立,于是原方程组的通积分为 122C y A C B x C B A =--- 222C z B A C y A C B =--- 。

评注:借助于方程组的首次积分求解方程组的方法称为首次积分法。要得到通积分需要求得n 个独立的首次积分,n 为组成方程组的方程个数。用雅可比矩阵的秩来验证首次积分的独立性。 例7-2 求方程组 () () ???????-+--=-+-=11d 222 2y x y x dt dy y x x y dt x 的通解。 解 由原方程组可得 )1)((2222-++-=+y x y x dt dy y dt dx x 即 dt y x y x y x d )1)((2)(2 2 2 2 2 2 -++-=+ 这个方程关于变量t 和2 2 y x +是可以分离的,因此易求得它的通积分为 122 2221),,(C e y x y x t y x t =+-+=Φ 这是原方程组的一个首次积分。 再次利用方程组,得到 )(22y x dt dx y dt dy x +-=-, 即有 1arctan -=?? ? ?? x y dt d 由此得到原方程组的另一个首次积分 2arctan ),,(C t x y t y x =+=ψ 。 由于,雅可比矩阵为 ()( ) ???? ? ?????? ?++-++=????????? ????ψ??ψ ??Φ??Φ ?=?ψΦ?2222 222 222 2222),(),(y x x y x y e y x y e y x x y x y x y x t t ,

偏微分方程上机实验报告.doc

上机实验2:五点差分格式法 偏微分方程(Matlab )实验报告 ——五点差分格式法 一、 实验题目 设G 是形如下图的十字形域,由五个相等的单位正方形组成,用五点差分格式求下列边值问题的数值解: 22 2 21,u u G x y ??+=-???于u=0,于G 二、 实验原理 取定沿X 轴和Y 轴方向的步长1h 和2h ,() 12 22 1 2 h h h =+,作两族与坐 标轴平行的直线:x=i 1h ,y=j 2h ,,0,1,2,i j =±± 若(,i j x y )为正则内点,沿x,y 方向分别用二阶中心差商代替 xx yy u u 和则得 1,1,,1,1 2 212 22[ ]i j ij i j i j ij i j ij u u u u u u f h h +-+--+-+-+ = 特别取正方形网格:12h h h ==,则原差分方程可简化为 2 1,,11,,11()44 ij i j i j i j i j ij h u u u u u f --++-+++= 三、 实验程序 1)function uxy = EllIni2Uxl(x,y) format long ;

uxy = 0; 2)function uxy = EllIni2Uxr(x,y) format long; uxy = y*(2-y); 3)function uxy = EllIni2Uyl(x,y) format long; uxy = 0; 4)function uxy = EllIni2Uyr(x,y) format long; if x < 1 uxy = x; else uxy = 2 - x; end 5)function u = peEllip5(nx,minx,maxx,ny,miny,maxy) format long; hx = (maxx-minx)/(nx-1); hy = (maxy-miny)/(ny-1); u0 = zeros(nx,ny); for j=1:ny u0(j,1) = EllIni2Uxl(minx,miny+(j-1)*hy); u0(j,nx) = EllIni2Uxr(maxx,miny+(j-1)*hy); end for j=1:nx u0(1,j) = EllIni2Uyl(minx+(j-1)*hx,miny); u0(ny,j) = EllIni2Uyr(minx+(j-1)*hx,maxy); end A = -4*eye((nx-2)*(ny-2),(nx-2)*(ny-2)); b = ones((nx-2)*(ny-2),1).*(-1); for i=1:(nx-2)*(ny-2) if mod(i,nx-2) == 1 if i==1 A(1,2) = 1; A(1,nx-1) = 1; b(1) = - u0(1,2) - u0(2,1); else if i == (ny-3)*(nx-2)+1 A(i,i+1) = 1; A(i,i-nx+2) = 1;

偏微分方程数值解法

《偏微分方程数值解法》 课程设计 题目:六点对称差分格式解热传导方程的初边值 问题 姓名:王晓霜 学院:理学院 专业:信息与计算科学 班级:0911012 学号:091101218 指导老师:翟方曼 2012年12月14日

一、题目 用六点对称差分格式计算如下热传导方程的初边值问题 222122,01,01(,0),01 (0,),(1,),01x t t u u x t t x u x e x u t e u t e t +???=<<<≤?????=≤≤??==≤≤??? 已知其精确解为 2(,)x t u x t e += 二、理论 1.考虑的问题 考虑一维模型热传导方程 (1.1) )(22x f x u a t u +??=??,T t ≤<0 其中a 为常数。)(x f 是给定的连续函数。(1.1)的定解问题分两类: 第一,初值问题(Cauchy 问题):求足够光滑的函数()t x u ,,满足方程(1.1)和初始条件: (1.2) ()()x x u ?=0,, ∞<<∞-x 第二,初边值问题(也称混合问题):求足够光滑的函数()t x u ,,满足方程(1.1)和初始条件: ()13.1 ()()x x u ?=0,, l x l <<- 及边值条件 ()23.1 ()()0,,0==t l u t u , T t ≤≤0 假定()x f 和()x ?在相应的区域光滑,并且于()0,0,()0,l 两点满足相容条件,则上述问题有唯一的充分光滑的解。 现在考虑边值问题(1.1),(1.3)的差分逼近 取 N l h = 为空间步长,M T =τ为时间步长,其中N ,M 是自然数, jh x x j ==, ()N j ,,1,0Λ=; τk y y k ==, ()M k ,,1,0Λ=

偏微分中心差分格式实验报告(含matlab程序)

二阶常微分方程的中心差分求解 学校:中国石油大学(华东)理学院 姓名:张道德 一、 实验目的 1、 构造二阶常微分边值问题: 22,(),(), d u Lu qu f a x b dx u a u b αβ?=-+=<

11122 222222333222122112 100121012010012 00N N N u f q h h u f q h h h u f q h h h q u f h h ---???? ??+-???? ??? ???? ???????-+-? ?????? ???????????=-+? ?????? ???????????-???? ????????-+????? ?? ????? 可以看出系数矩阵为三对角矩阵,而对于系数矩阵为三对角矩阵的方程组可以用“追赶法”求解,则可以得出二阶常微分方程问题的数值解。 四、 举例求解 我们选取的二阶常微分方程边值问题为: 2 22242,01 (0)1,(1), x d u Lu x u e x dx u u e ?=-+=-<

偏微分方程数值及matlab实验报告(8)

偏微分方程数值及m a t l a b 实验报告(8) -标准化文件发布号:(9456-EUATWK-MWUB-WUNN-INNUL-DDQTY-KII

偏微分方程数值实验报告八 实验题目:利用有限差分法求解 . 0)1(,0)1(),()()(==-=+''-u u x f x u x u 真解为 ) 1()(22 x e x u x -=-实现算法:对于两点边值问题 , )(,)(,,d 22βα==∈=-b u a u l x f dx u (1) 其中),(b a l =f b a ),(<为],[b a l =上的连续函数,βα,为给定常数. 其相应的有限差分法的算法如下: 1.对求解区域做网格剖分,得到计算网格.在这里我们对区间l 均匀剖分n 段,每个剖分单元的剖分步长记为n a b h -= .2.对微分方程中的各阶导数进行差分离散,得到差分方程.运用的离散方法有:方法一:用待定系数和泰勒展开进行离散 )()()()(d ) (d 11112 2++--++≈i i i i i i i i x u x u x u x x u ααα方法二:利用差商逼近导数 2 112 2) ()(2)()(d )(d h x u x u x u x x u i i i i i -++-≈(2) 将(2)带入(1)可以得到 )()() ()(2)(2 11u R x f h x u x u x u i i i i i +=+-- -+, 其中)(u R i 为无穷小量,这时我们丢弃)(u R i ,则有在i x 处满足的计算公式: 1,...,1)() ()(2)(2 11-==+-- -+n i x f h x u x u x u i i i i ,(3) 3.根据边界条件,进行边界处理.由(1)可得 β α==n u u ,0(4) 称(3)(4)为逼近(1)的差分方程,并称相应的数值解向量1-n U 为差分解,i u 为)(i x u 的近似值.4.最后求解线性代数方程组,得到数值解向量1 -n U .

微分方程数值解实验报告

微分方程数值解法课程设计报告 班级:_______ 姓名: ___ 学号:__________ 成绩: 2017年 6月 21 日

摘要 自然界与工程技术中的很多现象,可以归结为微分方程定解问题。其中,常微分方程求解是微分方程的重要基础内容。但是,对于许多的微分方程,往往很难得到甚至不存在精确的解析表达式,这时候,数值解提供了一个很好的解决思路。,针对于此,本文对常微分方程数值解法进行了简单研究,主要讨论了一些常用的数值解法,如欧拉法、改进的欧拉法、Runge—Kutta方法、Adams法以及椭圆型方程、抛物型方程的有限差分方法等,通过具体的算例,结合MATLAB求解画图,初步给出了一般常微分方程数值解法的求解过程。同时,通过对各种方法的误差分析,让大家对各种方法的特点和适用范围有一个直观的感受。 关键词:微分方程数值解、MATLAB

目录 摘要 (2) 目录 (3) 第一章常微分方程数值解法的基本思想与原理 (4) 1.1 常微分方程数值解法的基本思路 (4) 1.2用matlab编写源程序 (4) 1.3 常微分方程数值解法应用举例及结果 (5) 第二章常系数扩散方程的经典差分格式的基本思想与原理 (6) 2.1 常系数扩散方程的经典差分格式的基本思路 (6) 2.2 用matlab编写源程序 (7)

2.3 常系数扩散方程的经典差分格式的应用举例及结果 (8) 第三章椭圆型方程的五点差分格式的基本思想与原理 (10) 3.1 椭圆型方程的五点差分格式的基本思路 (10) 3.2 用matlab编写源程序 (10) 3.3 椭圆型方程的五点差分格式的应用举例及结果 (12) 第四章总结 (12) 参考文献 (12)

一阶偏微分方程初步

第六章 一阶偏微分方程初步 6.2 一阶常微分方程的首次积分 1.求下列方程的首次积分及通积分。 (1)???????==z y dx dz y z dx dy 22 (上式为(1)式,下式为(2)式) 解:(2)式除以(1)式可得 3 3 z y dy dz = 即033=-dy y dz z 积分可得: 14 4c y z =- 从而求得一个首次积分:44y z -=φ 其次,由???==dx y zdz dx z ydy 22 (上式为(3)式,下式为(4)式) 将(3)和(4)相加可得:()dx z y zdz ydy 22+=+ 即dx z y zdz ydy 22222=++ 积分可得到又一个首次积分x e z y 22 2+=ψ。于是微分方程组的通积分为 ?????=+=-2222144c e z y c y z x (2)???????-=-=y x x dt dy y x y dt dx (上式为(1)式,下式为(2)式) 解:(1)式除以(2)式,可得,x y dy dx =即ydy xdx =

积分后得122c y x =-。从而求得一个首次积分2 2y x -=φ (2)-(1)可得1)(=-dt x y d 即dt x y d =-)( 积分可得2c t x y =--。从而有得到一个首次积分t x y --=ψ 于是微分方程组的通积分为 ? ??=--=-2122c t x y c y x (3)xy dz xz dy yz dx == 解:由xz dy yz dx =,可得0=-ydy xdx 。积分可得:122c y x =- 从而求得一个首次积分22y x -=φ 由xy dz xz dy =,可得0=-zdz ydy ,积分可得222c z y =-。于是又得到一个首次积分: 22z y -=ψ 于是微分方程组的通积分为 ???=-=-2 22122c z y c y x (4)x y dz z x dy y z dx -=-=- 解:利用合比定理可得: x y dz z x dy y z dx -=-=-=()0z y x d ++ 由此可得()0=++z y x d ,于是得到一个首次积分: z y x ++=φ 另外,我们有()()()() 2222222 22z y x d x y z zdz z x y ydy y z x xdx ++=-=-=- 于是得到另一个首次积分222z y x ++=ψ。于是微分方程的通积分为: ???=++=++2 2221c z y x c z y x

相关文档