我的以下代码中的泊松边界值问题的近似解的错误在哪里？

问题描述

我已经解决了以下

的边值问题

它也在下面的代码Matlab中定义，但是我的代码不起作用，我的意思是我没有系统的近似解决方案，我想知道代码中的问题是什么，或者版本是什么我无法编译使用的matlab版本，谢谢

我使用的方法的说明：基于对装配矩阵和移动性矩阵的研究，我使用了有限元方法或我们所谓的 Galerkin方法通过满足边界条件的加权函数对系统进行变异，然后我对元素进行了积分（对范围为-1,1 [的元素矩阵进行了积分），我得到了四个元素矩阵。有关该方法的更多信息，我使用了pleas check { {3}}（第6,7,8页）

注意，我在编译代码时遇到的错误是

““ MatElt2Nd”的当前使用与其在第7行中的先前使用或定义不一致”

代码

function [U] = EquaDiff2(n)
% ----------------------------------
% -d²u/dx² + 6*u   = (-4*x^2-6)exp(x^2)
% u(-1) = 0 u(1)= 0
%----------------------------------
function [Ke,Fe] = MatElt2Nd(x1,x2) % déclaration de la fonction,function of computing matrix and elementary matrix (assembly matrix)
% ----------------------------------
x = ]-1:2/n:1['; % modification d1 of bound d’intégration
K = zeros(n+1 ) ;
F = zeros(n+1,1) ;
for i = 1:n
j = i+1;
t = [i j];
x1 = x(i);
x2 = x(j);
[Ke,x2);
K(t,t) = K(t,t) + Ke;
F(t) = F(t) + Fe;
end;
K(1,:) = [];
K(:,1) = [];
F(1) = [];
U = K\F;
U = [0.0;U];
t = 0:0.01:1;
return
%-------------------------------------------
% calculation of matrix  Ke and vector Fe
%-------------------------------------------
function [Ke,Fe] = MatElt2Nd0(x1,x2) % NEWly named nested function is introduced

Ke1 = 1/(x2-x1)*[ 1 -1     % no modification done
                  -1 1 ] ; % essentiellement que les matrices
Ke2 =(x2-x1)* [ 2 1        % élémentaires
               1 2 ] ;
N = [(x-x2)/(x1-x2) (x-x1)/(x2-x1)] % function of form
Fe =simple( int(N' * (-4*x^2-6)*exp(x^2),x,x1,x2) ) % vecteur Fe ;
Ke = Ke1 + 6*Ke2 ;

return

编辑我已经有了一个通用代码，但是我无法对通用代码进行更改来解决我的系统，有帮助吗？

常规代码

% au'(x)+bu"(x)=0 for 0<=x<=d
% BC: u(0)=0 and u(d)=h
%==============================================================
% ======Example======
% Finding an approximate solution to the following BVP using 4 elements of
% equal length.
% u'(x)-u"(x)=0 : 0<=x<=1
% BC: u(0)=0 and u(1)=1
% Solution:
% >> Galerkin(4,1,-1,1)
% ==============================================================

% The output of this program is 
% 1- The approximate solution (plotted in blue)
% 2- The exact solution (plotted in red)
% 3- The percentage error (plotted in magenta)

%=======================Program Begin==========================
function Galerkin(ne1,a,b,d,h) % Declare function
clc % Clear workspace

% Define the Coefficients of the  exact solution
% The Exact solution is : u(x)=C1+C2*exp(-ax/b)
% where C2=h/(exp(-a*d/b)-1)and C1=-C2
C2=h/((exp(-a*d/b))-1);
C1=-C2;

% Define element length
le = d/ne1; 

% Define x matrix
x = zeros (ne1+1,1); % 
for i=2:ne1 +1
  x(i,1) = x(i-1,1)+le;
end

% K1 matrix corresponding to the diffusion term (u"(x))
K1 = (b/le) * [1,-1;-1,1]
% K2 matrix corresponding to the convection term (u'(x))
K2 =  a*[-1/2 1/2;-1/2 1/2]

% Element stiffness Matrix
Ke = K1+K2

% Global stiffness matrix 
%********************Begin Assembly***************************
k = zeros(ne1+1);
for i=1:ne1+1
  for j=1:ne1 +1
      if (i==j)
          if(i==1)
              k(i,j)=Ke(1,1);
          elseif(i==ne1+1)
              k(i,j)=Ke(2,2);
          else
              k(i,1)+Ke(2,2);
          end
      elseif(i==j+1)
          k(i,2);
      elseif(j==i+1)
          k(i,1);
      else
          k(i,j)=0;
      end
  end
end
%********************End Assembly*****************************


%The Global f Matrix
f = zeros(ne1+1,1);
%BC apply u(0) = 0
f(1,1) = 0;
%BC apply u(d) = h
f(ne1+1,1) = h;
% Display the Global stifness matrix before striking row
K_Global=k
%Striking first row (u1=0)
k(1,1) = 1;
for i=2:ne1+1
  k(1,i) = 0;
  k(ne1+1,i) = 0;
end
k(ne1+1,ne1+1) = 1;
% Display the solvable stifness matrix 
K_strike=k
%solving the result and finding the displacement matrix,{u}
u=inv(k)*f
hold on

% ======Calculating Approximate Solution and plotting============
syms X
U_sym=sym(zeros(ne1,1));
dU_sym=sym(zeros(ne1,1));

for i=1:ne1
    N1x=1-((X-x(i))/le);
    N2x=(X-x(i))/le;
    U_X=(u(i)*N1x)+(u(i+1)*N2x);
    U_sym(i)=U_X;  
    dU_sym(i)=diff(U_sym(i));
    subplot(3,1)
    hold on
    ezplot(U_sym(i),[x(i) x(i+1)])
    subplot(3,2)
    hold on
    % du/dx approximate
    ezplot(dU_sym(i),[x(i) x(i+1)])
    
end

解决方法

暂无找到可以解决该程序问题的有效方法，小编努力寻找整理中！

如果你已经找到好的解决方法，欢迎将解决方案带上本链接一起发送给小编。

小编邮箱:dio#foxmail.com (将#修改为@）