I am trying to write a program (based on finite difference method) able to predict the Temperature within the skin of a human being when a laser is pointing toward the skin. The energy associated with the laser decays exponentially with the depth.

I am having troubles to find and/or solve the PDE at the boundary. The surface of the skin is submitted to convection, radiation and a heat flux from the laser energy.

My problem is when I write my energy balance on the first control volume (to solve for the temperature at the surface, the boundary) I end up with a hyperbolic equation (1st order derivative in time, and first order derivative in space). The PDE for the rest of the domain (within the skin) are parabolic as expected.

Should my PDE at the boundary be hyperbolic? If yes, do I use a first order upwind scheme to solve it? If yes, to express the first order difference in space what temperature do I use (T at the surface and the first node within the domain temperature, or the temperature of a node OUT of the domain???).

Thanks

If it is the conduction equation, upwind has no meaning. That has to do with advection (or convection terms) - energy carried by moving fluid. Those terms might be generated if you're considering energy due to blood flow within the skin.

Thanks for your concern,

Let me show you my PDE for the boundary:

Energy stored=conduction (only going out, nothing going in, because at the boundary what is going in is convection)+convection+radiation+energy from laser

so I have something like:

dT/dt=a dT/dx + b (Tinfinity-T)+c

with a, b and c coefficient. My concern is that this end up being a hyperbolic equation cause the conduction is a first order derivative. Could you tell me how I am supposed to deal with the conduction term without being in conflict with the physics of the problem.. Tanks PS : Once again unlike within the domain, at the boundary the conduction term is a first order derivative.

Go back to the elemental energy balance to derive your equation. The time derivitive starts as the rate of change of energy stored in a small volume,

density x volume x c_p x dT/dt = balance on energy crossing boundaries of the volume.

For a flat surface, the volume is S x h, the thickness of the volume. The transfer across S is - cond x S x (T(h) - T_s) / h + b (T_inf - T_s) + c.

Now, when you consider the limit as h -> 0 the right hand side comes out as you've written it. The left hand side however, is

S x h x density x c_p x dT/dt, which goes to 0 in the limit as h -> 0!

That is: in the limit there is no energy storage at the surface and the explicit time dependence vanishes from the equation.

If I'm thinking about this incorrectly, I hope some other folks will jump in.

Thank you very much It is what I did a few weeks ago (neglected the storage energy as the control volume was small). I thought it was wrong because for some reason when I change my space increment, the temperature changes drastically. However including the storage energy on the first control volume doesn t solve the problem (besides, with this term included I have huge time step increment limitation as my size step increment has to be extreamly small).

Maybe my space increment is not small enough which could explain the differences obtained when I change this space increment