如何求解包含三对角矩阵的方程组？MATLAB

Question

我不知道这是否是典型的行为，但我用后向差分法解决了一个有限差分问题。

我用适当的对角线项填充了一个稀疏矩阵(沿着中心对角线，上面和下面的一个)，我试图使用MATLAB的内置方法(B=A\x)来解决这个问题，似乎MATLAB就是弄错了。

此外，如果使用inv()并使用三对角矩阵的逆，则得到正确的解。

为什么这种行为会这样？

更多信息：

http://pastebin.com/AbuEW6CR (值被标出以便更容易读懂)刚度矩阵K：

1   0   0   0
-0.009  1.018   -0.009  0
0   -0.009  1.018   -0.009
0   0   0   1

D值：

0
15.55
15.55
86.73

内置输出：

-1.78595556155136e-05
0.00196073713853244
0.00196073713853244
0.0108149483252210

使用inv(K)的输出：

0
15.42
16.19
86.73

手工输出：

0
15.28
16.18
85.16

代码

nx = 21; %number of spatial steps
nt = 501; %number of time steps (varies between 501 and 4001)
p = alpha * dt / dx^2; %arbitrary constant
a = [0 -p*ones(1,nx-2) 0]'; %diagonal below central diagonal
b = (1+2*p)*ones(nx,1); %central diagonal
c = [1 -p*ones(1,nx-2) 1]'; %diagonal above central diagonal
d = zeros(nx, 1); %rhs values
% Variables a,b,c,d are used for the manual tridiagonal method for
% comparison with MATLAB's built-in functions. The variables represent
% diagonals and the rhs of the matrix    
% The equation is K*U(n+1)=U(N)

U = zeros(nx,nt);
% Setting initial conditions
U(:,[1 2]) = (60-32)*5/9;
K = sparse(nx,nx);
% Indices of the sparse matrix which correspond to the diagonal
diagonal = 1:nx+1:nx*nx; 

% Populating diagonals
K(diagonal) =1+2*p;
K(diagonal(2:end)-1) =-p;
K(diagonal(1:end-1)+1) =-p;

% Applying dirichlet condition at final spatial step, the temperature is
% derived from a table for predefined values during the calculation
K(end,end-1:end)=[0 1];

% Applying boundary conditions at first spatial step
K(1,1:2) = [1 0];


% Populating rhs values and applying boundary conditions, d=U(n)
d(ivec) = U(ivec,n);
d(nx) = R; %From table
d(1) = 0;

U(:,n+1) = tdm(a,b,c,d); % Manual solver, gives correct answer

U(:,n+1) = d\K; % Built-in solver, gives wrong answer

Answer 1

以下一行：

U(:,n+1) = d\K;

应该是

U(:,n+1) = K\d;

我错误地把它们搞错了，没有注意到，它显然改变了数学表达式，从而产生了错误的答案。