EGPoroFullSatTimeDerivative

Kernel = biot_coefficient*d(volumetric_strain)/dt + (1/biot_modulus)*d(porepressure)/dt. This is the time-derivative for poromechanics for a single-phase, fully-saturated fluid with constant bulk modulus

Description

TODO

Example Input Syntax

[AuxKernels<<<{"href": "../../syntax/AuxKernels/index.html"}>>>]
  [stress_yy]
    type = RankTwoAux<<<{"description": "Access a component of a RankTwoTensor", "href": "../auxkernels/RankTwoAux.html"}>>>
    rank_two_tensor<<<{"description": "The rank two material tensor name"}>>> = stress
    variable<<<{"description": "The name of the variable that this object applies to"}>>> = stress_yy
    index_i<<<{"description": "The index i of ij for the tensor to output (0, 1, 2)"}>>> = 1
    index_j<<<{"description": "The index j of ij for the tensor to output (0, 1, 2)"}>>> = 1
  []
  [tot_force]
    type = ParsedAux<<<{"description": "Sets a field variable value to the evaluation of a parsed expression.", "href": "../auxkernels/ParsedAux.html"}>>>
    args = 'stress_yy EGpressure'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
    variable<<<{"description": "The name of the variable that this object applies to"}>>> = tot_force
    function = '-stress_yy+0.6*EGpressure'
  []
  [EGpressure]
    type = ParsedAux<<<{"description": "Sets a field variable value to the evaluation of a parsed expression.", "href": "../auxkernels/ParsedAux.html"}>>>
    args = 'CGpressure DGpressure'
    variable<<<{"description": "The name of the variable that this object applies to"}>>> = EGpressure
    function = 'CGpressure + DGpressure'
  []
[]

Input Parameters

displacementsThe displacements appropriate for the simulation geometry and coordinate system
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:The displacements appropriate for the simulation geometry and coordinate system
vCoupled variable
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:Coupled variable
variableThe name of the variable that this residual object operates on
C++ Type:NonlinearVariableName
Unit:(no unit assumed)
Controllable:No
Description:The name of the variable that this residual object operates on

Required Parameters

blockThe list of blocks (ids or names) that this object will be applied
C++ Type:std::vector<SubdomainName>
Controllable:No
Description:The list of blocks (ids or names) that this object will be applied
lumpingFalseTrue for mass matrix lumping, false otherwise
Default:False
C++ Type:bool
Controllable:No
Description:True for mass matrix lumping, false otherwise
matrix_onlyFalseWhether this object is only doing assembly to matrices (no vectors)
Default:False
C++ Type:bool
Controllable:No
Description:Whether this object is only doing assembly to matrices (no vectors)

Optional Parameters

absolute_value_vector_tagsThe tags for the vectors this residual object should fill with the absolute value of the residual contribution
C++ Type:std::vector<TagName>
Controllable:No
Description:The tags for the vectors this residual object should fill with the absolute value of the residual contribution
extra_matrix_tagsThe extra tags for the matrices this Kernel should fill
C++ Type:std::vector<TagName>
Controllable:No
Description:The extra tags for the matrices this Kernel should fill
extra_vector_tagsThe extra tags for the vectors this Kernel should fill
C++ Type:std::vector<TagName>
Controllable:No
Description:The extra tags for the vectors this Kernel should fill
matrix_tagssystem timeThe tag for the matrices this Kernel should fill
Default:system time
C++ Type:MultiMooseEnum
Options:nontime, system, time
Controllable:No
Description:The tag for the matrices this Kernel should fill
vector_tagstimeThe tag for the vectors this Kernel should fill
Default:time
C++ Type:MultiMooseEnum
Options:nontime, time
Controllable:No
Description:The tag for the vectors this Kernel should fill

Contribution To Tagged Field Data Parameters

