PFLOTRAN Geomechanics help

309 views
Skip to first unread message

Sarah Reid

unread,
Dec 15, 2020, 9:46:01 AM12/15/20
to pflotra...@googlegroups.com
Good morning,

I'm an MSc student in the Department of Geoscience at the University of Calgary, and I'm using PFLOTRAN for my thesis project. My project involves modelling the formation of a groundwater discharge feature characterized by a thin and fragile crust over liquified mud. I've successfully used the subsurface deck to model fluid flow and pressure buildup, and next I'm trying to couple it to geomechanics, but I'm having some difficulty.

I've written my own MATLAB scripts to create the input files (unstructured grid, regions, mapping file) for the geomechanics, and after comparing them to the Terzaghi geomechanics example included with PFLOTRAN, I'm pretty sure they're correct. However, I'm getting an error when PFLOTRAN goes to read my unstructured mesh (see attached). From what I can tell, the unstructured grid is in the correct format, with the correct number of cells and vertices. I've searched in the google user group for any similar error, but no luck. Do you have any suggestion for what I might be doing wrong?

Thank you!

Sarah


Sarah Reid
She/Her
MSc Candidate
Department of Geoscience
University of Calgary
2500 University Drive NW, Calgary, Alberta, CANADA T2N 1N4

Screen Shot 2020-11-26 at 2.37.21 PM.png
terzaghi_sr.in

Hammond, Glenn E

unread,
Dec 15, 2020, 10:40:34 AM12/15/20
to pflotra...@googlegroups.com

Sarah,

 

Please send the complete input deck.  If large, please send a link (e.g. Google drive) where we may download the files.

 

Glenn

 

From: pflotra...@googlegroups.com <pflotra...@googlegroups.com> On Behalf Of Sarah Reid
Sent: Tuesday, December 15, 2020 6:46 AM
To: pflotra...@googlegroups.com
Subject: [pflotran-users: 6081] PFLOTRAN Geomechanics help

 

Check twice before you click! This email originated from outside PNNL.

 

--
You received this message because you are subscribed to the Google Groups "pflotran-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to pflotran-user...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/pflotran-users/YQXPR01MB343045783240E6DBBEAF700AA0C60%40YQXPR01MB3430.CANPRD01.PROD.OUTLOOK.COM.

Sarah Reid

unread,
Dec 21, 2020, 11:56:36 AM12/21/20
to pflotran-users
Hi Glenn,

Here is a link to a google drive with the files. They may be called from somewhere different in the input file. 
https://drive.google.com/drive/folders/11ttM9VuozptBX5sVY209bAQyvFFkXg5N?usp=sharing

Thank you!
Sarah

Hammond, Glenn E

unread,
Jan 5, 2021, 8:21:13 PM1/5/21
to pflotra...@googlegroups.com

Sarah,

 

Something about your ASCII mesh file is corrupted.  Try opening it in ‘nano’ on Linux.  On my Ubuntu machine, it appears binary, while with vi or emacs it appears normal.  But PFLOTRAN clearly did not like it.  When I opened the file in gedit and re-saved as a different ASCII format, the mesh was read properly.  See the following link for the updated version of the mesh file.

 

https://drive.google.com/file/d/1bKGSqtt3lWj1a45dBnHOEQIxp9TwO6pO/view?usp=sharing

 

Glenn

 

Sarah Reid

unread,
Jan 8, 2021, 12:47:56 PM1/8/21
to pflotra...@googlegroups.com
Thank you!

Sarah

You received this message because you are subscribed to a topic in the Google Groups "pflotran-users" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/pflotran-users/AvItrR7ZftU/unsubscribe.
To unsubscribe from this group and all its topics, send an email to pflotran-user...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/pflotran-users/PH0PR09MB8603E7C8E2BF08CE2B08FDBE9AD00%40PH0PR09MB8603.namprd09.prod.outlook.com.

Sarah Reid

unread,
Jan 15, 2021, 2:36:44 PM1/15/21
to pflotran-users
Hello again,

I successfully recreated the terzaghi example using mesh/region/mapping files that I created myself, but I'm missing something when I try to increase the mesh to more than 1 element in the x and y directions. I get an error "argument out of range". Are the vertices being numbered incorrectly in the mesh or mapping file? My input deck is in the attached link. 


