(整理)微分方程数值解
- 格式:doc
- 大小:56.00 KB
- 文档页数:8
第七章常微分方程数值解本章介绍求解微分方程数值解的基本思想和方法.●常微分方程含有自变量、未知函数和它的一阶导数和高阶导数的方程.它是描述运动、变化规律的重要数学方法之一,分为两类:1.初值问题,即给出未知函数及导数在初始点的值2.边值问题,即给出未知函数及(或)它的某些导数在区间两个端点的值●常微分方程初值问题00(,),, ||(),,y f x y a x b y y a y y R '=<≤<∞⎧⎨=∈⎩ 其中f 为,x y 的已知函数,0y 为给定的初值. 这里仅讨论一阶标量微分方程初值问题的数值解法.而高阶微分方程通常可化为一阶微分方程组来研究.数值解法寻求微分方程初值问题之解()y x 在一系列离散点012N a x x x x b =<<<<=上的近似值:,,,210y y yN y , 的方法.{}n y : 问题的数值解数值解所满足的离散方程统称为差分格式. 步长:1i i i h x x -=- ,一般取定步长i h h =● 初值问题的适定性(其解是否唯一存在)记(带形)区域:[,]a b R ⨯为G ,即G =[,]a b R ⨯. 设:f G R →为连续映射,若存在常数0L >使得不等式1212|(,)(,)|||f x y f x y L y y -≤- 对一切12(,),(,)x y x y G ∈都成立,则称(,)f x y 在G 上关于y 满足Lipschitz 条件,而式中的常数L 称为Lipschitz 常数.定理 初值问题00(,),, ||(),,y f x y a x b y y a y y R '=<≤<∞⎧⎨=∈⎩,当(,)f x y 在G 上连续,且关于y 满足Lipschitz 条件,则其解存在且唯一.§ 7.1 Euler 方法● Euler 公式将初值问题⎩⎨⎧∈=∞<≤<=',,)(||,),,(00R y y a y y b x a y x f y 的求解区间],[b a N 等分,分点:),,2,1,0(,N n nh a x n =+= ,其中N ab h -=将),(y x f y ='写成等价的积分方程形式:⎰++=+hx x d y f x y h x y τττ))(,()()(,在上式中令n x x =,并用左矩形公式计算右端积分,得到n n n n n R x y x hf x y h x y ++=+))(,()()(, (*)⎰+-=1))(,())(,(n nx x n n n x y x hf dx x y x f R ——余项将(*)中的余项n R 截去,可得))(,()()(1n n n n x y x hf x y x y +≈+则有)(n x y 的近似值n y 的递推公式1,,1,0),,(1-=+=+N n y x hf y y n n n n (*1)----------Euler 公式n R 表示当)(n n x y y =为精确值时,利用(*1)即Euler 公式计算)(1+n x y 时的误差.n n n y x y -=)(ε 为Euler 方法的局部截断误差.例1 用Euler 公式解初值问题⎪⎩⎪⎨⎧=≤<-='1)0(10,2y x y x y y解:取1.0=h ,Euler 公式的具体形式为)2(),(1nnn n n n n n y x y h y y x hf y y -+=+=+其中)10,,1,0,1.0 ===n n nh x n ( 已知10=y ,则有1.11.01)2(0001=+=-+=y x y h y y191818.1)1.12.01.1(1.01.1)2(11112=-+=-+=y x y h y y… … … 依次计算可得109876543,,,,,,,y y y y y y y y其部分结果见下表可见Euler 方法的计算结果精度不太高。
微分方程数值解法微分方程数值解法微分方程数值解法【1】摘要:本文结合数例详细阐述了最基本的解决常微分方程初值问题的数值法,即Euler方法、改进Euler法,并进行了对比,总结了它们各自的优点和缺点,为我们深入探究微分方程的其他解法打下了坚实的基础。
关键词:常微分方程数值解法 Euler方法改进Euler法1、Euler方法由微分方程的相关概念可知,初值问题的解就是一条过点的积分曲线,并且在该曲线上任一点处的切线斜率等于函数的值。
根据数值解法的基本思想,我们取等距节点,其中h为步长,在点处,以为斜率作直线交直线于点。
如果步长比较小,那么所作直线与曲线的偏差不会太大,所以可用的近似值,即:,再从点出发,以为斜率作直线,作为的近似值,即:重复上述步骤,就能逐步求出准确解在各节点处的近似值。
一般地,若为的近似值,则过点以为斜率的直线为:从而的近似值为:此公式就是Euler公式。
因为Euler方法的思想是用折线近似代替曲线,所以Euler方法又称Euler折线法。
Euler方法是初值问题数值解中最简单的一种方法,由于它的精度不高,当步数增多时,由于误差的积累,用Euler方法作出的折线可能会越来越偏离曲线。
举例说明:解: ,精确解为:1.2 -0.96 -1 0.041.4 -0.84 -0.933 0.9331.6 -0.64 -0.8 0.161.8 -0.36 -0.6 0.242.0 0 -0.333 0.332.2 0.44 0 0.44通过上表可以比较明显地看出误差随着计算在积累。
2、改进Euler法方法构造在常微分方程初值问题 ,对其从到进行定积分得:用梯形公式将右端的定积分进行近似计算得:用和来分别代替和得计算格式:这就是改进的Euler法。
解:解得:由于 ,是线形函数可以从隐式格式中解出问题的精确解是误差0.2 2.421403 2.422222 0.000813 0.021400.4 2.891825 2.893827 0.00200 0.051830.6 3.422119 3.425789 0.00367 0.094112.0 10.38906 10.43878 0.04872 1.1973通过比较上表的第四列与第五列就能非常明显看出改进Euler方法精度比Euler方法精度高。
微分方程的数值解法微分方程是描述自然界中众多现象和规律的重要数学工具。
然而,许多微分方程是很难或者无法直接求解的,因此需要使用数值解法来近似求解。
本文将介绍几种常见的微分方程数值解法。
1. 欧拉方法欧拉方法是最简单的数值解法之一。
它将微分方程转化为差分方程,通过计算离散点上的导数来逼近原方程的解。
欧拉方法的基本思想是利用当前点的导数值来估计下一个点的函数值。
具体步骤如下:首先,将自变量区间等分为一系列的小区间。
然后,根据微分方程的初始条件,在起始点确定初始函数值。
接下来,根据导数的定义,计算每个小区间上函数值的斜率。
最后,根据初始函数值和斜率,递推计算得到每个小区间上的函数值。
2. 龙格-库塔方法龙格-库塔方法是一种常用的高阶精度数值解法。
它通过进行多次逼近和修正来提高近似解的准确性。
相比于欧拉方法,龙格-库塔方法在同样的步长下可以获得更精确的解。
具体步骤如下:首先,确定在每个小区间上的步长。
然后,根据微分方程的初始条件,在起始点确定初始函数值。
接下来,根据当前点的导数值,使用权重系数计算多个中间点的函数值。
最后,根据所有中间点的函数值,计算出当前点的函数值。
3. 改进欧拉方法(改进的欧拉-克罗默法)改进欧拉方法是一种中阶精度数值解法,介于欧拉方法和龙格-库塔方法之间。
它通过使用两公式递推来提高精度,并减少计算量。
改进欧拉方法相对于欧拉方法而言,增加了一个估计项,从而减小了局部截断误差。
具体步骤如下:首先,确定在每个小区间上的步长。
然后,根据微分方程的初始条件,在起始点确定初始函数值。
接下来,利用欧拉方法计算出中间点的函数值。
最后,利用中间点的函数值和斜率,计算出当前点的函数值。
总结:微分方程的数值解法为我们研究和解决实际问题提供了有力的工具。
本文介绍了欧拉方法、龙格-库塔方法和改进欧拉方法这几种常见的数值解法。
选择合适的数值解法取决于微分方程的性质以及对解的精确性要求。
在实际应用中,我们应该根据具体情况选择最合适的数值解法,并注意控制步长以尽可能减小误差。
微分方程组的数值求解方法微分方程组数值求解方法微分方程组是数学中非常重要的一个分支,它描述了许多自然界和社会生活中的现象,例如电路的运行、天体的运行、生命体的生长等等。
我们需要对微分方程组进行求解,才能够得到它们的解析解,从而更好地理解和应用它们。
然而,大多数微分方程组不可能用解析法求解,因此,我们需要采用数值方法来求解微分方程组。
常见的微分方程组数值求解方法包括欧拉法、龙格库塔法和变步长法等。
下面,我们将逐一介绍它们的基本原理和优缺点。
一、欧拉法欧拉法是微分方程组数值求解方法中最简单的一种。
它的基本思想是将微分方程组中的各个变量离散化,然后根据微分方程组的导数计算每一步的值。
具体来讲,欧拉法的数值求解公式为:\begin{aligned} &x_{n+1}=x_n+hf_n(x_n,y_n,z_n),\\&y_{n+1}=y_n+hf_n(x_n,y_n,z_n),\\&z_{n+1}=z_n+hf_n(x_n,y_n,z_n), \end{aligned}其中,$x(t)$,$y(t)$,$z(t)$是微分方程组的解,$f_n(x_n,y_n,z_n)$是微分方程组导数在点$(x_n,y_n,z_n)$处的值,$h$为时间步长。
欧拉法的优点是简单易懂,方便实现,缺点是误差较大,计算不够精确。
因此,在实际应用中,往往需要采用更加精确的数值方法。
二、龙格库塔法龙格库塔法是微分方程组数值求解方法中比较常用的一种。
它的基本思想是通过多次计算微分方程组中的导数,以获得更加精确的数值解。
具体来讲,龙格库塔法的求解公式为:\begin{aligned}&k_{1x}=hf_n(x_n,y_n,z_n),k_{1y}=hf_n(x_n,y_n,z_n),k_{1z}=hf_n (x_n,y_n,z_n),\\&k_{2x}=hf_n(x_n+\frac{h}{2},y_n+\frac{k_{1y}}{2},z_n+\frac{k_ {1z}}{2}),k_{2y}=hf_n(x_n+\frac{h}{2},y_n+\frac{k_{1y}}{2},z_n+ \frac{k_{1z}}{2}),k_{2z}=hf_n(x_n+\frac{h}{2},y_n+\frac{k_{1y}}{ 2},z_n+\frac{k_{1z}}{2}),\\&k_{3x}=hf_n(x_n+\frac{h}{2},y_n+\frac{k_{2y}}{2},z_n+\frac{k_ {2z}}{2}),k_{3y}=hf_n(x_n+\frac{h}{2},y_n+\frac{k_{2y}}{2},z_n+ \frac{k_{2z}}{2}),k_{3z}=hf_n(x_n+\frac{h}{2},y_n+\frac{k_{2y}}{ 2},z_n+\frac{k_{2z}}{2}),\\&k_{4x}=hf_n(x_n+h,y_n+k_{3y},z_n+k_{3z}),k_{4y}=hf_n(x_n+h,y_n+k_{3y},z_n+k_{3z}),k_{4z}=hf_n(x_n+h,y_n+k_{3y},z_n+k_{3 z}),\\&x_{n+1}=x_n+\frac{k_{1x}}{6}+\frac{k_{2x}}{3}+\frac{k_{3x}}{ 3}+\frac{k_{4x}}{6},\\&y_{n+1}=y_n+\frac{k_{1y}}{6}+\frac{k_{2y}}{3}+\frac{k_{3y}}{ 3}+\frac{k_{4y}}{6},\\&z_{n+1}=z_n+\frac{k_{1z}}{6}+\frac{k_{2z}}{3}+\frac{k_{3z}}{ 3}+\frac{k_{4z}}{6}, \end{aligned}其中,$k_{1x}$,$k_{1y}$,$k_{1z}$,$k_{2x}$,$k_{2y}$,$k_{2z}$,$k_{3x}$,$k_{3y}$,$k_{3z}$,$k_{4x}$,$k_{4y}$,$k_{4z}$是微分方程组中导数的值。
微分方程数值解使用数值方法求解微分方程微分方程是描述自然现象中变化的数学模型,它是数学和科学研究中的重要工具。
然而,许多微分方程并没有精确的解析解,因此需要使用数值方法来近似求解。
本文将介绍一些常用的数值方法来求解微分方程,包括欧拉方法、改进的欧拉方法和龙格-库塔方法。
一、欧拉方法欧拉方法是最简单、最基础的数值方法之一。
它基于微分方程解的定义,通过离散化自变量和因变量来逼近解析解。
假设我们要求解的微分方程为dy/dx = f(x, y),初始条件为y(x0) = y0。
将自变量x分割成若干个小区间,步长为h,得到x0, x1, x2, ..., xn。
根据微分方程的定义,我们可以得到递推公式 yn+1 = yn + h*f(xn, yn)。
用代码表示即为:```def euler_method(f, x0, y0, h, n):x = [x0]y = [y0]for i in range(n):xn = x[i]yn = y[i]fn = f(xn, yn)xn1 = xn + hyn1 = yn + h*fnx.append(xn1)y.append(yn1)return x, y```二、改进的欧拉方法欧拉方法存在着局部截断误差,即在每个小区间上的误差。
改进的欧拉方法是对欧拉方法的改进,可以减小截断误差。
它的递推公式为yn+1 = yn + h*(f(xn, yn) + f(xn+1, yn+1))/2。
用代码表示即为:```def improved_euler_method(f, x0, y0, h, n):x = [x0]y = [y0]for i in range(n):xn = x[i]yn = y[i]fn = f(xn, yn)xn1 = xn + hyn1 = yn + h*(fn + f(xn1, yn + h*fn))/2x.append(xn1)y.append(yn1)return x, y```三、龙格-库塔方法龙格-库塔方法是一种更加精确的数值方法,它通过计算多个递推式的加权平均值来逼近解析解。
微分方程的数值解法微分方程(Differential Equation)是描述自然界中变化的现象的重要工具,具有广泛的应用范围。
对于一般的微分方程,往往很难找到解析解,这时候就需要使用数值解法来近似求解微分方程。
本文将介绍几种常见的微分方程数值解法及其原理。
一、欧拉方法(Euler's Method)欧拉方法是最基本也是最容易理解的数值解法之一。
它的基本思想是将微分方程转化为差分方程,通过给定的初始条件,在离散的点上逐步计算出函数的近似值。
对于一阶常微分方程dy/dx = f(x, y),利用欧拉方法可以得到近似解:y_n+1 = y_n + h * f(x_n, y_n)其中,h是步长,x_n和y_n是已知点的坐标。
欧拉方法的优点在于简单易懂,但是由于是一阶方法,误差较大,对于复杂的微分方程可能不够准确。
二、改进的欧拉方法(Improved Euler's Method)改进的欧拉方法又称为改进的欧拉-柯西方法,是对欧拉方法的一种改进。
它通过在每一步计算中利用两个不同点的斜率来更准确地逼近函数的值。
对于一阶常微分方程dy/dx = f(x, y),改进的欧拉方法的迭代公式为:y_n+1 = y_n + (h/2) * [f(x_n, y_n) + f(x_n+1, y_n + h * f(x_n, y_n))]相较于欧拉方法,改进的欧拉方法具有更高的精度,在同样的步长下得到的结果更接近真实解。
三、四阶龙格-库塔方法(Fourth-Order Runge-Kutta Method)四阶龙格-库塔方法是一种更高阶的数值解法,通过计算多个点的斜率进行加权平均,得到更为准确的解。
对于一阶常微分方程dy/dx = f(x, y),四阶龙格-库塔方法的迭代公式为:k1 = h * f(x_n, y_n)k2 = h * f(x_n + h/2, y_n + k1/2)k3 = h * f(x_n + h/2, y_n + k2/2)k4 = h * f(x_n + h, y_n + k3)y_n+1 = y_n + (k1 + 2k2 + 2k3 + k4)/6四阶龙格-库塔方法是数值解法中精度最高的方法之一,它的计算复杂度较高,但是能够提供更为准确的结果。
常微分方程的数值解法常微分方程是研究变量的变化率与其当前状态之间的关系的数学分支。
它在物理、工程、经济等领域有着广泛的应用。
解常微分方程的精确解往往十分困难甚至不可得,因此数值解法在实际问题中起到了重要的作用。
本文将介绍常见的常微分方程的数值解法,并比较其优缺点。
1. 欧拉方法欧拉方法是最简单的数值解法之一。
它基于近似替代的思想,将微分方程中的导数用差商近似表示。
具体步骤如下:(1)确定初始条件,即问题的初值。
(2)选择相应的步长h。
(3)根据微分方程的定义使用近似来计算下一个点的值。
欧拉方法的计算简单,但是由于误差累积,精度较低。
2. 改进欧拉方法为了提高欧拉方法的精度,改进欧拉方法应运而生。
改进欧拉方法通过使用两个点的斜率的平均值来计算下一个点的值。
具体步骤如下:(1)确定初始条件,即问题的初值。
(2)选择相应的步长h。
(3)根据微分方程的定义使用近似来计算下一个点的值。
改进欧拉方法相较于欧拉方法而言,精度更高。
3. 龙格-库塔法龙格-库塔法(Runge-Kutta)是常微分方程数值解法中最常用的方法之一。
它通过迭代逼近精确解,并在每一步中计算出多个斜率的加权平均值。
具体步骤如下:(1)确定初始条件,即问题的初值。
(2)选择相应的步长h。
(3)计算各阶导数的导数值。
(4)根据权重系数计算下一个点的值。
与欧拉方法和改进欧拉方法相比,龙格-库塔法的精度更高,但计算量也更大。
4. 亚当斯法亚当斯法(Adams)是一种多步法,它利用之前的解来近似下一个点的值。
具体步骤如下:(1)确定初始条件,即问题的初值。
(2)选择相应的步长h。
(3)通过隐式或显式的方式计算下一个点的值。
亚当斯法可以提高精度,并且比龙格-库塔法更加高效。
5. 多步法和多级法除了亚当斯法,还有其他的多步法和多级法可以用于解常微分方程。
多步法通过利用多个点的值来逼近解,从而提高精度。
而多级法则将步长进行分割,分别计算每个子问题的解,再进行组合得到整体解。
微分方程数值解法微分方程是数学中的重要概念,它描述了物理系统中变量之间的关系。
解微分方程是许多科学领域中常见的问题,其中又可以分为解析解和数值解两种方法。
本文将重点介绍微分方程的数值解法,并详细讨论其中的常用方法和应用。
一、微分方程的数值解法概述微分方程的解析解往往较为复杂,难以直接求解。
在实际问题中,我们通常利用计算机进行数值计算,以获得方程的数值解。
数值解法的基本思想是将微分方程转化为一组离散的数值问题,通过逼近连续函数来获得数值解。
二、常见的数值解法1. 欧拉法欧拉法是最基础的数值解法之一,其核心思想是将微分方程转化为差分方程,通过逼近连续函数来获得数值解。
欧拉法的基本形式为:yn+1 = yn + h·f(xn, yn)其中,yn表示第n个时间步的数值解,h为时间步长,f为微分方程右端的函数。
欧拉法的精度较低,但计算简单,适用于初步估计或简单系统的求解。
2. 改进的欧拉法(Heun法)改进的欧拉法(Heun法)是对欧拉法的改进,其关键在于求解下一个时间步的近似值时,利用了两个斜率的平均值。
Heun法的基本形式为:yn+1 = yn + (h/2)·(k1 + k2)k1 = f(xn, yn),k2 = f(xn+h, yn+h·k1)Heun法较欧拉法的精度更高,但计算量较大。
3. 龙格-库塔法(RK方法)龙格-库塔法是一类常用的数值解法,包含了多个不同阶数的方法。
其中,最常用的是经典四阶龙格-库塔法(RK4法),其基本形式为:k1 = f(xn, yn)k2 = f(xn + h/2, yn + (h/2)·k1)k3 = f(xn + h/2, yn + (h/2)·k2)k4 = f(xn + h, yn + h·k3)yn+1 = yn + (h/6)·(k1 + 2k2 + 2k3 + k4)RK4法实现较为复杂,但精度较高,适用于解决大多数常微分方程问题。
微分方程数值解差分法微分方程是自然科学和工程技术中广泛使用的工具,它们描述了许多物理过程的动力学行为。
对于复杂的微分方程,解析解往往很难或者不可能得到。
此时我们需要数值解差分法来解决问题。
一、微分方程数值解的方法1.分裂法分裂法是将一个复杂的微分方程分解为多个简单的方程。
例如,将一个偏微分方程分解成几个常微分方程,从而可以方便地使用数值方法计算解。
2.有限差分法有限差分法是一种常见的微分方程数值计算方法。
它将一维或多维的连续函数离散为一系列离散点,然后使用差分方程近似微分方程,最后用迭代法计算数值解。
3.有限元法有限元法是一种广泛使用的数值计算方法,它可以用于求解各种类型的微分方程。
该方法将求解区域分割成多个小区域,然后对每个小区域进行离散化和近似处理。
二、数值解差分法数值解差分法是微分方程数值解的基本方法之一。
它是一种基于差分方程的离散化方法,可以对微分方程进行近似,并将微分方程转化为一个差分方程。
数值解的差分法可以分为前向差分、后向差分和中心差分三种方法。
1.前向差分法前向差分法使用前一时间步的值,计算当前时间步的值。
它的近似误差随着时间步长的增大而增大。
前向差分的公式如下:y_i+1 = y_i + hf_i(x_i,y_i)其中,h是时间步长,f_i是微分方程的左侧。
2.后向差分法后向差分法使用后一时间步的值,计算当前时间步的值。
它的近似误差随着时间步长的增大而减小。
后向差分的公式如下:y_i+1=y_i + hf_i(x_i+1,y_i+1)3.中心差分法中心差分法使用前一时间步和后一时间步的值,计算当前时间步的值。
它的近似误差随着时间步长的增大而增大。
中心差分的公式如下:y_i+1=y_i + 1/2hf_i(x_i,y_i) + 1/2hf_i(x_i+1,y_i+1)三、差分法的优缺点差分法作为微分方程数值解的一种基本方法,具有以下优缺点:1.优点(1)简单易实现:差分法的实现很简单,只需要计算微分方程的离散值和靠近值即可。
微分方程的解析与数值解法微分方程既是数学分析的重要分支,也是许多学科领域的基础。
在实际问题的求解中,我们常常需要寻找微分方程的解析解或者数值解。
本文将围绕微分方程的解析和数值解法展开讨论。
一、微分方程的解析解解析解指的是通过代数计算得到的方程的解。
对于某些简单的微分方程,我们可以通过分离变量、变量代换等方法得到解析解。
下面以一阶线性常微分方程为例,讨论解的求解过程。
考虑一阶线性常微分方程形式如下:$$\frac{dy}{dx} + P(x)y = Q(x)$$其中,$P(x)$和$Q(x)$为已知函数。
我们可以通过以下步骤求解该微分方程:1. 将方程改写为标准形式:$\frac{dy}{dx} + P(x)y - Q(x) = 0$2. 求解齐次线性微分方程:$\frac{dy}{dx} + P(x)y = 0$。
记其解为$y_h$,即$y_h = Ce^{-\int P(x)dx}$,其中$C$为常数。
3. 利用常数变易法,假设原方程的解为$y = u(x)y_h$,其中$u(x)$为待定函数。
4. 将$y = u(x)y_h$代入原方程,得到关于$u(x)$的方程。
5. 求解$u(x)$的方程,得到$u(x)$的表达式。
6. 将$u(x)$代入$y = u(x)y_h$,得到原方程的解析解。
上述过程就是一阶线性常微分方程求解的一般步骤。
对于其他类型的微分方程,也有相应的解析解求解方法。
但并非所有微分方程都存在解析解。
二、微分方程的数值解法对于一些复杂的微分方程,无法找到解析解,此时我们需要借助数值方法求解。
常见的数值解法包括欧拉法、改进的欧拉法、四阶龙格-库塔法等。
1. 欧拉法欧拉法是一种较为简单的数值解法,其基本思想是通过离散化微分方程,将微分方程转化为差分方程。
具体步骤如下:将求解区间$[a, b]$等分成$n$个小段,步长为$h = \frac{b-a}{n}$。
利用微分方程的导数定义,将微分方程转化为差分方程,即$y_{i+1} = y_i + h \cdot f(x_i, y_i)$,其中$f(x, y)$为微分方程右端的函数。
实验1
题目:
两点边值问题的差分求解
目的与要求:
掌握中心差分格式的程序实现
掌握分析算法误差的方法
实验内容:
(i) 分别在步长h=1/20,1/40,1/80,1/160情形下用中心差分格式计算齐次两点边值问题-u''=f,u(0)=u(1)=0。
其中f(x) = 100*exp(-10*x),精确解为u(x) = 1 - (1-exp(-10))*x - exp(-10*x)
(ii) 给出差分解近似精确解在无穷范数和离散L2范数下的误差阶。
% 实验报告提交期限 2011.10.10
提示:关于(ii)必须给出类似下面的表格
实验2
题目:
五点差分格式
目的与要求:
掌握五点差分格式的程序实现
掌握分析算法误差的方法
实验内容:
(i) 分别在步长h=1/16,1/32,1/64,1/128情形下用五点差分格式计算二维椭圆问题-delta(u) = f,求解区域为[0,1]*[0,1],边值条件为Dirichlet 齐次边值条件。
建议五点差分格式的求解使用Gauss-Seidel 迭代法
(ii) 给出系数矩阵的非零元图像,使用函数spy(A)
(iii) 分析差分解近似精确解在无穷范数和L2范数下的误差阶。
% 实验报告提交期限 2010.11.20
提示:
关于(i )要先写出2个排序函数并测试
关于(iii )必须给出类似下面的表格
实验3
题目:最简差分格式
目的与要求:
掌握向前向后差分格式的程序实现
掌握分析算法误差的方法
实验内容:
(i) 分别用向前和向后差分格式计算一维抛物问题
22
u u f t x ∂∂=+∂∂ 求解区域为[0,1]*[0,T],其中T=1,精确解和右端项分别为
(,)sin()t u x t e x π-=
2(,)(1)sin()t f x t e x ππ-=-
定解条件为齐次边值条件及初值条件
(,0)sin()u x x π=
(ii) 在空间步长h=1/10, 1/20, 1/40,网比r 固定为1/2情形下,分析差分解近似精确解在无穷范数和L2范数下关于空间步长和时间步长的误差阶。
%实验一
clc;clear all;
digits(32);
f=inline('100*exp(-10*x)');
u=inline('1-(1-exp(-10))*x-exp(-10*x)');
W=[20 40 80 160];
iLabel={'步长为1/20' '步长为1/40' '步长为1/80' '步长为1/160'};
%四种步长下的差值结果
for l=1:4
%N为该次分割区间数n为插值节点个数
N=W(l);
h=1/N;
n=N-1;
a=-1*ones(1,n);
b=2*ones(1,n);
c=a;
x=a;U=a;
for k=1:n
F(k)=f(k*h)*h^2;
U(k)=u(k*h);
end
%求解插值节点的值
for i=1:n-1
a(i+1)=a(i+1)/b(i);
b(i+1)=b(i+1)-a(i+1)*c(i);
F(i+1)=F(i+1)-a(i+1)*F(i);
end
x(n)=F(n)/b(n);
for k=n-1:-1:1
x(k)=(F(k)-c(k)*x(k+1))/b(k);
end
%寻找最大误差
max=abs(U(1)-x(1));
for i=2:n
if max<abs(U(i)-x(i));
max= abs(U(i)-x(i));
end
end
max
L2=sqrt((U-x)*(U-x)'*h) %绘制拟合曲线
t=1/N:h:(1-1/N);
subplot(2,2,l);
plot(t,x,'r*',t,U,'g');
title(iLabel{l});
end
%%向后差分实验
clc;clear all;
format long;
r=1/2; shuzu=[10,20,40];
for bu=1:3
N=shuzu(bu)-1;
h=1/(N+1);
t=h^2*r;
M=1/t;
V1=zeros(1,N+2);
x=0:h:1;
V=sin(pi*x);
U=exp(-1)*sin(pi*x);
%
for k=1:M
a=-r*ones(1,N-1);
b=(1+2*r)*ones(1,N);
f=(pi^2-1)*exp(-k*t)*sin(pi*x);
F=f*t+V;F=F(2:N+1);
X=zeros(N,N);
X=diag(b)+diag(a,1)+diag(a,-1);
V1=inv(X)*F';V(2:N+1)=V1;
end
e=V(2:N+1)-U(2:N+1);
%%寻找最大误差
mmax=abs(e(1));
for i=1:N
mmax= max(abs(e(i)),mmax);
end
mmax;
%%L2误差
L2=e*e';
L2=sqrt(L2*h);
%求误差的阶Morder表示最大误差阶,Lorder 表示L2误差阶
if bu>1
Morder=log2(zmax/mmax)
Lorder=log2(zL2/L2)
end zmax=mmax zL2=L2
end。