How to Construct a Matrix for Diffusion PDE in MATLAB?

henxan
Messages
46
Reaction score
2
I am going to do a numerical simulation of diffusion in matlab. The diffusion coefficient is concentration dependant, and i use an array operation to calculate D(x), so it is known.

Based on Fick's second equation:

$$
\frac{\partial C}{\partial t} = \frac{\partial}{\partial x} D \frac{\partial C}{\partial x}
$$

I am going to use the Forward time centered scheme for a variable D:

$$
u_i^{n+1} = u_i^{n} + \frac{\Delta t}{2\Delta x^2}\left( \left( D_{i+1}^n + D_{i}^n \right) \left(u_{i+1}^n -u_i^n \right) - \left( D_{i}^n + D_{i-1}^n \right) \left(u_{i}^n -u_{i-1}^n \right)\right)
$$

How can i construct a matrix which I can put into matlab? I can do the following with the rightmost expression:

$$
\left( \left( D_{i+1}^n + D_{i}^n \right) \left(u_{i+1}^n -u_i^n \right) - \left( D_{i}^n + D_{i-1}^n \right) \left(u_{i}^n -u_{i-1}^n \right)\right) = \begin{pmatrix} D_{i-1}^{n} & D_{i}^{n} & D_{i+1}^{n} \end{pmatrix} \begin{pmatrix} 1 & -1 & 0\\ 1 & -2 & 1\\ 0 & -1 & 1 \end{pmatrix} \begin{pmatrix} u_{i-1}^{n} \\ u_{i}^{n} \\ u_{i+1}^{n} \end{pmatrix}
$$

But how can I apply this to the entire row n+1 and include boundary conditions?
 
Last edited:
Physics news on Phys.org
What I maybe should have written is, "How do I expand the following matrix expressions to calculate:
<br /> \Delta u_i^{n+1} = \begin{pmatrix} D_{i-1}^{n} &amp; D_{i}^{n} &amp; D_{i+1}^{n} \end{pmatrix} \begin{pmatrix} 1 &amp; -1 &amp; 0\\ 1 &amp; -2 &amp; 1\\ 0 &amp; -1 &amp; 1 \end{pmatrix} \begin{pmatrix} u_{i-1}^{n} \\ u_{i}^{n} \\ u_{i+1}^{n} \end{pmatrix}<br />
for all i \in [2, .. N-1] where ##N## is the size of the array ##u \& D##"
 
Thread 'Need help understanding this figure on energy levels'
This figure is from "Introduction to Quantum Mechanics" by Griffiths (3rd edition). It is available to download. It is from page 142. I am hoping the usual people on this site will give me a hand understanding what is going on in the figure. After the equation (4.50) it says "It is customary to introduce the principal quantum number, ##n##, which simply orders the allowed energies, starting with 1 for the ground state. (see the figure)" I still don't understand the figure :( Here is...
Thread 'Understanding how to "tack on" the time wiggle factor'
The last problem I posted on QM made it into advanced homework help, that is why I am putting it here. I am sorry for any hassle imposed on the moderators by myself. Part (a) is quite easy. We get $$\sigma_1 = 2\lambda, \mathbf{v}_1 = \begin{pmatrix} 0 \\ 0 \\ 1 \end{pmatrix} \sigma_2 = \lambda, \mathbf{v}_2 = \begin{pmatrix} 1/\sqrt{2} \\ 1/\sqrt{2} \\ 0 \end{pmatrix} \sigma_3 = -\lambda, \mathbf{v}_3 = \begin{pmatrix} 1/\sqrt{2} \\ -1/\sqrt{2} \\ 0 \end{pmatrix} $$ There are two ways...
Back
Top