control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Controllable:No
Description:Adds user-defined labels for accessing object parameters via control logic.
diag_save_inThe name of auxiliary variables to save this Kernel's diagonal Jacobian contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
C++ Type:std::vector<AuxVariableName>
Unit:(no unit assumed)
Controllable:No
Description:The name of auxiliary variables to save this Kernel's diagonal Jacobian contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Controllable:Yes
Description:Set the enabled status of the MooseObject.
implicitTrueDetermines whether this object is calculated using an implicit or explicit form
Default:True
C++ Type:bool
Controllable:No
Description:Determines whether this object is calculated using an implicit or explicit form
save_inThe name of auxiliary variables to save this Kernel's residual contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
C++ Type:std::vector<AuxVariableName>
Unit:(no unit assumed)
Controllable:No
Description:The name of auxiliary variables to save this Kernel's residual contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
seed0The seed for the master random number generator
Default:0
C++ Type:unsigned int
Controllable:No
Description:The seed for the master random number generator
use_displaced_meshFalseWhether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Default:False
C++ Type:bool
Controllable:No
Description:Whether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.

Advanced Parameters

prop_getter_suffixAn optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:An optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
use_interpolated_stateFalseFor the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.
Default:False
C++ Type:bool
Controllable:No
Description:For the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.

Material Property Retrieval Parameters

Input Files

(test/tests/enrichedGalerkin/EG_Terzaghi.i)
(test/tests/enrichedGalerkin/EG_Tides.i)
(test/tests/enrichedGalerkin/EG_Mandel.i)

(test/tests/enrichedGalerkin/EG_Mandel.i)

# Mandel's problem of consolodation of a drained medium
#
# A sample is in plane strain.
# -a <= x <= a
# -b <= y <= b
# It is squashed with constant force by impermeable, frictionless plattens on its top and bottom surfaces (at y=+/-b)
# Fluid is allowed to leak out from its sides (at x=+/-a)
# The porepressure within the sample is monitored.
#
# As is common in the literature, this is simulated by
# considering the quarter-sample, 0<=x<=a and 0<=y<=b, with
# impermeable, roller BCs at x=0 and y=0 and y=b.
# Porepressure is fixed at zero on x=a.
# Porepressure and displacement are initialised to zero.
# Then the top (y=b) is moved downwards with prescribed velocity,
# so that the total force that is inducing this downwards velocity
# is fixed.  The velocity is worked out by solving Mandel's problem
# analytically, and the total force is monitored in the simulation
# to check that it indeed remains constant.
#
# Here are the problem's parameters, and their values:
# Soil width.  a = 1
# Soil height.  b = 0.1
# Soil's Lame lambda.  la = 0.5
# Soil's Lame mu, which is also the Soil's shear modulus.  mu = G = 0.75
# Soil bulk modulus.  K = la + 2*mu/3 = 1
# Drained Poisson ratio.  nu = (3K - 2G)/(6K + 2G) = 0.2
# Soil bulk compliance.  1/K = 1
# Fluid bulk modulus.  Kf = 8
# Fluid bulk compliance.  1/Kf = 0.125
# Soil initial porosity.  phi0 = 0.1
# Biot coefficient.  alpha = 0.6
# Biot modulus.  M = 1/(phi0/Kf + (alpha - phi0)(1 - alpha)/K) = 4.705882
# Undrained bulk modulus. Ku = K + alpha^2*M = 2.694118
# Undrained Poisson ratio.  nuu = (3Ku - 2G)/(6Ku + 2G) = 0.372627
# Skempton coefficient.  B = alpha*M/Ku = 1.048035
# Fluid mobility (soil permeability/fluid viscosity).  k = 1.5
# Consolidation coefficient.  c = 2*k*B^2*G*(1-nu)*(1+nuu)^2/9/(1-nuu)/(nuu-nu) = 3.821656
# Normal stress on top.  F = 1
#
# The solution for porepressure and displacements is given in
# AHD Cheng and E Detournay "A direct boundary element method for plane strain poroelasticity" International Journal of Numerical and Analytical Methods in Geomechanics 12 (1988) 551-572.
# The solution involves complicated infinite series, so I shall not write it here

[Mesh]
  type = GeneratedMesh
  dim = 3
  nx = 10
  ny = 1
  nz = 1
  xmin = 0
  xmax = 1
  ymin = 0
  ymax = 0.1
  zmin = 0
  zmax = 1
[]

[GlobalParams]
  displacements = 'disp_x disp_y disp_z'
  block = 0
[]

[Variables]
  [disp_x]
  []
  [disp_y]
  []
  [disp_z]
  []
  [CGpressure]
  []
  [./DGpressure]
    family = MONOMIAL
    order = CONSTANT
  [../]
