0

(已编辑,我已更改代码)嗯,我有复合方程,我需要在 matlab 中求解并找到结果。我尝试了不同的技术,但都失败了。方程是:

u(j-1)-2(u(j))+u(j+1)= -4*h^2*pi^2 * sin(2*pi*xj)

在哪里

n=100

j=1 到 n

xj=jh

h=1/n

u(0)==u(n)==0

我需要解方程并绘制结果。这样我就可以将结果与确切的解决方案进行比较。这是我到目前为止写的代码...

function c= prob1()
n=100;

c=(0);   % variable to store all results
u = linspace(1,n-1);
    for k=3:90
    jay=k;
    h=1/k;
    syms xj  
    eqn6 = u(jay-1) -2*u(jay)+u(jay+1)==-4*(h^2)*(pi^2)*sin(2*pi*xj);
    A = solve(eqn6, xj); % solving the equation with xj as unknown
      if(~(A==0))
      c=vertcat(c,A);  % just filtering out the results with 0 output
      end
    end
end

现在我在 A 中得到这样的答案“ (625*asin(1/9877545463176224))/3927 ”。我无法绘制。

4

1 回答 1

1

通过将数学转换为 MATLAB 语言来建立方程组Au = b,如下所示:

n = 100;
h = 1/n;

j = 0:n; % include zero for the boundary condition
xj = j*h;

% set up right hand side
b = (-4*h^2*pi^2*sin(2*pi*xj))';
% overwrite the first and last value of b with the right hand side of the boundary conditions:
b([1 end]) = 0;

% This is the interesting part:
% set up A: 1*u(j-1) -2*u(j) + 1*u(j+1) and include the boundary conditions
main = [1; -2*ones(n-1,1); 1];
upper = [0; ones(n-1,1)];
lower = [ones(n-1,1); 0];
A = gallery('tridiag', lower, main, upper);

A如果您不明白为什么会这样,我建议您根据and写出至少 j = 0、n/2 和 n 的b方程式,并将它们与您的方程式进行比较。

现在,我们准备解决系统问题。系统很小,所以我使用反斜杠运算符(这是一种直接方法),但您也可以选择迭代方法,如bicgstab, gmres, qmr:

u = A\b;

绘制结果u

plot(xj,u)
于 2014-10-11T10:12:52.577 回答