0

Steady-state simulation of a low-type geothermal system

T1.sim

Dear TOUGH2 Users

Good morning/afternoon

This is my first post to this forum and I’d like to receive guidance and help to perform a steady-state natural simulation model of a low-type geothermal system.

Recently I am using TOUGH2 to simulate the steady-state condition of a low-type geothermal system by understanding the temperature and pressure distribution beneath the study area.

This study area contains in total five geological layers with different physical characteristics, plus “Air” layer plus a “heat source” layer and a “fault” layer. All these layers/materials have physical characteristics mentioned in this Table:

 

Materials

Geology

Density (kg/m3)

Porosity (%)

Permeability (m2)

Wet Heat Conductivity (W/m K)

Specific Heat (J/kg K)

Material 1

Quaternary rocks

1650

0.1

1.0E-15

2.5

500

Material 2

Upper Tertiary rocks

2500

0.3

1.0E-13

3.44

500

Material 3

Lower Tertiary rocks

2700

0.3

1.0E-13

3.44

500

Material 4

Upper Cretaceous rocks

2600

0.3

1.0E-13

3.44

500

Material 5

Basement rocks

2800

0.3

1.0E-13

3.60

500

Material 6

Fault

2400

0.3

1.0E-13

3.50

500

Material 7

Heat source

2800

0.3

1.0E-13

3.80

500

Material 8

Air

1225

0.003

1.0E-13

2.50

500

 

Our questions are:

1)      In this model, I want to create initial conditions using a hydrostatic pressure of 1.013 105 Pa and a geothermal gradient of 3 oC/100m for all the model except “Air” layer and an average air temperature of 28 oC. I have difficulties in including such information in the initial conditions (Function). The model is running from top to -8000m.

2)      We included a fault in this model. We increased the vertical (Z) permeability to permit the vertical flow of hot waters. However we dont see any upflow in this model, may be due to initial conditions or other parameters.

3)      We included a heat source at the bottom of the model, at some specific blocks.

4)      We run the model for 1 million years but couldn’t get good results (the temperature is low at the botton and high at the top of the model). May be problem due to initial conditions, or boundary conditions.

I attach here the TOUGH2 input file.

Thank you,

Saibi

19 replies

