Radiative/Convective Boundary Conditions for Heat Equation

[sharpybox]

Hi everyone,
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$\mu)$$T^{t+1}_{i,j,k}$ - μ($T^{t+1}_{i+1,j,k}$ + $T^{t+1}_{i-1,j,k}$ + $T^{t+1}_{i,j+1,k}$+$T^{t+1}_{i,j-1,k}$+$T^{t+1}_{i,j,k+1}$ + $T^{t+1}_{i,j,k-1}$) = (1-6$\mu$)$T^{t}_{i,j,k}$ - μ($T^{t}_{i+1,j,k}$ + $T^{t}_{i-1,j,k}$ + $T^{t}_{i,j+1,k}$+$T^{t}_{i,j-1,k}$+$T^{t}_{i,j,k+1}$ + $T^{t}_{i,j,k-1}$)

Where T represents the temperature field and $\mu$ = $\frac{tα}{2h^{2}}$ (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$\frac{\partial T}{\partial x}$ = hc(T-$T_{a}$) + $\epsilon$$\sigma$($T^{4}$-$T^{4}_{a}$)
(hc = convective heat transfer coefficient, k = thermal conductivity, ε = emissivity, σ = Stefan Boltzmann constant, T = node temperature and $t_{a}$ is the ambient temperature.)

Applying a central difference approximation to the derivative at node $T_{0,j,k}$ yields:
-k$\frac{T^{t}_{1,j,k} - T^{t}_{-1,j,k}}{2h}$ = hc(T-$T_{a}$) + $\epsilon$$\sigma$($T^{4}$-$T^{4}_{a}$)

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

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 $T_{i-1,j,k}$ node for this boundary? Is the method the same when considering the x = maximum boundary when it is the $T_{i+1,j,k}$ 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

 

[hunt_mat]

It is possible to obtain the value of the derivative on the boundary if this is helpful to you.