[]

[BCs]
  [roller_xmin]
    type = DirichletBC
    variable = disp_x
    value = 0
    boundary = 'left'
  []
  [roller_ymin]
    type = DirichletBC
    variable = disp_y
    value = 0
    boundary = 'bottom'
  []
  [plane_strain]
    type = DirichletBC
    variable = disp_z
    value = 0
    boundary = 'back front'
  []
  [./xmax_drained]
    type = DGFunctionDiffusionDirichletBC
    variable = CGpressure
    function = 0
    boundary = right
    sigma = 10
    diff = 1.5
    epsilon = 0
    # v = DGpressure
  [../]
  [./xmax_drained2]
    type = DGFunctionDiffusionDirichletBC
    variable = DGpressure
    function = 0
    boundary = right
    sigma = 10
    diff = 1.5
    epsilon = 0
    # v = DGpressure
  [../]
  [top_velocity]
    type = FunctionDirichletBC
    variable = disp_y
    function = top_velocity
    boundary = top
  []
[]

[Functions]
  [top_velocity]
    type = PiecewiseLinear
    x = '0 0.002 0.006   0.014   0.03    0.046   0.062   0.078   0.094   0.11    0.126   0.142   0.158   0.174   0.19 0.206 0.222 0.238 0.254 0.27 0.286 0.302 0.318 0.334 0.35 0.366 0.382 0.398 0.414 0.43 0.446 0.462 0.478 0.494 0.51 0.526 0.542 0.558 0.574 0.59 0.606 0.622 0.638 0.654 0.67 0.686 0.702'
    y = '-0.041824842    -0.042730269    -0.043412712    -0.04428867     -0.045509181    -0.04645965     -0.047268246 -0.047974749      -0.048597109     -0.0491467  -0.049632388     -0.050061697      -0.050441198     -0.050776675     -0.051073238      -0.0513354 -0.051567152      -0.051772022     -0.051953128 -0.052113227 -0.052254754 -0.052379865 -0.052490464 -0.052588233 -0.052674662 -0.052751065 -0.052818606 -0.052878312 -0.052931093 -0.052977751 -0.053018997 -0.053055459 -0.053087691 -0.053116185 -0.053141373 -0.05316364 -0.053183324 -0.053200724 -0.053216106 -0.053229704 -0.053241725 -0.053252351 -0.053261745 -0.053270049 -0.053277389 -0.053283879 -0.053289615'
  []
[]

[AuxVariables]
  [stress_yy]
    order = CONSTANT
    family = MONOMIAL
  []
  [tot_force]
    order = CONSTANT
    family = MONOMIAL
  []
  [EGpressure]
    family = MONOMIAL
    order = FIRST
  []
[]

[AuxKernels]
  [stress_yy]
    type = RankTwoAux
    rank_two_tensor = stress
    variable = stress_yy
    index_i = 1
    index_j = 1
  []
  [tot_force]
    type = ParsedAux
    args = 'stress_yy EGpressure'
    execute_on = timestep_end
    variable = tot_force
    function = '-stress_yy+0.6*EGpressure'
  []
  [EGpressure]
    type = ParsedAux
    args = 'CGpressure DGpressure'
    variable = EGpressure
    function = 'CGpressure + DGpressure'
  []
[]

[Kernels]
  [./grad_stress_x]
    type = StressDivergenceTensors
    variable = disp_x
    component = 0
  [../]
  [./grad_stress_y]
    type = StressDivergenceTensors
    variable = disp_y
    component = 1
  [../]
  [./grad_stress_z]
    type = StressDivergenceTensors
    variable = disp_z
    component = 2
  [../]
    [./poro_x]
    type = PoroMechanicsCoupling
    variable = disp_x
    component = 0
    porepressure = CGpressure
  [../]
  [./poro_y]
    type = PoroMechanicsCoupling
    variable = disp_y
    component = 1
    porepressure = CGpressure
  [../]
  [./poro_z]
    type = PoroMechanicsCoupling
    variable = disp_z
    component = 2
    porepressure = CGpressure
  [../]
  [./poro_xDG]
    type = PoroMechanicsCoupling
    variable = disp_x
    component = 0
    porepressure = DGpressure
  [../]
  [./poro_yDG]
    type = PoroMechanicsCoupling
    variable = disp_y
    component = 1
    porepressure = DGpressure
  [../]
  [./poro_zDG]
    type = PoroMechanicsCoupling
    variable = disp_z
    component = 2
    porepressure = DGpressure
  [../]
  [./poro_timederiv]
    type = PoroFullSatTimeDerivative
    variable = CGpressure
  [../]
  [./poro_timederiv2]
    type = PoroFullSatTimeDerivative
    variable = DGpressure
  [../]
  [./darcy_flow]
    type = CoefDiffusion
    variable = CGpressure
    coef = 1.5
  [../]
  [./time_derivEG1]
    type = EGPoroFullSatTimeDerivative
    variable = CGpressure
    v = DGpressure
  [../]
  [./time_derivEG2]
    type = EGPoroFullSatTimeDerivative
    variable = DGpressure
    v = CGpressure
  [../]
