solvers_groundwaterFoam

groundwaterFoam solver

Description

Transient solver for groundwater flow in porous media (Richards' equation with isotropic permeability volScalarField)

The non-linearities of relative permeability and capillary pressure are handled at each timestep by successive non linear methods:

1. Picard's iterations first to converge under a user-defined Picard tolerance
2. Newton's iterations at each time step to converge under a user-defined Newton tolerance

Each algorithm can be skipped by setting tolerance equal to 1. If the convergence is not reached in maxIter iterations, the current time iteration starts again with a smaller timestep.

Time step is managed using an estimation of the time truncation error (see below).

Seepage can be activated by specifying the name of the patch patchDEM (representing the top of the watershed). At given time iteration, if a neighbour cell is saturated and cell flux < 0 : solver sets h=0 at this cell and computes exfiltration outflow.

Configuration files

constant/transportProperties :

thetamin  thetamin [0 0 0 0 0 0 0] 0.102; // minimal saturation
thetamax thetamax [0 0 0 0 0 0 0] 0.368; // maximal saturation

Ss Ss [0 -1 0 0 0 0 0] 0; // specific storage

patchDEM nameOfThePatch; // patch used for seepage options

phase.theta
{
    rho     rho [1 -3 0 0 0 0 0] 1e3; // density of the fluid
    mu      mu [1 -1 -1 0 0 0 0] 1e-3; // dynamic viscosity of the fluid
}

relativePermeabilityModel  VanGenuchten;

capillarityModel    VanGenuchten;

VanGenuchtenCoeffs // parameters of the Van Genuchten kr/pc model
{
    m       0.5;
    alpha 3.35;
}

system/fvSolution :

solvers
{
   "h" // for Picard's solver
    {
        solver          PCG;
        preconditioner  DIC;
        tolerance       1e-12;
        relTol          0;
    }
    "deltah" // for Newton's solver
    {
        solver          PBiCG;
        preconditioner  DILU;
        tolerance       1e-12;
        relTol          0;
    }

}

Picard
{
    tolerance 0.001; // tolerance for the Picard's algorithm
    maxIter 10; // maximum iterations for the Picard's algorithm
}

Newton
{
    tolerance 1e-06; // tolerance for the Newton's algorithm
    maxIter 10; // maximum iterations for the Newton's algorithm
}

system/controlDict

adjustTimeStep  yes;

truncationError   0.005; // theta variation instruction for computing the time step

eventTimeTracking false; // to force the solver to  explicitly compute event time solutions (instead of time linear interpolations)

CSVoutput       true; // active the waterMassBalance.csv output

eventFileOutput outputList.dat; // to specify the writing time outputs (replaces usual write() function of OpenFOAM)

Required fields

• 0/h : The potential field
• 0/Utheta : The velocity field
• constant/g : gravity field
• constant/K : permeability field

Optional fields

• 0/thetamin and 0/thetamax : spatialized minimal and maximal saturation (replace thetamin and thetamax in transportProperties)
• 0/m and 0/alpha : spatialized Van Genuchten parameters (replace m and alpha in transportProperties)

Timestep managing

The computation of timestep for next iteration is directly computed using truncation error related to the time scheme used. Due to the Newton's method time dicretization, only the Euler scheme is available with:

deltaT = Foam::pow(2 x truncationError x Hmax[speciesi]/dH2dT2max[speciesi],1./3.)

where dH2dT2max is the maximal value of the second order time derivative and Hmax the value of hwater in this cell.

Return to porousMultiphaseFoam main page