MATLAB: Facing problem to solve convection-diffusion equation

convection diffusion equationcrank-nicolson methodfinite difference method

I want to solve simple diffusion equation MATLAB.
 (del_C/del_t)=(D^2)(del^2_C/del_x^2)-U(del_C/del_x)
I.C: C(x,0)=C0; B.C: C(0,t)=C0, C(L,t)=C1
My program:
numr =50;         %number of grid points in x direction
numt = 10000;       %number of time steps 
L = 0.01;           %thickness of a flat plate (m)
dx = R/(numr - 1);  %grid size in x direction
dt = 0.0005;
D = 1.0e-05;        %species diffusivity 
C0 = 1.0;           %species initial concentration
Cinf = 0.1;         %species surface concentration 
x = 0:dx:L;         %vector of x positions
C = zeros(numr,numt);  %initialize the concentration C(x,t) matrix
%specify initial conditions along the radius(C=C0,x=0 to L,t=0)
t(1) = 0;                 
for i=1:numr
  C(i,1)=C0; 
end
% uniform boundary condition at the centre(C=C0,x=0,t=0 to tf)constant mass
% source term
 for j = 1:numt
     C(1,j)=C0; 
 end
% boundary condition at particle surface (C=Cinf,x=L,t>0)
% specis concentration is fixed at the boundary
for j=2:numt
  C(numr,j)=Cinf; 
end
%Stability criteria
%dx<(2*D)/u
%dt<(dx^2)/(2*D)
u=0.5;
s=D*dt/(2*dx^2);
p=u*dt/(4*dx);
%iteration based on forward differencing scheme 
for j=2:numt
   t(j) = t(j-1) + dt;
   for i=2:numr-1
C(i,j+1)=C(i,j)+(s-p)*(C(i+1,j+1)+C(i+1,j))-(2.0*s)*(C(i,j+1)+C(i,j))+(s+p)*(C(i-1,j+1)+C(i-1,j));
end
end
plot(C,x)
xlabel('x direction')
ylabel('C')

Best Answer

Your line
C(i,j+1)=C(i,j)+(((s-p)*C(i+1,j+1)+(s-p)*C(i+1,j)+((-2*s)*C(i,j+1))+((-2*s)*C(i,j))+((s+p)*C(i-1,j+1))+((s+p)*C(i-1,j));
is missing two closing brackets somewhere
C(i,j+1)=C(i,j)+(((s-p)*C(i+1,j+1)+(s-p)*C(i+1,j)+((-2*s)*C(i,j+1))+((-2*s)*C(i,j))+((s+p)*C(i-1,j+1))+((s+p)*C(i-1,j));
 1     0  1   0 123   2  3       2 3   2  3     2 34    3  4     32 34    3  4   32 34   3  4       32 34   3  4     32
so at the end of the expression you have 2 unclosed brackets.
Related Question