I looked at a similar question to find out how to discretize the heat equation in 2 spacial dimensions at the boundary between two materials. I am required to use explicit method (forward-time-centered-space) to solve. I believe I have arrived at a reasonable approximation with the following:

\begin{align}\frac{T^{k+1}_{i,j}-T^k_{i,j}}{\Delta t}&=\frac1{h^2}\left[\left(\alpha^k_{i+\frac12,j}-\alpha^k_{i-\frac12,j}\right)\left(T^k_{i+\frac12,j}-T^k_{i-\frac12,j}\right)+\left(\alpha^k_{i,j+\frac12}-\alpha^k_{i,j-\frac12}\right)\left(T^k_{i,j+\frac12}-T^k_{i,j-\frac12}\right)\right.\\ &\quad\left.+\alpha^k_{i,j}\left(T^k_{i-1,j}-2T^k_{i,j}+T^k_{i+1,j}\right)+\alpha^k_{i,j}\left(T^k_{i,j-1}-2T^k_{i,j}+T^k_{i,j+1}\right)\right] \end{align}

where any index $\pm\frac12$ represents the average of the value at current index and the one right in front of or behind, i.e.: $$\frac12(T_{i\pm1}+T_{i}).$$

$k$ is the time index, $i,j$ are the spacial indices, $h$ is a grid size, and $\Delta t$ is a time step.

This was determined through the generalized heat equation in $x$ and $y$ dimensions $$\frac{\partial T}{\partial t}=\nabla\cdot\left(\alpha\nabla T\right).$$

I have two questions:

  1. Is the discretized form of the heat equation above correct? If not, how do I fix it?
  2. Is it possible to use a nondimensional time step $\tau = \frac{\alpha\Delta t}{h^2}$ for this problem even though the $\alpha$ is changing across materials? If the nondimensional time step is changing, how would I then extract the actual solution in space and time after it is finished solving?
  • $\begingroup$ For #2, I'd suspect that you'd use the maximum, as it's done in the CFL condition. $\endgroup$ – Kyle Kanos Apr 13 '18 at 1:22
  • $\begingroup$ Should I then use a constant time step $\tau$ ? And then scale it by whichever $\alpha$ range I'm in? $\endgroup$ – AAC Apr 13 '18 at 2:58
  • $\begingroup$ At the discontinuity between the materials, you have $$k_1\left(\frac{\partial T}{\partial x}\right)_1=k_2\left(\frac{\partial T}{\partial x}\right)_2$$where k is the thermal conductivity, and the subscripts refer to the two sides of the interface. $\endgroup$ – Chet Miller Apr 13 '18 at 14:56

I discuss in that answer that the FDM for the 1D diffusion equation with variable coefficient ends up being, $$ \frac{T^{n+1}_i-T^n_i}{dt}=\frac{1}{dx^2}\left[\alpha_{i+\frac{1}{2}}\left(T^n_{i+1}-T^n_{i}\right)-\alpha_{i-\frac{1}{2}}\left(T^n_i-T^n_{i-1}\right)\right] $$ with $\alpha_{i+\frac12}=\left(\alpha_{i+1}+\alpha_{i}\right)/2$. Extending this 2D should be direct: just add $j$ subscript to all terms, then add a similar $\left[\cdots\right]/dy^2$ term. So it seems you've got a few extra terms, at least at an initial, cursory look. I also discuss the Crank-Nicolson version in this answer.

For the timestep, $dt$, you need to satisfy the Courant-Freidrich-Lewy condition (CFL) for diffusion equations: $$ CFL\leq\max\left(\alpha\right)\cdot\frac{dt}{\left(dx^2+dy^2\right)} $$ where, typically, $CFL=1/2$. So there will be one time step and its value depends on only the maximum value of the (variable) diffusion coefficient.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.