有限体积法习题
对方程
-d/dx(p*du/dx)+q*u=f(x),
不妨使p=-1,q=0，f(x)=6x；a=0，b=1，u(a)=0,u’(b)=3.
则易知正解为u(x)=x^3.
为得到有限体积法收敛阶估计，由于已知正解，只需对一h值用有限体积元法求出其数值解uh，估计r1=|uh-u|1；然后分别对h2=h/2,h3=h2/2,h4=h3/2,…用有限体积元法求出其数值解uhi估计ri=|uhi-u|1，i=2,3,4,…；
则log2(ri/r(i-1))即为对收敛阶的估计。
编程如下：
函数：
function [u]=f(x)
u=x^3;
主程序：
r=zeros(1,10);%取十个h得到的误差储存空间
p=zeros(1,10);%误差阶估计
for i=1:10
a=0;%左端点
b=1;%右端点
n=2^(i+1);%分划数
h=(b-a)/n;%分划
x=zeros(1,n+1);%整分点
x_b=zeros(1,n+1);%半整分点
u=zeros(1,n+1);%整分点函数值
u_b=zeros(1,n+1);%半整分点函数值
x(0)=a;
for j=1:n
x(j)=x(j-1)+h;
x_b(j-1)=x(j)-h/2;
end
function [k]=f(x)
k=(x-x(i))/h;
end
因篇幅问题不能全部显示，请点此查看更多更全内容