null
    • Staff Scientist
    • Christine_Doughty
    • 4 yrs ago
    • Reported - view

    Dear Saibi,

    I downloaded your zip file and extracted what should be the input file, but when I tried to read it I found it is not an ASCII (i.e. plain text) file, which is what a TOUGH2 input file should be.  Are you using TOUGH2 through an interface that helps you create the input file?  If you can post the actual input file, I will look at it.

    Christine

      • Hakim
      • 4 yrs ago
      • Reported - view

      Christine Doughty Dear Christine, thank you very much for your reply and help. Actually I am using Petrasim, which has TOUGH2. I attach also an Excel file (Al Ain). I use to run my models in the past (15 years ago) using Hydrotherm, then I moved to TOUGH2. Recently, I purchased Petrasim and try to develop same model in Petrasim (TOUGH2). Thank you again, Regards, Hakim

      AL AIN
    • Staff Scientist
    • Christine_Doughty
    • 4 yrs ago
    • Reported - view

    Dear Hakim,

    I looked at your Excel file, and as far as I can tell, it all looks reasonable.  However, to have a chance of helping you, I need to look at the actual TOUGH input file.  I think Petrasim has an option for you to save your TOUGH input file as an ASCII file or possibly as several ASCII files.  See if you can do this.  Christine

      • Hakim
      • 4 yrs ago
      • Reported - view

      Christine Doughty Dear Christine, thank you a lot for your help. I checked in Petrasim and it accepts to save under .sim file type only. I attached a screen shot of my folder from Petrasim. I attach also T1.dat

      files T1.dat
    • Keurfon_Luu
    • 4 yrs ago
    • Reported - view

    There is something off with your initial pressure (see attached figure). Your initial temperature does not seem to follow the temperature gradient you wanted (28 - 0.03 * z which means 268 at the bottom of the model).

    I never used Petrasim so I can't tell how to solve your issue in Petrasim.

      • Hakim
      • 4 yrs ago
      • Reported - view

      Keurfon Luu Thank you very much. Can you tell me which software/code did you use?

      • Keurfon_Luu
      • 4 yrs ago
      • Reported - view

      Hakim Hi Hakim,

      I am using a code of mine (https://github.com/keurfonluu/toughio).

      I simply read your input file and then exported it to an unstructured VTK file and visualized it with ParaView (https://www.paraview.org/).

      Here is the script to generate such file from your input file: 

      import numpy
import toughio

parameters = toughio.read_input("T1.dat")
points = numpy.array([v["center"] for v in parameters["elements"].values()])
incons = numpy.array([v["values"] for v in parameters["initial_conditions"].values()])

mesh = toughio.meshmaker.triangulate(points)
mesh.add_point_data("pressure", incons[:, 0])
mesh.add_point_data("temperature", incons[:, 2])
mesh.write("T1.vtu")
    • Staff Scientist
    • Christine_Doughty
    • 4 yrs ago
    • Reported - view

    Dear Hakim,

    The file you attached T1.dat is the correct format for a TOUGH input file.  However, since I haven't used Petrasim for a really long time, I cannot be sure this is the input file that Petrasim is using, but let's assume that it is.  I see your model has 24 layers, with 900 grid blocks per layer.  It is built from the bottom up, so the first elements shown in the ELEME block are the deepest and the last ones represent the air.  I agree with Keurfon that there is something wrong with your initial conditions (INCON block).  For a normal hydrostatic pressure gradient, I would expect the pressure at the bottom of the model (8000 m depth) to be about 800 bars, but yours is 97 bars.  Actually, I am not sure TOUGH's normal equation of states go to pressures as high as 800 bars.  What EOS are you using?  Also, as Keurfon pointed out, your deep temperatures are far too low to represent 3C per 100 m.  The properties in your ROCKS block do not correspond to those shown in the table you posted to the forum.  In fact most permeabilities are so small (1E-18 to 1E-19 m2) that I can see why you might get no fluid flow.  Please check the units that Petrasim expects for permeabilities, to see if a misunderstanding is causing your rock permeabilities to be so low.  On the other hand, the air permeability is very big (1.25E-6 m2).  When you specify an air boundary condition, you should realize that TOUGH is not set up to actually model what is going on in the air - you are just holding T and P fixed there.  I don't think you should specify a permeability any bigger than 1E-12 m2 for air.  Also, the volumes of your air grid blocks are not really big enough to act like a constant T, constant P boundary.  They are only about  10 to 100 times bigger than other grid blocks in the problem, whereas they should be 1E20 to 1E30 times bigger. 

    One other thing, your input file (first line of the PARAM block) asks the code to take 9999 steps and create a complete printout of all gridblock conditions every step.  This is probably not necessary and is making your output file way too big.  I recommend asking for a printout every 9999 steps (so you will just get one if you run all 9999 steps) and using the TIMES block to get printouts at a few specified times.  Also, set MOP(1)=1 and MOP(7)=1 to get more helpful information about your simulation.

    I know TOUGH input parameters get set by Petrasim based on how you fill out a form, so it may not be easy for you to adjust them, but I think the Petrasim User's Guide does provide helpful information.  Perhaps you could contact the folks at Petrasim to make sure the interface is working the way you expect it to.  Or maybe there are forum readers who are experienced in Petrasim and are willing to help.

    Christine

      • Hakim
      • 4 yrs ago
      • Reported - view

      Christine Doughty Dear Christine, thank you so much for your valuable comments. We used EOS1. I agree that there is a problem in the input of initial conditions. I am a new user of Petrasim. Can you suggest other codes to perform such steady-simulations for T and P ?

    • Staff Scientist
    • Christine_Doughty
    • 4 yrs ago
    • Reported - view

    Dear Hakim,

    I ran a sample problem injecting into a closed volume with EOS1 and the pressure was able to go way above 800 bars without any obvious problems, so your hydrostatic profile for 8000 m depth should be okay. 

    I glanced through the Petrasim User's Manual, and I think one way for you to create your initial conditions would be to specify them by layer.  This will be kind of tedious as you have 24 layers, but if you are patient I think it will work.  I recommend you use Excel or something similar to calculate columns of z, P(z), and T(z) according to the hydrostatic and geothermal gradients you want.  Pick the z values for each row to correspond to the main part of your layer depths, ignoring the little region where the layers are not flat.  You might want to have a column that gives the name of the layer for each row too, to make it easy keep track of where you are for entering them in Petrasim.  Then have Excel open in one window and Petrasim open in another window, and in Petrasim go to Edit Properties and Specify by Layer and then copy the P and T values from Excel into Petrasim one layer at a time.  (See Petrasim User's Manual, Ch. 9 Boundary and Initial Conditions, sub-section Layer and Regional Initial Conditions).  

    When you run TOUGH, if you have elements with a huge volume (like the top layer), the P and T values you entered will not change - they become a fixed boundary condition.  You might want to hold the temperature at part of the bottom layer fixed also (at the locations where there is not a heat source).  To keep temperature fixed while not keeping pressure fixed, have a normal volume and a huge rock density.  Then it will be a closed boundary for fluid flow and a constant-temperature boundary.  But remember to not do this for your heat source elements, they must have normal volume and normal rock density.

    Does this sound reasonable?  When I run TOUGH directly (without the Petrasim interface), I usually write a utility program to create initial conditions cell by cell, where I read the z value of each cell from the ELEME block, and calculate T=T0 + gradT*z and P=P0 + gradP*z, then write an INCON block.   Since your grid layers are mostly flat, what I suggest you do by layer in Petrasim is essentially equivalent to this.  And as TOUGH runs, it will adjust the temperatures for the non-flat part of the grid.

    If any other Petrasim users have a better idea, I would be happy to hear it.

    Christine

      • Keurfon_Luu
      • 4 yrs ago
      • Reported - view

      Christine Doughty Also, for Dirichlet boundary conditions, nodal distances should be set to a really small value.

      Here is a Python script to modify your input file with Christine's suggestions: 

      import copy
import toughio

# Read input file
parameters = toughio.read_input("T1.dat")

# Loop over elements
air = set()
for label, element in parameters["elements"].items():
    # Get coordinates of current element
    x, y, z = element["center"]

    # Hydrostatic pressure
    parameters["initial_conditions"][label]["values"][0] = max(1.013e5, 1.013e5 - 9810.0 * (z + 1700.0))

    # Vertical thermal gradient
    parameters["initial_conditions"][label]["values"][2] = max(28.0, 28.0 - 0.03 * (z + 1700.0))

    # Fix pressure and temperature at top
    if element["material"] == "MAT0":
        element["volume"] *= 1.0e50
        air.add(label)

    # Fix temperature at bottom
    if element["material"] == "MAT 7" and z < -7800.0:
        element["material"] = "MAT8"

# Create new material at bottom with large density
parameters["rocks"]["MAT8"] = copy.deepcopy(parameters["rocks"]["MAT 7"])
parameters["rocks"]["MAT8"]["density"] = 1.0e4

# Loop over connections
for label, connection in parameters["connections"].items():
    # Set nodal distance to 1.0e-9 if element is air
    if label[:5] in air:
        connection["nodal_distances"][0] = 1.0e-9
    elif label[5:] in air:
        connection["nodal_distances"][1] = 1.0e-9

# Write modified input file
toughio.write_input("T1_modified.dat", parameters)

      You may find the modified input file below.

      T1_modified.dat
      • Hakim
      • 4 yrs ago
      • Reported - view

      Keurfon Luu Dear Keurfon Luu, Thank you so much for your code and helpful comments. Hakim

      • Hakim
      • 4 yrs ago
      • Reported - view

      Keurfon Luu Dear Keurfon, hope you are doing well during this covid19 situation. We successfuly installed Paraview and toughio. Now, is it possible to run the simulation in your code (toughio) or its just for showing results? in this case, I need to run it in tough2/petrasim. Please help.

      • Keurfon_Luu
      • 4 yrs ago
      • Reported - view

      Hakim toughio is just a Python pre- and post-processor for TOUGH. It can read/write both TOUGH input and output files. I guess that you can prepare the mesh and TOUGH input files with Petrasim, modify the input files with toughio if necessary, and run your simulation either directly using TOUGH or through Petrasim. For post-processing, Petrasim is probably easier to use.

      • Hakim
      • 4 yrs ago
      • Reported - view

      Keurfon Luu Thanks for your rapid reply. if you know how to open the modified file (T1_modified.dat) in Petrasim, please let me know.

      • JuDi
      • 4 yrs ago
      • Reported - view

      Hakim , I am not aware of a way to import a modified *.dat file into Petrasim. But you should be able to run the simulation by running Tough2 executables (you receive with Petrasim) via command line. Here are instructions on how to do that. I hope they are still valid.

      Julia

    • Staff Scientist
    • Christine_Doughty
    • 4 yrs ago
    • Reported - view

    Dear Hakim, 

    I might be confused, but it looks like you are using the two-waters option for EOS1, since your MULTI block is 2 3 2 6.  According to the User's Guide, INCON is expecting variables in the order (P, T, X2), where X2 is the mass fraction of the second water, but the INCON block you have is in the order (P, 0, T).  Could you be using EOS3?  If you really are using EOS1, I recommend getting started with the one-water version (MULTI 1 2 2 6), because there are things you need to be careful about with the two-water option.

    Christine

    P.S. Thanks Keurfon for the Python script and pointer on Dirichlet boundary conditions.

      • Hakim
      • 4 yrs ago
      • Reported - view

      Christine Doughty Dear Christine, thank you and sorry for this long absence and hope you are all doing well during this covid-19 situation. I asked Keurfon about the way to run T1_modified file in his code or do I need to go back to Petrasim. Thanks

    • Staff Scientist
    • Christine_Doughty
    • 4 yrs ago
    • Reported - view

    Dear Hakim, You should use the Command Line option in Petrasim, where you specify the name of the input file you want to use (in your case, the T1_modified file).  Check "Chapter 14 Command Line Execution of TOUGH2" in the Petrasim Manual. Christine

Login to reply

Content aside

  • 4 yrs agoLast active
  • 19Replies
  • 468Views
  • 6 Following