偏微分方程学习..docxVIP

1、椭圆形方程差分方法应用带微分余项的Taylor公式可得将上两式相加，得使用标准的二阶中心差分格式，可得 我们可以分别应用符号四阶紧致逼近算子到方程1中的二阶导和。离散的二维泊松方程可以表示为如下形式和用来表示在理查德森外推过程中忽略掉的误差，表示六阶精度。应用符号算子，设等于0，略去高阶项，方程可以写成如果我们设置网格比为，则方程可以写成系数是， ， ， .当网格比为1时，差分格式为如下形式用矩阵可表示为由于该差分格式是九点紧致格式，四阶收敛，所以误差可以表示为这里表示在步长为h的情况下的误差，A为精确解，表示在步长为h的情况下的数值解。由此可以看出，当步长降低一倍时，误差比值为16，这个可以有效的检验计算结果的准确性。数值算例该问题的精确解为.表1给出了不同步长取值下的误差，误差比值及计算时间。表1 计算误差表h最大误差误差比值运行时间0.14.8242e-5－0.0665280.053.0074e-616.04110.0457960.0251.8818e-715.98160.0534910.01251.1759e-816.00240.154051误差比值接近16，因此可以看出计算结果是比较准确的。图1、图2、图3分别为步长0.05时的精确解曲面，数值解曲面和误差曲面。图1 精确解曲面图2 步长为0.05时的数值解曲面图3 步长为0.05时的误差曲面MATLAB程序function [ error ] = jiudian( h )%UNTITLED Summary of this function goes here% Detailed explanation goes heretic;n = 1/h;for i=1:n+1 %boundary condition u(1,i) = sin(pi*h*(i-1)); u(n+1,i) = exp(1)*sin(pi*h*(i-1)); u(i,1) = 0; u(i,n+1) = 0;end%coefficient matrix AB = sparse(20*eye(n-1));D = sparse(-4*eye(n-1));for i = 1:n-2 B(i,i+1) = -4; B(i+1,i) = -4;endfor i = 1:n-2 D(i,i+1) = -1; D(i+1,i) = -1;enda = cell(n-1);for i = 1:n-1for j = 1:n-1 a{i,j} = sparse(zeros(n-1));endendfor i = 1:n-1 a(i,i) = {B};endfor i = 1:n-2 a{i,i+1} = D; a{i+1,i} = D;end%Cell Array to MatrixA = cell2mat(a);A = sparse(A);% right vector Kfor i = 1:n+1for j = 1:n+1 f(i,j) = (pi*pi-1)*exp((i-1)*h)*sin(pi*(j-1)*h);endendfor i = 2:nfor j = 2:n g(i-1,j-1) = 8*f(i,j) + f(i,j+1) + f(i,j-1) +f(i+1,j) +f(i-1,j);endendk = zeros(n-1);for i = 1:n-1 k(1,i) = 4*u(1,i+1) + u(1,i) + u(1,i+2) + k(1,i); k(i,1) = 4*u(i+1,1) + u(i,1) + u(i+2,1) + k(i,1); k(n-1,i) = 4*u(n+1,i+1) + u(n+1,i) + u(n+1,i+2) + k(n-1,i); k(i,n-1) = 4*u(i+1,n+1) + u(i,n+1) + u(i+2,n+1) + k(i,n-1);endK = h*h/2*g +k;K = reshape(K, (n-1)*(n-1), 1);%solve XX = A\K;for i = 2:nfor j = 2:n u(j,i) = X((n-1)*(i-2) + (j-1));endend%plot X, ex = 0:h:1;y = 0:h:1;subplot(2,2,1); surf(x,y,u);xlabel(x);ylabel(y);zlabel(u);title(差分解)for i = 1:n+1for j = 1:n+1