I am modeling the temperature of the groundwater using heat equation. I have Dirichlet BC at the top but at the bottom I have constant temperature equal 12 degrees C (see attached pic). It is look like at the bottom we have zero flux. When I specify the Dirichlet BC I have a "knee" right before "neutral zone" which is not the case. So far I have 2 ideas how to implement it:

  1. I was thinking maybe the Robyn BC will help here?
  2. Maybe I need to implement somehow the [-hh * (T - T_inf)/k] as gradient there. When we have not 12C we will have the flux and when we have 12 then we will have 0 gradient.

what do you think? is it correct way of tackling this problem?

attached image)

  • $\begingroup$ This stackexchange isn't really the right place to ask about the correct mathematical statement of the problem. You need to know that before you try computing anything. $\endgroup$ – David Ketcheson Oct 7 '15 at 7:36
  • $\begingroup$ I'm voting to close this question as off-topic because it is about the physics or engineering of the problem, not about its computational solution. $\endgroup$ – David Ketcheson Oct 7 '15 at 7:37

You cannot (in general) both enforce zero flux and constant value at the same boundary.

Additionally, you don't need a Robin boundary condition. The boundary condition depends on what you wish to enforce. At the top, a Dirichlet condition seems natural. You have already suggested two options at the bottom: Dirichlet (constant, 12 degrees as you said) or Neumann (zero flux, as you said).

In this problem, if the bottom boundary condition is a great enough distance from the top there will be only a very small difference between a zero flux Neumann condition and a Dirichlet condition at the mean value of the annual surface temperature. How far is "a great enough distance" depends on the skin depth as mentioned by @GeoMatt22. If you are seeing a bend just above the "neutral zone" with the Dirichlet condition it is probably best to switch to Neumann.


Do you have variable thermal diffusivity with depth?

If not, then if I understand you problem correctly, there is an analytical solution (see here).

In either case, assuming I understand your problem correctly, then typically the lower boundary condition is zero heat-flux "at infinity" (i.e. Neumann). Here "infinity" means at $\frac{z_0-z}{\lambda}\gg1$, where $\lambda$ is the "skin depth" (see the above link). With this condition, the temperature will end up as the time-average of your surface-forcing.

Edit: Your comment suggests you have a piecewise-constant diffusivity, dry soil vs. saturated soil (ignoring the vadose zone). Across a step change in diffusivity, the temperature $T$ and heat flux $q=-KT_z$ will be continuous, but the temperature gradient $Tz$ will be discontinuous. Is the kink in your results at the interface (i.e. water table)?

Numerically, piecewise constant conductivities would be handled in a finite volume method by assigning the interface-conductivity as the harmonic average of the two adjacent cell-conductivities when computing the flux. (These notes give an overview.)

  • $\begingroup$ yes, I have. I have water at the bottom and air at the top of the soil. $\endgroup$ – Igor Markelov Oct 6 '15 at 14:00

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.