下载文档原格式
毕业论文（设计）内容介绍论文（设计）题目微分方程数值解选题时间2010.11.20完成时间2011.3.9论文（设计）字数33840关键词微分方程，周期解，边值问题论文（设计）题目的来源、理论和实践意义：微分方程是数学科学联系实际问题的主要桥梁之一，它是含有未知函数及其导数的方程。
常微分方程的求解是现代科学研究和工程技术中经常遇到的实际问题，然而，从实际问趣中建立出来的微分方程往往具有非常复杂的形式，有些解析式难以计算，有些则根本不能用解析式来表达，所以利用数值解法‘叫求解实际问题就显得非常重要。
论文（设计）的主要内容及创新点：如果未知函数的自变量是一个，称为常微分方程;自变量多于一个，称为偏微分方程。
在科学研究和工程计算中碰到的许多微分方程，根本不存在解析解，或者求解析解的代价很大，求解过程过于复杂，在这种情况下，我们只能借助于数值计算来求方程的数值解。
附：论文（设计）本人签名：年月日目录中文摘要 (5)第一章常微分方程的解 (6)第一节常微分方程的基本概念 (6)第二节常微分方程的12步骤 (10)第三节偏导数的方程 (14)第二章递增方程的应用 (17)第一节递增数列 (17)第二节数列的极限 (20)第三章与积分有关的数列的极限问题 (24)第一节积分的应用 (24)第二节单调定性的松弛法 (26)第三节松弛算法法的证明 (33)第四章简单的单步法及基本概念 (36)第一节解初值问题的梯形法 (36)第二节左矩形公式 (39)第三节隐式Euler方法 (40)第四节预估–校正Euler方法 (42)参考文献 (45)摘要：常微分方程的形成和发展是与力学、天文学、物理学及其他自然科学技术的发展互相促进和互相推动的。
本文第一章讲述了常微分方程的发展历史，第二章介绍了一系列常微分方程的周期解和边值问题，说明了其研究现状，第三章举例说明了其在生态学和军事上的应用。
无论在数学研究还是在自然科学以及其他应用科学，常微分方程都显现出其重要的理论和应用价值。
第1页共18页微分方程数值解实验指导书李光云数学与计算科学学院O 一0十月第2页共18页一、概述本课程实验指导书是根据陆金甫，关治编著的《偏微分方程数值解法（第二版）》编写的。
通过上机实验，可帮助学生理解偏微分方程数值计算方法的基本思想，巩固学生所学过的算法，同时让学生进一步学习和掌握Matlab软件在数值计算方面中的应用。
二、实验环境本书选择的实验环境是计算机以及软件Matlab（版本6.5以上）一套。
三、实验课时安排共7次实验，16课时，试验五为综合设计性试验4课时，其它实验2个课时。
四、实验要求上机完成实验指导书中所规定的内容，自行按实验指导书要求完成程序设计和调试，并提交每次实验的实验报告，附带算法程序清单和算法输出结果。
五、实验考核要求上机完成试验内容，并提交一份算法程序清单和数值结果。
第3页共18页实验一 双曲型方程的迎风格式、实验目的1、 掌握双曲型方程的迎风格式的算法思想 ，格式稳定的条件;2、 培养编程与上机调试能力。
、实验课时：2个课时 三、实验准备1、 了解偏微分方程的分类，以及双曲型方程的特点；2、 熟悉求解偏微分方程有限差分法的基本原理，双曲型方程迎风格式的推导和 特点，熟悉用迎风格式求解双曲型方程的算法步骤；3、 熟悉Matlab 软件的基本命令及其数值计算中的基本命令和函数。
4、 了解双曲型方程迎风格式稳定的条件。
四、实验内容用双曲型方程的迎风格式n nU j -U j 」0, hn 1 nnnU j —U jU j1—U j c a 0,Th其中，• ,h 分别为时间步长和空间步长。
求解初值问题賦”,ts,I u(x,0 )=u°(x ),其中,p!,x",U。
八 0,x 0.五、基本思想及主要步骤考虑对流方程:U :ua 0,x . ,t 0* *-t:x其a 为给定常数。
其迎风格式的基本思想是在方程中关于空间的偏导数用在特征线n 1 nU j Uja 0 a<0.计算〔0,1〕区间的u 值至t n 二0.5,并画出解的图形。
微分方程数值解法实验报告班级：姓名：学号：日期：一、实验目的1、熟悉微分方程（组）数值解的Euler算法，改进的Euler算法和Runge-Kutta算法，利用matlab软件实现微分方程数值解法来求解具体试题；2、虽然求解常微分方程有各种各样的解析解，但解析方法只能用来求解一些特殊类型的方程，通常它们无法求出解析解，而需要数值方法来近似求解。
因此产生了常微分方程初值问题的数值计算方法，常微分方程数值解法是通过计算机便捷的求解近似值。
二、基本理论及背景1、在数值求解常微分方程中，主要有有限差分计算和有限元计算两大类方法，其中在有限差分计算方法中有一类方法称为龙格－库塔(Runge_Kutta)方法。
四阶的龙格－库塔方法为最佳的计算格式。
2、参考三中的代码，分别用Euler算法，改进的Euler算法和Runge-Kutta 算法实现微分方程（组）的数值求解，完成下列题目：三、算法设计及实现1、算法设计，通过Euler算法，改进的Euler算法和Runge-Kutta三种算法来实现微分方程（组）的数值求解；2、程序文件及功能清单：(1) Euler Method:function [x,y]=EulerDSolve(f,ab,y0,h)x=(ab(1):h:ab(2))';n=length(y0);y=zeros(length(x),n);y(1,:)=y0';for k=2:length(x)y(k,:)=y(k-1,:)+h*feval(f,x(k-1),y(k-1,:)')';end；(2) Improved Euler Method:function [x,y]=MEulerDSolve(f,ab,y0,h)x=(ab(1):h:ab(2))';n=length(x);y=zeros(n,length(y0));y(1,:)=y0';for k=2:nyp=y(k-1,:)+h*feval(f,x(k-1),y(k-1,:));yc=y(k-1,:)+h*feval(f,x(k),yp);y(k,:)=(yp+yc)/2;end(3) Runge-Kutta Method:function [TOut,YOut]=Runge_Kutta(f,ab,y0,h)TOut=(ab(1):h:ab(2))';n=length(TOut);YOut=zeros(n,length(y0));YOut(1,:)=y0';for k=2:nx=TOut(k-1); y=YOut(k-1,:)';K1=feval(f,x,y);K2=feval(f,x+h/2,y+K1*h/2);K3=feval(f,x+h/2,y+K2*h/2);K4=feval(f,x+h,y+K3*h);YOut(k,:)=(y+(K1+2*K2+2*K3+K4)*h/6)';end四、实验步骤1、打开MATLAB软件，新建 *.m文件，在m文件的窗口中编辑Euler算法的函数程序，另建一m文件，编辑自己改进的Euler算法的函数程序，再新建一m文件，在窗口中编辑Runge-Kutta算法的函数程序，并全部保存在指定的文件夹下；2、将MATLAB软件的工作页面的工具栏下的目标文件指向指定的文件夹；3、分别调用上述三种算法的函数，实现微分方程（组）的数值求解完成给定的实验题目；4、输出结果和初步分析说明（见附页）。
微分方程数值解法实验报告2篇微分方程数值解法实验报告（一）在实际科学与工程问题中，我们经常会遇到微分方程的求解。
然而，许多微分方程往往没有解析解，这就需要我们利用数值方法来获得近似解。
本篇实验报告将介绍两种常见的微分方程数值解法：欧拉方法和改进的欧拉方法。
一、欧拉方法欧拉方法是最简单的微分方程数值解法之一。
其基本原理为离散化微分方程，将微分方程中的导数用差商代替。
设要求解的微分方程为dy/dx = f(x, y)，步长为h，则可用以下公式进行递推计算：y_{n+1} = y_n +hf(x_n, y_n)二、算法实现为了对欧拉方法进行数值实验，我们以一阶线性常微分方程为例：dy/dx = x - y, y(0) = 1步骤如下：（1）选择合适的步长h和求解区间[a, b]，这里我们取h=0.1，[a, b] = [0, 1]；（2）初始化y_0 = 1；（3）利用欧拉方法递推计算y_{n+1} = y_n + 0.1(x_n - y_n)；（4）重复步骤（3），直到x_n = 1。
三、实验结果与讨论我们通过上述步骤得到了在[0, 1]上的近似解y(x)。
下图展示了欧拉方法求解的结果。
从图中可以看出，欧拉方法得到的近似解与精确解有一定的偏差。
这是因为欧拉方法只是通过递推计算得到的近似解，并没有考虑到更高阶的误差。
所以在需要高精度解时，欧拉方法并不理想。
四、改进的欧拉方法针对欧拉方法的不足，我们可以考虑使用改进的欧拉方法（也称为改进的欧拉-柯西方法）。
它是通过利用前后两个步长欧拉方法得到的结果作为差商的中间项，从而提高了求解精度。
一阶线性常微分方程的改进欧拉方法可以表示为：y_{n+1} = y_n + h(f(x_n, y_n) + f(x_{n+1}, y_n + hf(x_n,y_n)))/2五、算法实现与结果展示将改进的欧拉方法应用于前述的一阶线性常微分方程，我们同样选择h=0.1，[a, b] = [0, 1]。
微分方程数值解实验报告实验目的：掌握微分方程数值解的基本方法，能够利用计算机编程求解微分方程。
实验原理：微分方程是自然科学与工程技术中常见的数学模型，它描述了变量之间的关系及其随时间、空间的变化规律。
解微分方程是研究和应用微分方程的基础，但有很多微分方程无法找到解析解，只能通过数值方法进行求解。
本实验采用欧拉方法和改进的欧拉方法求解微分方程的初值问题：$$\begin{cases}\frac{dy}{dt}=f(t,y)\\y(t_0)=y_0\end{cases}$$其中，$f(t,y)$是给定的函数，$y(t_0)=y_0$是已知的初值条件。
欧拉方法是最基本的数值解法，其步骤如下：1.给定$t_0$和$y_0$2.计算$t_{i+1}=t_i+h$，其中$h$是步长3. 计算$y_{i+1}=y_i+hf(t_i,y_i)$4.重复步骤2、3直到达到终止条件改进的欧拉方法是对欧拉方法进行改进，通过利用函数$y(t)$在$t+\frac{1}{2}h$处的斜率来更准确地估计$y_{i+1}$，其步骤如下：1.给定$t_0$和$y_0$2.计算$t_{i+1}=t_i+h$，其中$h$是步长3. 计算$y_*=y_i+\frac{1}{2}hf(t_i,y_i)$4. 计算$y_{i+1}=y_i+hf(t_i+\frac{1}{2}h,y_*)$5.重复步骤2、3、4直到达到终止条件实验步骤：1.编写程序实现欧拉方法和改进的欧拉方法2.给定微分方程和初值条件3.设置步长和终止条件4.利用欧拉方法和改进的欧拉方法求解微分方程5.比较不同步长下的数值解与解析解的误差6.绘制误差-步长曲线，分析数值解的精度和收敛性实验结果：以一阶常微分方程$y'=3ty+t$为例，给定初值$y(0)=1$，取步长$h=0.1$进行数值求解。
利用欧拉方法求解微分方程得到的数值解如下：\begin{array}{cccc}t & y_{\text{exact}} & y_{\text{Euler}} & \text{误差} \\ \hline0.0&1.000&1.000&0.000\\0.1&1.035&1.030&0.005\\0.2&1.104&1.108&0.004\\0.3&1.212&1.217&0.005\\0.4&1.360&1.364&0.004\\0.5&1.554&1.559&0.005\\0.6&1.805&1.810&0.005\\0.7&2.131&2.136&0.005\\0.8&2.554&2.560&0.006\\0.9&3.102&3.107&0.006\\1.0&3.807&3.812&0.005\\\end{array}利用改进的欧拉方法求解微分方程得到的数值解如下：\begin{array}{cccc}t & y_{\text{exact}} & y_{\text{Improved Euler}} & \text{误差} \\\hline0.0&1.000&1.000&0.000\\0.1&1.035&1.035&0.000\\0.2&1.104&1.103&0.001\\0.3&1.212&1.211&0.001\\0.4&1.360&1.358&0.002\\0.5&1.554&1.552&0.002\\0.6&1.805&1.802&0.003\\0.7&2.131&2.126&0.005\\0.8&2.554&2.545&0.009\\0.9&3.102&3.086&0.015\\1.0&3.807&3.774&0.032\\\end{array}误差-步长曲线如下：实验结论：通过对比欧拉方法和改进的欧拉方法的数值解与解析解的误差，可以发现改进的欧拉方法具有更高的精度和收敛性。
,nn a 均不为零，,)nn a ，并将)D D -+ ，令1x B x =,,0,1,2,n k =(0),)Tn x 为初始向量。
,)nn a ，）可以写成()D L -2(B D =为迭代矩阵构成的迭代法,,0,1,2,n k =112()((1))()2/(11)cos D L R D D L bu hϖϖϖϖϖμπ---+--=+-=解法的迭代次数以及差分解(,)(1/64,1/128)h u x y h =与精确解()(,)sin sin x y u x y e x y πππ+=的精度。
2.在MATLAB a R 2014软件上编写代码；3.运行得出结论。
【实验过程】（实验步骤、记录、数据、分析）1.分析问题。
2.在MATLAB a R 2014上编写代码，调试。
3.运行代码，发现问题后改进代码。
4.分别输入n 值为64和128，得出结论。
【实验结论】（结果）n=64时: 1.Jacobi 迭代 迭代次数 i =78112.Guass-Seidel 迭代 迭代次数 i =42113.SOR 超松弛法 迭代次数 i =121n=128时：1.Jacobi迭代迭代次数i =200012.Guass-Seidel迭代迭代次数i =156603.SOR超松弛法迭代次数i =242【实验小结】（收获体会）通过本次实验，掌握了求解Poisson方程第一边值问题的五点差分格式，比较了三种解法的迭代次数以及差分解与精确解的精度，更加熟悉的掌握了利用MATLAB求解数学问题的方法。
附录1：源程序附录2：实验报告填写说明1．实验项目名称：要求与实验教学大纲一致。
2．实验目的：目的要明确，要抓住重点，符合实验教学大纲要求。
3．实验原理：简要说明本实验项目所涉及的理论知识。
4．实验环境：实验用的软、硬件环境。
5．实验方案（思路、步骤和方法等）：这是实验报告极其重要的内容。
概括整个实验过程。
对于验证性实验，要写明依据何种原理、操作方法进行实验，要写明需要经过哪几个步骤来实现其操作。
微分方程数值解 的教学研究与实践杨韧,杨光崇,谢海英(成都信息工程学院数学学院,四川成都,610225)收稿日期:2008-12-08;修改日期:2009-06-18.作者简介:杨韧(1957-),女,江苏大丰人,学士,副教授,研究方向为微分方程和数值计算,E _mail:ryan g@ ;杨光崇(1963-),男,四川泸州人,硕士,教授,研究方向为微分方程数值解,E _mail:gcyang@ ;谢海英(1970-),女,四川遂宁人,硕士,讲师,研究方向为优化与决策,E _mail:hyxie@ .摘要 结合近几年我院信息与计算科学专业微分方程数值解课程教学改革的研究与实践,对微分方程数值解课程的教学内容、课程体系、教学方法和教学手段进行了深入的探讨.实践表明,通过运用多元化的现代教学手段并加强教学实践环节培养了学生综合应用知识解决问题的能力,教学效果良好.关键词 微分方程数值解;教学研究;教学实践.中图分类号 G420;O241.81 引言微分方程数值解是我院信息与计算科学本科专业的专业方向课,它是一门具有较强的实际背景、专门研究科学计算的课程.我们针对学院定位,构建适合本专业的课程内容体系,强调理论联系实际,重视数值计算方法在实际问题中的应用.经过多年的探索、实践和发展,我院的微分方程数值解课程形成了自己的课程内容体系,注重理论联系实际,重视数学建模的思想方法,坚持课堂教学与应用训练紧密结合,使学生在应用中扩充和贯通相关知识,给学生以更多的自主思考空间,激发学生的学习欲望,鼓励学生参与分析讨论,着重掌握解决具体问题的实际方法和计算原理.2 教学内容根据教育部课程教学指导委员会颁发的信息与计算科学专业规范对微分方程数值解课程的基本要求和我院本科学生的实际情况,对课程内容作了一定的取舍并有所侧重.经过调整,教学内容趋于全面,且突出了重点,尤其是增加了学生应用理论知识解决问题的实践教学环节部分.微分方程数值解课程的理论性很强,为了化简难点,便于学生能够正确理解教材内容,微分方程数值解课程的实验课设置了常微分方程初值问题的数值解、热传导方程初边值问题的差分模拟、Laplace 方程Dirichlet 边值问题的差分模拟等3个验证性实验和用龙格 库塔方法求解实际问题、用古典隐式差分格式求解抛物型方程初边值问题、五点差分格式求解波动方程混合问题等3个综合性实验,验证性实验的计算结果验证数值格式的合理性;综合性实验解决实际问题,训练由物理问题建立数学模型的能力,学习应用微分方程数值解的数值方法编制程序的动手能力,训练分析、归纳总结问题的综合能力.3 教学改革与教学研究微分方程数值解这门课理论性强,公式推导烦琐、枯燥,计算量大.这就需要采用适合学生的教学方法进行教学,在教学方法与教学手段改革方面,我院微分方程数值解课程坚持教师的主导作用与学生主体作用相结合,充分调动学生学习的主动性与自觉性,突出个性,因材施教,将教学改革与课程建设建立在现代教育技术的平台上,形成比较一致的教学思想和指导方针.3.1 设计思想微分方程数值解是我院信息与计算科学专业学生必修的核心专业课之一.针对新世纪科技人才对大学生数学素质的要求,通过对一般本科院校信息与计算科学专业课程教学现状的分析,我们提出了信息与计算科学专业微分方程数值解的课程教学内容与体系结构改革的设计思想.(1)适应当前高等教育从精英教育过渡到大众教育的需要,针对一般本科院校的教学实际,选择适当的教学定位.(2)根据各门课程教学内容之间的内在本质联系,整合课程体系结构.(3)优化教学内容,合理安排理论体系.在保证教学内容的科学性与系统性的前提下逐步培养学生综合运用所学知识分析、解决实际问题的能力.(4)加强实践性教学,培养学生用科学计算的方法来抽象和处理实际问题的能力.将数学软件MA TLAB 引入课堂教学,展示计算机解决实际问题的动态过程和手段,培养学生利用计算机处理实际问题的能力,解决现有教材和黑板教学不能实现的非常复杂的科学计算问题,使教学内容直观、明了、生动.124高等数学研究ST U DIES IN COL L EG E M A T H EM AT ICSVo l.13,N o.1Jan.,20103.2重点、难点的处理方法(1)微分方程数值解内容抽象,学生难于理解在保证教材内容科学性的前提下,安排由浅入深的内容次序以及简捷、直观的理论体系.课程始终贯以连续问题离散化的基本思想,注重与先行课程的衔接,力图达到与相关学科的相互渗透和利用.采用数值分析的基本技巧详细分析微分方程的离散过程,比如用差商代替微商、数值积分、数值微分、Taylor展式等基本方法推导各种差分格式,使理论的叙述和推导过程简化.在分析差分格式稳定性时,分别用高等代数中的矩阵理论和数学分析中的Fourier级数来进行讨论.根据一般本科院校教学的实际需要,我们适当地调整了部分理论的推导证明.例如,对于同一类微分方程方程的多种差分格式,通过详细分析推导一种差分格式的基本理论(稳定性、收敛性、误差估计等),衍生到其它的差分格式,深入浅出,启发思维,重点关注分析、比较异同,让学生参与讨论,调动学生学习积极性,尝试达到教与学融为一体的教学成效.针对教材理论性强的特点,结合各章节内容增设部分例题,利用计算机演示计算过程并绘制出解的图形,充分发挥计算机计算速度快、容量大,图形直观等优势,解决用笔难以实现的复杂的科学计算问题,使学生能够更加直观地理解教材内容.(2)计算问题复杂,学生难于掌握微分方程数值解的计算问题复杂,学生难于掌握,为此,我们加强了计算机编程能力的训练,通过算法分析和编制程序,由计算机来实现问题的求解.同时,用多种计算方法来解决同一种计算问题,研究其共性和差异,有利于学生灵活掌握各种算法.3.3创新与特点针对微分方程数值解这门课具有 理论性强,公式推导烦琐、枯燥,计算量大 的特点以及学生实际情况,灵活应用多种教学方法.(1)在教学中着重于对基本概念、基本理论和思想方法的讲解.对于引入的实际问题,采用启发式教学,通过师生讨论,建立模型,由老师展示计算机实现计算过程,最后分析结果,给学生以更多的自主思考空间,激发学生的学习欲望,鼓励学生参与分析讨论,着重掌握解决具体问题的实际方法和计算原理.(2)采用 课堂(多媒体+计算机动态演示+黑板) 实验室 的多元化教学模式.课堂教学以讲授为主,提问、讨论等多种方法进行,注重启迪学生思维,根据不同的教学内容,有意识地尝试不同的教学方式,将多种不同的教学形式进行优化组合,充分调动学生的主观能动性和思维的积极性,培养创新意识和创新能力以及自我更新知识的能力.(3)在课程教学中溶入数学建模思想,培养学生建模的能力.对理工科学生而言从实际问题中提取相关信息建立数学模型尤为重要,微分方程数值解教学将提高学生建模能力作为重要的教学目标,在教学过程中注重培养学生的数学意识和运用数学的能力.(4)在微分方程数值解课程中开设实验课,理论联系实际,让学生体会高维方程组的多种求解方法以及程序设计思想.并用适当的实际问题激发学生对数值方法理论内容的进一步实验论证和应用的兴趣,让学生感受到理论分析的重要性.使用计算机软件模拟克服和补充理论教学中的不足,更形象更生动的去验证和表现相关理论.特别是数值解的模拟中,让学生非常容易理解和体会理论解和数值解之间的联系和差别.(5)结合课程内容,将微分方程数值解的相关知识应用于毕业设计,培养学生综合应用所学知识解决问题的能力.多名学生应用所学微分方程数值解的相关知识在老师的指导下,完成!流体力学边界层理论中Falker Skan方程的数值解∀、!流体力学边界层理论中混合对流方程的数值解∀、!二阶常微分方程的一个差分格式∀、!利用MAT LAB软件实现偏微分方程的可视化###波动方程∀、!傅立叶变换的可视化及其应用∀等多篇毕业设计论文并取得较好的效果.4结束语通过微分方程数值解的教学研究与实践,我院的微分方程数值解在课程建设与教学改革研究方面已初步走向科学化与规范化的进程.该课程在教育思想转变、教学内容、课程体系、教学方法和教学手段的改革等方面取得了实际成效,基本实现教学过程立体化.但课程教学立体化尚需不断总结,不断完善;教学方法与教学手段改革尚需进一步深化.在今后的教学实践中根据我院本专业的定位,我们还需进一步调整和优化相应的课程内容体系;推进教学方法、教学手段现代化,提高教、学双方的效率;注重计算机实现科学计算与计算机绘图能力的培养,使我院微分方程数值解的教学更具特色.参考文献[1]张宏伟.注重培养研究能力的!微分方程数值解法∀课程教学研究与实践[J].大学数学.2006,22(6):4- 6.[2]教育部数学与统计学教学指导委员会数学类专业教学指导分委员会.信息与计算科学专业教学规范[J].大学数学.2003,19(6):6-11.[3]戴嘉尊.微分方程数值解法[M].南京:东南大学出版社,2002.[4]p utational M ethods in Engineering BoundaryV alue P roblems[M].A cademic Press,N ew Y ork,1979.125第13卷第1期杨韧,杨光崇,谢海英: 微分方程数值解 的教学研究与实践。
微分方程数值解法实验报告姓名： 班级： 学号：一：问题描述求解边值问题：()2(sin cos cos sin (0,1)(0,1)0,(,)x y u e x y x y G u x y G ππππππ+⎧⎫∆=+⎪⎪∈=⨯⎨⎬⎪⎪=∈∂⎩⎭（x,y) 其精确解为)sin()sin(),()(y x e y x u y x πππ+=问题一：取步长h=k=1/64,1/128,作五点差分格式，用Jacobi 迭代法，Gauss_Seidel 迭代法，SOR 迭代法（w=1.45）。
求解差分方程，以前后两次重合到小数点后四位的迭代值作为解的近似值，比较三种解法的迭代次数以及差分解)128/1,64/1)(,(=h y x u h 与精确解的精度。
问题二：取步长h=k=1/64,1/128,作五点差分格式，用单参数和双参数PR 法解差分方程，近似到小数点后四位。
与SOR 法比较精度和迭代步数。
问题三：取步长h=k=1/64,1/128,作五点差分格式，用共轭梯度法和预处理共轭梯度法解差分方程，近似到小数点后四位。
与SOR法与PR 法比较精度和迭代步数。
二．实验目的：分别使用五点差分法（Jacobi 迭代，Gauss_Seidel 迭代，SOR迭代），PR 交替隐式差分法（单参数，双参数），共轭梯度法，预共轭梯度法分别求椭圆方程的数值解。
三．实验原理：(1) Jacobi 迭代法设线性方程组(1)的系数矩阵A 可逆且主对角元素均不为零,令 并将A 分解成(2) 从而(1)可写成令其中. (3) 以为迭代矩阵的迭代法(公式)(4)称为雅可比(Jacobi)迭代法(公式),用向量的分量来表示,(4)为(5) 其中为初始向量. (2) Guass-Seidel 迭代法由雅可比迭代公式可知,在迭代的每一步计算过程中是用的全部分量来计算的所有分量,显然在计算第i 个分量时,已经b Ax =nn a ,...,a ,a 2211()nn a ,...,a ,a diag D 2211=()D D A A +-=()b x A D Dx +-=11f x B x +=b D f ,A D I B 1111--=-=1B ()()111f x B x k k +=+⎩⎨⎧[],...,,k ,n ,...,i x a b a x n ij j )k (j j i i ii )k (i 21021111==∑-=≠=+()()()()()Tn x ,...x ,x x 002010=()k x ()1+k x ()1+k i x计算出的最新分量没有被利用，从直观上看,最新计算出的分量可能比旧的分量要好些.因此，对这些最新计算出来的第次近似的分量加以利用,就得到所谓解方程组的高斯—塞德（Gauss-Seidel ）迭代法.把矩阵A 分解成(6)其中,分别为的主对角元除外的下三角和上三角部分,于是,方程组(1)便可以写成 即其中(7)以为迭代矩阵构成的迭代法(公式)(8)称为高斯—塞德尔迭代法(公式),用 量表示的形式为(3) SOR 迭代(4) 交替方向迭代法（PR 法）迭代格式为：()()1111+-+k i k x ,...,x 1+k ()1+k x ()1+k j x U L D A --=()nn a ,...,a ,a diag D 2211=U ,L --A ()b Ux x L D +=-22f x B x +=()()b L D f ,U L D B 1212---=-=2B ()()221f x B x k k +=+⎩⎨⎧[],...,,k ,n ,,i x a x a b a x i j n i j )k (j ij )k (j ij i ii )k (i 21021111111==∑∑--=-=+=++Λ))1(()(1D R L D T ωωω-+-=-b )(1--=L D d ωωhu πμωcos )11/(22opt =-+=2121,,1,1,1,,122L L L L u u u L u u u j i j i j i j i j i j i +==+-=+-+-+-对于单参数PR 法，对于多参数，(5) 共轭梯度法 算法步骤如下： ［预置步］任意，计算，并令取：指定算法终止常数，置，进入主步；［主步］ （１）如果，终止算法，输出；否则下行；（２）计算：（３）计算：（４）置，转入（１）．(6) 预共轭梯度法b uL I uL I b u L I uL I k k k k k k k k k k ττττττ+-=++-=++++211122211)()()()(hh optπτsin 22=2sin a ....2,1)11(421k 221h k a h k πρρτ==+-=--其中［预置步］任意，计算，并令取：指定算法终止常数，置，进入主步；［主步］（１）计算：，（２）如果，转入（3）．否则，终止算法，输出计算结果（３）计算：（４）置，转入（1）注：在算法［主步］中，引入变量，及，可以简化计算。
微分方程数值解实验报告微分方程数值解实验报告一、引言微分方程是数学中一类重要的方程，广泛应用于各个科学领域。
在实际问题中，往往难以得到微分方程的解析解，因此需要借助数值方法来求解。
本实验旨在通过数值解法，探索微分方程的数值解及其应用。
二、数值解法介绍常用的微分方程数值解法有欧拉法、改进欧拉法、四阶龙格-库塔法等。
在本实验中，我们将采用改进欧拉法进行数值解的求取。
改进欧拉法是一种一阶的显式迭代法，其基本思想是将微分方程的导数用差商来近似表示，并通过迭代逼近真实解。
具体迭代公式如下：\[y_{n+1} = y_n + h \cdot f(x_n + \frac{h}{2}, y_n + \frac{h}{2} \cdot f(x_n, y_n))\]其中，\(y_n\)表示第n步的近似解，\(h\)表示步长，\(f(x_n, y_n)\)表示微分方程的导数。
三、实验步骤1. 确定微分方程及初始条件在本实验中，我们选择经典的一阶常微分方程：\[y' = -2xy\]并给定初始条件\(y(0) = 1\)。
2. 设置步长和迭代次数为了得到较为准确的数值解，我们需要合理选择步长和迭代次数。
在本实验中，我们将步长设置为0.1，迭代次数为10。
3. 迭代计算数值解根据改进欧拉法的迭代公式，我们可以通过编写计算程序来求解微分方程的数值解。
具体计算过程如下：- 初始化：设定初始条件\(y_0 = 1\)，并给定步长\(h = 0.1\)。
- 迭代计算：使用改进欧拉法的迭代公式，通过循环计算得到\(y_1, y_2, ...,y_{10}\)。
- 输出结果：将计算得到的数值解输出，并进行可视化展示。
四、实验结果与分析通过以上步骤，我们得到了微分方程的数值解\(y_1, y_2, ..., y_{10}\)。
将这些数值解进行可视化展示，可以更直观地观察到解的变化趋势。
根据实验结果，我们可以发现随着迭代次数的增加，数值解逐渐逼近了真实解。
实验：微分方程数值解1、实验目的了解微分方程初值问题的数学模型，掌握使用MATLAB 求解常微分方程的数值方法，学会归纳不同参数和仿真之间的规律。
2、实验的基础MATLAB 求常微分方程初值问题00(,)()y f x y y x y '=⎧⎨=⎩数值方法是先创建函数文件,用以描述微分方程右端二元函数,然后用ode23()求出数值解。
常微分方程组初值问题一阶常微分方程组初值问题数值求解方法[T,y] = ode23(' F ',Tspan,y 0)其中, F 是函数文件,表示微分方程右端函数Tspan = [t 0Tfinal] ——求解区域;y 0——初始条件注: 函数F(t,y) 必须返回列向量.数值解y 的每一行对应于列向量T 中的每一行数据00(,)()dy f t y t t dt y t y ⎧=≥⎪⎨⎪=⎩例1马尔萨斯模型,以1994 年我国人口为12亿为初值，求解常微分方程N (t )表示人口数量,取人口变化率r =0.015,微分方程ode23('fun1',[1994,2020],12)[T,N]=ode23('fun1',[1994,2020],12)199019952000200520102015202012141618命令窗口→function z=fun1(t,N)z=0.015*N; 编辑窗口→0.015dN N dt =(1994)12N =3、微分方程数值解实验：捕食者与被捕食者问题一、实验内容假设有足够多的青草供野兔享用，而狐狸仅以野兔为食物。
x为野兔数量，y表示狐狸数量。
假定在没有狐狸的情况下，野兔增长率为100%.如果没有野兔，狐狸将被饿死，死亡率为100%。
狐狸和野兔相互作用的关系是，狐狸的存在使野兔受到威胁，且狐狸越多野兔增长受到阻碍越大，设野兔数量的负增长系数为0.015.而野兔的存在又为狐狸提供食物，设狐狸数量的增长与野兔的数量成正比，比例系数为0.001.数学模型为计算x (t )，y (t )当t ∈[0，20]时的数据。
微分方程数值解法实验报告实验题目：数值解微分方程的实验研究引言：微分方程是描述自然现象、科学问题和工程问题中变量之间的关系的重要数学工具。
然而，大部分微分方程很难找到解析解，因此需要使用数值方法来近似求解。
本实验旨在研究数值解微分方程的方法和工具，并通过实验验证其有效性和准确性。
实验步骤：1.了解微分方程的基本概念和求解方法，包括欧拉法、改进的欧拉法和龙格-库塔法。
2. 配置实验环境，准备实验所需的工具和软件，如Python编程语言和相关数值计算库。
3.选择一种微分方程进行研究和求解，可以是一阶、二阶或更高阶的微分方程。
4.使用欧拉法、改进的欧拉法和龙格-库塔法分别求解选定的微分方程，并比较其结果的准确性和稳定性。
5.计算数值解与解析解之间的误差，并进行误差分析和讨论。
6.对比不同数值解法的性能，包括计算时间和计算精度。
7.结果展示和总结，根据实验结果对数值解方法进行评价和选取。
实验结果：以一阶线性常微分方程为例，我们选择经典的“衰减振荡”问题进行实验研究。
该问题的微分方程形式为：dy/dt = -λy其中，λ为正实数。
我们首先使用Python编程语言实现了欧拉法、改进的欧拉法和龙格-库塔法。
进一步，我们选择了λ=0.5和初始条件y(0)=1，使用这三种数值解法求解该微分方程，并比较结果的准确性。
通过对比数值解和解析解可以发现，在短时间内，欧拉法、改进的欧拉法和龙格-库塔法的结果与解析解非常接近。
但随着时间的增加，欧拉法的结果开始偏离解析解，而改进的欧拉法和龙格-库塔法仍然能够提供准确的近似解。
这是因为欧拉法采用线性逼近的方式，误差随着步长的增加而累积，而改进的欧拉法和龙格-库塔法采用更高阶的逼近方式，可以减小误差。
为了更直观地比较不同方法的性能，我们还计算了它们的计算时间。
实验结果显示，欧拉法计算时间最短，而龙格-库塔法计算时间最长。
这表明在计算时间要求较高的情况下，可以选择欧拉法作为数值解方法。
微分方程数值解实验报告 姓名 丁建伟 学号 200708020211 日期2010.12.30 实验项目 五点差分格式 指导教师徐强 一、上机实验的问题和要求（需求分析）：实验内容：(I) 分别在步长h=1/16,1/32,1/64,1/128情形下用五点差分格式计算⎩⎨⎧Ω∂=Ω=-∆on u in f u ;0; 其中Ω=（0，1）×（0，1）； 2=f πsin πxsin πy ；精确解为 u=sin πxsin πy ；(I)给出系数矩阵A 的图像。
(II) 用G-S 迭代法求解AV=F ，并分析无穷范数和L2范数下的误差阶。
二、程序设计的基本思想，原理和算法描述：基本思想及原理：三、主要程序代码或命令：#include<stdio.h>#include<math.h>#define MAX 200#define M 3.1415926void main(){int n,i,j;float V[MAX][MAX]={{0}};float F[MAX][MAX],U[MAX][MAX],e[MAX][MAX];float h,x,y,err,temp;printf("请输入分割数n的值:");scanf("%d",&n);h=1/float(n+1);for(i=1,x=h;i<=n;i++,x+=h)for(j=1,y=h;j<=n;j++,y+=h){F[i][j]=2*M*M*sin(M*x)*sin(M*y);U[i][j]=sin(M*x)*sin(M*y);}err=1.0;while(err>1.0e-6){err=0.0;for(j=1;j<=n;j++)for(i=1;i<=n;i++){temp=0.25*(V[i+1][j]+V[i-1][j]+V[i][j-1]+V[i][j+1])+(h*h/4)*F[i][j]; err=(temp-V[i][j])*(temp-V[i][j])+err;V[i][j]=temp;}err=sqrt(err);}float e1=0,e2=0;for(i=1;i<=n;i++)for(j=1;j<=n;j++){e[i][j]=fabs(V[i][j]-U[i][j]);e2+=h*e[i][j]*e[i][j];if(e1<e[i][j])e1=e[i][j];}e2=sqrt(e2);printf("请输出在无穷范数和二范数下的误差:e1=%f和e2=%f\n",e1,e2);for(i=1;i<=n;i++)for(j=1;j<=n;j++)printf("请输出V[%d][%d]的值:%f\n",i,j,V[i][j]);}四、调试和运行程序过程中产生的问题及采取的措施：1.编译时出错，若想在主函数和被调用函数都使用一些变量，必须把这些变量设为全局变量。
实验4 常微分方程数值解化工系毕啸天2010011811 【实验目的】1. 练习数值微分的计算；2. 掌握用MATLAB 软件求微分方程初值问题数值解的方法；3. 通过实例学习用微分方程模型解决简化的实际问题；4. 了解欧拉方法和龙格-库塔方法的基本思想和计算公式，及稳定性等概念。
【实验内容】题目3小型火箭初始重量为1400kg，其中包括1080kg 燃料。
火箭竖直向上发射时燃料燃烧率为18kg/s，由此产生32000N 的推力，火箭引擎在燃料用尽时关闭。
设火箭上升时空气阻力正比于速度的平方，比例系数为0.4kg/m，求引擎关闭瞬间火箭的高度、速度、加速度，及火箭到达最高点时的高度和加速度，并画出高度、速度、加速度随时间变化的图形。
3.1 燃料燃烧过程物理模型分析设火箭质量为m，高度为h，速度为v，加速度为a，火箭推力为F，重力加速度为g，阻力为f。
1.由火箭上总共携带燃料1080kg，燃料燃烧率为18kg/s，可知火箭上升时间t=60s时，燃料全部烧尽。
2.由阻力正比于速度的平方，比例系数0.4kg/m，可知阻力表达式为f=0.4v2。
3.由于燃料燃烧，火箭的质量是时间的函数，易知m(t)=m0-18t4.5.根据牛顿第二运动定律，有。
代入数据有解出由以上5条分析，我们得到了一个常微分方程组：初值条件为：v0=0,h(0)=0,t60s.3.2 程序代码根据常微分方程组的初值问题，在MATLAB中计算数值解。
记x(1) = h，x(2)= v，x =(x(1), x(2))T首先编写M文件function dx = Rocket(t,x)dx=[x(2);(32000-0.4*x(2)^2)/(1400-18*t)-9.8];%以向量形式表示微分方程endts=0:60 %终点时间为60s，步长定义为1即可x0=[0,0];[t,x]=ode45(@Rocket,ts,x0);[t,x]plot(t,x(:,1)),grid,title('图1.高度-时间')xlabel('t/s')ylabel('h/m')pauseplot(t,x(:,2)),grid,title('图2.速度-时间')xlabel('t/s')ylabel('v/(m/s)')pausea=(32000-0.4*x(:,2).^2)./(1400-18*t)-9.8plot(t,a),grid,title('图3.加速度-时间')xlabel('t/s')ylabel('a/(m/s2)')[t,x(:,1),x(:,2),a]由MATLAB计算得到从开始上升到关闭引擎瞬间的情况如下表：由数据可以看出，引擎关闭瞬间，火箭的高度为h=12190m，速度v=267.26m/s，加速度a=0.91701m/s2。
微分方程数值解实验报告专业 信息与计算科学 班级 0703 姓名 老武 学号 2007060*** 协作队员 实验日期2010 年 10 月成绩评定 教师签名 批改日期 题目 热传导方程的有限差分数值模拟一、 问题提出一根长为L 的均匀导热细杆，侧面绝热，内部无热源。
其热传导系数为k ，比热为c ，线密度为ρ。
求细杆内温度变化的规律。
二、 模型建立设杆长方向为x 轴，考虑杆上从x 到x+△x 的一段。
其质量为△m=ρ△x ，热容量为c △m 。
设杆中的热流沿x 轴正向，热流强度为q(x,t)，热量为Q(x,t),温度分布为 u(x,t)。
△x 内细杆吸收热量的来源只有热传导(无热源)。
由热传导的Fourier 定律,有),(),(t x ku t x q x -= (1)由能量守恒定律，在△t 内细杆[x,x+ △x]上的能量有Q u m c ∆=∆∆即t t x x q t x q u x c ∆∆+-=∆∆)],(),([ρ于是有),(),(t x q t x u c x t -=ρ (2)结合(1)和(2)得xx t u a u 2= (3)其中 ρc k a /2=三、 求解方法使用古典显格式：)2(111n m n m n m n m n m U U U U U -+++-+=τ其中 22/h k a =τ (k 和h 分别为时间与空间方向的步长，取k=0.005，h=0.1使得2/1/2≤h k )取 12=a L=1，细杆各处的初始温度为 )sin(x π，两端截面上的温度为0。
Matlab 程序如下：%古典显格式clc;x0=0:0.1:1t=0:0.005:0.1;n=length(t)Un=sin(pi*x0)for i=1:nun=[];u=[];for r=1:11u1=exp(-pi^2*t(i))*sin(pi*x0(r));u=[u u1];endufor j=2:10Un1=Un(j)+0.5*(Un(j+1)-2*Un(j)+Un(j-1));un=[un Un1];endun=[0 un 0]e=abs(u-Un)Un=un;End四、 输出结果表一：数值模拟结果 n m U n t m x 0 0.005 0.010 0.015 0.020 … 0.995 (1e-004) 1.000(1e-004) 0.0000 0.1000 0.2000 0.3000 0.4000 0.5000 0.6000 0.7000 0.8000 0.9000 1.0000 0.0000 0.3090 0.5878 0.8090 0.9511 1.0000 0.9511 0.8090 0.5878 0.3090 0.0000 0.0000 0.2939 0.5590 0.7694 0.9045 0.9511 0.9045 0.7694 0.5590 0.2939 0.0000 0.0000 0.2795 0.5317 0.7318 0.8602 0.9045 0.8602 0.7318 0.5317 0.2795 0.0000 0.0000 0.2658 0.5056 0.6959 0.8181 0.8602 0.8181 0.6959 0.5056 0.2658 0.0000 0.0000 0.2528 0.4809 0.6619 0.7781 0.8181 0.7781 0.6619 0.4809 0.2528 0.0000 … … … … … … … … … … … 0.0000 0.1422 0.2706 0.3724 0.4378 0.4603 0.4378 0.3724 0.2706 0.1422 0.00000.0000 0.1353 0.2573 0.3542 0.4164 0.4378 0.4164 0.3542 0.2573 0.1353 0.0000表二：真实值 n m U _n t m x 00.005 0.010 0.015 0.020 … 0.995 (1e-004) 1.000 (1e-004) 0.0000 0.1000 0.2000 0.3000 0.4000 0.5000 0.6000 0.7000 0.8000 0.9000 1.0000 0.0000 0.3090 0.5878 0.8090 0.9511 1.0000 0.9511 0.8090 0.5878 0.3090 0.0000 0.0000 0.2941 0.5595 0.7701 0.9053 0.9518 0.9053 0.7701 0.5595 0.2941 0.0000 0.0000 0.2800 0.5325 0.7330 0.8617 0.9060 0.8617 0.7330 0.5325 0.2800 0.0000 0.0000 0.2665 0.5069 0.6977 0.8202 0.8624 0.8202 0.6977 0.5069 0.2665 0.0000 0.0000 0.2537 0.4825 0.6641 0.7807 0.8209 0.7807 0.6641 0.4825 0.2537 0.0000 … … … … … … … … … … … 0.0000 0.1679 0.3194 0.4396 0.5168 0.5434 0.5168 0.4396 0.3194 0.1679 0.00000.0000 0.1598 0.3040 0.4184 0.4919 0.5172 0.4919 0.4184 0.3040 0.1598 0.0000误差1e = 0 0 0 0 0 0 0 0 0 02e ＝1.0e-003 *0 0.2451 0.4663 0.6418 0.7545 0.7933 0.7545 0.6418 0.4663 0.2451 0.00003e ＝ 0 0.0005 0.0009 0.0012 0.0014 0.00150.0014 0.0012 0.0009 0.0005 0.0000e＝0 0.0007 0.0013 0.0017 0.0020 0.0022 40.0020 0.0017 0.0013 0.0007 0.0000e＝0 0.0008 0.0016 0.0022 0.0026 0.0027 50.0026 0.0022 0.0016 0.0008 0.0000……….e＝1.0e-005 *2000 0.2567 0.4883 0.6721 0.7901 0.83080.7901 0.6721 0.4883 0.2567 0.0000e＝1.0e-005 *2010 0.2455 0.4670 0.6427 0.7555 0.79440.7555 0.6427 0.4670 0.2455 0.0000五、结果分析在同一个时间下，细杆内的温度分布为：细杆内中间的温度最高，往两边逐渐下降到0，并且温度值关于x=0.5这条直线对称。
微分方程数值解实验指导（四）实验四: 二阶双曲型方程差分格式的程序实现1. 实验题目 就05.0,31=∆=∆t x π分别对15,2==a a ，用显式格式求如下周期自由振动问题的数值解：⎪⎪⎪⎩⎪⎪⎪⎨⎧==∂∂=>∞<<∞-=∂∂-∂∂).,2(),0(,sin )0,(,0)0,(,0,,022222t u t u x a x t u x u t x x u a t u π其中a 为常数。
观察数值解的收敛情况,分析解释计算结果。
2. 实验要求按照给定的差分格式编程实现求出数值解；比较数值解与精确解的误差；结合格式的相容性、稳定性和收敛性条件简单分析计算结果。
要求在实验课上算出数值结果；按要求格式写出实验报告；下次实验课前交本次实验的实验报告。
3. 实验原理与实验过程第一步: 对求解区域作网格剖分。
由于所求问题是自由振动的周期边值问题，我们可以在一个空间变量周期内求出该问题的数值解。
记求解区域G={(x,t)|0<x<2pi,0<t<T},分割G 为一些均匀的矩形，其边与坐标轴平行，取空间步长为J h x /2π==∆,时间步长为M T t /==∆τ，网格比为h a r /τ=，其中J ，M 为正整数，用两族平行直线J j jh x x j ...,,2,1,===平行于t 轴作竖直线, M n n t t n ...,,2,1,===τ平行于X 轴作水平线, 将矩形区域}0,20|),{(T t x t x G ≤≤≤≤=π分割成矩形网格。
第二步: 差分法的目的是：求周期自由振动问题的解),(t x u 在节点),(n j t x 的近似值nju （M n J j ...,,2,1,0,,,2,1,0== ）。
为此，需要构造逼近微分方程定解问题的差分格式。
采用显格式，J j h u u u a u u u n j n j n j n jnj n j ...,,2,1,0,222112211=+-=+--+-+τ（1） 初始条件用如下差分格式逼近)()(11100j jj j j x u u x u ϕτϕ=-=-（2）在（1）中取0=n ，与（2）联立解出)()()1()]()([2102101021j j j j jx x r x x r u τϕϕϕϕ+-++=+- 此处h a r /τ=表示网比。
微分方程数值解实验指导书李光云数学与计算科学学院二O一O 十月一、概述本课程实验指导书是根据陆金甫，关治编著的《偏微分方程数值解法(第二版)》编写的。
通过上机实验，可帮助学生理解偏微分方程数值计算方法的基本思想，巩固学生所学过的算法，同时让学生进一步学习和掌握Matlab软件在数值计算方面中的应用。
二、实验环境本书选择的实验环境是计算机以及软件Matlab(版本6.5以上)一套。
三、实验课时安排共7次实验，16课时，试验五为综合设计性试验4课时，其它实验2个课时。
四、实验要求上机完成实验指导书中所规定的内容，自行按实验指导书要求完成程序设计和调试，并提交每次实验的实验报告，附带算法程序清单和算法输出结果。
五、实验考核要求上机完成试验内容，并提交一份算法程序清单和数值结果。
实验一 双曲型方程的迎风格式一、实验目的1、掌握双曲型方程的迎风格式的算法思想,格式稳定的条件；2、培养编程与上机调试能力。
二、实验课时：2个课时 三、实验准备1、了解偏微分方程的分类，以及双曲型方程的特点；2、熟悉求解偏微分方程有限差分法的基本原理，双曲型方程迎风格式的推导和特点，熟悉用迎风格式求解双曲型方程的算法步骤；3、熟悉Matlab 软件的基本命令及其数值计算中的基本命令和函数。
4、了解双曲型方程迎风格式稳定的条件。
四、实验内容用双曲型方程的迎风格式11110, 00, <0n n n nj jj j n n n njjj ju u u u a a h uuu u aa hττ+-++--+=>--+=其中，,h τ分别为时间步长和空间步长。
求解初值问题[]()()00,,0,,,0,u ux t T t xu x u x ∂∂⎧+=∈∈⎪∂∂⎨⎪=⎩其中()01,0,0,0.x u x x ≤⎧=⎨>⎩取10.01,2h λ==.计算[]0,1区间的u 值至0.5n t =,并画出解的图形。
五、基本思想及主要步骤考虑对流方程0,,0u u a x t t x∂∂+=∈>∂∂ 其a 为给定常数。
其迎风格式的基本思想是在方程中关于空间的偏导数用在特征线方向一侧的单边差商代替，得到对流方程的迎风格式11110, 00, <0n n n nj jj j n nn njjj ju u u u a a h uuu u aa hττ+-++--+=>--+=用Fourier 方法讨论格式的稳定性可知，以上两个格式都是条件稳定的，而且都在1,a hτλλ≤=时稳定。
程序设计步骤如下：(1)首先由a 的正负确定用哪个格式，再把格式改写为显示格式； (2)利用空间步长定义一个零向量用于存储u 值； (3)将初值离散存于u 向量内；(4)利用差分格式沿着t 坐标逐层计算u 值，每计算一次覆盖存储于u 向量内，直到算到0.5n t =层，然后输出结果；(5)利用得到的u 向量画出在0.5n t =时u 的函数图像。
六、实验提示：此算法可以编写M 程序，使得该程序可以实现对不同的a ，时间步长和空间步长求解，由于不用符号的a 用到的格式不同，迎风格式在条件1,a hτλλ≤=时稳定，故可以在程序里加一个判断部分，当不稳定时及时报警跳出。
七、实验结果记录和要求1、实验结果输出，画出解的图形；2、把得到的解和初值问题的解析解比较，观察比较结果；3、递交实验报告。
八、思考：若考虑对流方程的差分格式11110, <00, >0n n n nj jj j n n n njjj ju u u u a a h uuu u aa hττ+-++--+=--+=是否也能得到差分格式稳定？实验二 常系数扩散方程的经典差分格式一、实验目的1、了解抛物型方程的经典差分格式,格式稳定的条件；2、掌握常系数扩散方程初边值问题的加权隐式格式的求解方法；3、培养编程与上机调试能力。
二、实验课时：2个课时 三、实验准备1、了解偏微分方程的分类方法，抛物型方程的特点；2、了解常见的抛物型方程的差分格式；3、熟悉抛物型方程的加权隐式格式及其推导过程；4、熟悉常系数抛物型方程初边值问题的初值和边界离散方法。
四、实验内容考虑常系数扩散方程的初边值问题()()()22,01,0,,0sin ,01,0,1,0,0.u ux t t x u x x x u t u t t π⎧∂∂=<<>⎪∂∂⎪⎪=≤≤⎨⎪==≥⎪⎪⎩其解析解为()2,sin ,01,0.t u x t e x x t ππ-=≤≤≥用加权隐式格式近似求解()11111111222210n n n n n n n n j jj j j j j j u u u u u u u u h h θθτ----+-+-⎡⎤--+-+-++=⎢⎥⎢⎥⎣⎦其中01θ≤≤，取()110,,0,1,...,,10j J h x jh j J τ====为时间步长，2h τλ=为网格比，对不同的时间步长(11,,1,2,4,842λ=)，计算当10,1,2θ=时初边值问题的解()0.4,0.4u ，并且与精确解比较，分析比较结果。
五、基本思想及主要步骤用有限差分法解常系数扩散方程22u ua t x∂∂=∂∂有加权隐式差分格式()11111111222210n n n n n n n n j jj j j j j j u u u u u u u u a h h θθτ----+-+-⎡⎤--+-+-++=⎢⎥⎢⎥⎣⎦其中01θ≤≤，当12θ=时为Crank-Nicolson 格式，当1θ=时为向后差分格式，当0θ=时为向前差分格式。
加权隐式格式稳定的条件是12,121, 1.2a λθθθ≤≤-≤≤1，当0<2无限制当加权隐式格式是两层隐式格式，用第n 层计算第n+1层节点值的时候，要解线性方程组。
实验步骤如下：(1)输入.θλ，确定加权隐式格式的参数；(2)定义向量v ，把初边值条件离散，得到0,0,1,...,j u j J =的值存入向量v ； (3)利用差分格式由第n 层计算第n+1层，建立相应线性方程组，求解并且存入向量v ；(4)计算到0.4t =，输出()0.4,0.4u 。
六、实验提示：求解程序可编写M 文件，以便改变.θλ，得到不用参数下()0.4,0.4u 的值，然后比较与真实值的误差，判断格式的稳定性。
七、实验结果记录和要求1、输出不同参数下求得的节点的近似值；2、比较近似值和真实值的误差，印证格式稳定条件；3、实验结束后，递交实验报告。
八、思考：抛物型方程和双曲型方程的初边值问题的边界处理有什么不同，用隐式格式计算与显示格式相比，有什么优劣。
实验三 椭圆型方程的五点格式一、实验目的1、掌握poisson 方程定解问题的五点格式的求解方法；2、培养编程与上机调试能力。
二、实验课时：2个课时 三、实验准备1、了解偏微分方程的分类方法，椭圆型方程的特点；2、了解常见的poisson 方程的差分格式；3、熟悉poisson 方程的五点格式及其推导过程。
四、实验内容给定如下Laplace 方程（poisson 方程的特殊情况）的定解问题{}()()()()0,017,010,0,17,,00,,10100.u x y u y u y u x u x ∆=Ω=<<<<====利用椭圆型方程的五点格式，计算该问题的近似解，并且画出近似解的图形。
五、基本思想及主要步骤对Laplace 方程的第一边值问题()()()()22220,,,,,,u uu x y Dx y u x y x y x y Dα∂∂∆=+=∈∂∂=∈∂ 利用taylor 展开可得逼近它的五点差分格式的差分逼近()()1, 1.,1.122220,,,,i j ij i ji j ij i j h i j hij ij i j hu u u u u u u x y D h k u x y D α+-+--+-+=+=∈=∈∂其中,h k 分别为x 轴和y 轴步长，边界条件可以由()(),,u x y x y α=离散可得，当(,)x y D ∈∂时有(),ij i j x y αα=。
注意五点格式计算节点是由边界的已知节点，计算内部节点，计算时需要联立大型方程组，该方程组可以用迭代法求解。
主要步骤：(1)首先取定,h k ，对求解区域划分网格，按照网格定义矩阵v1，使得矩阵里面每个元素对应求解区域中的每个节点；(2)由边界条件定义v1矩阵中边界元素的值，其余元素定义为零；(3)定义与v1同型的零矩阵v2；(4)用五点格式公式通过矩阵v1迭代计算矩阵v2，迭代精度为0.1；(5)画图。
六、实验提示：由v1迭代计算矩阵v2时，v2元素通过v1元素利用五点格式计算，当v1v2对应元素误差在控制范围之内时，终止迭代，否则继续，直到满足误差条件。
七、实验结果记录和要求1、输出所有节点的计算结果，并且画图；2、实验结束后，递交实验报告。
八、思考：如果求解区域是一般区域，如何划分网格，处理边界条件。
实验四 pde 工具箱在热传导，静电学等方面的应用一、实验目的1、全面了解pde 工具箱对其他偏微分方程的求解方法；2、培养编程与上机调试能力。
二、实验课时：2个课时 三、实验准备1、了解静电学,热传导等方面的偏微分方程；2、熟悉pde 工具箱的参数的具体意义；3、熟悉Matlab 软件中pde 工具箱在求解其他偏微分方程定解问题的使用方法。
四、实验内容(1)一块有矩形裂纹的金属块,左侧被加热到100度，右侧以恒定速率降低到周围空气温度，所有其他边界独立，则有如下定解问题：0,,100,(Dirichlet ),10,(Neumann ),0,(Neumann )udu u tu un un∂-∆=∈Ω∂=∂=-∂∂=∂右侧条件左侧条件其它边界条件 设起始时刻0t 时刻金属块的温度为零度。
金属块长度为1.6，宽度为1，裂纹长度为0.8，宽度为0.1，指定起始时间为0，希望研究开始5s 内的热散发问题。
(2)一个内部充满空气的方框，内边长为0.2m ，外边长为0.5m 。
内边界上电势为1000v ，外边界上电势为0v ，域内没有电荷，现在求方框中的电势。
此问题即求方程0v ∇=边界条件为Dirichlet 条件，内边界上1000v =，外边界上0v =。
五、主要方法及提示这些问题都可归为求解偏微分方程的定解问题，我们用matlab 中的pde 工具箱求解此类问题的近似解。
对问题1，首先打开图形用户界面，画CSG 模型：首先画矩形（R1）四角坐标为[]0.5,0.5,0.5,0.5,x =--[]y 0.8,0.8,0.8,0.8=--，后画矩形（R2）四角坐标为[]=--，选择Options选项，打y0.4,0.4,0.4,0.4x=--[]0.05,0.05,0.05,0.05,开Grid Spacing对话框，在-0.05和0.05处输入x轴的附加短线，以帮助画出代表裂纹的矩形，然后显示网格和“Snap-to-gird”风格，金属块的CSG模型现在可以用公式R1-R2来表示。
