0%

Legendre 配置法及其算例

谱配置方法也称配点法,选取节点为 Legendre 高斯节点称为 Legendre 配置方法。这里介绍两点边值问题的 Legendre 配置方法,然后 MATLAB 编程算出数值解以及误差。

两点边值问题

考虑两点边值问题:

Legendre Collocation 方法

在区间 $[-1,1]$ 上选取配置点 $\left\{x_{j}\right\}_{j=0}^{N}$ 为 Legendre-Gauss-Lobatto 节点,其中 $x_{0}=-1$,$x_{N}=1$. 以上方程的配置法 (\label{equcm})

逼近解展开

其中 $h_j$ 是 Lagrange 基多项式(也称为节点基函数),即 $h_{j} \in P_{N}$,$h_{j}(x_{k})=\delta_{k j}$,令 $D=(d_{kj}:=h_j^{\prime}(x_k))_{k,j=1,2,\cdots,N.}$. 记 $w_j=u_{N}(x_j)$, 将 $u_{N}(x)$ 展开为 $u_{N}(x)=\sum_{j=0}^{N}w_j h_j(x)$, 可知

将以上公式代入方程 (equcm) 可得

线性系统可以简化为

数值例子

对于前面提到的两点边值问题,取 $\nu=1$,$\mu=1, \rho=1$,如果取精确解为 $\sin (k \pi x)$,则右端项为 $f=\mu k^2 \pi^2 \sin(k\pi x)+\nu k\pi \cos(k\pi x)+\rho \sin(k\pi x)$.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
% LegenCM2.m
% Legendre-collocation method for the model equation:
% -mu*u''(x)+nu*u'(x)+rho*u(x)=f(x), x in (-1,1),
% boundary condition: u(-1)=u(1)=0;
% exact solution: u=sin(kw*pi*x);
% RHS: f(x)=mu*kw^2*pi^2*sin(kw*pi*x)+nu*kw*pi*cos(kw*pi*x)+rho*sin(kw*pi*x)
% Rmk: Use routines lepoly(); legslb(); legslbdm();
clear all; close all;
kw=10;
nu=1;
mu=1;
rho=1;
Nvec=32:2:76;
% Initialization for error
L2_Err=[]; Max_Err=[];
% Loop for various modes N to calculate numerical errors
for N=Nvec
[xv,wv]=legslb(N); % Legendre-Gauss-Lobatto points and weights
u=sin(kw*pi*xv); % exact solution
f=mu*kw*kw*pi^2*sin(kw*pi*xv)+nu*kw*pi*cos(kw*pi*xv)+rho*sin(kw*pi*xv); % RHS
% Solve the collocation system
D1=legslbdiff(N,xv); % 1st order differentiation matrix
%D1=legslbdm(N); % 1st order differentiation matrix
D2=D1*D1; % 2nd order differentiation matrix
% Compositing the coefficient matrix
D=-mu*D2(2:N-1,2:N-1)+nu*D1(2:N-1,2:N-1)+rho*eye(N-2);
b=f(2:N-1); % RHS
un=D\b;
un=[0;un;0];

L2_error=sqrt(((un-u).^2)'*wv); % L^2 error
Max_error=norm(abs(un-u),inf); % maximum pointwise error
L2_Err=[L2_Err;L2_error];
Max_Err=[Max_Err;Max_error];
end
% Plot L^2 and maximum pointwise error
plot(Nvec,log10(L2_Err),'bo-','MarkerFaceColor','w','LineWidth',1)
hold on
plot(Nvec,log10(Max_Err),'rd-','MarkerFaceColor','w','LineWidth',1)
grid on
legend('L^2 error','L^{\infty} error','location','NorthEast')
set(gca,'fontsize',12)
xlabel('N','fontsize', 14), ylabel('log_{10}Error','fontsize',14)
% title('Convergence of Legendre-collocation method','fontsize',12)

% sets axis tick and axis limits
xticks(30:10:80)
yticks(-15:3:0)
xlim([30 80])
ylim([-15 0])

% print -dpng -r600 LegenCM2.png
% print -depsc2 LegenCM2.eps

LegenCM2.png

参考书籍

  • Spectral Methods: Algorithms, Analysis and Applications