ipeidg1.m

来自「various matlab programs to slove various」· M 代码 · 共 46 行

function [u,x]=ipeidg1(I,n,a,u0,bc)
nx = n(1); 
h = (I(2)-I(1))/nx; 
x = [I(1):h:I(2)]; 
t = I(3); 
um = (eval(u0))';
u = []; 
xx=[];
for i = 1:nx+1
   u = [u, um(i), um(i)];   
   xx = [xx, x(i), x(i)];   
end
u = u'; 
nt = n(2); 
k = (I(4)-I(3))/nt;
lambda = k/h; 
e = ones(2*nx+2,1);
B = spdiags([1/6*e,1/3*e,1/6*e],-1:1,2*nx+2,2*nx+2);
dd = 1/3+0.5*a*lambda; 
du = 1/6+0.5*a*lambda;
dl = 1/6-0.5*a*lambda; 
A=sparse([]);
A(1,1) = dd; 
A(1,2) = du; 
A(2,1) = dl; 
A(2,2) = dd;
for i=3:2:2*nx+2
   A(i,i-1)   =-a*lambda;
   A(i,i)     = dd;   
   A(i,i+1)   = du;
   A(i+1,i)   = dl;   
   A(i+1,i+1) = A(i,i);
end
[L,U]=lu(A);
for t = I(3)+k:k:I(4)
   f = B*u;
   if a > 0
      f(1) = a*bc(1)+f(1);
   elseif a <= 0
      f(nx) = a*bc(2)+f(nx);
   end
   y = L \ f;
   u = U \ y;
end
x = xx;
return