Thank you!
Sarah

Peter Alt-Epping

unread,
Jan 19, 2021, 3:20:59 PM1/19/21
to pflotra...@googlegroups.com
Hi all,

I am running a simulation in TH mode using a database with temperature
points that are different from the ones in the default database
hanford.dat. In the latest version of PFLOTRAN this seems to be a
problem, because I get the following error message:

ERROR: The number of temperature and/or the database temperatures do not
match the hardwired values in PFLOTRAN for Debye Hueckel. Please use 8
temperatures at 0, 25, 60, 100, 150, 200, 250 and 300 C.

Did anything change recently in the way the D-H parameters are
calculated for non-isothermal problems?

Thanks,

Peter




--
______________________________________
Dr. Peter Alt-Epping
University of Bern
Institute of Geological Sciences
Baltzerstrasse 3
3012 Bern
Switzerland

T: 0041 31 631 4531
F: 0041 31 631 4843

Peter Lichtner

unread,
Jan 19, 2021, 7:20:52 PM1/19/21
to pflotran-users
Hi Peter: currently the Debye-Hueckel parameters are hardwired in the code so until this is refactored to a more general implementation you will need to use the values for the temperatures as listed.
There have been no recent changes in implementation of the D-H parameters.
…Peter
> --
> You received this message because you are subscribed to the Google Groups "pflotran-users" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to pflotran-user...@googlegroups.com.
> To view this discussion on the web visit https://groups.google.com/d/msgid/pflotran-users/d4d2c8fc-6fcc-8b77-efb8-d524d917fedc%40geo.unibe.ch.

Peter Alt-Epping

unread,
Jan 20, 2021, 4:33:01 AM1/20/21
to pflotra...@googlegroups.com
Hi Peter,

I thought the D-H parameters are interpolated for intermediate
temperatures (see reaction_database.F90) ...

Thanks,

Peter


   if (option%reference_temperature <= 0.01d0) then
     reaction%debyeA = 0.4939d0
     reaction%debyeB = 0.3253d0
     reaction%debyeBdot = 0.0374d0
   else if (option%reference_temperature > 0.d0 .and. &
            option%reference_temperature <= 25.d0) then
     temp_low = 0.d0
     temp_high = 25.d0
     call Interpolate(temp_high,temp_low,option%reference_temperature, &
                      0.5114d0,0.4939d0,reaction%debyeA)
     call Interpolate(temp_high,temp_low,option%reference_temperature, &
                      0.3288d0,0.3253d0,reaction%debyeB)
     call Interpolate(temp_high,temp_low,option%reference_temperature, &
                      0.0410d0,0.0374d0,reaction%debyeBdot)
   else if (option%reference_temperature > 25.d0 .and. &
            option%reference_temperature <= 60.d0) then
     temp_low = 25.d0
     temp_high = 60.d0
     call Interpolate(temp_high,temp_low,option%reference_temperature, &
                      0.5465d0,0.5114d0,reaction%debyeA)
     call Interpolate(temp_high,temp_low,option%reference_temperature, &
                      0.3346d0,0.3288d0,reaction%debyeB)
     call Interpolate(temp_high,temp_low,option%reference_temperature, &
                      0.0440d0,0.0410d0,reaction%debyeBdot)
   else if (option%reference_temperature > 60.d0 .and. &
            option%reference_temperature <= 100.d0) then
     temp_low = 60.d0
     temp_high = 100.d0

....

Hammond, Glenn E

unread,
Jan 20, 2021, 12:17:45 PM1/20/21
to pflotra...@googlegroups.com
I believe the issue was that the D-H parameters need to be aligned with the log Ks temperature-wise. Otherwise, there is a double interpolation effect for D-H.

Glenn

> -----Original Message-----
> From: pflotra...@googlegroups.com <pflotran-
> us...@googlegroups.com> On Behalf Of Peter Alt-Epping
> Sent: Wednesday, January 20, 2021 1:33 AM
> To: pflotra...@googlegroups.com
> Subject: Re: [pflotran-users: 6106] T-dependent D-H parameters
>
> Check twice before you click! This email originated from outside PNNL.
>
>
> > ...Peter
> https://gcc02.safelinks.protection.outlook.com/?url=https%3A%2F%2Fgroup
> s.google.com%2Fd%2Fmsgid%2Fpflotran-users%2Fd4d2c8fc-6fcc-8b77-efb8-
> d524d917fedc%2540geo.unibe.ch&amp;data=04%7C01%7Cglenn.hammond
> %40pnnl.gov%7C69c199ac6b6e48136de608d8bd266757%7Cd6faa5f90ae24033
> 8c0130048a38deeb%7C0%7C0%7C637467319938635633%7CUnknown%7CTW
> FpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJX
> VCI6Mn0%3D%7C3000&amp;sdata=qQy2oSszRjvW1jCo4bo4QFhvLn1cQYsCS
> Crb%2FV7k8Yw%3D&amp;reserved=0.
>
> --
> ______________________________________
> Dr. Peter Alt-Epping
> University of Bern
> Institute of Geological Sciences
> Baltzerstrasse 3
> 3012 Bern
> Switzerland
>
> T: 0041 31 631 4531
> F: 0041 31 631 4843
>
> --
> You received this message because you are subscribed to the Google Groups
> "pflotran-users" group.
> To unsubscribe from this group and stop receiving emails from it, send an
> email to pflotran-user...@googlegroups.com.
> To view this discussion on the web visit
> https://gcc02.safelinks.protection.outlook.com/?url=https%3A%2F%2Fgroup
> s.google.com%2Fd%2Fmsgid%2Fpflotran-users%2F9d0aa02e-7b82-c9e1-
> dcf0-
> e9d3637fe53b%2540geo.unibe.ch&amp;data=04%7C01%7Cglenn.hammond
> %40pnnl.gov%7C69c199ac6b6e48136de608d8bd266757%7Cd6faa5f90ae24033
> 8c0130048a38deeb%7C0%7C0%7C637467319938645588%7CUnknown%7CTW
> FpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJX
> VCI6Mn0%3D%7C3000&amp;sdata=%2BlLqJ9Vr4ufcpqYg5tIqJScMaGIiAmaTA
> DNLq7T6vaI%3D&amp;reserved=0.

Peter Alt-Epping

unread,
Jan 20, 2021, 1:28:29 PM1/20/21
to pflotra...@googlegroups.com
Hi Glenn,

If the reference temperatures in the database are different from the
hardwired ones, do I get the correct D-H parameters if I change temp_low
and temp_high in the interpolation in reaction_database.F90 according to
the intervals between these new reference temperatures -  and derive the
D-H parameters at these new reference temperatures via an external
interpolation (e.g. in Excel)?

Thanks,

Peter

Hammond, Glenn E

unread,
Jan 21, 2021, 11:24:15 AM1/21/21
to pflotra...@googlegroups.com
> -----Original Message-----
> From: pflotra...@googlegroups.com <pflotran-
> us...@googlegroups.com> On Behalf Of Peter Alt-Epping
> Sent: Wednesday, January 20, 2021 10:28 AM
> To: pflotra...@googlegroups.com
> Subject: Re: [pflotran-users: 6108] T-dependent D-H parameters
>
> Hi Glenn,
>
> If the reference temperatures in the database are different from the
> hardwired ones, do I get the correct D-H parameters if I change temp_low
> and temp_high in the interpolation in reaction_database.F90 according to the
> intervals between these new reference temperatures - 

Yes, this should work.

> and derive the D-H
> parameters at these new reference temperatures via an external
> interpolation (e.g. in Excel)?

I believe that this will work, but I am not 100% confident about whether it is technically correct (is it still double interpolation?).

Ideally, we have temperature-dependent D-H parameters associated with and read from each database. We could then refactor the database setup to read these parameters. We could accomplish this by adding a DEBYE_HUCKEL block at the top of the database.

Glenn
> s.google.com%2Fd%2Fmsgid%2Fpflotran-users%2Fb7ed651d-630c-43c6-
> f989-
> 7d1694075aa9%2540geo.unibe.ch&amp;data=04%7C01%7Cglenn.hammond
> %40pnnl.gov%7C707a2d8c81ad471e0c1108d8bd713963%7Cd6faa5f90ae24033
> 8c0130048a38deeb%7C0%7C0%7C637467641291254228%7CUnknown%7CTW
> FpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJX
> VCI6Mn0%3D%7C1000&amp;sdata=ifxFUjL92b8ecc%2BZhO5%2Bzvs6NUA%2
> B3XiWyr7TWiJKp3Y%3D&amp;reserved=0.

