Re: [DuMuX] using el2p with fluid flow
Dear Timo,
Let me introduce myself in this talk.
Edscott has been able to successfully develop a program based on Dumux 2.12 for
the problem of low salinity water injection in the last year.
In these moments we have a geomechanics project where we are even going to
model the effects of fractures. We start with the linear elasticity model
(elastic) and until now everything has worked as specified in the
documentation. But with the problems "el1p2c" and "el2p" we have serious
difficulties. I do not know if you could be so kind to help us as it seems that
they only work for the particular test cases.
In the problem el1p2c we would like it to be simplify to "el1p", i.e., a single
component, and we would like to implement a Mandel-type problem (poroelastic
problem). We have the same problem already implemented in COMSOL and in FEniCS
and it gives the correct results. I wonder if I send you the data you or
anybody else from your team could try to implement it?
Thank you very much in advance.
Regards,
Martin
--
__
Dr. Martín A. Díaz Viera
Instituto Mexicano del Petróleo
Gerencia de Ingeniería de Recuperación Adicional
Edificio 6, Cubículo 112
Eje Central Lázaro Cárdenas Norte 152,
San Bartolo Atepehuacan,
Gustavo A. Madero,
Ciudad de México, C.P. 07730
Tel.: (5255) 9175-6473
E-mail: mdi...@imp.mx, mdiaz...@yahoo.com.mx
WEB Page: http://www.esmg-mx.org/
El viernes, 22 de febrero de 2019 15:25:16 GMT-6, Timo Koch
escribió:
Hi Edscott,
Martin Beck implemented exactly that in his PhD thesis in Dumux 2.12 and the
code should be somewhere (check dumux-pub) or ask Martin.
Also this should be much easier to do in Dumux 3.0, as the structure was
developed with something like sequential coupling in mind.
Timo
Am 22.02.2019 um 21:58 schrieb Ed Scott Wilson Garcia :
#yiv1963932208 #yiv1963932208 -- _filtered #yiv1963932208
{font-family:Calibri;panose-1:2 15 5 2 2 2 4 3 2 4;}#yiv1963932208
#yiv1963932208 p.yiv1963932208MsoNormal, #yiv1963932208
li.yiv1963932208MsoNormal, #yiv1963932208 div.yiv1963932208MsoNormal
{margin:0cm;margin-bottom:.0001pt;font-size:12.0pt;font-family:New
serif;}#yiv1963932208 a:link, #yiv1963932208 span.yiv1963932208MsoHyperlink
{color:blue;text-decoration:underline;}#yiv1963932208 a:visited, #yiv1963932208
span.yiv1963932208MsoHyperlinkFollowed
{color:purple;text-decoration:underline;}#yiv1963932208
span.yiv1963932208EstiloCorreo17
{font-family:sans-serif;color:#1F497D;}#yiv1963932208
.yiv1963932208MsoChpDefault {font-size:10.0pt;} _filtered #yiv1963932208
{margin:70.85pt 3.0cm 70.85pt 3.0cm;}#yiv1963932208
div.yiv1963932208WordSection1 {}#yiv1963932208
I’ve kind of given up on this approach. What I want to try next is to create a
multiprocess run, with an elastic model coupled sequentially with a 2p model.
The first process will solve for the deformations (preferably with the el2p
pdelab formula). This should give me a pressure distribution in the grid. Then
use this pressure distribution to feed a 2p model for a single iteration with
an outflow condition. This would have to spin until pressure converges. Once
that happens, go back to the elastic and start again.
Basically all I want to do for the moment is to emulate a Mandel problem, where
we have a constant vertical pressure on top, no flow on the bottom and on the
left, no vertical deformation on the left side and no horizontal deformation on
the bottom. The right side is left open both the deformation and fluid flow.
Fluid will be squished out.
De: Dumux [mailto:dumux-boun...@listserv.uni-stuttgart.de]En nombre de Timo Koch
Enviado el: viernes, 22 de febrero de 2019 02:47 p. m.
Para: DuMuX User Forum
Asunto: Re: [DuMuX] using el2p with fluid flow
Am 21.02.2019 um 23:50 schrieb Ed Scott Wilson Garcia :
I'm trying to use the el2p example program to construct a program with
inflow/outflow. In order to test the feasibility, after simplifying to 2
dimensions instead of 3, I have used the data from a working problem with the
2p model. On the first test, set the boundary conditions allDirichlet for uxIdx
and uyIdx and equal to zero. This would be equivalent to null elasticity and
the model should reduce to the 2p model.
Hi Edscott,
why would it reduce to 2p? You can still have deformations inside the domain.
Formulation was changed from the default to pNsW. And the fluidsystem was
changed to the oil/water fluidsystem used a working 2p problem.
The boundary conditions are set to:
void boundaryTypesAtPos(BoundaryTypes , const GlobalPosition&
globalPos) const
{
values.setDirichlet(uxIdx);
values.setDirichlet(uyIdx);
Scalar top = this->bBoxMax()[1];
Scalar bottom = this->bBoxMin()[1];
if (globalPos[1] > top - eps_)
{
values.setDirichlet(pressureIdx); // nonwetting pressure
Re: [DuMuX] using el2p with fluid flow
Hi Edscott,
Martin Beck implemented exactly that in his PhD thesis in Dumux 2.12 and the
code should be somewhere (check dumux-pub) or ask Martin.
Also this should be much easier to do in Dumux 3.0, as the structure was
developed with something like sequential coupling in mind.
Timo
> Am 22.02.2019 um 21:58 schrieb Ed Scott Wilson Garcia :
>
> I’ve kind of given up on this approach. What I want to try next is to create
> a multiprocess run, with an elastic model coupled sequentially with a 2p
> model. The first process will solve for the deformations (preferably with the
> el2p pdelab formula). This should give me a pressure distribution in the
> grid. Then use this pressure distribution to feed a 2p model for a single
> iteration with an outflow condition. This would have to spin until pressure
> converges. Once that happens, go back to the elastic and start again.
>
> Basically all I want to do for the moment is to emulate a Mandel problem,
> where we have a constant vertical pressure on top, no flow on the bottom and
> on the left, no vertical deformation on the left side and no horizontal
> deformation on the bottom. The right side is left open both the deformation
> and fluid flow. Fluid will be squished out.
>
>
>
>
> De: Dumux [mailto:dumux-boun...@listserv.uni-stuttgart.de] En nombre de Timo
> Koch
> Enviado el: viernes, 22 de febrero de 2019 02:47 p. m.
> Para: DuMuX User Forum
> Asunto: Re: [DuMuX] using el2p with fluid flow
>
>
> Am 21.02.2019 um 23:50 schrieb Ed Scott Wilson Garcia :
>
>
> I'm trying to use the el2p example program to construct a program with
> inflow/outflow. In order to test the feasibility, after simplifying to 2
> dimensions instead of 3, I have used the data from a working problem with the
> 2p model. On the first test, set the boundary conditions allDirichlet for
> uxIdx and uyIdx and equal to zero. This would be equivalent to null
> elasticity and the model should reduce to the 2p model.
> Hi Edscott,
>
> why would it reduce to 2p? You can still have deformations inside the domain.
>
>
> Formulation was changed from the default to pNsW. And the fluidsystem was
> changed to the oil/water fluidsystem used a working 2p problem.
>
> The boundary conditions are set to:
>
> void boundaryTypesAtPos(BoundaryTypes , const GlobalPosition&
> globalPos) const
> {
> values.setDirichlet(uxIdx);
> values.setDirichlet(uyIdx);
>
> Scalar top = this->bBoxMax()[1];
> Scalar bottom = this->bBoxMin()[1];
> if (globalPos[1] > top - eps_)
> {
>values.setDirichlet(pressureIdx); // nonwetting pressure
>values.setOutflow(Indices::contiNEqIdx);
> }
> else // Neumann for the remaining boundaries
>values.setAllNeumann();
> }
>
>
> void neumannAtPos(PrimaryVariables , const GlobalPosition&
> globalPos) const
> {
> if (globalPos[1] < eps_)
> {
> auto wflow = -1.37e-3; // water flow at bottom
> values[Indices::contiWEqIdx] = wflow;
> values[Indices::contiNEqIdx] = 0; // oil flow at bottom
> } else {
> // no-flow on the remaining Neumann-boundaries.
> values[Indices::contiWEqIdx] = 0;
> values[Indices::contiNEqIdx] = 0;
> }
>
> }
>
> void dirichletAtPos(PrimaryVariables , const GlobalPosition&
> globalPos) const
> {
> values = 0.0;
> Scalar top = this->bBoxMax()[1];
> Scalar bottom = this->bBoxMin()[1];
> Scalar right = this->bBoxMin()[0];
> if (globalPos[1] > top - eps_){
> values[Indices::pressureIdx] = 1.72e7;
> }
> }
>
> Since "void initialAtPos(PrimaryVariables , const
> GlobalPosition )" is ignored, the initial pressure is set
> withinitial pressure is set with:
>
> void initializePressure()
> {
> for(const auto& vertex : vertices(gridView_))
> {
> int vIdxGlobal = this->vertexMapper().index(vertex);
> GlobalPosition globalPos = vertex.geometry().corner(0);
>
> // initial approximate pressure distribution at start of
> initialization run
> pInit_[vIdxGlobal] = -1.72e7;
> }
> }
>
> Gravity is disabled so depthBore should not have any effect.
>
> Expected result: Should reproduce the 2p flow problem.
>
>
> why?
>
>
>
>
> Result: Segment violation crash. This is the gdb output following compilation
> with -O0.
>
> Looks like the fluid state is uninitialized. Did you debug and check why it
> is uninitialized? Is the error in Dumux code (does it occur in an unmodified
> setup)? You seem to use a custom fluid system with a member function
> “saturation()”?! It’s not clear to me why any fluid system should have such a
> function...
>
>
>
> (gdb) run
> Starting program:
>
Re: [DuMuX] using el2p with fluid flow
I’ve kind of given up on this approach. What I want to try next is to create a
multiprocess run, with an elastic model coupled sequentially with a 2p model.
The first process will solve for the deformations (preferably with the el2p
pdelab formula). This should give me a pressure distribution in the grid. Then
use this pressure distribution to feed a 2p model for a single iteration with
an outflow condition. This would have to spin until pressure converges. Once
that happens, go back to the elastic and start again.
Basically all I want to do for the moment is to emulate a Mandel problem, where
we have a constant vertical pressure on top, no flow on the bottom and on the
left, no vertical deformation on the left side and no horizontal deformation on
the bottom. The right side is left open both the deformation and fluid flow.
Fluid will be squished out.
De: Dumux [mailto:dumux-boun...@listserv.uni-stuttgart.de] En nombre de Timo
Koch
Enviado el: viernes, 22 de febrero de 2019 02:47 p. m.
Para: DuMuX User Forum
Asunto: Re: [DuMuX] using el2p with fluid flow
Am 21.02.2019 um 23:50 schrieb Ed Scott Wilson Garcia
mailto:edsc...@imp.mx>>:
I'm trying to use the el2p example program to construct a program with
inflow/outflow. In order to test the feasibility, after simplifying to 2
dimensions instead of 3, I have used the data from a working problem with the
2p model. On the first test, set the boundary conditions allDirichlet for uxIdx
and uyIdx and equal to zero. This would be equivalent to null elasticity and
the model should reduce to the 2p model.
Hi Edscott,
why would it reduce to 2p? You can still have deformations inside the domain.
Formulation was changed from the default to pNsW. And the fluidsystem was
changed to the oil/water fluidsystem used a working 2p problem.
The boundary conditions are set to:
void boundaryTypesAtPos(BoundaryTypes , const GlobalPosition&
globalPos) const
{
values.setDirichlet(uxIdx);
values.setDirichlet(uyIdx);
Scalar top = this->bBoxMax()[1];
Scalar bottom = this->bBoxMin()[1];
if (globalPos[1] > top - eps_)
{
values.setDirichlet(pressureIdx); // nonwetting pressure
values.setOutflow(Indices::contiNEqIdx);
}
else // Neumann for the remaining boundaries
values.setAllNeumann();
}
void neumannAtPos(PrimaryVariables , const GlobalPosition&
globalPos) const
{
if (globalPos[1] < eps_)
{
auto wflow = -1.37e-3; // water flow at bottom
values[Indices::contiWEqIdx] = wflow;
values[Indices::contiNEqIdx] = 0; // oil flow at bottom
} else {
// no-flow on the remaining Neumann-boundaries.
values[Indices::contiWEqIdx] = 0;
values[Indices::contiNEqIdx] = 0;
}
}
void dirichletAtPos(PrimaryVariables , const GlobalPosition&
globalPos) const
{
values = 0.0;
Scalar top = this->bBoxMax()[1];
Scalar bottom = this->bBoxMin()[1];
Scalar right = this->bBoxMin()[0];
if (globalPos[1] > top - eps_){
values[Indices::pressureIdx] = 1.72e7;
}
}
Since "void initialAtPos(PrimaryVariables , const GlobalPosition
)" is ignored, the initial pressure is set withinitial pressure is
set with:
void initializePressure()
{
for(const auto& vertex : vertices(gridView_))
{
int vIdxGlobal = this->vertexMapper().index(vertex);
GlobalPosition globalPos = vertex.geometry().corner(0);
// initial approximate pressure distribution at start of
initialization run
pInit_[vIdxGlobal] = -1.72e7;
}
}
Gravity is disabled so depthBore should not have any effect.
Expected result: Should reproduce the 2p flow problem.
why?
Result: Segment violation crash. This is the gdb output following compilation
with -O0.
Looks like the fluid state is uninitialized. Did you debug and check why it is
uninitialized? Is the error in Dumux code (does it occur in an unmodified
setup)? You seem to use a custom fluid system with a member function
“saturation()”?! It’s not clear to me why any fluid system should have such a
function...
(gdb) run
Starting program:
/home/FreeBSD12/home/edscott/geomecanica/projects/el2p/build-cmake/src/el2p
-ParameterFile ./pd1-el2p.input
[Thread debugging using libthread_db enabled]
Using host libthread_db library "/usr/lib/libthread_db.so.1".
[Detaching after fork from child process 26597]
[New Thread 0x747fc700 (LWP 26605)]
[New Thread 0x73eeb700 (LWP 26609)]
Don't panic... !
El2P_TestProblem: Initializing the fluid system for the el2p model
myOilBrine-nt: Using non tabulated H20 properties.
episode set to: 1
Initializing problem 'el2p'
problem init()
el2pmodel calls: initialPressSat
InitialPressSat: setPressure function called
numDofs = 396
Writing result file for
[DuMuX] typetag and property visualization app (2.12)
One of the main hurdles I have had with Dumux code is figuring out where properties and typetags are defined, because this sends me opening and closing files and by the time I 've reached what I was looking for three more questions have popped up in my mind. The moment I just got plain tired of this scheme, I decided to create a parser which would gobble up all the included files and create a structure xml file of where, what and how things were constructed. This is the script "parse7.pl". I suppose there could be an xml viewer out there to view the results, but I found it much easier to construct my own app specific to the created xml. It works for me. If it works for someone else, that's great. You can download the code from github.com from https://github.com/edscott/dumux/tree/master/structure There are some screenshots there as well. Any bugreports or feature enhancement suggestions are welcome. This is the content of the README file: To build, just do a "make" Requirements, gtk3 (developer) and perl. Once built you will have two programs, the parser (parse7.pl) and the visualizer (structure). The parser's job is to create an xml file which will detail the way template headers are included, properties which are defined and typetags and their inheritance relationship. It is executed with "perl parser7.pl source.cc [problemTypeTag]" These three items, "Files", "properties" and "Typetags" are the top level elements in a gtk treeview. As the "Files" treeview is expanded, you can see the order and from where each template header is included. See Snapshot-1. When the "Properties" treeview is expanded, you can see all the properties which are defined in the included template headers, their value and location. Snapshot-2. When the "Typetags" treeview is expanded (Snapshot-3), you can see all the defined typetags. If the optional problemTypeTag is specified, then the xml will use it to create an inheritance chain and you can see the main type tag with an "emblem-important" icon. First level inheritance have a "emblem-default" icon, which by default is a green checkmark. Second level has a solid background checkmark, third level a simple checkmark and all the rest a blue dot. Whatever does not have an icon is of no use to the problemTypeTag specified. This has helped me speed up looking for what properties are defined and their values, as well as to get a general idea of the structure of a Dumux aplication. If I continue with Dumux, I will try to get a Dumux-3.0 parser running, but at the moment the future is unclear. And as usual, the code is distributed with GNU version 3 license, hoping that it be of use to others, but with no warranty nor suitability of use. feb-21-201 Enjoy. ___ Dumux mailing list Dumux@listserv.uni-stuttgart.de https://listserv.uni-stuttgart.de/mailman/listinfo/dumux
Re: [DuMuX] el2p initialization run
Hi Edscott,
I’m assuming you use version 2.12. Can you be more specific? What are these
boundary condition? Why wouldn’t it be possible to solve it with the el2p model?
What does it mean “Newton does not converge on the first step”?
Timo
> Am 22.02.2019 um 20:57 schrieb Ed Scott Wilson Garcia :
>
> I am studying th el2p test example to see if it is possible to use it without
> modifications to the geomechanics/el2p/model in order to validate the code
> with a well known published example of the Mandel problem (Sangnimnuan 2018).
> My question is whether the el2p model can be used without further
> modifications for these boundary conditions.
>
> The code for the initialization has no problem whatsoever when modifying the
> fixed parameters (Young Modulus, Poisson ratio, porosity, permeability rock
> density).
> But when I modify the boundary conditions from Dirichlet to Neumann to allow
> an external force to be applied on the top boundary and expansion on the left
> boundary, the Newton method does not converge on the first step.
>
> IOW:
>
>void boundaryTypesAtPos(BoundaryTypes , const GlobalPosition&
> globalPos) const
>{
>values.setAllNeumann();
>
>// The solid displacement normal to the lateral boundaries is fixed.
>if(globalPos[0] < eps_ || globalPos[0] > this->bBoxMax()[0]-eps_)
>{
>values.setDirichlet(uxIdx);
>if(initializationRun_ == false)
>{
>values.setDirichlet(pressureIdx);
>values.setDirichlet(saturationIdx);
>}
>}
>// The solid displacement normal to the lateral boundaries is fixed.
>// Only for 3 d
>if (dim == 3) {
>if(globalPos[1] < eps_ || globalPos[1] > this->bBoxMax()[1]-eps_)
>{
>values.setDirichlet(uyIdx);
>if(initializationRun_ == false)
>{
>values.setDirichlet(pressureIdx);
>values.setDirichlet(saturationIdx);
>}
>}}
>
>// Lower boundary closed for brine and CO2 flux, uz is fixed.
>if(globalPos[dimWorld-1] < eps_)
>{
>if(dim == 2)
>values.setDirichlet(uyIdx);
>if(dim == 3)
>values.setDirichlet(uzIdx);
>}
> #define SANGNIMNUAN
> #ifdef SANGNIMNUAN
>// right side
>if(globalPos[0] > this->bBoxMax()[0]-eps_) {
>values.setNeumann(uxIdx);
>values.setNeumann(uyIdx);
>}
>// top side
>if(globalPos[1] > this->bBoxMax()[1]-eps_){
>values.setNeumann(uxIdx);
>values.setNeumann(uyIdx);
>}
> #endif
>
>// for the initialization run the pressure and saturation
>// values are only given at the top boundary.
>if(globalPos[dimWorld-1] > this->bBoxMax()[dimWorld-1]-eps_)
>{
>values.setDirichlet(pressureIdx);
>values.setDirichlet(saturationIdx);
>}
>}
>
>void neumannAtPos(PrimaryVariables , const GlobalPosition&
> globalPos) const
>{
>values = 0;
> #ifdef SANGNIMNUAN
>// right side
>if(globalPos[0] > this->bBoxMax()[0]-eps_) {
>// no force
>}
>// top side
>if(globalPos[1] > this->bBoxMax()[1]-eps_){
>// assuming area of 1 and load of 616 psi. P=F/A
>values[uyIdx] = 4.2472e+6; // Newtons
>// no force in ux
>
>}
> #endif
>if(globalPos[1] < eps_ || globalPos[1] > this->bBoxMax()[1]-eps_){
>// assuming area of 1 and load of 616 psi. P=F/A
>//values[uyIdx] = 4.2472e+6; // Newtons
>}
>}
>
> I am assuming that the el2p model, as is, will not suffice to run a simple
> Mandel problem without some deep modifications. Yes?
>
> best regards,
>
> Edscott
>
>
>
>
>
> ___
> Dumux mailing list
> Dumux@listserv.uni-stuttgart.de
> https://listserv.uni-stuttgart.de/mailman/listinfo/dumux
___
Dumux mailing list
Dumux@listserv.uni-stuttgart.de
https://listserv.uni-stuttgart.de/mailman/listinfo/dumux
[DuMuX] el2p initialization run
I am studying th el2p test example to see if it is possible to use it without
modifications to the geomechanics/el2p/model in order to validate the code with
a well known published example of the Mandel problem (Sangnimnuan 2018). My
question is whether the el2p model can be used without further modifications
for these boundary conditions.
The code for the initialization has no problem whatsoever when modifying the
fixed parameters (Young Modulus, Poisson ratio, porosity, permeability rock
density).
But when I modify the boundary conditions from Dirichlet to Neumann to allow an
external force to be applied on the top boundary and expansion on the left
boundary, the Newton method does not converge on the first step.
IOW:
void boundaryTypesAtPos(BoundaryTypes , const GlobalPosition&
globalPos) const
{
values.setAllNeumann();
// The solid displacement normal to the lateral boundaries is fixed.
if(globalPos[0] < eps_ || globalPos[0] > this->bBoxMax()[0]-eps_)
{
values.setDirichlet(uxIdx);
if(initializationRun_ == false)
{
values.setDirichlet(pressureIdx);
values.setDirichlet(saturationIdx);
}
}
// The solid displacement normal to the lateral boundaries is fixed.
// Only for 3 d
if (dim == 3) {
if(globalPos[1] < eps_ || globalPos[1] > this->bBoxMax()[1]-eps_)
{
values.setDirichlet(uyIdx);
if(initializationRun_ == false)
{
values.setDirichlet(pressureIdx);
values.setDirichlet(saturationIdx);
}
}}
// Lower boundary closed for brine and CO2 flux, uz is fixed.
if(globalPos[dimWorld-1] < eps_)
{
if(dim == 2)
values.setDirichlet(uyIdx);
if(dim == 3)
values.setDirichlet(uzIdx);
}
#define SANGNIMNUAN
#ifdef SANGNIMNUAN
// right side
if(globalPos[0] > this->bBoxMax()[0]-eps_) {
values.setNeumann(uxIdx);
values.setNeumann(uyIdx);
}
// top side
if(globalPos[1] > this->bBoxMax()[1]-eps_){
values.setNeumann(uxIdx);
values.setNeumann(uyIdx);
}
#endif
// for the initialization run the pressure and saturation
// values are only given at the top boundary.
if(globalPos[dimWorld-1] > this->bBoxMax()[dimWorld-1]-eps_)
{
values.setDirichlet(pressureIdx);
values.setDirichlet(saturationIdx);
}
}
void neumannAtPos(PrimaryVariables , const GlobalPosition&
globalPos) const
{
values = 0;
#ifdef SANGNIMNUAN
// right side
if(globalPos[0] > this->bBoxMax()[0]-eps_) {
// no force
}
// top side
if(globalPos[1] > this->bBoxMax()[1]-eps_){
// assuming area of 1 and load of 616 psi. P=F/A
values[uyIdx] = 4.2472e+6; // Newtons
// no force in ux
}
#endif
if(globalPos[1] < eps_ || globalPos[1] > this->bBoxMax()[1]-eps_){
// assuming area of 1 and load of 616 psi. P=F/A
//values[uyIdx] = 4.2472e+6; // Newtons
}
}
I am assuming that the el2p model, as is, will not suffice to run a simple
Mandel problem without some deep modifications. Yes?
best regards,
Edscott
___
Dumux mailing list
Dumux@listserv.uni-stuttgart.de
https://listserv.uni-stuttgart.de/mailman/listinfo/dumux
