Skip to main content

STDBE511 example 1

  • Both fluid and porous medium are considered
  • L = 0.1 m
  • dx = 0.1/50
  • K0 = E-10 m2
  • Da = E-8
  • V_TOP = 0.01
  • dt = dx / V_TOP = 0.2
  • StabOption=1

Mesh:

In this example we test following:

  • [[STDBE511_#SetSTDBE511Param]]
  • [[LinSolver_#SetLinSolverParam]]
  • [[Domain_#Initiate]]

[TOC]

!!! note "Some simulation parameter"

#define _K0_ 1.0E-10
#define _PHI_ 0.5_DFP
#define _DT_ 0.2
#define _VTOP_ 0.01_DFP
#define _MU_FLUID_ 0.001
#define _RHO_FLUID_ 1000.0
#define _TIMESTEPS_ maxI4B
#define _MAX_ITER_ 20
#define _tDBC4P_ 0
#define _REF_PRESSURE_NODE_ 90
#define _STAB_PARAM_OPTION_ 1
#define _MATERTIAL_INTERFACES_ [7]
#define _SAVE_STATE_FREQ_ 100
#define _WRITE_FREQ_ 1
varNameComment
_K0_permeability of the soil
_PHI_porosity of the soil
_DT_time step size
_VTOP_velocity at top
_MU_FLUID_dynamic viscosity of fluid
_RHO_FLUID_mass-density of fluid
_TIMESTEPS_number of time steps
_MAX_ITER_maximum number of Newton iteration
_tDBC4P_total Dirichlet BC for pressure
_REF_PRESSURE_NODE_reference node number for pressure
_STAB_PARAM_OPTION_stabilization option
_MATERTIAL_INTERFACES_mesh-id of fluid-soil interfaces
_SAVE_STATE_FREQ_save-state frequency
_WRITE_FREQ_write-frequency

!!! note "Include modules"

PROGRAM main
USE easifemBase
USE easifemClasses
USE easifemMaterials
USE easifemKernels
USE STDBE511_Class
IMPLICIT NONE

Declare variables

  TYPE(STDBE511_) :: obj
TYPE(ParameterList_) :: param
TYPE(HDF5File_) :: domainFile
TYPE( HDF5File_ ) :: outfile
TYPE( HDF5File_ ) :: statefile
TYPE(MeshSelection_) :: region
TYPE(Domain_), TARGET :: dom
CLASS(DirichletBC_), POINTER :: dbc => NULL()
CHARACTER(LEN=*), PARAMETER :: &
& domainFileName = "./mesh.h5"
CHARACTER( LEN = * ), PARAMETER :: &
& outfilePrefix="./output"
CHARACTER( LEN = * ), PARAMETER :: &
& statePrefix="./state"
INTEGER( I4B ), PARAMETER :: CoordinateSystem = KERNEL_2D
INTEGER( I4B ), PARAMETER :: SaveStateFrequency = _SAVE_STATE_FREQ_
INTEGER( I4B ), PARAMETER :: writeFrequency = _WRITE_FREQ_
INTEGER( I4B ) :: backupNumber
var-namecomment
objInstance of [[STDBE511_]] kernel
paramparameters instance of [[ParameterList_]]
domainFiledomain file [[HDF5File_]]
outfileoutput file [[HDF5File_]]
statefilestate file [[HDF5File_]]
regionmesh region [[MeshSelection_]]
domdomain [[Domain_]]
dbcinstance of dirichlet boundary condition [[NeumannBC_]]
domainFileNamedomain file name
outfilePrefixoutput file prefix
statePrefixstate file prefix
coordinateSystemcoordinate system
saveStateFrequencysaving frequency for state files
writeFrequencywrite frequency for output files
backupNumberinteral variable

!!! note "Linear solver parameters"

  INTEGER(I4B), PARAMETER :: LinSolver_name = LIS_GMRES
INTEGER(I4B), PARAMETER :: KrylovSubspaceSize = 50
INTEGER(I4B), PARAMETER :: maxIter_LinSolver = -1
REAL(DFP), PARAMETER :: rtol_LinSolver=1.0D-6
REAL(DFP), PARAMETER :: atol_LinSolver=1.0D-10
INTEGER(I4B), PARAMETER :: preconditionOption = NO_PRECONDITION
INTEGER(I4B), PARAMETER :: convergenceIn_LinSolver = convergenceInRes
INTEGER(I4B), PARAMETER :: convergenceType_LinSolver = relativeConvergence
LOGICAL( LGT ), PARAMETER :: relativeToRHS = .FALSE.
var-namecomment
LinSolver_namename of linear solver
KrylovSubspaceSizea parameter for linsol
maxIter_LinSolvermax number of iterations in linsol
rtol_LinSolverrelative tolerance for linsol
atol_LinSolverabsolute tolerance for linsol
preconditionOptionprecondition option
convergenceIn_LinSolverconvergence criteria for linsol
convergenceType_LinSolverconvergence criteria for linsol: relative or absolute.
relativeToRHSif convergence is relative then relative to RHS or not.

Currently, NO_PRECONDITON means diagonal preconditioning will be selected.

!!! note "Iteration parameters"

  INTEGER(I4B), PARAMETER :: nnt = 2
INTEGER( I4B ), PARAMETER :: stabParamOption = _STAB_PARAM_OPTION_
INTEGER( I4B ), PARAMETER :: maxIter= _MAX_ITER_
INTEGER(I4B), PARAMETER :: TotalTimeSteps = _TIMESTEPS_
REAL(DFP), PARAMETER :: startTime=0.0_DFP
REAL(DFP), PARAMETER :: dt= _DT_
REAL(DFP), PARAMETER :: endTime=dt*TotalTimeSteps
REAL( DFP ), PARAMETER :: atol_Pressure=1.0D-10
REAL( DFP ), PARAMETER :: rtol_Pressure=1.0D-10
REAL( DFP ), PARAMETER :: atol_Velocity=1.0D-6
REAL( DFP ), PARAMETER :: rtol_Velocity=1.0D-6
REAL( DFP ), PARAMETER :: tol_steadyState = 1.0E-8
LOGICAL( LGT ), PARAMETER :: resetIteration = .FALSE.
LOGICAL( LGT ), PARAMETER :: resetTimeStep = .FALSE.
var-namecomment
nntnumber of nodes in time
stabParamOptionstabilization parameter option
maxItermaximum number of iteration
TotalTimeStepstotal number of time steps for the simulation
startTimestart time of simulation
dttime step size of simulation
endTimeend time of simulation
atol_Pressureabsolute tolerance for checking convergence in pressure
rtol_Pressurerelative tolerance for checking convergence in pressure
atol_Velocityabsolute tolerance for checking convergence in velocity
rtol_Velocityrelative tolerance for checking convergence in velocity
tol_steadyStatetolerance for checking steady state
resetIterationreset variables to zero after each iteration
resetTimeStepreset variables to zero after each time step

!!! note "Materials and boundary conditions"

  INTEGER(I4B), PARAMETER :: tPorousMaterials = 2
INTEGER(I4B), PARAMETER :: tFluidMaterials = 1
INTEGER(I4B), PARAMETER :: tDirichletBCForVelocity = 3
INTEGER(I4B), PARAMETER :: tDirichletBCForPressure = _tDBC4P_
INTEGER( I4B ), PARAMETER :: refPressureNode = _REF_PRESSURE_NODE_
INTEGER( I4B ), PARAMETER :: materialInterfaces(1) = _MATERTIAL_INTERFACES_
REAL(DFP), PARAMETER :: porosity= _PHI_
REAL(DFP), PARAMETER :: permeability = _K0_
REAL(DFP), PARAMETER :: massdensity_soil=1700.0
REAL(DFP), PARAMETER :: massdensity_fluid = _RHO_FLUID_
REAL(DFP), PARAMETER :: dynamicViscosity = _MU_FLUID_
REAL(DFP), PARAMETER :: V_TOP= _VTOP_
var-namecomment
tPorousMaterialstotal porous media
tFluidMaterialstotal fluid material
tDirichletBCForVelocitytotal dirichlet boundary condition for velocity
tDirichletBCForPressuretotal dirichlet boundary condition for pressure
refPressureNodereference node number for pressure
porosityporosity of soil
massdensity_soilmass density of soil
permeabilitypermeability of soil
massdensity_fluidmass density of fluid
dynamicViscositydynamic viscosity of fluid
V_TOPvelocity at the top

!!! note "Internal variables"

  LOGICAL( LGT ) :: convg
INTEGER( I4B ) :: iter
INTEGER( I4B ) :: its
REAL( DFP ) :: t1, t2, t2t1

SetSTDBE511Param

  • [[STDBE511_#SetSTDBE511Param]]
  CALL FPL_INIT(); CALL param%Initiate()
!!
CALL SetSTDBE511Param( &
& param=param, &
& materialInterfaces=materialInterfaces, &
& stabParamOption=stabParamOption, &
& engine="NATIVE_SERIAL", &
& nnt=nnt, &
& startTime=startTime, &
& endTime=endTime, &
& dt=dt, &
& coordinateSystem=coordinateSystem, &
& maxIter=maxIter, &
& atoleranceForPressure=atol_pressure, &
& atoleranceForVelocity=atol_velocity, &
& rtoleranceForPressure=rtol_pressure, &
& rtoleranceForVelocity=rtol_velocity, &
& toleranceForSteadyState=tol_steadyState, &
& tPorousMaterials=tPorousMaterials, &
& tFluidMaterials=tFluidMaterials, &
& tDirichletBCForVelocity=tDirichletBCForVelocity, &
& tDirichletBCForPressure=tDirichletBCForPressure, &
& domainFile=domainFileName)

SetLinSolverParam

  CALL SetLinSolverParam( &
& param=param, &
& solverName=LinSolver_name,&
& preconditionOption=preconditionOption, &
& convergenceIn=convergenceIn_LinSolver, &
& convergenceType=convergenceType_LinSolver, &
& maxIter=maxIter_LinSolver, &
& relativeToRHS=relativeToRHS, &
& KrylovSubspaceSize=KrylovSubspaceSize, &
& rtol=rtol_LinSolver, &
& atol=atol_LinSolver )

Domain/Initiate

Domain for pressure and velocity field

  CALL domainFile%Initiate( filename=domainFileName, MODE="READ")
CALL domainFile%Open()
CALL dom%Initiate(domainFile, "")
CALL domainFile%Deallocate()

Debug

  CALL param%print()

!!! example "result"

[0] Key = LinSolver/preconditionOption, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = 0
[3] Key = LinSolver/engine, Data Type = DLCA, Dimensions = 0, Bytes = 13, Value = NATIVE_SERIAL
[3] Key = LinSolver/maxIter, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = -1
[4] Key = STDBE511/nnt, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = 2
[5] Key = STDBE511/isCommonDomain, Data Type = L, Dimensions = 0, Bytes = 4, Value = T
[5] Key = STDBE511/baseInterpolationForVelocity, Data Type = DLCA, Dimensions = 0, Bytes = 21, Value = LagrangeInterpolation
[5] Key = LinSolver/atol, Data Type = R8P, Dimensions = 0, Bytes = 8, Value = 0.100000000000000E-009
[9] Key = STDBE511/baseContinuityForPressure, Data Type = DLCA, Dimensions = 0, Bytes = 2, Value = H1
[14] Key = STDBE511/tFluidMaterials, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = 1
[14] Key = STDBE511/quadratureTypeForVelocity, Data Type = DLCA, Dimensions = 0, Bytes = 13, Value = GaussLegendre
[15] Key = STDBE511/baseInterpolationForPressure, Data Type = DLCA, Dimensions = 0, Bytes = 21, Value = LagrangeInterpolation
[22] Key = STDBE511/timeDependency, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = 2
[22] Key = LinSolver/rtol, Data Type = R8P, Dimensions = 0, Bytes = 8, Value = 0.100000000000000E-005
[23] Key = LinSolver/convergenceIn, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = 1
[24] Key = STDBE511/stabParamOption, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = 1
[24] Key = STDBE511/quadratureTypeForPressure, Data Type = DLCA, Dimensions = 0, Bytes = 13, Value = GaussLegendre
[25] Key = STDBE511/startTime, Data Type = R8P, Dimensions = 0, Bytes = 8, Value = 0.000000000000000E+000
[25] Key = LinSolver/solverName, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = 9
[26] Key = STDBE511/domainFileForVelocity, Data Type = DLCA, Dimensions = 0, Bytes = 9, Value = ./mesh.h5
[27] Key = STDBE511/tMaterialInterfaces, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = 1
[36] Key = STDBE511/domainFileForPressure, Data Type = DLCA, Dimensions = 0, Bytes = 9, Value = ./mesh.h5
[43] Key = STDBE511/materialInterfaces, Data Type = I4P, Dimensions = 1, Bytes = 4, Value = 7
[44] Key = STDBE511/baseContinuityForSpace, Data Type = DLCA, Dimensions = 0, Bytes = 2, Value = H1
[50] Key = STDBE511/baseInterpolationForSpace, Data Type = DLCA, Dimensions = 0, Bytes = 21, Value = LagrangeInterpolation
[51] Key = STDBE511/baseContinuityForTime, Data Type = DLCA, Dimensions = 0, Bytes = 2, Value = H1
[51] Key = STDBE511/toleranceForSteadyState, Data Type = R8P, Dimensions = 0, Bytes = 8, Value = 0.999999993922529E-008
[53] Key = LinSolver/KrylovSubspaceSize, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = 50
[57] Key = STDBE511/baseInterpolationForTime, Data Type = DLCA, Dimensions = 0, Bytes = 21, Value = LagrangeInterpolation
[58] Key = LinSolver/convergenceType, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = 1
[59] Key = STDBE511/quadratureTypeForSpace, Data Type = DLCA, Dimensions = 0, Bytes = 13, Value = GaussLegendre
[62] Key = STDBE511/tPorousMaterials, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = 2
[64] Key = STDBE511/atoleranceForVelocity, Data Type = R8P, Dimensions = 0, Bytes = 8, Value = 0.100000000000000E-005
[65] Key = LinSolver/relativeToRHS, Data Type = L, Dimensions = 0, Bytes = 4, Value = F
[66] Key = STDBE511/quadratureTypeForTime, Data Type = DLCA, Dimensions = 0, Bytes = 13, Value = GaussLegendre
[74] Key = STDBE511/atoleranceForPressure, Data Type = R8P, Dimensions = 0, Bytes = 8, Value = 0.100000000000000E-009
[77] Key = STDBE511/coordinateSystem, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = 3
[78] Key = STDBE511/endTime, Data Type = R8P, Dimensions = 0, Bytes = 8, Value = 0.429496735800000E+009
[79] Key = STDBE511/tDirichletBCForVelocity, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = 3
[81] Key = STDBE511/rtoleranceForVelocity, Data Type = R8P, Dimensions = 0, Bytes = 8, Value = 0.100000000000000E-005
[81] Key = STDBE511/tNeumannBCForVelocity, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = 0
[84] Key = STDBE511/dt, Data Type = R8P, Dimensions = 0, Bytes = 8, Value = 0.200000002980232E+000
[85] Key = STDBE511/name, Data Type = DLCA, Dimensions = 0, Bytes = 8, Value = STDBE511
[89] Key = STDBE511/tDirichletBCForPressure, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = 0
[91] Key = STDBE511/rtoleranceForPressure, Data Type = R8P, Dimensions = 0, Bytes = 8, Value = 0.100000000000000E-009
[91] Key = STDBE511/tNeumannBCForPressure, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = 0
[98] Key = STDBE511/engine, Data Type = DLCA, Dimensions = 0, Bytes = 13, Value = NATIVE_SERIAL
[98] Key = STDBE511/maxIter, Data Type = I4P, Dimensions = 0, Bytes = 4, Value = 20
[99] Key = STDBE511/baseContinuityForVelocity, Data Type = DLCA, Dimensions = 0, Bytes = 2, Value = H1

Next

  • [[STDBE511_test_2#STDBE511/Initiate]]
END PROGRAM main