$$ \frac {\partial C}{dt} = D \frac {\partial^2 C}{dx^2} - w \frac {\partial C}{dx} - \lambda C + P $$

Discretise the equation

$$\frac {C^{t+1}_{x} - C^{t}_{x}}{\Delta t} = \frac {D}{2} \Bigg(\frac {C^{t+1}_{x+1} - 2C^{t+1}_{x} + C^{t+1}_{x-1}}{\Delta x^2} + \frac {C^{t}_{x+1} - 2C^{t}_{x} + C^{t}_{x-1}}{\Delta x^2} \Bigg) -\frac {w}{2} \Bigg( \frac {C^{t+1}_{x+1} - C^{t+1}_{x-1}}{2 \Delta x} + \frac {C^{t}_{x+1} - C^{t}_{x-1}}{2 \Delta x} \Bigg) - \frac{\lambda}{2} \Big( C^{t+1}_{x} + C^t_{x} \Big) + P^t_x$$$$C^{t+1}_{x} = C^{t}_{x} + \frac {D \Delta t}{2 \Delta x^2} \Bigg( C^{t+1}_{x+1} - 2C^{t+1}_{x} + C^{t+1}_{x-1} + C^{t}_{x+1} - 2C^{t}_{x} + C^{t}_{x-1} \Bigg) - \frac {w \Delta t}{4\Delta x} \Bigg(C^{t+1}_{x+1} - C^{t+1}_{x-1} + C^{t}_{x+1} - C^{t}_{x-1} \Bigg) - \frac{\lambda \Delta t}{2} \Big( C^{t+1}_{x} + C^t_{x} \Big) + \Delta t P^t_x$$

if $\frac {w \Delta t}{4\Delta x} = \sigma$, $\frac {D \Delta t}{2 \Delta x^2} = \rho$ and $\frac{\lambda \Delta t}{2} = \mu$

$$C^{t+1}_{x} = C^{t}_{x} + \rho C^{t+1}_{x+1} - 2\rho C^{t+1}_{x} + \rho C^{t+1}_{x-1} + \rho C^{t}_{x+1} - 2\rho C^{t}_{x} + \rho C^{t}_{x-1} - \sigma C^{t+1}_{x+1} + \sigma C^{t+1}_{x-1} - \sigma C^{t}_{x+1} + \sigma C^{t}_{x-1} - \mu C^{t+1}_{x} - \mu C^t_{x} + \Delta t P^t_x$$

and rearranging,

$$-(\sigma + \rho)C^{t+1}_{x-1} + (1 + 2\rho + \mu)C^{t+1}_{x} + (\sigma - \rho)C^{t+1}_{x+1} = (\sigma + \rho) C^{t}_{x-1} + (1 - 2\rho - \mu)C^{t}_{x} + (\rho - \sigma) C^{t}_{x+1} + \Delta t P^t_x$$

gives the discretised equation for $x > 1, x < L$.

$$-(\sigma + \rho)C^{t+1}_{-1} + (1 + 2\rho + \mu)C^{t+1}_{0} + (\sigma - \rho)C^{t+1}_{1} = (\sigma + \rho) C^{t}_{-1} + (1 - 2\rho - \mu)C^{t}_{0} + (\rho - \sigma) C^{t}_{1} + \Delta t P^t_x$$$$-(\sigma + \rho)C^{t+1}_{0} + (1 + 2\rho + \mu)C^{t+1}_{0} + (\sigma - \rho)C^{t+1}_{1} = (\sigma + \rho)C^{t}_{0} + (1 - 2\rho - \mu)C^{t}_{0} + (\rho - \sigma) C^{t}_{1} + \Delta t P^t_x$$$$(1 - \sigma + \rho + \mu)C^{t+1}_{0} + (\sigma - \rho)C^{t+1}_{1} = (1 + \sigma - \rho - \mu)C^{t}_{0} + (\rho - \sigma) C^{t}_{1} + \Delta t P^t_x$$$$-(\sigma + \rho)C^{t+1}_{L-2} + (1 + 2\rho + \mu)C^{t+1}_{L-1} + (\sigma - \rho)C^{t+1}_{L} = (\sigma + \rho) C^{t}_{L-2} + (1 - 2\rho - \mu)C^{t}_{L-1} + (\rho - \sigma) C^{t}_{L} + \Delta t P^t_x$$$$-(\sigma + \rho)C^{t+1}_{L-2} + (1 + 2\rho + \mu)C^{t+1}_{L-1} + (\sigma - \rho)C^{t+1}_{L-1} = (\sigma + \rho)C^{t}_{L-2} + (1 - 2\rho - \mu)C^{t}_{L-1} + (\rho - \sigma) C^{t}_{L-1} + \Delta t P^t_x$$$$-(\sigma + \rho)C^{t+1}_{L-2} + (1 + \sigma + \rho +\mu)C^{t+1}_{L-1} = (\sigma + \rho)C^{t}_{L-2} + (1 - \sigma - \rho - \mu)C^{t}_{L-1} + \Delta t P^t_x$$

In [ ]: