Minus pressure in THM process

Hi all,
I have some problems in a THM process simulation. I use a 2D model to simulate the injection process in EGS. However, I found from the first result vtu file that the pressure distribution contains a large area of large minus values like the second picture. In fact, I set the distribution of pressure as an approximate hydrostatic condition initially, like the first picture. I attach my prj file here. Can anyone tell me what’s the reason for my abnormal inversion pressure result?

THM_test_complex.zip (3.3 KB)

From the second figure it looks like you have 0 Dirichlet boundary conditions for pressure all around the domain. Try to specify boundary conditions which are compatible with the initial conditions.

Hi, @dmitri.naumov
I’m sorry that I made a mistake in my post. the first picture shows the initial condition and the second one is the first result vtu file. I have corrected the mistake now.
I add Dirichlet boundary conditions at the four boundary lines in my model, but the results are just like the abnormal pictures…

I understood it correctly, the first figure is for initial conditions, second the first time step.
I would think the Dirichlet 0 for pressure on all sides of the domain is not correct, it should be pressure gradient on left and right, 7e7 on the bottom, and some reasonable value on the top, depending on what model you are trying to solve.

the two result figures are just the results with the boundary conditions set as you have described in the reply… when I first run the prj file I have just set them. For the top set the atmospheric pressure 1.013e5pa, for the bottom is 7e7pa, for the left and right sides it is gradient. You can see the details in the two following figures. However, it seems no help. In my model the coordinate Y is minus because the real area is underground, so you can see that my gradient pressure value expression is 1.013e5-1e4*y. Maybe there is something here? I’m not sure.

Your gravity is showing upwards. Should be 0 -9.81

Hi, I have corrected the gravity, but the result doesn’t change…

do you use an initial effective stress that matches the initial pressure distribution, the boundary conditions and gravity?

EDIT: otherwise your initial state is not in equilibrium and you will get sudden deformation which messes up your pore pressure field due to the HM coupling.

Hi, @Thomas_Nagel
The initial pressure distribution is obtained by a HT process initialization that I have run using ogs6 before…so the pressure distribution is not so regular that perhaps I can’t figure out the correct initial stress… Could you give me some advice?

There is an option to compensate intial non-equilibrium forces on the residual level. We have an improvement in the pipeline that is more suitable for your problem. @wenqing.wang can you advise, pls?

You can try to use compensate_non_equilibrium_initial_residuum as

            <process ref="THERMO_HYDRO_MECHANICS">

However, the variable distribution figure show that only pressure has problem. You can try to apply the same pressure conditions on the top, left and right boundaries as that for HT for the initial condition.

I see. I will try it later. Thank you!

Thank you! I will try @wenqing.wang 's advice.

Sorry to reply to you after such a long time for I had been busy with my graduation the past two months. I have tried @wenqing.wang 's code to compensate initial non-equilibrium forces but the results still have some problems. First, I cancel the pressure source term and nodal temperature boundary condition, where no additional forces are added during the whole process. with the compensate initial non-equilibrium code the variables stay constant during the simulation time. Then I added only the source term to run the project, it is confusing that the variables are still constant. Next I add both source term and nodal temperature boundary condition. Besides, following wang’s advice I add pressure boundary condition for four lines.the displacement seems well but the pressure,simga and temperature are still unreasonable. I put some distributions of these variables here. Could you please help me to figure out the problems?


source terms and Neumann conditions enter each flux balance and are indeed “nulled” when applying the compensation option. The easiest is that you assign a curve to each of them so that in time step 0 they are still 0, and are activated only after the first time step. This way, your source terms and Neumann boundary conditions are not “compensated away” and everything should work. Switching to Dirichlet conditions (fixed T or p) is physically very different, so you’d be solving a different problem.

Thanks for your advice and I have solved the problem. When I check the variables’ distribution, I notice some errors, too. First, the pressure around the injection and production well seems too large. I think it is due to the pressure source term values. For it’s a 2D model, I assume an injection distributed over a 300m section of open well, so I set the inject rate as 3.44e-5m³/s.is that right? Or I want to know what’s the standard to set an inject rate for a 2D model?
Second, the first result file shows that the compressive stress distribution has a sudden increase during the first time step, which I think is not reasonable because I haven’t add source term and nodal temperature boundary condition in this step. The first figure shows the distribution of sigma_xx in initial time and the second shows that in first time step. Could you explain to me that how this happens and what to do to solve it?

I canceled the source term and nodal temperature boundary condition in the whole process to run the project again, the sudden change of the sigma happened again in the first time step and it remained unchanged later. Doesn’t the compensation code consider the stress compensation?