1.Introduction
1.5.1 Finite-Difference Versus Spectral-Collocation P16
初学小白,书16页数值实现,有限差分法和谱配点法的对比。
生疏了,很不精简,有什么奇怪的地方,就忽略它吧!
function p16
clear;clc;
u0=inline('log(2+sin(x))','x');
du0=inline('(cos(x))./(2+sin(x))','x');%点除点乘啊啊啊,不然代入就是输出矩阵啊啊
% algebraic convergence for FEM and FDM
figure('Position',[0 0 1000 500] ), subplot(1,2,1);
for Nvec =10:4:100;
errfd =fdm(u0,du0,Nvec);%有限差分
errsp= spm(u0,du0,Nvec);
loglog(Nvec, errfd, 'r+', 'markersize', 8), hold on%双对数坐标变换
loglog(Nvec, errsp, 'bo', 'markersize', 8), hold on
end
title('Convergence of 4th-order FD &Fourier-collection methed');
grid on, xlabel N, ylabel error
for Nvec =10:0.1:100;
loglog(Nvec, (2.*pi./Nvec).^4, 'b-','markersize', 14),hold on%对数坐标系中h^{4}
end
text(100, 1e-5, 'h^{4}', 'fontsize', 18);
% spectral convergence for spectral method
subplot(1,2,2);%1行2列第二图
for Nvec =2:4:50;%图二横坐标
errsp= spm(u0,du0,Nvec);%Fourier 谱方法
[Nvec errsp];
semilogy(Nvec, errsp, 'r+', 'markersize', 8)
hold on;%y坐标系为对数坐标系
end
grid on, xlabel N, ylabel error
title('Convergence of Fourier-collection method');
end
function [u,du,x]= EvalSolEx(uinit,duin,N)
h=2*pi./N; x=(0:h:2*pi-h)';%xi
u=feval(uinit, x);%已知数据代入函数句柄,uinit即u0
du=feval(duin,x);
end
function errsp = spm(uinit,du0,N)%Fourier spectral
%uinit原函数,Nvec取点数
[u,du,x] = EvalSolEx(uinit,du0,N);
h=2*pi./N;
H=zeros(N,N);
for k=0:1:N-1
for j=0:1:N-1
if k~=j
H(j+1,k+1)=(1/2).*((-1).^(k+j))*(cot((j-k).*pi./N));
else
H(j+1,k+1)=0;
end
end
end
w=H*u;
errsp= norm(du-w, 'inf');%返回矩阵的无穷范数,即最大一行的和
end
% algebraic convergence for f1
function err= fdm(uinit,duin,N)%差分
errfd =[];
A=diag(repmat([0], 1, N))+diag(repmat([8],1, N-1), 1)+diag(repmat([-1], 1, N-2), 2)+diag(repmat([-8], 1, N-1), -1)+diag(repmat([1], 1, N-2), -2);
A(1,N-1)=1;A(2,N)=1;
A(1,N)=-8;
A(N-1,1)=-1;A(N,2)=-1;
A(N,1)=8;
[u,du,x] = EvalSolEx(uinit,duin, N);
h=2.*pi./N;
w=(1./(12.*h))*A*u;
err = norm(du-w, 'inf');
errfd = [errfd err];
end