acroucher / PyTOUGH

A Python library for automating TOUGH2 simulations of subsurface fluid and heat flow
GNU Lesser General Public License v3.0
96 stars 36 forks source link

Use mulgrid() to read existing TOUGH input file #21

Closed jxalfoudqt closed 5 years ago

jxalfoudqt commented 5 years ago

Hi Adrian,

How can I visualize the grid of radial model example by using 'mulgrid()'? There is 'dat.grid' but no 'geo' as specified in Hydrostatics 3D model. I think this method can be used for grid visualization of all existing TOUGH input files.

Thanks. Tao

acroucher commented 5 years ago

It depends whether it is a 1-D radial model or a 2-D radial model (i.e. r-z cylindrical coordinates).

For a 1-D radial model, you can do the visualisation with simple line plots vs. radius, using e.g. matplotlib. There is an example of this in the PyTOUGH wiki.

For a 2-D radial model in r-z coordinates, you can create a dummy 2-D rectangular mulgrid object representing the r-z geometry, purely for visualisation purposes. Give it the same radial and vertical block sizes that you pass to t2grid.radial() to create your TOUGH2 grid, and a dummy (arbitrary) block size in the y-direction. Then you can do visualisation of your results with e.g. mulgrid.slice_plot().

jxalfoudqt commented 5 years ago

Many thanks Adrian. Actually, I take your idea to generate a mulgrid object and write it to a .vtu file for visualization.

chenmingzhang commented 4 years ago

just a comment for others who may face similar situations https://tough.forumbee.com/t/18hx69

jousheinfo commented 2 years ago

Hi @acroucher, one question please. I have one radial grid very similar to the example in the link shared by @chenmingzhang, but I would like to know if it is possible to make a 3D plot for visualization. I have been reading other questions related to this one, but I could not find any useful so far. I would very much appreciate your help!

acroucher commented 2 years ago

hi @jousheinfo , just to clarify, are you wanting to make a 3D plot of a radial grid? I generally don't do that, because a radial grid is really only 2D (r, z), so it can be visualized better using 2D plots, as explained above.

If you are wanting to make a 3D plot of a 3D grid, then you can do that in PyTOUGH by using the write_vtk() method of a mulgrid object (or a t2listing object if you are wanting to plot results).

jousheinfo commented 2 years ago

@acroucher thank you for your response. I got another question, please. I have the following radial grid made in PYTOUGH:

from t2data import * from math import log10

mesh parameters:

well_radius = (6 / 2) * 0.0254 reservoir_radius = 1000. well_depth = 2500 thickness = 50 num_blocks = int((well_depth / thickness) + 1)

layer_thickness = np.repeat(thickness,num_blocks-1)

dr = np.logspace(log10(well_radius), log10(reservoir_radius), num_blocks) dz = layer_thickness

create TOUGH2 input data:

dat = t2data() dat.title = 'model_injection' dat.grid = t2grid().radial(dr, layer_thickness, origin = [0.,0.,0.])

Even tough I defined the max desired radius 1000 m, in the input file, I can see the last block has a radius of about 10 times bigger (5.288e+3*2 m). Why is this happenning? The desired number of cells is 50, and the thickness of each cell is constant equal to 50m.

radial_grid

I do not really understand how pytough computes the logarithmic increment and how can I find the increment factor?

I am trying to replicate the mesh file with meshmaker, but I'm struggling. I would very much appreciate your help. Thank you in advance!

acroucher commented 2 years ago

Note that the dr values are radial block thicknesses, not radii of the centres. Your last dr value is 1000 which means the last block will have a radial thickness of 1000 m. The outer radius of the grid will be the sum of the dr values which in your case is 5788.5 m at the moment.

Also, I'm not sure how many layers you were meaning to create, but if you just want 50 blocks (i.e. one layer) then your layer_thickness value should have just one element. You don't need to specify thicknesses for each column. At the moment you are creating a grid with 50 * 51 = 2550 blocks.