# Solving 1D convection using method of lines

I'm interested in solving the following 1D-advection equation using method of lines.

$$\frac{\partial C}{\partial t} = - v\frac{\partial C}{\partial x}$$ The spatial domain has been discretized into N nodes.Using forward difference for discretizing the spatial domain gives,

$$\frac{dC}{dt} = - v\frac{C_{i+1} - C_i}{\Delta x}$$

At the Nth node, Neumann boundary condition is applied

$$(\frac{dC}{dt})_N = 0$$

Implies, $$C_{N+1} = C_N$$

However, the following system doesn't converge. I am not sure if the above procedure is correct.

function convection()
global v;
N= 5;
C0 = [5 0 0 0 0];
v = 10000;
delx = 6;
tspan = [0 10];
[t C]  = ode15s(@(t,s) fun(t,s), tspan , C0);
plot(t,C)
function dC= fun(t,C)
for i = 1:nnode-1
dC(i,1) = -v*(C(i+1) - C(i))/delx;
end
dC(nnode,1) = 0;
end
end


I'd like to know whether boundary condition has to be applied on both ends. I have used Neumann BC only at the terminal node.

I think you are using a downwind- instead of an upwind finite difference. This leads to your code imposing a boundary condition where it is not allowed. The solution to your convection equation is basically (ignoring the left BC for the moment)

$$C(x,t) = C_0(x - v t)$$

where $$C_0$$ is your initial value.

Thus, if $$v > 0$$, it is a rightward travelling signal. The value at the right end of your domain, in your case wherever your node $$x_N$$ is located, is completely determined from your initial value and left boundary condition. Imposing a Neumann BC there overdetermines the system and will lead to instability.

Second, your finite difference is biased in the wrong way. For positive $$v$$ (as in your example) it should be

$$\frac{d C}{dt} = -v \frac{C_{i} - C_{i-1}}{\Delta x}$$

Changing this will also fix the issue with your boundary condition at $$x_N$$, since the upwind discretization won't need it.

• Thanks a lot. I could impose Dirichlet BC at the left boundary by setting dC(1,1) = 0. I'd like to confirm whether Robin BC reduces to Dirichlet BC when there is no diffusion. Nov 19 '19 at 2:31
• No, Robin BC is typically a combination of Dirichlet and Neumann BC. That is, something like u(0) = 2 u'(0), for example. I don't think diffusion has any effect here. Nov 19 '19 at 7:41
• Sorry if I wasn't clear. When both convection and diffusion are involved, constant flux can be set at the left boundary using Robin BC, $J_L = vC_0 - DdC/dx$. However, when diffusive effects aren't involved I'm not sure how constant flux can be imposed at the left boundary. , Often constant flux BC is imposed using $J_L = vC_{at ghost node = -1}$ .This implies $vC_{-1} = vC_{0}$ in the absence of diffusion.Therefore, $C_{-1} =C_0$; dC/dt( at i=0) = 0. I am not sure if this is the right way to proceed if one has to implement a conatnt flux BC in the absence of diffusion. Nov 19 '19 at 10:25
• Talking about fluxes suggests to me that you are now considering finite volumes, which are quite different to finite difference. There, you can impose a BC via a flux at the boundary. You will find more on this in lots of numerical maths textbooks. A personal favourite is the one by Dale Durran, but other tastes may vary. Nov 19 '19 at 11:53