[]
[DGKernels]
  [./DG]
    type = DGDiffusion
    variable = DGpressure
    sigma = 10
    diff = 1.5
    epsilon = 0
  [../]
  [./EG1]
    type = EGDiffusion
    variable = DGpressure
    sigma = 10
    diff = 1.5
    epsilon = 0
    v = CGpressure
  [../]
  [./EG2]
    type = EGDiffusion
    variable = CGpressure
    v = DGpressure
    sigma = 10
    diff = 1.5
    epsilon = 0
  [../]
[]

[Materials]
  [elasticity_tensor]
    type = ComputeElasticityTensor
    C_ijkl = '0.5 0.75'
    # bulk modulus is lambda + 2*mu/3 = 0.5 + 2*0.75/3 = 1
    fill_method = symmetric_isotropic
  []
  [./strain]
    type = ComputeSmallStrain
    displacements = 'disp_x disp_y disp_z'
  [../]
  [./stress]
    type = ComputeLinearElasticStress
  [../]
  [./poro_material]
    type = PoroFullSatMaterial
    porosity0 = 0.1
    biot_coefficient = 0.6
    solid_bulk_compliance = 1
    fluid_bulk_compliance = 0.125
    constant_porosity = true
    porepressure = EGpressure
  [../]
[]

[Postprocessors]
  [./p0]
    type = PointValue
    outputs = csv
    point = '0.0 0 0'
    variable = EGpressure
  [../]
  [./p1]
    type = PointValue
    outputs = csv
    point = '0.1 0 0'
    variable = EGpressure
  [../]
  [./p2]
    type = PointValue
    outputs = csv
    point = '0.2 0 0'
    variable = EGpressure
  [../]
  [./p3]
    type = PointValue
    outputs = csv
    point = '0.3 0 0'
    variable = EGpressure
  [../]
  [./p4]
    type = PointValue
    outputs = csv
    point = '0.4 0 0'
    variable = EGpressure
  [../]
  [./p5]
    type = PointValue
    outputs = csv
    point = '0.5 0 0'
    variable = EGpressure
  [../]
  [./p6]
    type = PointValue
    outputs = csv
    point = '0.6 0 0'
    variable = EGpressure
  [../]
  [./p7]
    type = PointValue
    outputs = csv
    point = '0.7 0 0'
    variable = EGpressure
  [../]
  [./p8]
    type = PointValue
    outputs = csv
    point = '0.8 0 0'
    variable = EGpressure
  [../]
  [./p9]
    type = PointValue
    outputs = csv
    point = '0.9 0 0'
    variable = EGpressure
  [../]
  [./p99]
    type = PointValue
    outputs = csv
    point = '1 0 0'
    variable = EGpressure
  [../]
  [./xdisp]
    type = PointValue
    outputs = csv
    point = '1 0.1 0'
    variable = disp_x
  [../]
  [./ydisp]
    type = PointValue
    outputs = csv
    point = '1 0.1 0'
    variable = disp_y
  [../]
  [./total_downwards_force]
     type = ElementAverageValue
     outputs = csv
     variable = tot_force
  [../]

  [./dt]
    type = FunctionValuePostprocessor
    outputs = console
    function = if(0.15*t<0.01,0.15*t,0.01)
  [../]
[]

