[haakon]

Hi,

I am trying to solve a temperature transport problem in a steady flow. The case in question is groundwater flow, where hot water (303 K) is injected through some wells, and water is extracted at some other wells (cold wells).

The flow is governed by the laplace equation (laplacian(p) = 0), and then the volume fluxes q is q = - K & grad(p) (volume fluxes = average velocity). We assume that the water and soil phases is in thermal equilibrium at all times. For details regarding this, please see the USGS Hydrotherm manual: ftp://brrftp.cr.usgs.gov/pub/klkipp/...HT_Doc-3.1.pdf pages 2-2 to 2-4.

My intention was to solve this problem in OpenFOAM, and I created a solver for this purpose with basis in the "scalarTransportFoam" solver. The flow/pressure equation is solved once at the beginning, and then this vector field is used to solve the transient thermal equation:

(
    (porosity*rhow*Cpw + (1-porosity)*rhos*Cps)*fvm::ddt(T) 
  - fvm::laplacian(Gamma, T)
  + rhow*Cpw*fvm::div(phi, T)
);

• porosity - average soil porosity
• rhow - density of water
• rhos - density of soil
• Cpw - specific heat of water
• Cps - specific heat of soil
• Gamma - average thermal condiutivity of both phases combined
• phi - fvc::interpolate(q) & mesh.Sf();

The problem is that more or less independently on how I make my mesh and witch numerical schemes I use, the resulting temperature field diverges. Please see attached screenshot for the field after a while.

The mesh is 3D, and I have tried a lot of different resolutions, adding of prism layers around the wells, making the wells square etc. The attached screenshot is a slice through the middle of the domain. The wells are created by setting a fixed pressure on the well surfaces, and then adjusting this to the desired volume flow. As previously mentioned, this works great. The BC for the temperature is a fixed temperature of 303 kelvin at the inlet wells (left) and zeroGradient at the outlet wells (right).

What I find most strange is the extremely low-temperature regions around the hot wells, almost like a flower...

I have figured out that this must be a problem with the convective heat transport (the diffusion is not very important here). Can any of you help me understand why this fails so miserably?

BTW: I have also tried to set the wells as sources and sinks directly in the equation with SuSp-terms instead of resolving the mesh around the wells, and again, the flow was nice and smooth, while the temperature went bananas...

[Kumudu]

Originally Posted by haakon

Hi,

I am trying to solve a temperature transport problem in a steady flow. The case in question is groundwater flow, where hot water (303 K) is injected through some wells, and water is extracted at some other wells (cold wells).

The flow is governed by the laplace equation (laplacian(p) = 0), and then the volume fluxes q is q = - K & grad(p) (volume fluxes = average velocity). We assume that the water and soil phases is in thermal equilibrium at all times. For details regarding this, please see the USGS Hydrotherm manual: ftp://brrftp.cr.usgs.gov/pub/klkipp/...HT_Doc-3.1.pdf pages 2-2 to 2-4.

My intention was to solve this problem in OpenFOAM, and I created a solver for this purpose with basis in the "scalarTransportFoam" solver. The flow/pressure equation is solved once at the beginning, and then this vector field is used to solve the transient thermal equation:

(
    (porosity*rhow*Cpw + (1-porosity)*rhos*Cps)*fvm::ddt(T) 
  - fvm::laplacian(Gamma, T)
  + rhow*Cpw*fvm::div(phi, T)
);

• porosity - average soil porosity
• rhow - density of water
• rhos - density of soil
• Cpw - specific heat of water
• Cps - specific heat of soil
• Gamma - average thermal condiutivity of both phases combined
• phi - fvc::interpolate(q) & mesh.Sf();

The problem is that more or less independently on how I make my mesh and witch numerical schemes I use, the resulting temperature field diverges. Please see attached screenshot for the field after a while.

The mesh is 3D, and I have tried a lot of different resolutions, adding of prism layers around the wells, making the wells square etc. The attached screenshot is a slice through the middle of the domain. The wells are created by setting a fixed pressure on the well surfaces, and then adjusting this to the desired volume flow. As previously mentioned, this works great. The BC for the temperature is a fixed temperature of 303 kelvin at the inlet wells (left) and zeroGradient at the outlet wells (right).

What I find most strange is the extremely low-temperature regions around the hot wells, almost like a flower...

I have figured out that this must be a problem with the convective heat transport (the diffusion is not very important here). Can any of you help me understand why this fails so miserably?

BTW: I have also tried to set the wells as sources and sinks directly in the equation with SuSp-terms instead of resolving the mesh around the wells, and again, the flow was nice and smooth, while the temperature went bananas...

Hi,

I am interested in your problem. Please tell me how you add these source terms to the laplacian equation. I am really stuck in my thesis , and can you send me the all the files you changed in the LaplacinFoam.

It is really, helpful if you help me in this regards
Kumudu

[haakon]

Unfortunately I cannot help you. I never found out this problem, and I hadn't much time to work on it either, so I gave it up. I might have the files laying around in my archives somewhere, but I do not think that I will share them at the moment as I never finished this project.

When it comes to the modification of solvers and equations, I know that there is a lot of useful resources around on the Internet, and I seriously doubt that I can write anything better here.

Good luck with your thesis.