$$ \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$.

Discretise the boundaries

$$D \frac {\partial C}{dx} - wC = 0 $$$$D \frac {C^{t}_{x+1} - C^{t}_{x-1}}{2 \Delta x} - w C^{t}_{x} = 0$$$$D \frac {C^{t}_{1} - C^{t}_{-1}}{2 \Delta x} - w C^{t}_{0} = 0$$$$ \frac {D}{2 \Delta x}(C^{t}_{1} - C^{t}_{-1}) - w C^{t}_{0} = 0$$$$ \frac {D}{2 \Delta x} C^{t}_{1} - \frac {D}{2 \Delta x} C^{t}_{-1} - w C^{t}_{0} = 0$$$$C^{t}_{-1} = C^{t}_{1} - \frac {2 w \Delta x}{D} C^{t}_{0}$$$$-(\sigma + \rho) \Bigg(C^{t+1}_{1} - \frac {2 w \Delta x}{D} C^{t+1}_{0} \Bigg) + (1 + 2\rho + \mu)C^{t+1}_{0} + (\sigma - \rho)C^{t+1}_{1} = (\sigma + \rho) \Bigg(C^{t}_{1} - \frac {2 w \Delta x}{D} C^{t}_{0} \Bigg) + (1 - 2\rho - \mu)C^{t}_{0} + (\rho - \sigma) C^{t}_{1} + \Delta t P^t_x$$$$-(\sigma + \rho)C^{t+1}_{1} + (\sigma + \rho) \frac {2 w \Delta x}{D} C^{t+1}_{0} + (1 + 2\rho + \mu)C^{t+1}_{0} + (\sigma - \rho)C^{t+1}_{1} = (\sigma + \rho) C^{t}_{1} - (\sigma + \rho) \frac {2 w \Delta x}{D} C^{t}_{0} + (1 - 2\rho - \mu)C^{t}_{0} + (\rho - \sigma) C^{t}_{1} + \Delta t P^t_x$$$$\Bigg((\sigma + \rho) \frac {2 w \Delta x}{D} + (1 + 2\rho + \mu) \Bigg)C^{t+1}_{0} - 2 \rho C^{t+1}_{1} = \Bigg((1 - 2\rho - \mu) - (\sigma + \rho) \frac {2 w \Delta x}{D} \Bigg) C^{t}_{0} + 2 \rho C^{t}_{1} + \Delta t P^t_x$$$$(\sigma + \rho) \frac {2 w \Delta x}{D}$$$$\frac {w \Delta t}{4\Delta x} = \sigma, \frac {D \Delta t}{2 \Delta x^2} = \rho$$$$(\frac {w \Delta t}{4\Delta x} + \frac {D \Delta t}{2 \Delta x^2}) \frac {2 w \Delta x}{D}$$$$\frac {2 w \Delta x w \Delta t}{4\Delta x D} + \frac {2 w \Delta x D \Delta t}{2 \Delta x^2 D}$$$$(\sigma + \rho) \frac {2 w \Delta x}{D} = \frac {w^2 \Delta t}{2 D} + \frac {w \Delta t}{\Delta x}$$$$(\sigma + \rho) \frac {2 w \Delta x}{D} = \frac {w^2 \Delta t}{2 D} + 4 \sigma$$$$\Bigg(\frac {w^2 \Delta t}{2 D} + 4 \sigma + (1 + 2\rho + \mu) \Bigg)C^{t+1}_{0} - 2 \rho C^{t+1}_{1} = \Bigg((1 - 2\rho - \mu) - \frac {w^2 \Delta t}{2 D} - 4 \sigma \Bigg) C^{t}_{0} + 2 \rho 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$$