[Executioner]
  type = Transient
  solve_type =  NEWTON
  start_time = 0
  end_time = 0.7

  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
  petsc_options_value = 'lu       superlu_dist'

  [TimeStepper]
    type = PostprocessorDT
    postprocessor = dt
    dt = 0.001
  []
[]

[Outputs]
  execute_on = 'timestep_end'
  file_base = EG_Mandel
  exodus = false
  perf_graph = false
  [csv]
    interval = 3
    type = CSV
  []
[]

(test/tests/enrichedGalerkin/EG_Terzaghi.i)

# Terzaghi's problem of consolodation of a drained medium
#
# A saturated soil sample sits in a bath of water.
# It is constrained on its sides, and bottom.
# Its sides and bottom are also impermeable.
# Initially it is unstressed.
# A normal stress, q, is applied to the soil's top.
# The soil then slowly compresses as water is squeezed
# out from the sample from its top (the top BC for
# the porepressure is porepressure = 0).
#
# See, for example.  Section 2.2 of the online manuscript
# Arnold Verruijt "Theory and Problems of Poroelasticity" Delft University of Technology 2013
# but note that the "sigma" in that paper is the negative
# of the stress in TensorMechanics
#
# Here are the problem's parameters, and their values:
# Soil height.  h = 10
# Soil's Lame lambda.  la = 2
# Soil's Lame mu, which is also the Soil's shear modulus.  mu = 3
# Soil bulk modulus.  K = la + 2*mu/3 = 4
# Soil confined compressibility.  m = 1/(K + 4mu/3) = 0.125
# Soil bulk compliance.  1/K = 0.25
# Fluid bulk modulus.  Kf = 8
# Fluid bulk compliance.  1/Kf = 0.125
# Fluid mobility (soil permeability/fluid viscosity).  k = 1.5
# Soil initial porosity.  phi0 = 0.1
# Biot coefficient.  alpha = 0.6
# Soil initial storativity, which is the reciprocal of the initial Biot modulus.  S = phi0/Kf + (alpha - phi0)(1 - alpha)/K = 0.0625
# Consolidation coefficient.  c = k/(S + alpha^2 m) = 13.95348837
# Normal stress on top.  q = 1
# Initial porepressure, resulting from instantaneous application of q, assuming corresponding instantaneous increase of porepressure (Note that this is calculated by MOOSE: we only need it for the analytical solution).  p0 = alpha*m*q/(S + alpha^2 m) = 0.69767442
# Initial vertical displacement (down is positive), resulting from instantaneous application of q (Note this is calculated by MOOSE: we only need it for the analytical solution).  uz0 = q*m*h*S/(S + alpha^2 m)
# Final vertical displacement (down in positive) (Note this is calculated by MOOSE: we only need it for the analytical solution).  uzinf = q*m*h
#
# The solution for porepressure is
# P = 4*p0/\pi \sum_{k=1}^{\infty} \frac{(-1)^{k-1}}{2k-1} \cos ((2k-1)\pi z/(2h)) \exp(-(2k-1)^2 \pi^2 ct/(4 h^2))
# This series converges very slowly for ct/h^2 small, so in that domain
# P = p0 erf( (1-(z/h))/(2 \sqrt(ct/h^2)) )
#
# The degree of consolidation is defined as
# U = (uz - uz0)/(uzinf - uz0)
# where uz0 and uzinf are defined above, and
# uz = the vertical displacement of the top (down is positive)
# U = 1 - (8/\pi^2)\sum_{k=1}^{\infty} \frac{1}{(2k-1)^2} \exp(-(2k-1)^2 \pi^2 ct/(4 h^2))

[Mesh]
  type = GeneratedMesh
  dim = 3
  nx = 1
  ny = 1
  nz = 10
  xmin = -1
  xmax = 1
  ymin = -1
  ymax = 1
  zmin = 0
  zmax = 10
[]

[GlobalParams]
  displacements = 'disp_x disp_y disp_z'
  block = 0
[]

[Variables]
  [disp_x]
  []
  [disp_y]
  []
  [disp_z]
  []
  [CGpressure]
  []
  [./DGpressure]
    family = MONOMIAL
    order = CONSTANT
  [../]
[]

