下载文档原格式
一维抛物型偏微分方程初边值问题求解摘要：一、引言二、一维抛物型偏微分方程1.定义与性质2.初边值问题三、求解方法1.紧差分格式2.追赶法3.有限元算法四、Matlab程序实现1.紧差分格式程序2.追赶法程序五、结论与展望正文：一、引言在数学、物理等领域，偏微分方程是一类重要的方程。
其中，一维抛物型偏微分方程在科学研究和实际应用中具有广泛的意义。
本文将探讨一维抛物型偏微分方程的初边值问题的求解方法，并介绍相应的Matlab程序实现。
二、一维抛物型偏微分方程1.定义与性质一维抛物型偏微分方程是指具有如下形式的方程：u_t = a * u_xx其中，u(x, t) 表示未知函数，t 表示时间，x 表示空间坐标，a 为常数。
2.初边值问题初边值问题是指在给定的初始条件和边界条件下求解偏微分方程的问题。
在一维抛物型偏微分方程中，初边值问题可以表示为：u(x, 0) = u_0(x)u(x, t) = u_t(x, t) 在边界x=0，x=L上三、求解方法1.紧差分格式紧差分格式是一种求解偏微分方程的方法，其精度为O(h^(1/2) * Δt)，无条件稳定。
在这种方法中，我们首先需要建立离散的网格系统，然后通过数值积分求解离散化的偏微分方程。
2.追赶法追赶法是一种求解线性方程组的方法，也可以用于求解初边值问题。
在这种方法中，我们首先需要将偏微分方程转化为线性方程组，然后使用追赶法求解线性方程组。
3.有限元算法有限元算法是一种基于变分原理的求解方法，可以将偏微分方程问题转化为求解有限元空间的线性方程组。
这种方法在求解一维抛物型偏微分方程时具有较高的精度和可靠性。
引言偏微分方程定解问题有着广泛的应用背景。
人们用偏微分方程来描述、解释或者预见各种自然现象，并用于科学和工程技术的各个领域fll。
然而，对于广大应用工作者来说，从偏微分方程模型出发，使用有限元法或有限差分法求解都要耗费很大的工作量，才能得到数值解。
现在，MATLAB PDEToolbox已实现对于空间二维问题高速、准确的求解过程。
偏微分方程如果一个微分方程中出现的未知函数只含一个自变量，这个方程叫做常微分方程，也简称微分方程；如果一个微分方程中出现多元函数的偏导数，或者说如果未知函数和几个变量有关，而且方程中出现未知函数对几个变量的导数，那么这种微分方程就是偏微分方程。
常用的方法有变分法和有限差分法。
变分法是把定解问题转化成变分问题，再求变分问题的近似解；有限差分法是把定解问题转化成代数方程，然后用计算机进行计算；还有一种更有意义的模拟法，它用另一个物理的问题实验研究来代替所研究某个物理问题的定解。
虽然物理现象本质不同，但是抽象地表示在数学上是同一个定解问题，如研究某个不规则形状的物体里的稳定温度分布问题，由于求解比较困难，可作相应的静电场或稳恒电流场实验研究，测定场中各处的电势，从而也解决了所研究的稳定温度场中的温度分布问题。
随着物理科学所研究的现象在广度和深度两方面的扩展，偏微分方程的应用范围更广泛。
从数学自身的角度看，偏微分方程的求解促使数学在函数论、变分法、级数展开、常微分方程、代数、微分几何等各方面进行发展。
从这个角度说，偏微分方程变成了数学的中心。
一、MATLAB方法简介及应用1.1 MATLAB简介MATLAB是美国MathWorks公司出品的商业数学软件，用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言和交互式环境，主要包括MATLAB和Simulink两大部分。
1.2 Matlab主要功能数值分析数值和符号计算工程与科学绘图控制系统的设计与仿真数字图像处理数字信号处理通讯系统设计与仿真财务与金融工程1.3 优势特点1) 高效的数值计算及符号计算功能，能使用户从繁杂的数学运算分析中解脱出来；2) 具有完备的图形处理功能，实现计算结果和编程的可视化；3) 友好的用户界面及接近数学表达式的自然化语言，使学者易于学习和掌握；4) 功能丰富的应用工具箱(如信号处理工具箱、通信工具箱等) ，为用户提供了大量方便实用的处理工具。
引言偏微分方程定解问题有着广泛的应用背景。
人们用偏微分方程来描述、解释或者预见各种自然现象，并用于科学和工程技术的各个领域fll。
然而，对于广大应用工作者来说，从偏微分方程模型出发，使用有限元法或有限差分法求解都要耗费很大的工作量，才能得到数值解。
现在，MATLAB PDEToolbox已实现对于空间二维问题高速、准确的求解过程。
偏微分方程如果一个微分方程中出现的未知函数只含一个自变量，这个方程叫做常微分方程，也简称微分方程；如果一个微分方程中出现多元函数的偏导数，或者说如果未知函数和几个变量有关，而且方程中出现未知函数对几个变量的导数，那么这种微分方程就是偏微分方程。
常用的方法有变分法和有限差分法。
变分法是把定解问题转化成变分问题，再求变分问题的近似解；有限差分法是把定解问题转化成代数方程，然后用计算机进行计算；还有一种更有意义的模拟法，它用另一个物理的问题实验研究来代替所研究某个物理问题的定解。
虽然物理现象本质不同，但是抽象地表示在数学上是同一个定解问题，如研究某个不规则形状的物体里的稳定温度分布问题，由于求解比较困难，可作相应的静电场或稳恒电流场实验研究，测定场中各处的电势，从而也解决了所研究的稳定温度场中的温度分布问题。
随着物理科学所研究的现象在广度和深度两方面的扩展，偏微分方程的应用范围更广泛。
从数学自身的角度看，偏微分方程的求解促使数学在函数论、变分法、级数展开、常微分方程、代数、微分几何等各方面进行发展。
从这个角度说，偏微分方程变成了数学的中心。
一、MATLAB方法简介及应用1.1 MATLAB简介MATLAB是美国MathWorks公司出品的商业数学软件，用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言和交互式环境，主要包括MATLAB和Simulink两大部分。
1.2 Matlab主要功能数值分析数值和符号计算工程与科学绘图控制系统的设计与仿真数字图像处理数字信号处理通讯系统设计与仿真财务与金融工程1.3 优势特点1) 高效的数值计算及符号计算功能，能使用户从繁杂的数学运算分析中解脱出来；2) 具有完备的图形处理功能，实现计算结果和编程的可视化；3) 友好的用户界面及接近数学表达式的自然化语言，使学者易于学习和掌握；4) 功能丰富的应用工具箱(如信号处理工具箱、通信工具箱等) ，为用户提供了大量方便实用的处理工具。
偏微分方程数值解法上机报告（一）一、实验题目：用Ritz-Galerkin 方法求解边值问题2u '',01(0)0,(1)1u x x u u ⎧-+=<<⎨==⎩的第n 次近似()n u x ，基函数()sin(),1,2,...,i x i x i n ϕπ==.二、实验目的：通过本次上机实验，理解求解初值问题的变分问题的最重要的近似解法——Ritz-Galerkin 方法，以便为学习有限元法打好基础。
此外，要熟悉用Matlab 解决数学问题的基本编程方法，提高运用计算机解决问题的能力。
三、实验代码：n=5;syms x;for i=1:np(i)=sin(i*pi*x);q(i)=-i^2*pi^2*sin(i*pi*x);endfor i=1:nb(i)=2*int(p(i),0,1);for j=1:nA(i,j)=int((-q(j)+p(j))*p(i),0,1);endendt=inv(A)*b'四、运行结果：t=2251799813685248/3059521645650671/pi281474976710656/9481460623939047/pi281474976710656/43582901062631895/pi五、总结：通过本次上机，我了解了Ritz-Galerkin 方程 n j j p f cj p i p a n i i ,...,2,1)),(,())(),((1==∑=，明白了用Ritz-Galerkin 方法解决边值问题的变分问题的基本原理，并接近一步提高自己的编程动手能力，受益匪浅。
偏微分方程数值解法上机报告（二）一、 实验题目：用线性元求下列边值问题的数值解2''2sin ,0142(0)0,'(1)0y y x x y y ππ⎧-+=<<⎪⎨⎪==⎩二、 实验目的：通过本次上机，熟悉和掌握用Galerkin 法观点出发导出的求解处置问题数值解的线性有限元法。
MATLAB中的偏微分方程数值解法偏微分方程（Partial Differential Equations，PDEs）是数学中的重要概念，广泛应用于物理学、工程学、经济学等领域。
解决偏微分方程的精确解往往非常困难，因此数值方法成为求解这类问题的有效途径。
而在MATLAB中，有丰富的数值解法可供选择。
本文将介绍MATLAB中几种常见的偏微分方程数值解法，并通过具体案例加深对其应用的理解。
一、有限差分法（Finite Difference Method）有限差分法是最为经典和常用的偏微分方程数值解法之一。
它将偏微分方程的导数转化为差分方程，通过离散化空间和时间上的变量，将连续问题转化为离散问题。
在MATLAB中，使用有限差分法可以比较容易地实现对偏微分方程的数值求解。
例如，考虑一维热传导方程（Heat Equation）：∂u/∂t = k * ∂²u/∂x²其中，u为温度分布随时间和空间的变化，k为热传导系数。
假设初始条件为一段长度为L的棒子上的温度分布，边界条件可以是固定温度、热交换等。
有限差分法可以将空间离散化为N个节点，时间离散化为M个时刻。
我们可以使用中心差分近似来计算二阶空间导数，从而得到以下差分方程：u(i,j+1) = u(i,j) + Δt * (k * (u(i+1,j) - 2 * u(i,j) + u(i-1,j))/Δx²)其中，i表示空间节点，j表示时间步。
Δt和Δx分别为时间和空间步长。
通过逐步迭代更新节点的温度值，我们可以得到整个时间范围内的温度分布。
而MATLAB提供的矩阵计算功能，可以大大简化有限差分法的实现过程。
二、有限元法（Finite Element Method）有限元法是另一种常用的偏微分方程数值解法，特点是适用于复杂的几何形状和边界条件。
它将求解区域离散化为多个小单元，通过构建并求解代数方程组来逼近连续问题。
在MATLAB中，我们可以使用Partial Differential Equation Toolbox提供的函数进行有限元法求解。
偏微分方程的MATLAB求解精讲©MA TLAB求解微分/偏微分方程，一直是一个头大的问题，两个字，“难过”，由于MA TLAB对LaTeX的支持有限，所有方程必须化成MA TLAB可接受的标准形式，不支持像其他三个数学软件那样直接傻瓜式输入，这个真把人给累坏了！不抱怨了，还是言归正传，回归我们今天的主体吧！MA TLAB提供了两种方法解决PDE问题，一是pdepe()函数，它可以求解一般的PDEs，据用较大的通用性，但只支持命令行形式调用。
二是PDE工具箱，可以求解特殊PDE问题，PDEtool有较大的局限性，比如只能求解二阶PDE问题，并且不能解决偏微分方程组，但是它提供了GUI界面，从繁杂的编程中解脱出来了，同时还可以通过File->Save As直接生成M代码一、一般偏微分方程组(PDEs)的MA TLAB求解 (3)1、pdepe函数说明 (3)2、实例讲解 (4)二、PDEtool求解特殊PDE问题 (6)1、典型偏微分方程的描述 (6)（1）椭圆型 (6)（2）抛物线型 (6)（3）双曲线型 (6)（4）特征值型 (7)2、偏微分方程边界条件的描述 (8)（1）Dirichlet条件 (8)（2）Neumann条件 (8)3、求解实例 (9)一、一般偏微分方程组(PDEs)的MATLAB 求解1、pdepe 函数说明MA TLAB 语言提供了pdepe()函数，可以直接求解一般偏微分方程(组)，它的调用格式为sol=pdepe(m,@pdefun,@pdeic,@pdebc,x,t)【输入参数】@pdefun ：是PDE 的问题描述函数，它必须换成下面的标准形式(,,)[(,,,)](,,,)()m m u u u uc x t x x f x t u s x t u x t x x x−∂∂∂∂∂=+∂∂∂∂∂式1 这样，PDE 就可以编写下面的入口函数 [c,f,s]=pdefun(x,t,u,du)m,x,t 就是对应于(式1)中相关参数，du 是u 的一阶导数，由给定的输入变量即可表示出出c,f,s 这三个函数@pdebc ：是PDE 的边界条件描述函数，必须先化为下面的形式(,,)(,,).*(,,,)0up x t u q x t u f x t u x∂+=∂ 于是边值条件可以编写下面函数描述为 [pa,qa,pb,qb]=pdebc(x,t,u,du)其中a 表示下边界，b 表示下边界@pdeic ：是PDE 的初值条件，必须化为下面的形式00(,)u x t u =我们使用下面的简单的函数来描述为 u0=pdeic(x)m,x,t ：就是对应于(式1)中相关参数【输出参数】sol ：是一个三维数组，sol(:,:,i)表示u i 的解，换句话说u k 对应x(i)和t(j)时的解为sol(i,j,k)通过sol ，我们可以使用pdeval()直接计算某个点的函数值2、实例讲解试求解下面的偏微分2111222221220.024()0.17()u u F u u t xu u F u u tx ∂∂=−− ∂∂ ∂∂ =−− ∂∂ 其中， 5.7311.46()x x F x e e −=−，且满足初始条件12(,0)1,(,0)0u x u x ==及边界条件1221(0,)0,(0,)0,(1,)1,(1,)0u ut u t u t t x x∂∂====∂∂【解】（1）对照给出的偏微分方程，根据标注形式，则原方程可以改写为111222120.024()1.*1()0.17u u F u u x u u F u u t t x ∂−−∂∂∂=+ ∂−∂∂∂可见m=0，且1122120.024()1,,1()0.17u F u u x c f s u F u u x ∂−− ∂===∂−∂%% 目标PDE 函数function [c,f,s]=pdefun (x,t,u,du) c=[1;1];f=[0.024*du(1);0.17*du(2)]; temp=u(1)-u(2);s=[-1;1].*(exp(5.73*temp)-exp(-11.46*temp));（2）边界条件改写为12011010.*.*00000u f f u −+=+=下边界上边界%% 边界条件函数function [pa,qa,pb,qb]=pdebc(xa,ua,xb,ub,t) %a 表示下边界，b 表示上边界 pa=[0;ua(2)];qa=[1;0]; pb=[ub(1)-1;0]; qb=[0;1];（3）初值条件改写为1210u u =%% 初值条件函数function u0=pdeic(x) u0=[1;0];（4）最后编写主调函数 clcx=0:0.05:1; t=0:0.05:2; m=0;sol=pdepe(m,@pdefun,@pdeic,@pdebc,x,t);figure('numbertitle','off','name','PDE Demo ——by Matlabsky') subplot(211)surf(x,t,sol(:,:,1)) title('The Solution of u_1') xlabel('X') ylabel('T') zlabel('U') subplot(212)surf(x,t,sol(:,:,2)) title('The Solution of u_2') xlabel('X') ylabel('T') zlabel('U')二、PDEtool 求解特殊PDE 问题MATLAB 的偏微分工具箱(PDE toolbox)可以比较规范的求解各种常见的二阶偏微分方程，但是惋惜的是只能求解特殊二阶的PDE 问题，并且不支持偏微分方程组！PDE toolbox 支持命令行形式求解PDE 问题，但是要记住那些命令以及调用形式真的很累人，还好MATLAB 提供了GUI 可视交互界面pdetool ，在pdetool 中可以很方便的求解一个PDE 问题，并且可以帮我们直接生成M 代码(File->Save As)。
文章标题：深入探讨 Matlab 中求解偏微分方程的方法和应用一、引言在现代科学和工程中，偏微分方程是一种重要的数学工具，用于描述各种自然现象和物理过程，如热传导、流体力学、电磁场等。
Matlab 是一个用于科学计算和工程应用的强大工具，提供了丰富的数值计算和数据可视化功能，其中包括求解偏微分方程的工具箱，本文将深入探讨在Matlab中求解偏微分方程的方法和应用。
二、基本概念偏微分方程（Partial Differential Equation, PDE）是关于多个变量的函数及其偏导数的方程。
在物理学和工程学中，PDE广泛应用于描述空间变量和时间变量之间的关系。
在Matlab中，求解PDE通常涉及到确定PDE类型、边界条件、初始条件和求解方法等步骤。
三、求解方法1. 有限差分法（Finite Difference Method）有限差分法是求解PDE的常用数值方法之一，它将PDE转化为差分方程组，并通过迭代求解得到数值解。
在Matlab中，可以使用pdepe 函数来求解具有一维、二维或三维空间变量的PDE，该函数可以直接处理边界条件和初始条件。
2. 有限元法（Finite Element Method）有限元法是另一种常用的数值方法，它将求解区域离散化为有限数量的单元，并通过单元之间的插值来逼近PDE的解。
Matlab提供了pdenonlin函数来求解非线性PDE，该函数支持各种复杂的几何形状和非线性材料参数。
3. 特征线法（Method of Characteristics）特征线法适用于一维双曲型PDE的求解，该方法基于特征线方程的性质来构造数值解。
在Matlab中，可以使用pdegplot函数来展示特征线，并通过构造特征线网格来求解PDE。
四、实际应用1. 热传导方程的求解假设我们需要求解一个长条形的材料中的热传导方程，可以通过在Matlab中定义边界条件和初始条件，然后使用pdepe函数来求解得到温度分布和热流线。
一、Maple V 系统Maple V是由Waterloo大学开发的数学系统软件，它不但具有精确的数值处理功能，而且具有无以伦比的符号计算功能。
Maple V的符号计算能力还是MathCAD和MATLAB等软件的符号处理的核心。
Maple提供了2000余种数学函数，涉及范围包括：普通数学、高等数学、线性代数、数论、离散数学、图形学。
它还提供了一套内置的编程语言，用户可以开发自己的应用程序，而且Maple自身的2000多种函数，基本上是用此语言开发的。
Maple采用字符行输入方式，输入时需要按照规定的格式输入，虽然与一般常见的数学格式不同，但灵活方便，也很容易理解。
输出则可以选择字符方式和图形方式，产生的图形结果可以很方便地剪贴到Windows应用程序内。
二、MATLAB 系统MATLAB原是矩阵实验室（Matrix Laboratory）在70年代用来提供Linpack和Eispack 软件包的接口程序，采用C语言编写。
从80年代出现3.0的DOS版本，逐渐成为科技计算、视图交互系统和程序语言。
MATLAB可以运行在十几个操作平台上，比较常见的有基于Wi ndows 9X/NT、OS/2、Macintosh、Sun、Unix、Linux等平台的系统。
MATLAB程序主要由主程序和各种工具包组成，其中主程序包含数百个内部核心函数，工具包则包括复杂系统仿真、信号处理工具包、系统识别工具包、优化工具包、神经网络工具包、控制系统工具包、μ分析和综合工具包、样条工具包、符号数学工具包、图像处理工具包、统计工具包等。
而且5.x版本还包含一套几十个的PDF文件，从MATLAB的使用入门到其他专题应用均有详细的介绍。
MATLAB是数值计算的先锋，它以矩阵作为基本数据单位，在应用线性代数、数理统计、自动控制、数字信号处理、动态系统仿真方面已经成为首选工具，同时也是科研工作人员和大学生、研究生进行科学研究的得力工具。
一维抛物线偏微分方程数值解法（2）上一篇文章请参看一维抛物线偏微分方程数值解法（1）解一维抛物线型方程(理论书籍可以参看孙志忠：偏微分方程数值解法）Ut-Uxx=0, 0<x<1,0<t<=1（Ut-aUxx=f(x,t),a>0)U(x,0)=e^x, 0<=x<=1,U(0,t)=e^t,U(1,t)=e^(1+t), 0<t<=1精确解为：U(x,t)=e^(x+t);Matlab程序：（此为向后差分法）function [u p e x t]=pwxywxh(h1,h2,m,n)%欧拉向后差分法解一维抛物线型偏微分方程%此程序用的是追赶法解线性方程组%h1为空间步长，h2为时间步长%m,n分别为空间，时间网格数%p为精确解，u为数值解,e为误差x=(0:m)*h1+0;t=(0:n)*h2+0;for(i=1:n+1)for(j=1:m+1)f(i,j)=0;endendfor(i=1:n+1)u(i,1)=exp(t(i));u(i,m+1)=exp(1+t(i));endfor(i=1:m+1)u(1,i)=exp(x(i));endr=h2/(h1*h1);for(i=2:n+1) %外循环，先固定每一时间层，每一时间层上解一线性方程组%a(1)=0;b(1)=1+2*r;c(1)=-r;d(1)=u(i-1,2)+h2*f(i,2)+r*u(i,1);for(k=2:m-2)a(k)=-r;b(k)=1+2*r;c(k)=-r;d(k)=u(i-1,k+1)+h2*f(i,k+1);%输入部分系数矩阵，为0的矩阵元素不输入%enda(m-1)=-r;b(m-1)=1+2*r;d(m-1)=u(i-1,m)+h2*f(i,m)+r*u(i,m+1);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);endu(i,m)=d(m-1)/b(m-1); %回代过程%for(k=m-2:-1:1)u(i,k+1)=(d(k)-c(k)*u(i,k+2))/b(k);endendfor(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为误差endend[u p e x t]=pwxywxh(0.1,0.005,10,200);surf(x,t,e);xlabel('x');ylabel('t');zlabel('e');>> title('误差曲面');plot(t,e)误差较之前的欧拉向前差分格式增长了两倍[u p e x t]=pwxywxh(0.1,0.05,10,20); plot(t,e)[u p e x t]=pwxywxh(0.01,0.05,100,20); plot(t,e)[u p e x t]=pwxywxh(0.01,0.005,100,200);plot(x,e)[u p e x t]=pwxywxh(0.005,0.005,200,200); plot(x,e)X=1时，出现了误差？？? 不是边界条件吗？不能理解这方法还是比前一种方法误差大呀不过可以随便改变时间、空间步长。
matlab求解偏微分方程
Matlab求解偏微分方程的步骤:
1、首先，定义偏微分方程，并确定微分方程的种类；
2、然后，选择Matlab解决方案，所有内置微分方程求解器都支持基于初始值的手算方案；
3、接着，指定偏微分方程的解决参数，如函数、初始值、区间、边界
条件和终止条件；
4、之后，启动Matlab微分方程求解器，以计算偏微分方程的解决结果，如需要则可以绘制曲线图；
5、最后，检查偏微分方程的解决结果是否准确，可以利用MATLAB
自带的代数系统软件Maple来检查数值结果。
总体来说，使用Matlab求解偏微分方程非常容易，用户可以根据实际
情况，快速地完成偏微分方程的解决。
Matlab提供了一系列灵活的解
决方案，可以满足日常研究工作的所有需求。
另外，Matlab的可视化
绘图，可以帮助用户更好地理解偏微分方程的结果。
一维抛物型偏微分方程初边值问题求解【原创实用版】目录一、引言二、一维抛物型偏微分方程的基本概念1.抛物型偏微分方程的定义2.一维抛物型偏微分方程的特点三、一维抛物型偏微分方程初边值问题的求解方法1.紧差分法2.追赶法3.有限元算法四、各种方法的适用范围和优缺点比较1.紧差分法的适用范围和优缺点2.追赶法的适用范围和优缺点3.有限元算法的适用范围和优缺点五、结论正文一、引言一维抛物型偏微分方程在物理、化学、生物等领域有着广泛的应用，例如热传导方程、扩散方程等。
求解一维抛物型偏微分方程初边值问题，对于理解现实世界中的各种现象具有重要意义。
本文将介绍一维抛物型偏微分方程初边值问题的求解方法，包括紧差分法、追赶法和有限元算法，并对这些方法的适用范围和优缺点进行比较。
二、一维抛物型偏微分方程的基本概念1.抛物型偏微分方程的定义抛物型偏微分方程是指描述抛物型函数的偏微分方程，其一般形式为：u_t = au_xx + bu_x + cu其中，u 表示函数值，t 表示时间，x 表示空间坐标，a、b、c 为常数。
2.一维抛物型偏微分方程的特点一维抛物型偏微分方程的特点是其系数矩阵 A 为对称矩阵，且其特征值均为实数。
这使得一维抛物型偏微分方程的求解较为简单。
三、一维抛物型偏微分方程初边值问题的求解方法1.紧差分法紧差分法是一种常用的求解一维抛物型偏微分方程初边值问题的数值方法。
该方法通过离散化方程，将连续的偏微分方程转化为离散的代数方程组，然后使用迭代法求解该代数方程组。
紧差分法的精度为O(h^12h^24)，无条件差分稳定。
2.追赶法追赶法是另一种求解一维抛物型偏微分方程初边值问题的数值方法。
该方法通过解线性方程组，得到数值解。
追赶法的优点是稳定性较好，适用于较大时间步长和空间步长的情况。
3.有限元算法有限元算法是一种求解一维抛物型偏微分方程初边值问题的数值方法，其基本思想是将整个计算域进行网格划分，然后在每个网格节点上求解微分方程。
一维抛物线偏微分方程数值解法（3）上一篇参看一维抛物线偏微分方程数值解法（2）（附图及matlab程序）解一维抛物线型方程(理论书籍可以参看孙志忠：偏微分方程数值解法）Ut-Uxx=0, 0<x<1,0<t<=1（Ut-aUxx=f(x,t),a>0)U(x,0)=e^x, 0<=x<=1,U(0,t)=e^t,U(1,t)=e^(1+t), 0<t<=1精确解为：U(x,t)=e^(x+t);此种方法精度为o(h1^2+h2^2)一：用追赶法解线性方程组（还可以用迭代法解）Matlab程序function [u p e x t]=CN(h1,h2,m,n)%Crank-Nicolson格式差分法解一维抛物线型偏微分方程%此程序用的是追赶法解线性方程组%h1为空间步长，h2为时间步长%m,n分别为空间，时间网格数%p为精确解，u为数值解,e为误差x=(0:m)*h1+0; x0=(0:m)*h1; %定义x0,t0是为了f(x,t)~=0的情况%t=(0:n)*h2+0; t0=(0:n)*h2+1/2*h2;syms f;for(i=1:n+1)for(j=1:m+1)f(i,j)=0; %f(i,j)=f(x0(j),t0(i))==0%endendfor(i=1:n+1)u(i,1)=exp(t(i));u(i,m+1)=exp(1+t(i));endfor(i=1:m+1)u(1,i)=exp(x(i));endr=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的矩阵元素不输入%enda(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);endu(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);endendfor(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为误差endend[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('误差曲面')plot(x,e)plot(t,e)误差较向前欧拉法减小一半但是运行时间较长，约39秒,而前两次运行只需l秒左右；[u p e x t]=CN(0.01,0.01,100,100);运行需三分钟左右,误差比前次提高五倍,运算量也提高五倍[u p e x t]=CN(0.1,0.1,10,10);surf(x,t,e) 运行需要2秒；精度还是挺高的；[u p e x t]=CN(0.1,0.2,10,5);surf(x,t,e)误差还可以接受此种方法精度高，计算量较大二：用迭代法解线性方程组：Matlab程序如下：function [u e p x t k]=CN1(h1,h2,m,n,kmax,ep) % 解抛物线型一维方程 C-N格式（Ut-aUxx=f(x,t),a>0) %用g-s(高斯-赛德尔)迭代法解%kmax为最大迭代次数%m,n为x,t方向的网格数，例如（2-0）/0.01=200;%e为误差，p为精确解syms temp;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));endfor(i=1:m+1)u(1,i)=exp(x(i));endfor(i=1:n+1)for(j=1:m+1)f(i,j)=0;endenda=zeros(n,m-1);r=h2/(h1*h1); %此处r=a*h2/(h1*h1）；a=1for(k=1:kmax)for(i=1:n)for(j=2:m)temp=((r/2*u(i,j-1)+(1-r)*u(i,j)+r/2*u(i,...j+1)+h2*f(i,j)+r/2*u(i+1,j-1)+r/2*u(i+1,j+1))/(1+r));a(i+1,j)=(temp-u(i+1,j))*(temp-u(i+1,j));u(i+1,j)=temp;%此处注意是u(i+1,j),,而不是u(i+1,j+1)% endenda(i+1,j)=sqrt(a(i+1,j));if(k>kmax)break;endif(max(max(a))<ep)break;endendfor(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));endend[u e p x t k]=CN1(0.1,0.005,10,200,10000,1e-10);运行速度：1秒迭代次数k =81surf(x,t,e)第二幅图为三角追赶法解方程作出的图，两者几乎一样；由于迭代法速度很快，所以可以将区间分得更小[u e p x t k]=CN1(0.01,0.01,100,100,10000,1e-12);surf(x,t,e);shading interp; k=6903。
一维抛物线偏微分方程数值解法（3）
上一篇参看一维抛物线偏微分方程数值解法（2）（附图及matlab程序）
解一维抛物线型方程(理论书籍可以参看孙志忠：偏微分方程数值解法）
Ut-Uxx=0, 0<x<1,0<t<=1（Ut-aUxx=f(x,t),a>0)
U(x,0)=e^x, 0<=x<=1,
U(0,t)=e^t,U(1,t)=e^(1+t), 0<t<=1
精确解为：U(x,t)=e^(x+t);
此种方法精度为o(h1^2+h2^2)
一：用追赶法解线性方程组（还可以用迭代法解）
Matlab程序
function [u p e x t]=CN(h1,h2,m,n)
%Crank-Nicolson格式差分法解一维抛物线型偏微分方程
%此程序用的是追赶法解线性方程组
%h1为空间步长，h2为时间步长
%m,n分别为空间，时间网格数
%p为精确解，u为数值解,e为误差
x=(0:m)*h1+0; x0=(0:m)*h1; %定义x0,t0是为了f(x,t)~=0的情况%
t=(0:n)*h2+0; t0=(0:n)*h2+1/2*h2;
syms f;
for(i=1:n+1)
for(j=1:m+1)
f(i,j)=0; %f(i,j)=f(x0(j),t0(i))==0%
end
end
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
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('误差曲面')
plot(x,e) plot(t,e)
误差较向前欧拉法减小一半
但是运行时间较长，约39秒,而前两次运行只需l秒左右；
[u p e x t]=CN(0.01,0.01,100,100);运行需三分钟左右,误差比前次提高五倍,运算量也提高五倍
[u p e x t]=CN(0.1,0.1,10,10);surf(x,t,e) 运行需要2秒；精度还是挺高的；[u p e x t]=CN(0.1,0.2,10,5);surf(x,t,e)
误差还可以接受
此种方法精度高，计算量较大
二：用迭代法解线性方程组：
Matlab程序如下：
function [u e p x t k]=CN1(h1,h2,m,n,kmax,ep) % 解抛物线型一维方程 C-N格式（Ut-aUxx=f(x,t),a>0) %用g-s(高斯-赛德尔)迭代法解
%kmax为最大迭代次数
%m,n为x,t方向的网格数，例如（2-0）/0.01=200;
%e为误差，p为精确解
syms temp;
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
a=zeros(n,m-1);
r=h2/(h1*h1); %此处r=a*h2/(h1*h1）；a=1
for(k=1:kmax)
for(i=1:n)
for(j=2:m)
temp=((r/2*u(i,j-1)+(1-r)*u(i,j)+r/2*u(i,...
j+1)+h2*f(i,j)+r/2*u(i+1,j-1)+r/2*u(i+1,j+1))/(1+r));
a(i+1,j)=(temp-u(i+1,j))*(temp-u(i+1,j));
u(i+1,j)=temp;%此处注意是u(i+1,j),,而不是u(i+1,j+1)% end
end
a(i+1,j)=sqrt(a(i+1,j));
if(k>kmax)
break;
end
if(max(max(a))<ep)
break;
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
[u e p x t k]=CN1(0.1,0.005,10,200,10000,1e-10);运行速度：1秒
迭代次数k =
81
surf(x,t,e)
第二幅图为三角追赶法解方程作出的图，两者几乎一样；
由于迭代法速度很快，所以可以将区间分得更小
[u e p x t k]=CN1(0.01,0.01,100,100,10000,1e-12);surf(x,t,e);shading interp; k=6903。