Peter Lichtner

unread,
Jan 21, 2021, 11:27:46 AM1/21/21
to pflotran-users
Another possibility is to calculate the D-H parameters directly at the desired temperatures. I can provide a fortran code for that. Might be useful as a check of the interpolation.
…Peter

Peter Alt-Epping

unread,
Jan 21, 2021, 11:42:23 AM1/21/21
to pflotra...@googlegroups.com
Hi Glenn and Peter,

Thanks for your replies. Yes, I think it would be useful to have the D-H
parameters in the database. Databases such as SOLTHERM have that as well.

Peter, I would appreciate very much if you could send me the code. I
will use it to check that the interpolation is correct.

Thanks,

Peter

benjami...@gmail.com

unread,
Jan 22, 2021, 11:03:16 AM1/22/21
to pflotran-users
Is the dielectric constant calculated in the IF97 or any EoS in PFLOTRAN?  If so, all you need to calculate Debye-Huckel A and B is density, temperature, and dielectric constant.  We have the FORTRAN code to calculate this, as well as Bdot, implemented in DBCreate and it could easily be absorbed by PFLOTRAN

DHparams.jpg

Hammond, Glenn E

unread,
Jan 24, 2021, 2:38:49 AM1/24/21
to pflotra...@googlegroups.com

Sarah,

 

I believe that the issue is that your subsurface structured grid has 200 grid cells 2x2x50 and your unstructured grid has 400 grid cells.  Geomechanics sets up a mapping from the structured to unstructured, but it cannot map beyond 200 cells.  The mapping assigns a value of -1 to the cells that is cannot map (unstructured geomech grid cell ids > 200) and that is the cause of the array bounds issue.  Please see if this is resolved by using an unstructured grid with 200 cells or a structured grid with 400 cells.

 

Glenn

 

Message has been deleted

Sarah Reid

unread,
Feb 4, 2021, 2:12:06 PM2/4/21
to pflotran-users
Hi Glenn,

In the geomech examples included with pflotran, the unstructured grids have more elements than the structured grid too, so I don't understand how to make it work with mine too. 

I'm trying to couple the unstructured geomech grid and the structured flow grid by placing the vertex of an unstructured element at the centre of the structured element, so there will always be more unstructured elements than structured elements. Is PFLOTRAN the right tool I should be using?

Thanks,
Sarah



Hammond, Glenn E

unread,
Feb 4, 2021, 2:21:58 PM2/4/21
to pflotra...@googlegroups.com

Satish,

 

Can you answer Sarah’s questions?

 

Sarah,

 

If you do not receive a response, I would choose a different tool.

 

Glenn

 

Swapnil Kar

unread,
Mar 18, 2021, 6:59:13 AM3/18/21
to pflotran-users
Respected Sir,
                           I am facing the following issue while trying to run the terzaghi's problem under the geomechanics example problem.I would be thankful to you if you can help me in this regard. I have downloaded the example problems from this bitbucket link :

Screenshot from 2021-03-18 16-19-36.pngScreenshot from 2021-03-18 16-20-30.png
Thanking You,
Swapnil Kar

Hammond, Glenn E

unread,
Mar 22, 2021, 11:10:53 AM3/22/21
to pflotra...@googlegroups.com

Swapnil,

 

Please send the path to the example input file that you are attempting to run.  I will forewarn you that support for the geomechanics capability in PFLOTRAN has been very slow.  I can try to help, but I did not implement the capability.


Glenn

 

From: pflotra...@googlegroups.com <pflotra...@googlegroups.com> On Behalf Of Swapnil Kar
Sent: Thursday, March 18, 2021 3:59 AM
To: pflotran-users <pflotra...@googlegroups.com>
Subject: Re: [pflotran-users: 6230] PFLOTRAN Geomechanics help

 

Respected Sir,

                           I am facing the following issue while trying to run the terzaghi's problem under the geomechanics example problem.I would be thankful to you if you can help me in this regard. I have downloaded the example problems from this bitbucket link :

 

Reply all
Reply to author
Forward
0 new messages