[BCs]
  [./confinex]
    type = DirichletBC
    variable = disp_x
    value = 0
    boundary = 'left right'
  [../]
  [./confiney]
    type = DirichletBC
    variable = disp_y
    value = 0
    boundary = 'bottom top'
  [../]
  [./basefixed]
    type = DirichletBC
    variable = disp_z
    value = 0
    boundary = back
  [../]
  [./topdrained]
    type = DGFunctionDiffusionDirichletBC
    variable = CGpressure
    function = 0
    boundary = front
    sigma = 10
    diff = 1.5
    epsilon = 0
  [../]
  [./topdrained2]
    type = DGFunctionDiffusionDirichletBC
    variable = DGpressure
    function = 0
    boundary = front
    sigma = 10
    diff = 1.5
    epsilon = 0
  [../]
  [./topload]
    type = NeumannBC
    variable = disp_z
    value = -1
    boundary = front
  [../]
[]

[AuxVariables]
  [EGpressure]
    family = MONOMIAL
    order = FIRST
  []
[]

[AuxKernels]
  [EGpressure]
    type = ParsedAux
    args = 'CGpressure DGpressure'
    variable = EGpressure
    function = 'CGpressure + DGpressure'
  []
[]

[Kernels]
  [./grad_stress_x]
    type = StressDivergenceTensors
    variable = disp_x
    component = 0
  [../]
  [./grad_stress_y]
    type = StressDivergenceTensors
    variable = disp_y
    component = 1
  [../]
  [./grad_stress_z]
    type = StressDivergenceTensors
    variable = disp_z
    component = 2
  [../]
    [./poro_x]
    type = PoroMechanicsCoupling
    variable = disp_x
    component = 0
    porepressure = CGpressure
  [../]
  [./poro_y]
    type = PoroMechanicsCoupling
    variable = disp_y
    component = 1
    porepressure = CGpressure
  [../]
  [./poro_z]
    type = PoroMechanicsCoupling
    variable = disp_z
    component = 2
    porepressure = CGpressure
  [../]
  [./poro_xDG]
    type = PoroMechanicsCoupling
    variable = disp_x
    component = 0
    porepressure = DGpressure
  [../]
  [./poro_yDG]
    type = PoroMechanicsCoupling
    variable = disp_y
    component = 1
    porepressure = DGpressure
  [../]
  [./poro_zDG]
    type = PoroMechanicsCoupling
    variable = disp_z
    component = 2
    porepressure = DGpressure
  [../]
  [./poro_timederiv]
    type = PoroFullSatTimeDerivative
    variable = CGpressure
  [../]
  [./poro_timederiv2]
    type = PoroFullSatTimeDerivative
    variable = DGpressure
  [../]
  [./darcy_flow]
    type = CoefDiffusion
    variable = CGpressure
    coef = 1.5
  [../]
  [./time_derivEG1]
    type = EGPoroFullSatTimeDerivative
    variable = CGpressure
    v = DGpressure
  [../]
  [./time_derivEG2]
    type = EGPoroFullSatTimeDerivative
    variable = DGpressure
    v = CGpressure
  [../]
[]
[DGKernels]
  [./DG]
    type = DGDiffusion
    variable = DGpressure
    sigma = 1
    diff = 1.5
    epsilon = 0
  [../]
  [./EG1]
    type = EGDiffusion
    variable = DGpressure
    sigma = 1
    diff = 1.5
    epsilon = 0
    v = CGpressure
  [../]
  [./EG2]
    type = EGDiffusion
    variable = CGpressure
    v = DGpressure
    sigma = 1
    diff = 1.5
    epsilon = 0
  [../]
[]

[Materials]
  [elasticity_tensor]
    type = ComputeElasticityTensor
    C_ijkl = '2 3'
    # bulk modulus is lambda + 2*mu/3 = 2 + 2*3/3 = 4
    fill_method = symmetric_isotropic
  []
  [./strain]
    type = ComputeSmallStrain
    displacements = 'disp_x disp_y disp_z'
  [../]
  [./stress]
    type = ComputeLinearElasticStress
  [../]
  [./poro_material]
    type = PoroFullSatMaterial
    porosity0 = 0.1
    biot_coefficient = 0.6
    solid_bulk_compliance = 0.25
    fluid_bulk_compliance = 0.125
    constant_porosity = true
    porepressure = EGpressure
  [../]
