Hi everyone,(adsbygoogle = window.adsbygoogle || []).push({});

I'm attempting to create a computer program to solve the transient 3d heat equation using the Crank Nicolson method.

I would like to model the boundaries of my domain as losing heat via convection and radiation due to the temperature difference between the boundary and the air in which the system I am modelling resides, but would like to check i have the correct method for incorporating these modes of heat transfer into my model.

At the moment for each internal node I have the following finite difference scheme:

(1+6[itex]\mu)[/itex][itex]T^{t+1}_{i,j,k}[/itex] - μ([itex]T^{t+1}_{i+1,j,k}[/itex] + [itex]T^{t+1}_{i-1,j,k}[/itex] + [itex]T^{t+1}_{i,j+1,k}[/itex]+[itex]T^{t+1}_{i,j-1,k}[/itex]+[itex]T^{t+1}_{i,j,k+1}[/itex] + [itex]T^{t+1}_{i,j,k-1}[/itex]) = (1-6[itex]\mu[/itex])[itex]T^{t}_{i,j,k}[/itex] - μ([itex]T^{t}_{i+1,j,k}[/itex] + [itex]T^{t}_{i-1,j,k}[/itex] + [itex]T^{t}_{i,j+1,k}[/itex]+[itex]T^{t}_{i,j-1,k}[/itex]+[itex]T^{t}_{i,j,k+1}[/itex] + [itex]T^{t}_{i,j,k-1}[/itex])

Where T represents the temperature field and [itex]\mu[/itex] = [itex]\frac{tα}{2h^{2}}[/itex] (t = time step, α = thermal diffusivity and h = step size in x/y/z dimensions).

As I understand it this type of boundary condition is a Neumann condition and can be represented by (assuming a 1d case along the x=0 boundary):

-k[itex]\frac{\partial T}{\partial x}[/itex] = hc(T-[itex]T_{a}[/itex]) + [itex]\epsilon[/itex][itex]\sigma[/itex]([itex]T^{4}[/itex]-[itex]T^{4}_{a}[/itex])

(hc = convective heat transfer coefficient, k = thermal conductivity, ε = emissivity, σ = Stefan Boltzmann constant, T = node temperature and [itex]t_{a}[/itex] is the ambient temperature.)

Applying a central difference approximation to the derivative at node [itex]T_{0,j,k}[/itex] yields:

-k[itex]\frac{T^{t}_{1,j,k} - T^{t}_{-1,j,k}}{2h}[/itex] = hc(T-[itex]T_{a}[/itex]) + [itex]\epsilon[/itex][itex]\sigma[/itex]([itex]T^{4}[/itex]-[itex]T^{4}_{a}[/itex])

[itex]T^{t}_{-1,j,k}[/itex] = [itex]\frac{2h(hc * (T^{t}_{0,j,k} - T_{a}) + \epsilon \sigma (T^{4t}_{0,j,k} - T^{4}_{a})}{-k}[/itex]

Am I correct in thinking that the statement above is then substituted into the original Crank Nicolson FD scheme quoted earlier in place of the [itex]T_{i-1,j,k}[/itex] node for this boundary? Is the method the same when considering the x = maximum boundary when it is the [itex]T_{i+1,j,k}[/itex] node that needs replacing?

And finally is it necessary to repeat the process at the t+1th time step as well as the time t?

Thanks for your help

**Physics Forums - The Fusion of Science and Community**

The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

# Radiative/Convective Boundary Conditions for Heat Equation

Loading...

Similar Threads - Radiative Convective Boundary | Date |
---|---|

A Pseudospectral method using FFT to solve diff equation | Sep 20, 2016 |

Boundary conditions for heat transfer in the pipe | Sep 2, 2015 |

Stability condition for solving convection equation by FDM | May 25, 2015 |

Differentiating Planck's radiation law | Dec 14, 2010 |

Planck's law of radiation? | Apr 11, 2004 |

**Physics Forums - The Fusion of Science and Community**