[]

[Postprocessors]
  [./p0]
    type = PointValue
    outputs = csv
    point = '0 0 0'
    variable = EGpressure
  [../]
  [./p1]
    type = PointValue
    outputs = csv
    point = '0 0 1'
    variable = EGpressure
  [../]
  [./p2]
    type = PointValue
    outputs = csv
    point = '0 0 2'
    variable = EGpressure
  [../]
  [./p3]
    type = PointValue
    outputs = csv
    point = '0 0 3'
    variable = EGpressure
  [../]
  [./p4]
    type = PointValue
    outputs = csv
    point = '0 0 4'
    variable = EGpressure
  [../]
  [./p5]
    type = PointValue
    outputs = csv
    point = '0 0 5'
    variable = EGpressure
  [../]
  [./p6]
    type = PointValue
    outputs = csv
    point = '0 0 6'
    variable = EGpressure
  [../]
  [./p7]
    type = PointValue
    outputs = csv
    point = '0 0 7'
    variable = EGpressure
  [../]
  [./p8]
    type = PointValue
    outputs = csv
    point = '0 0 8'
    variable = EGpressure
  [../]
  [./p9]
    type = PointValue
    outputs = csv
    point = '0 0 9'
    variable = EGpressure
  [../]
  [./p99]
    type = PointValue
    outputs = csv
    point = '0 0 10'
    variable = EGpressure
  [../]
  [./zdisp]
    type = PointValue
    outputs = csv
    point = '0 0 10'
    variable = disp_z
  [../]
  [./dt]
    type = FunctionValuePostprocessor
    outputs = console
    function = if(0.5*t<0.1,0.5*t,0.1)
  [../]
[]

[Executioner]
  type = Transient
  solve_type =  NEWTON
  start_time = 0
  end_time = 10

  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
  petsc_options_value = 'lu       superlu_dist'

  [./TimeStepper]
    type = PostprocessorDT
    postprocessor = dt
    dt = 0.0001
  [../]
[]

[Outputs]
  execute_on = 'timestep_end'
  file_base = EG_Terzaghi
  exodus = false
  perf_graph = false
  [./csv]
    type = CSV
  [../]
[]

(test/tests/enrichedGalerkin/EG_Tides.i)

# A confined aquifer is fully saturated with water
# Earth tides apply strain to the aquifer and the resulting porepressure changes are recorded
#
# To replicate standard poroelasticity exactly:
# (1) the PorousFlowBasicTHM Action is used;
# (2) multiply_by_density = false;
# (3) PorousFlowConstantBiotModulus is used
[Mesh]
  type = GeneratedMesh
  dim = 3
  nx = 1
  ny = 1
  nz = 1
  xmin = 0
  xmax = 1
  ymin = 0
  ymax = 1
  zmin = 0
  zmax = 1
[]

[GlobalParams]
  displacements = 'disp_x disp_y disp_z'
  block = 0
[]

[Variables]
  [disp_x]
  []
  [disp_y]
  []
  [disp_z]
  []
  [CGpressure]
  []
  [./DGpressure]
    family = MONOMIAL
    order = CONSTANT
  [../]
[]

[BCs]
  # zmin is called 'back'
  # zmax is called 'front'
  # ymin is called 'bottom'
  # ymax is called 'top'
  # xmin is called 'left'
  # xmax is called 'right'
  [strain_x]
    type = FunctionDirichletBC
    variable = disp_x
    function = earth_tide_x
    boundary = 'left right'
  []
  [strain_y]
    type = FunctionDirichletBC
    variable = disp_y
    function = earth_tide_y
    boundary = 'bottom top'
  []
  [strain_z]
    type = FunctionDirichletBC
    variable = disp_z
    function = earth_tide_z
    boundary = 'back front'
  []
[]

[AuxVariables]
  [EGpressure]
    family = MONOMIAL
    order = FIRST
  []
[]

[AuxKernels]
  [EGpressure]
    type = ParsedAux
    args = 'CGpressure DGpressure'
    variable = EGpressure
    function = 'CGpressure + DGpressure'
  []
[]

[Kernels]
  [./grad_stress_x]
    type = StressDivergenceTensors
    variable = disp_x
    component = 0
  [../]
  [./grad_stress_y]
    type = StressDivergenceTensors
    variable = disp_y
    component = 1
  [../]
  [./grad_stress_z]
    type = StressDivergenceTensors
    variable = disp_z
    component = 2
  [../]
    [./poro_x]
    type = PoroMechanicsCoupling
    variable = disp_x
    component = 0
    porepressure = CGpressure
  [../]
  [./poro_y]
    type = PoroMechanicsCoupling
    variable = disp_y
    component = 1
    porepressure = CGpressure
  [../]
  [./poro_z]
    type = PoroMechanicsCoupling
    variable = disp_z
    component = 2
    porepressure = CGpressure
  [../]
  [./poro_xDG]
    type = PoroMechanicsCoupling
    variable = disp_x
    component = 0
    porepressure = DGpressure
  [../]
  [./poro_yDG]
    type = PoroMechanicsCoupling
    variable = disp_y
    component = 1
    porepressure = DGpressure
  [../]
  [./poro_zDG]
    type = PoroMechanicsCoupling
    variable = disp_z
    component = 2
    porepressure = DGpressure
  [../]
  [./poro_timederiv]
    type = PoroFullSatTimeDerivative
    variable = CGpressure
  [../]
  [./poro_timederiv2]
    type = PoroFullSatTimeDerivative
    variable = DGpressure
  [../]
  [./darcy_flow]
    type = DarcyFluxPressure
    variable = CGpressure
    gravity = '0 0 0'
  [../]
  [./time_derivEG1]
    type = EGPoroFullSatTimeDerivative
    variable = CGpressure
    v = DGpressure
  [../]
  [./time_derivEG2]
    type = EGPoroFullSatTimeDerivative
    variable = DGpressure
    v = CGpressure
  [../]
[]
[DGKernels]
  [./DG]
    type = DGDiffusion
    variable = DGpressure
    sigma = 10
    diff = conductivity
    epsilon = 0
  [../]
  [./EG1]
    type = EGDiffusion
    variable = DGpressure
    sigma = 10
    diff = conductivity
    epsilon = 0
    v = CGpressure
  [../]
  [./EG2]
    type = EGDiffusion
    variable = CGpressure
    v = DGpressure
    sigma = 10
    diff = conductivity
    epsilon = 0
  [../]
[]

[Functions]
  [earth_tide_x]
    type = ParsedFunction
    value = 'x*1E-8*(5*cos(t*2*pi) + 2*cos((t-0.5)*2*pi) + 1*cos((t+0.3)*0.5*pi))'
  []
  [earth_tide_y]
    type = ParsedFunction
    value = 'y*1E-8*(7*cos(t*2*pi) + 4*cos((t-0.3)*2*pi) + 7*cos((t+0.6)*0.5*pi))'
  []
  [earth_tide_z]
    type = ParsedFunction
    value = 'z*1E-8*(7*cos((t-0.5)*2*pi) + 4*cos((t-0.8)*2*pi) + 7*cos((t+0.1)*4*pi))'
  []
[]

[Materials]
  [./elasticity_tensor]
    type = ComputeIsotropicElasticityTensor
    bulk_modulus = 10.0E9 # drained bulk modulus
    poissons_ratio = 0.25
  [../]
  [./strain]
    type = ComputeSmallStrain
  [../]
  [./stress]
    type = ComputeLinearElasticStress
  [../]
  [./poro_material]
    type = PoroFullSatMaterial
    porosity0 = 0.1
    biot_coefficient = 0.6
    solid_bulk_compliance = 1E-10
    fluid_bulk_compliance = 5E-10
    constant_porosity = true
    porepressure = EGpressure
  [../]
  [density]
    type = GenericConstantMaterial
    prop_names = density
    prop_values = 1000.0
  []
  [conductivity]
    type = GenericConstantMaterial
    prop_names = conductivity
    prop_values = 1E-9
  []
[]

[Postprocessors]
  [pp]
    type = PointValue
    point = '0.5 0.5 0.5'
    variable = EGpressure
  []
[]

[Executioner]
  type = Transient
  solve_type = Newton
  dt = 0.01
  end_time = 2

  petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
  petsc_options_value = 'lu       superlu_dist'
[]

[Outputs]
  file_base = EG_Tides
  console = true
  csv = true
[]