Density-driven convective mixing
PorousFlow may be used to study the convective mixing of fluids. In this page we discuss the density-driven convective mixing of CO into brine. CO exists initially in both the liquid and gas phase. As the CO diffuses into the brine, the gaseous CO disappears and the brine's density increases, which drives convective mixing of the two fluids and accelerates the dissolution of the gas phase (Emami-Meybodi et al., 2015). Input files for similar situations may be found here.
These examples should be run using multiple cores due to the computational expense. It is recommended that the single-phase model is run using at least four cores, and the two-phase model with at least eight cores.
Because this problem involves the disappearance of the gaseous phase, the simulation utilises persistent primary variables and a vapor-liquid flash.
[AuxVariables]
[xnacl]
initial_condition = 0.01
[]
[saturation_gas]
order = FIRST
family = MONOMIAL
[]
[xco2l]
order = FIRST
family = MONOMIAL
[]
[density_liquid]
order = FIRST
family = MONOMIAL
[]
[porosity]
order = CONSTANT
family = MONOMIAL
[]
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)Here pgas
is the gas-state pressure (the capillary pressure is zero, so this is also the brine porepressure) and zi
is the total mass fraction of the CO summed over the two phases. In this input file, the salt content of the brine is fixed at 0.01:
[xnacl]
initial_condition = 0.01
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)and the temperature at 45C:
[temperature]
type = PorousFlowTemperature
temperature = '45'
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)The porepressure is initialised to approximate hydrostatic conditions, and the total mass fraction of CO is zero except towards the top of the model domain where it is 0.2. This corresponds to a gas saturation of 0.29 (the gas is almost all CO and only about 0.2% brine) and the mass-fraction of CO in the liquid phase is 0.048. The porosity is slightly randomised by sampling from a uniform distribution between the minimum and maximum values specified in the RandomIC initial condition, and this promotes initiation of the beautiful "fingers" of CO-rich brine seen in the results. The random field seed can be changed in the input file, or on the command line using command line overriding. The porosity field is output in the Exodus file that is created, so can be visualized in the MOOSE GUI peacock or in an external visualization package.
[ICs]
[pressure]
type = FunctionIC
function = 10e6-9.81*1000*y
variable = pgas
[]
[zi]
type = BoundingBoxIC
variable = zi
x1 = 0
x2 = 2
y1 = 1.95
y2 = 2
inside = 0.2
outside = 0
[]
[porosity]
type = RandomIC
variable = porosity
min = 0.25
max = 0.275
seed = 0
[]
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)The equations of state for the liquid brine, the CO and the liquid brine-CO mixture are the high-precision versions supplied by the fluid_properties module:
[Modules]
[FluidProperties]
[co2sw]
type = CO2FluidProperties
[]
[co2]
type = TabulatedFluidProperties
fp = co2sw
[]
[brine]
type = BrineFluidProperties
[]
[]
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i) [fs]
type = PorousFlowBrineCO2
brine_fp = brine
co2_fp = co2
capillary_pressure = pc
(modules/porous_flow/examples/lava_lamp/2phase_convection.i) [brineco2]
type = PorousFlowFluidState
gas_porepressure = 'pgas'
z = 'zi'
temperature_unit = Celsius
xnacl = 'xnacl'
capillary_pressure = pc
fluid_state = fs
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)Quantities of interest are tracked using AuxVariables
populated using the following AuxKernels
[AuxKernels]
[saturation_gas]
type = PorousFlowPropertyAux
variable = saturation_gas
property = saturation
phase = 1
execute_on = 'timestep_end'
[]
[xco2l]
type = PorousFlowPropertyAux
variable = xco2l
property = mass_fraction
phase = 0
fluid_component = 1
execute_on = 'timestep_end'
[]
[density_liquid]
type = PorousFlowPropertyAux
variable = density_liquid
property = density
phase = 0
execute_on = 'timestep_end'
[]
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)The usual 2-phase Kernels are used in this example with zero molecular dispersion (but the diffusion coefficient is nonzero)
[Kernels]
[mass0]
type = PorousFlowMassTimeDerivative
fluid_component = 0
variable = pgas
[]
[flux0]
type = PorousFlowAdvectiveFlux
fluid_component = 0
variable = pgas
[]
[diff0]
type = PorousFlowDispersiveFlux
fluid_component = 0
variable = pgas
disp_long = '0 0'
disp_trans = '0 0'
[]
[mass1]
type = PorousFlowMassTimeDerivative
fluid_component = 1
variable = zi
[]
[flux1]
type = PorousFlowAdvectiveFlux
fluid_component = 1
variable = zi
[]
[diff1]
type = PorousFlowDispersiveFlux
fluid_component = 1
variable = zi
disp_long = '0 0'
disp_trans = '0 0'
[]
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i) [diffusivity]
type = PorousFlowDiffusivityConst
diffusion_coeff = '2e-9 2e-9 2e-9 2e-9'
tortuosity = '1 1'
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)Various postprocessors help with tracking the total quantity of CO in the two phases, and the flux of CO into the brine:
[Postprocessors]
[total_co2_in_gas]
type = PorousFlowFluidMass
phase = 1
fluid_component = 1
[]
[total_co2_in_liquid]
type = PorousFlowFluidMass
phase = 0
fluid_component = 1
[]
[numdofs]
type = NumDOFs
[]
[delta_xco2]
type = ChangeOverTimePostprocessor
postprocessor = total_co2_in_liquid
[]
[dt]
type = TimestepSize
[]
[flux]
type = FunctionValuePostprocessor
function = flux
[]
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)Mesh adaptivity is used to refine the mesh as the fingers form, based on the jump in the gradient of the total mass fraction of CO summed over the two phases:
[Adaptivity]
max_h_level = 2
marker = marker
initial_marker = initial
initial_steps = 2
[Indicators]
[indicator]
type = GradientJumpIndicator
variable = zi
[]
[]
[Markers]
[marker]
type = ErrorFractionMarker
indicator = indicator
refine = 0.8
[]
[initial]
type = BoxMarker
bottom_left = '0 1.95 0'
top_right = '2 2 0'
inside = REFINE
outside = DO_NOTHING
[]
[]
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)The evolution of the CO mass fraction is shown in the animation below

Figure 1: Left: Mass fraction of CO in the liquid phase. Right: Saturation of gas. The adaptive mesh is overlayed on each figure.
References
- Hamid Emami-Meybodi, Hassan Hassanzadeh, Christopher P Green, and Jonathan Ennis-King.
Convective dissolution of co$_2$ in saline aquifers: progress in modeling and experiments.
International Journal of Greenhouse Gas Control, 40:238–266, 2015.[BibTeX]
@article{emami, author = "Emami-Meybodi, Hamid and Hassanzadeh, Hassan and Green, Christopher P and Ennis-King, Jonathan", title = "Convective dissolution of CO$\_2$ in saline aquifers: Progress in modeling and experiments", journal = "International Journal of Greenhouse Gas Control", volume = "40", pages = "238--266", year = "2015", publisher = "Elsevier" }
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)
# Two phase density-driven convection of dissolved CO2 in brine
#
# Initially, the model has a gas phase at the top with a saturation of 0.29
# (which corresponds to an initial value of zi = 0.2).
# Diffusion of the dissolved CO2
# component from the saturated liquid to the unsaturated liquid below reduces the
# amount of CO2 in the gas phase. As the density of the CO2-saturated brine is greater
# than the unsaturated brine, a gravitational instability arises and density-driven
# convection of CO2-rich fingers descend into the unsaturated brine.
#
# The instability is seeded by a random perturbation to the porosity field.
# Mesh adaptivity is used to refine the mesh as the fingers form.
#
# Note: this model is computationally expensive, so should be run with multiple cores,
# preferably on a cluster.
[GlobalParams]
PorousFlowDictator = 'dictator'
gravity = '0 -9.81 0'
[]
[Adaptivity]
max_h_level = 2
marker = marker
initial_marker = initial
initial_steps = 2
[Indicators]
[indicator]
type = GradientJumpIndicator
variable = zi
[]
[]
[Markers]
[marker]
type = ErrorFractionMarker
indicator = indicator
refine = 0.8
[]
[initial]
type = BoxMarker
bottom_left = '0 1.95 0'
top_right = '2 2 0'
inside = REFINE
outside = DO_NOTHING
[]
[]
[]
[Mesh]
type = GeneratedMesh
dim = 2
ymax = 2
xmax = 2
ny = 40
nx = 40
bias_y = 0.95
[]
[Kernels]
[mass0]
type = PorousFlowMassTimeDerivative
fluid_component = 0
variable = pgas
[]
[flux0]
type = PorousFlowAdvectiveFlux
fluid_component = 0
variable = pgas
[]
[diff0]
type = PorousFlowDispersiveFlux
fluid_component = 0
variable = pgas
disp_long = '0 0'
disp_trans = '0 0'
[]
[mass1]
type = PorousFlowMassTimeDerivative
fluid_component = 1
variable = zi
[]
[flux1]
type = PorousFlowAdvectiveFlux
fluid_component = 1
variable = zi
[]
[diff1]
type = PorousFlowDispersiveFlux
fluid_component = 1
variable = zi
disp_long = '0 0'
disp_trans = '0 0'
[]
[]
[AuxVariables]
[xnacl]
initial_condition = 0.01
[]
[saturation_gas]
order = FIRST
family = MONOMIAL
[]
[xco2l]
order = FIRST
family = MONOMIAL
[]
[density_liquid]
order = FIRST
family = MONOMIAL
[]
[porosity]
order = CONSTANT
family = MONOMIAL
[]
[]
[AuxKernels]
[saturation_gas]
type = PorousFlowPropertyAux
variable = saturation_gas
property = saturation
phase = 1
execute_on = 'timestep_end'
[]
[xco2l]
type = PorousFlowPropertyAux
variable = xco2l
property = mass_fraction
phase = 0
fluid_component = 1
execute_on = 'timestep_end'
[]
[density_liquid]
type = PorousFlowPropertyAux
variable = density_liquid
property = density
phase = 0
execute_on = 'timestep_end'
[]
[]
[Variables]
[pgas]
[]
[zi]
scaling = 1e4
[]
[]
[ICs]
[pressure]
type = FunctionIC
function = 10e6-9.81*1000*y
variable = pgas
[]
[zi]
type = BoundingBoxIC
variable = zi
x1 = 0
x2 = 2
y1 = 1.95
y2 = 2
inside = 0.2
outside = 0
[]
[porosity]
type = RandomIC
variable = porosity
min = 0.25
max = 0.275
seed = 0
[]
[]
[UserObjects]
[dictator]
type = PorousFlowDictator
porous_flow_vars = 'pgas zi'
number_fluid_phases = 2
number_fluid_components = 2
[]
[pc]
type = PorousFlowCapillaryPressureConst
pc = 0
[]
[fs]
type = PorousFlowBrineCO2
brine_fp = brine
co2_fp = co2
capillary_pressure = pc
[]
[]
[Modules]
[FluidProperties]
[co2sw]
type = CO2FluidProperties
[]
[co2]
type = TabulatedFluidProperties
fp = co2sw
[]
[brine]
type = BrineFluidProperties
[]
[]
[]
[Materials]
[temperature]
type = PorousFlowTemperature
temperature = '45'
[]
[brineco2]
type = PorousFlowFluidState
gas_porepressure = 'pgas'
z = 'zi'
temperature_unit = Celsius
xnacl = 'xnacl'
capillary_pressure = pc
fluid_state = fs
[]
[porosity]
type = PorousFlowPorosityConst
porosity = porosity
[]
[permeability]
type = PorousFlowPermeabilityConst
permeability = '1e-11 0 0 0 1e-11 0 0 0 1e-11'
[]
[relperm_water]
type = PorousFlowRelativePermeabilityCorey
phase = 0
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[relperm_gas]
type = PorousFlowRelativePermeabilityCorey
phase = 1
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[diffusivity]
type = PorousFlowDiffusivityConst
diffusion_coeff = '2e-9 2e-9 2e-9 2e-9'
tortuosity = '1 1'
[]
[]
[Preconditioning]
active = basic
[mumps_is_best_for_parallel_jobs]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = ' lu mumps'
[]
[basic]
type = SMP
full = true
petsc_options_iname = '-ksp_type -pc_type -sub_pc_type -sub_pc_factor_shift_type -pc_asm_overlap'
petsc_options_value = 'gmres asm lu NONZERO 2 '
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
end_time = 1e6
nl_max_its = 25
l_max_its = 100
dtmax = 1e4
nl_abs_tol = 1e-6
[TimeStepper]
type = IterationAdaptiveDT
dt = 10
growth_factor = 2
cutback_factor = 0.5
[]
[]
[Functions]
[flux]
type = ParsedFunction
vals = 'delta_xco2 dt'
vars = 'dx dt'
value = 'dx/dt'
[]
[]
[Postprocessors]
[total_co2_in_gas]
type = PorousFlowFluidMass
phase = 1
fluid_component = 1
[]
[total_co2_in_liquid]
type = PorousFlowFluidMass
phase = 0
fluid_component = 1
[]
[numdofs]
type = NumDOFs
[]
[delta_xco2]
type = ChangeOverTimePostprocessor
postprocessor = total_co2_in_liquid
[]
[dt]
type = TimestepSize
[]
[flux]
type = FunctionValuePostprocessor
function = flux
[]
[]
[Outputs]
print_linear_residuals = false
perf_graph = true
exodus = true
csv = true
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)
# Two phase density-driven convection of dissolved CO2 in brine
#
# Initially, the model has a gas phase at the top with a saturation of 0.29
# (which corresponds to an initial value of zi = 0.2).
# Diffusion of the dissolved CO2
# component from the saturated liquid to the unsaturated liquid below reduces the
# amount of CO2 in the gas phase. As the density of the CO2-saturated brine is greater
# than the unsaturated brine, a gravitational instability arises and density-driven
# convection of CO2-rich fingers descend into the unsaturated brine.
#
# The instability is seeded by a random perturbation to the porosity field.
# Mesh adaptivity is used to refine the mesh as the fingers form.
#
# Note: this model is computationally expensive, so should be run with multiple cores,
# preferably on a cluster.
[GlobalParams]
PorousFlowDictator = 'dictator'
gravity = '0 -9.81 0'
[]
[Adaptivity]
max_h_level = 2
marker = marker
initial_marker = initial
initial_steps = 2
[Indicators]
[indicator]
type = GradientJumpIndicator
variable = zi
[]
[]
[Markers]
[marker]
type = ErrorFractionMarker
indicator = indicator
refine = 0.8
[]
[initial]
type = BoxMarker
bottom_left = '0 1.95 0'
top_right = '2 2 0'
inside = REFINE
outside = DO_NOTHING
[]
[]
[]
[Mesh]
type = GeneratedMesh
dim = 2
ymax = 2
xmax = 2
ny = 40
nx = 40
bias_y = 0.95
[]
[Kernels]
[mass0]
type = PorousFlowMassTimeDerivative
fluid_component = 0
variable = pgas
[]
[flux0]
type = PorousFlowAdvectiveFlux
fluid_component = 0
variable = pgas
[]
[diff0]
type = PorousFlowDispersiveFlux
fluid_component = 0
variable = pgas
disp_long = '0 0'
disp_trans = '0 0'
[]
[mass1]
type = PorousFlowMassTimeDerivative
fluid_component = 1
variable = zi
[]
[flux1]
type = PorousFlowAdvectiveFlux
fluid_component = 1
variable = zi
[]
[diff1]
type = PorousFlowDispersiveFlux
fluid_component = 1
variable = zi
disp_long = '0 0'
disp_trans = '0 0'
[]
[]
[AuxVariables]
[xnacl]
initial_condition = 0.01
[]
[saturation_gas]
order = FIRST
family = MONOMIAL
[]
[xco2l]
order = FIRST
family = MONOMIAL
[]
[density_liquid]
order = FIRST
family = MONOMIAL
[]
[porosity]
order = CONSTANT
family = MONOMIAL
[]
[]
[AuxKernels]
[saturation_gas]
type = PorousFlowPropertyAux
variable = saturation_gas
property = saturation
phase = 1
execute_on = 'timestep_end'
[]
[xco2l]
type = PorousFlowPropertyAux
variable = xco2l
property = mass_fraction
phase = 0
fluid_component = 1
execute_on = 'timestep_end'
[]
[density_liquid]
type = PorousFlowPropertyAux
variable = density_liquid
property = density
phase = 0
execute_on = 'timestep_end'
[]
[]
[Variables]
[pgas]
[]
[zi]
scaling = 1e4
[]
[]
[ICs]
[pressure]
type = FunctionIC
function = 10e6-9.81*1000*y
variable = pgas
[]
[zi]
type = BoundingBoxIC
variable = zi
x1 = 0
x2 = 2
y1 = 1.95
y2 = 2
inside = 0.2
outside = 0
[]
[porosity]
type = RandomIC
variable = porosity
min = 0.25
max = 0.275
seed = 0
[]
[]
[UserObjects]
[dictator]
type = PorousFlowDictator
porous_flow_vars = 'pgas zi'
number_fluid_phases = 2
number_fluid_components = 2
[]
[pc]
type = PorousFlowCapillaryPressureConst
pc = 0
[]
[fs]
type = PorousFlowBrineCO2
brine_fp = brine
co2_fp = co2
capillary_pressure = pc
[]
[]
[Modules]
[FluidProperties]
[co2sw]
type = CO2FluidProperties
[]
[co2]
type = TabulatedFluidProperties
fp = co2sw
[]
[brine]
type = BrineFluidProperties
[]
[]
[]
[Materials]
[temperature]
type = PorousFlowTemperature
temperature = '45'
[]
[brineco2]
type = PorousFlowFluidState
gas_porepressure = 'pgas'
z = 'zi'
temperature_unit = Celsius
xnacl = 'xnacl'
capillary_pressure = pc
fluid_state = fs
[]
[porosity]
type = PorousFlowPorosityConst
porosity = porosity
[]
[permeability]
type = PorousFlowPermeabilityConst
permeability = '1e-11 0 0 0 1e-11 0 0 0 1e-11'
[]
[relperm_water]
type = PorousFlowRelativePermeabilityCorey
phase = 0
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[relperm_gas]
type = PorousFlowRelativePermeabilityCorey
phase = 1
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[diffusivity]
type = PorousFlowDiffusivityConst
diffusion_coeff = '2e-9 2e-9 2e-9 2e-9'
tortuosity = '1 1'
[]
[]
[Preconditioning]
active = basic
[mumps_is_best_for_parallel_jobs]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = ' lu mumps'
[]
[basic]
type = SMP
full = true
petsc_options_iname = '-ksp_type -pc_type -sub_pc_type -sub_pc_factor_shift_type -pc_asm_overlap'
petsc_options_value = 'gmres asm lu NONZERO 2 '
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
end_time = 1e6
nl_max_its = 25
l_max_its = 100
dtmax = 1e4
nl_abs_tol = 1e-6
[TimeStepper]
type = IterationAdaptiveDT
dt = 10
growth_factor = 2
cutback_factor = 0.5
[]
[]
[Functions]
[flux]
type = ParsedFunction
vals = 'delta_xco2 dt'
vars = 'dx dt'
value = 'dx/dt'
[]
[]
[Postprocessors]
[total_co2_in_gas]
type = PorousFlowFluidMass
phase = 1
fluid_component = 1
[]
[total_co2_in_liquid]
type = PorousFlowFluidMass
phase = 0
fluid_component = 1
[]
[numdofs]
type = NumDOFs
[]
[delta_xco2]
type = ChangeOverTimePostprocessor
postprocessor = total_co2_in_liquid
[]
[dt]
type = TimestepSize
[]
[flux]
type = FunctionValuePostprocessor
function = flux
[]
[]
[Outputs]
print_linear_residuals = false
perf_graph = true
exodus = true
csv = true
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)
# Two phase density-driven convection of dissolved CO2 in brine
#
# Initially, the model has a gas phase at the top with a saturation of 0.29
# (which corresponds to an initial value of zi = 0.2).
# Diffusion of the dissolved CO2
# component from the saturated liquid to the unsaturated liquid below reduces the
# amount of CO2 in the gas phase. As the density of the CO2-saturated brine is greater
# than the unsaturated brine, a gravitational instability arises and density-driven
# convection of CO2-rich fingers descend into the unsaturated brine.
#
# The instability is seeded by a random perturbation to the porosity field.
# Mesh adaptivity is used to refine the mesh as the fingers form.
#
# Note: this model is computationally expensive, so should be run with multiple cores,
# preferably on a cluster.
[GlobalParams]
PorousFlowDictator = 'dictator'
gravity = '0 -9.81 0'
[]
[Adaptivity]
max_h_level = 2
marker = marker
initial_marker = initial
initial_steps = 2
[Indicators]
[indicator]
type = GradientJumpIndicator
variable = zi
[]
[]
[Markers]
[marker]
type = ErrorFractionMarker
indicator = indicator
refine = 0.8
[]
[initial]
type = BoxMarker
bottom_left = '0 1.95 0'
top_right = '2 2 0'
inside = REFINE
outside = DO_NOTHING
[]
[]
[]
[Mesh]
type = GeneratedMesh
dim = 2
ymax = 2
xmax = 2
ny = 40
nx = 40
bias_y = 0.95
[]
[Kernels]
[mass0]
type = PorousFlowMassTimeDerivative
fluid_component = 0
variable = pgas
[]
[flux0]
type = PorousFlowAdvectiveFlux
fluid_component = 0
variable = pgas
[]
[diff0]
type = PorousFlowDispersiveFlux
fluid_component = 0
variable = pgas
disp_long = '0 0'
disp_trans = '0 0'
[]
[mass1]
type = PorousFlowMassTimeDerivative
fluid_component = 1
variable = zi
[]
[flux1]
type = PorousFlowAdvectiveFlux
fluid_component = 1
variable = zi
[]
[diff1]
type = PorousFlowDispersiveFlux
fluid_component = 1
variable = zi
disp_long = '0 0'
disp_trans = '0 0'
[]
[]
[AuxVariables]
[xnacl]
initial_condition = 0.01
[]
[saturation_gas]
order = FIRST
family = MONOMIAL
[]
[xco2l]
order = FIRST
family = MONOMIAL
[]
[density_liquid]
order = FIRST
family = MONOMIAL
[]
[porosity]
order = CONSTANT
family = MONOMIAL
[]
[]
[AuxKernels]
[saturation_gas]
type = PorousFlowPropertyAux
variable = saturation_gas
property = saturation
phase = 1
execute_on = 'timestep_end'
[]
[xco2l]
type = PorousFlowPropertyAux
variable = xco2l
property = mass_fraction
phase = 0
fluid_component = 1
execute_on = 'timestep_end'
[]
[density_liquid]
type = PorousFlowPropertyAux
variable = density_liquid
property = density
phase = 0
execute_on = 'timestep_end'
[]
[]
[Variables]
[pgas]
[]
[zi]
scaling = 1e4
[]
[]
[ICs]
[pressure]
type = FunctionIC
function = 10e6-9.81*1000*y
variable = pgas
[]
[zi]
type = BoundingBoxIC
variable = zi
x1 = 0
x2 = 2
y1 = 1.95
y2 = 2
inside = 0.2
outside = 0
[]
[porosity]
type = RandomIC
variable = porosity
min = 0.25
max = 0.275
seed = 0
[]
[]
[UserObjects]
[dictator]
type = PorousFlowDictator
porous_flow_vars = 'pgas zi'
number_fluid_phases = 2
number_fluid_components = 2
[]
[pc]
type = PorousFlowCapillaryPressureConst
pc = 0
[]
[fs]
type = PorousFlowBrineCO2
brine_fp = brine
co2_fp = co2
capillary_pressure = pc
[]
[]
[Modules]
[FluidProperties]
[co2sw]
type = CO2FluidProperties
[]
[co2]
type = TabulatedFluidProperties
fp = co2sw
[]
[brine]
type = BrineFluidProperties
[]
[]
[]
[Materials]
[temperature]
type = PorousFlowTemperature
temperature = '45'
[]
[brineco2]
type = PorousFlowFluidState
gas_porepressure = 'pgas'
z = 'zi'
temperature_unit = Celsius
xnacl = 'xnacl'
capillary_pressure = pc
fluid_state = fs
[]
[porosity]
type = PorousFlowPorosityConst
porosity = porosity
[]
[permeability]
type = PorousFlowPermeabilityConst
permeability = '1e-11 0 0 0 1e-11 0 0 0 1e-11'
[]
[relperm_water]
type = PorousFlowRelativePermeabilityCorey
phase = 0
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[relperm_gas]
type = PorousFlowRelativePermeabilityCorey
phase = 1
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[diffusivity]
type = PorousFlowDiffusivityConst
diffusion_coeff = '2e-9 2e-9 2e-9 2e-9'
tortuosity = '1 1'
[]
[]
[Preconditioning]
active = basic
[mumps_is_best_for_parallel_jobs]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = ' lu mumps'
[]
[basic]
type = SMP
full = true
petsc_options_iname = '-ksp_type -pc_type -sub_pc_type -sub_pc_factor_shift_type -pc_asm_overlap'
petsc_options_value = 'gmres asm lu NONZERO 2 '
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
end_time = 1e6
nl_max_its = 25
l_max_its = 100
dtmax = 1e4
nl_abs_tol = 1e-6
[TimeStepper]
type = IterationAdaptiveDT
dt = 10
growth_factor = 2
cutback_factor = 0.5
[]
[]
[Functions]
[flux]
type = ParsedFunction
vals = 'delta_xco2 dt'
vars = 'dx dt'
value = 'dx/dt'
[]
[]
[Postprocessors]
[total_co2_in_gas]
type = PorousFlowFluidMass
phase = 1
fluid_component = 1
[]
[total_co2_in_liquid]
type = PorousFlowFluidMass
phase = 0
fluid_component = 1
[]
[numdofs]
type = NumDOFs
[]
[delta_xco2]
type = ChangeOverTimePostprocessor
postprocessor = total_co2_in_liquid
[]
[dt]
type = TimestepSize
[]
[flux]
type = FunctionValuePostprocessor
function = flux
[]
[]
[Outputs]
print_linear_residuals = false
perf_graph = true
exodus = true
csv = true
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)
# Two phase density-driven convection of dissolved CO2 in brine
#
# Initially, the model has a gas phase at the top with a saturation of 0.29
# (which corresponds to an initial value of zi = 0.2).
# Diffusion of the dissolved CO2
# component from the saturated liquid to the unsaturated liquid below reduces the
# amount of CO2 in the gas phase. As the density of the CO2-saturated brine is greater
# than the unsaturated brine, a gravitational instability arises and density-driven
# convection of CO2-rich fingers descend into the unsaturated brine.
#
# The instability is seeded by a random perturbation to the porosity field.
# Mesh adaptivity is used to refine the mesh as the fingers form.
#
# Note: this model is computationally expensive, so should be run with multiple cores,
# preferably on a cluster.
[GlobalParams]
PorousFlowDictator = 'dictator'
gravity = '0 -9.81 0'
[]
[Adaptivity]
max_h_level = 2
marker = marker
initial_marker = initial
initial_steps = 2
[Indicators]
[indicator]
type = GradientJumpIndicator
variable = zi
[]
[]
[Markers]
[marker]
type = ErrorFractionMarker
indicator = indicator
refine = 0.8
[]
[initial]
type = BoxMarker
bottom_left = '0 1.95 0'
top_right = '2 2 0'
inside = REFINE
outside = DO_NOTHING
[]
[]
[]
[Mesh]
type = GeneratedMesh
dim = 2
ymax = 2
xmax = 2
ny = 40
nx = 40
bias_y = 0.95
[]
[Kernels]
[mass0]
type = PorousFlowMassTimeDerivative
fluid_component = 0
variable = pgas
[]
[flux0]
type = PorousFlowAdvectiveFlux
fluid_component = 0
variable = pgas
[]
[diff0]
type = PorousFlowDispersiveFlux
fluid_component = 0
variable = pgas
disp_long = '0 0'
disp_trans = '0 0'
[]
[mass1]
type = PorousFlowMassTimeDerivative
fluid_component = 1
variable = zi
[]
[flux1]
type = PorousFlowAdvectiveFlux
fluid_component = 1
variable = zi
[]
[diff1]
type = PorousFlowDispersiveFlux
fluid_component = 1
variable = zi
disp_long = '0 0'
disp_trans = '0 0'
[]
[]
[AuxVariables]
[xnacl]
initial_condition = 0.01
[]
[saturation_gas]
order = FIRST
family = MONOMIAL
[]
[xco2l]
order = FIRST
family = MONOMIAL
[]
[density_liquid]
order = FIRST
family = MONOMIAL
[]
[porosity]
order = CONSTANT
family = MONOMIAL
[]
[]
[AuxKernels]
[saturation_gas]
type = PorousFlowPropertyAux
variable = saturation_gas
property = saturation
phase = 1
execute_on = 'timestep_end'
[]
[xco2l]
type = PorousFlowPropertyAux
variable = xco2l
property = mass_fraction
phase = 0
fluid_component = 1
execute_on = 'timestep_end'
[]
[density_liquid]
type = PorousFlowPropertyAux
variable = density_liquid
property = density
phase = 0
execute_on = 'timestep_end'
[]
[]
[Variables]
[pgas]
[]
[zi]
scaling = 1e4
[]
[]
[ICs]
[pressure]
type = FunctionIC
function = 10e6-9.81*1000*y
variable = pgas
[]
[zi]
type = BoundingBoxIC
variable = zi
x1 = 0
x2 = 2
y1 = 1.95
y2 = 2
inside = 0.2
outside = 0
[]
[porosity]
type = RandomIC
variable = porosity
min = 0.25
max = 0.275
seed = 0
[]
[]
[UserObjects]
[dictator]
type = PorousFlowDictator
porous_flow_vars = 'pgas zi'
number_fluid_phases = 2
number_fluid_components = 2
[]
[pc]
type = PorousFlowCapillaryPressureConst
pc = 0
[]
[fs]
type = PorousFlowBrineCO2
brine_fp = brine
co2_fp = co2
capillary_pressure = pc
[]
[]
[Modules]
[FluidProperties]
[co2sw]
type = CO2FluidProperties
[]
[co2]
type = TabulatedFluidProperties
fp = co2sw
[]
[brine]
type = BrineFluidProperties
[]
[]
[]
[Materials]
[temperature]
type = PorousFlowTemperature
temperature = '45'
[]
[brineco2]
type = PorousFlowFluidState
gas_porepressure = 'pgas'
z = 'zi'
temperature_unit = Celsius
xnacl = 'xnacl'
capillary_pressure = pc
fluid_state = fs
[]
[porosity]
type = PorousFlowPorosityConst
porosity = porosity
[]
[permeability]
type = PorousFlowPermeabilityConst
permeability = '1e-11 0 0 0 1e-11 0 0 0 1e-11'
[]
[relperm_water]
type = PorousFlowRelativePermeabilityCorey
phase = 0
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[relperm_gas]
type = PorousFlowRelativePermeabilityCorey
phase = 1
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[diffusivity]
type = PorousFlowDiffusivityConst
diffusion_coeff = '2e-9 2e-9 2e-9 2e-9'
tortuosity = '1 1'
[]
[]
[Preconditioning]
active = basic
[mumps_is_best_for_parallel_jobs]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = ' lu mumps'
[]
[basic]
type = SMP
full = true
petsc_options_iname = '-ksp_type -pc_type -sub_pc_type -sub_pc_factor_shift_type -pc_asm_overlap'
petsc_options_value = 'gmres asm lu NONZERO 2 '
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
end_time = 1e6
nl_max_its = 25
l_max_its = 100
dtmax = 1e4
nl_abs_tol = 1e-6
[TimeStepper]
type = IterationAdaptiveDT
dt = 10
growth_factor = 2
cutback_factor = 0.5
[]
[]
[Functions]
[flux]
type = ParsedFunction
vals = 'delta_xco2 dt'
vars = 'dx dt'
value = 'dx/dt'
[]
[]
[Postprocessors]
[total_co2_in_gas]
type = PorousFlowFluidMass
phase = 1
fluid_component = 1
[]
[total_co2_in_liquid]
type = PorousFlowFluidMass
phase = 0
fluid_component = 1
[]
[numdofs]
type = NumDOFs
[]
[delta_xco2]
type = ChangeOverTimePostprocessor
postprocessor = total_co2_in_liquid
[]
[dt]
type = TimestepSize
[]
[flux]
type = FunctionValuePostprocessor
function = flux
[]
[]
[Outputs]
print_linear_residuals = false
perf_graph = true
exodus = true
csv = true
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)
# Two phase density-driven convection of dissolved CO2 in brine
#
# Initially, the model has a gas phase at the top with a saturation of 0.29
# (which corresponds to an initial value of zi = 0.2).
# Diffusion of the dissolved CO2
# component from the saturated liquid to the unsaturated liquid below reduces the
# amount of CO2 in the gas phase. As the density of the CO2-saturated brine is greater
# than the unsaturated brine, a gravitational instability arises and density-driven
# convection of CO2-rich fingers descend into the unsaturated brine.
#
# The instability is seeded by a random perturbation to the porosity field.
# Mesh adaptivity is used to refine the mesh as the fingers form.
#
# Note: this model is computationally expensive, so should be run with multiple cores,
# preferably on a cluster.
[GlobalParams]
PorousFlowDictator = 'dictator'
gravity = '0 -9.81 0'
[]
[Adaptivity]
max_h_level = 2
marker = marker
initial_marker = initial
initial_steps = 2
[Indicators]
[indicator]
type = GradientJumpIndicator
variable = zi
[]
[]
[Markers]
[marker]
type = ErrorFractionMarker
indicator = indicator
refine = 0.8
[]
[initial]
type = BoxMarker
bottom_left = '0 1.95 0'
top_right = '2 2 0'
inside = REFINE
outside = DO_NOTHING
[]
[]
[]
[Mesh]
type = GeneratedMesh
dim = 2
ymax = 2
xmax = 2
ny = 40
nx = 40
bias_y = 0.95
[]
[Kernels]
[mass0]
type = PorousFlowMassTimeDerivative
fluid_component = 0
variable = pgas
[]
[flux0]
type = PorousFlowAdvectiveFlux
fluid_component = 0
variable = pgas
[]
[diff0]
type = PorousFlowDispersiveFlux
fluid_component = 0
variable = pgas
disp_long = '0 0'
disp_trans = '0 0'
[]
[mass1]
type = PorousFlowMassTimeDerivative
fluid_component = 1
variable = zi
[]
[flux1]
type = PorousFlowAdvectiveFlux
fluid_component = 1
variable = zi
[]
[diff1]
type = PorousFlowDispersiveFlux
fluid_component = 1
variable = zi
disp_long = '0 0'
disp_trans = '0 0'
[]
[]
[AuxVariables]
[xnacl]
initial_condition = 0.01
[]
[saturation_gas]
order = FIRST
family = MONOMIAL
[]
[xco2l]
order = FIRST
family = MONOMIAL
[]
[density_liquid]
order = FIRST
family = MONOMIAL
[]
[porosity]
order = CONSTANT
family = MONOMIAL
[]
[]
[AuxKernels]
[saturation_gas]
type = PorousFlowPropertyAux
variable = saturation_gas
property = saturation
phase = 1
execute_on = 'timestep_end'
[]
[xco2l]
type = PorousFlowPropertyAux
variable = xco2l
property = mass_fraction
phase = 0
fluid_component = 1
execute_on = 'timestep_end'
[]
[density_liquid]
type = PorousFlowPropertyAux
variable = density_liquid
property = density
phase = 0
execute_on = 'timestep_end'
[]
[]
[Variables]
[pgas]
[]
[zi]
scaling = 1e4
[]
[]
[ICs]
[pressure]
type = FunctionIC
function = 10e6-9.81*1000*y
variable = pgas
[]
[zi]
type = BoundingBoxIC
variable = zi
x1 = 0
x2 = 2
y1 = 1.95
y2 = 2
inside = 0.2
outside = 0
[]
[porosity]
type = RandomIC
variable = porosity
min = 0.25
max = 0.275
seed = 0
[]
[]
[UserObjects]
[dictator]
type = PorousFlowDictator
porous_flow_vars = 'pgas zi'
number_fluid_phases = 2
number_fluid_components = 2
[]
[pc]
type = PorousFlowCapillaryPressureConst
pc = 0
[]
[fs]
type = PorousFlowBrineCO2
brine_fp = brine
co2_fp = co2
capillary_pressure = pc
[]
[]
[Modules]
[FluidProperties]
[co2sw]
type = CO2FluidProperties
[]
[co2]
type = TabulatedFluidProperties
fp = co2sw
[]
[brine]
type = BrineFluidProperties
[]
[]
[]
[Materials]
[temperature]
type = PorousFlowTemperature
temperature = '45'
[]
[brineco2]
type = PorousFlowFluidState
gas_porepressure = 'pgas'
z = 'zi'
temperature_unit = Celsius
xnacl = 'xnacl'
capillary_pressure = pc
fluid_state = fs
[]
[porosity]
type = PorousFlowPorosityConst
porosity = porosity
[]
[permeability]
type = PorousFlowPermeabilityConst
permeability = '1e-11 0 0 0 1e-11 0 0 0 1e-11'
[]
[relperm_water]
type = PorousFlowRelativePermeabilityCorey
phase = 0
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[relperm_gas]
type = PorousFlowRelativePermeabilityCorey
phase = 1
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[diffusivity]
type = PorousFlowDiffusivityConst
diffusion_coeff = '2e-9 2e-9 2e-9 2e-9'
tortuosity = '1 1'
[]
[]
[Preconditioning]
active = basic
[mumps_is_best_for_parallel_jobs]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = ' lu mumps'
[]
[basic]
type = SMP
full = true
petsc_options_iname = '-ksp_type -pc_type -sub_pc_type -sub_pc_factor_shift_type -pc_asm_overlap'
petsc_options_value = 'gmres asm lu NONZERO 2 '
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
end_time = 1e6
nl_max_its = 25
l_max_its = 100
dtmax = 1e4
nl_abs_tol = 1e-6
[TimeStepper]
type = IterationAdaptiveDT
dt = 10
growth_factor = 2
cutback_factor = 0.5
[]
[]
[Functions]
[flux]
type = ParsedFunction
vals = 'delta_xco2 dt'
vars = 'dx dt'
value = 'dx/dt'
[]
[]
[Postprocessors]
[total_co2_in_gas]
type = PorousFlowFluidMass
phase = 1
fluid_component = 1
[]
[total_co2_in_liquid]
type = PorousFlowFluidMass
phase = 0
fluid_component = 1
[]
[numdofs]
type = NumDOFs
[]
[delta_xco2]
type = ChangeOverTimePostprocessor
postprocessor = total_co2_in_liquid
[]
[dt]
type = TimestepSize
[]
[flux]
type = FunctionValuePostprocessor
function = flux
[]
[]
[Outputs]
print_linear_residuals = false
perf_graph = true
exodus = true
csv = true
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)
# Two phase density-driven convection of dissolved CO2 in brine
#
# Initially, the model has a gas phase at the top with a saturation of 0.29
# (which corresponds to an initial value of zi = 0.2).
# Diffusion of the dissolved CO2
# component from the saturated liquid to the unsaturated liquid below reduces the
# amount of CO2 in the gas phase. As the density of the CO2-saturated brine is greater
# than the unsaturated brine, a gravitational instability arises and density-driven
# convection of CO2-rich fingers descend into the unsaturated brine.
#
# The instability is seeded by a random perturbation to the porosity field.
# Mesh adaptivity is used to refine the mesh as the fingers form.
#
# Note: this model is computationally expensive, so should be run with multiple cores,
# preferably on a cluster.
[GlobalParams]
PorousFlowDictator = 'dictator'
gravity = '0 -9.81 0'
[]
[Adaptivity]
max_h_level = 2
marker = marker
initial_marker = initial
initial_steps = 2
[Indicators]
[indicator]
type = GradientJumpIndicator
variable = zi
[]
[]
[Markers]
[marker]
type = ErrorFractionMarker
indicator = indicator
refine = 0.8
[]
[initial]
type = BoxMarker
bottom_left = '0 1.95 0'
top_right = '2 2 0'
inside = REFINE
outside = DO_NOTHING
[]
[]
[]
[Mesh]
type = GeneratedMesh
dim = 2
ymax = 2
xmax = 2
ny = 40
nx = 40
bias_y = 0.95
[]
[Kernels]
[mass0]
type = PorousFlowMassTimeDerivative
fluid_component = 0
variable = pgas
[]
[flux0]
type = PorousFlowAdvectiveFlux
fluid_component = 0
variable = pgas
[]
[diff0]
type = PorousFlowDispersiveFlux
fluid_component = 0
variable = pgas
disp_long = '0 0'
disp_trans = '0 0'
[]
[mass1]
type = PorousFlowMassTimeDerivative
fluid_component = 1
variable = zi
[]
[flux1]
type = PorousFlowAdvectiveFlux
fluid_component = 1
variable = zi
[]
[diff1]
type = PorousFlowDispersiveFlux
fluid_component = 1
variable = zi
disp_long = '0 0'
disp_trans = '0 0'
[]
[]
[AuxVariables]
[xnacl]
initial_condition = 0.01
[]
[saturation_gas]
order = FIRST
family = MONOMIAL
[]
[xco2l]
order = FIRST
family = MONOMIAL
[]
[density_liquid]
order = FIRST
family = MONOMIAL
[]
[porosity]
order = CONSTANT
family = MONOMIAL
[]
[]
[AuxKernels]
[saturation_gas]
type = PorousFlowPropertyAux
variable = saturation_gas
property = saturation
phase = 1
execute_on = 'timestep_end'
[]
[xco2l]
type = PorousFlowPropertyAux
variable = xco2l
property = mass_fraction
phase = 0
fluid_component = 1
execute_on = 'timestep_end'
[]
[density_liquid]
type = PorousFlowPropertyAux
variable = density_liquid
property = density
phase = 0
execute_on = 'timestep_end'
[]
[]
[Variables]
[pgas]
[]
[zi]
scaling = 1e4
[]
[]
[ICs]
[pressure]
type = FunctionIC
function = 10e6-9.81*1000*y
variable = pgas
[]
[zi]
type = BoundingBoxIC
variable = zi
x1 = 0
x2 = 2
y1 = 1.95
y2 = 2
inside = 0.2
outside = 0
[]
[porosity]
type = RandomIC
variable = porosity
min = 0.25
max = 0.275
seed = 0
[]
[]
[UserObjects]
[dictator]
type = PorousFlowDictator
porous_flow_vars = 'pgas zi'
number_fluid_phases = 2
number_fluid_components = 2
[]
[pc]
type = PorousFlowCapillaryPressureConst
pc = 0
[]
[fs]
type = PorousFlowBrineCO2
brine_fp = brine
co2_fp = co2
capillary_pressure = pc
[]
[]
[Modules]
[FluidProperties]
[co2sw]
type = CO2FluidProperties
[]
[co2]
type = TabulatedFluidProperties
fp = co2sw
[]
[brine]
type = BrineFluidProperties
[]
[]
[]
[Materials]
[temperature]
type = PorousFlowTemperature
temperature = '45'
[]
[brineco2]
type = PorousFlowFluidState
gas_porepressure = 'pgas'
z = 'zi'
temperature_unit = Celsius
xnacl = 'xnacl'
capillary_pressure = pc
fluid_state = fs
[]
[porosity]
type = PorousFlowPorosityConst
porosity = porosity
[]
[permeability]
type = PorousFlowPermeabilityConst
permeability = '1e-11 0 0 0 1e-11 0 0 0 1e-11'
[]
[relperm_water]
type = PorousFlowRelativePermeabilityCorey
phase = 0
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[relperm_gas]
type = PorousFlowRelativePermeabilityCorey
phase = 1
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[diffusivity]
type = PorousFlowDiffusivityConst
diffusion_coeff = '2e-9 2e-9 2e-9 2e-9'
tortuosity = '1 1'
[]
[]
[Preconditioning]
active = basic
[mumps_is_best_for_parallel_jobs]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = ' lu mumps'
[]
[basic]
type = SMP
full = true
petsc_options_iname = '-ksp_type -pc_type -sub_pc_type -sub_pc_factor_shift_type -pc_asm_overlap'
petsc_options_value = 'gmres asm lu NONZERO 2 '
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
end_time = 1e6
nl_max_its = 25
l_max_its = 100
dtmax = 1e4
nl_abs_tol = 1e-6
[TimeStepper]
type = IterationAdaptiveDT
dt = 10
growth_factor = 2
cutback_factor = 0.5
[]
[]
[Functions]
[flux]
type = ParsedFunction
vals = 'delta_xco2 dt'
vars = 'dx dt'
value = 'dx/dt'
[]
[]
[Postprocessors]
[total_co2_in_gas]
type = PorousFlowFluidMass
phase = 1
fluid_component = 1
[]
[total_co2_in_liquid]
type = PorousFlowFluidMass
phase = 0
fluid_component = 1
[]
[numdofs]
type = NumDOFs
[]
[delta_xco2]
type = ChangeOverTimePostprocessor
postprocessor = total_co2_in_liquid
[]
[dt]
type = TimestepSize
[]
[flux]
type = FunctionValuePostprocessor
function = flux
[]
[]
[Outputs]
print_linear_residuals = false
perf_graph = true
exodus = true
csv = true
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)
# Two phase density-driven convection of dissolved CO2 in brine
#
# Initially, the model has a gas phase at the top with a saturation of 0.29
# (which corresponds to an initial value of zi = 0.2).
# Diffusion of the dissolved CO2
# component from the saturated liquid to the unsaturated liquid below reduces the
# amount of CO2 in the gas phase. As the density of the CO2-saturated brine is greater
# than the unsaturated brine, a gravitational instability arises and density-driven
# convection of CO2-rich fingers descend into the unsaturated brine.
#
# The instability is seeded by a random perturbation to the porosity field.
# Mesh adaptivity is used to refine the mesh as the fingers form.
#
# Note: this model is computationally expensive, so should be run with multiple cores,
# preferably on a cluster.
[GlobalParams]
PorousFlowDictator = 'dictator'
gravity = '0 -9.81 0'
[]
[Adaptivity]
max_h_level = 2
marker = marker
initial_marker = initial
initial_steps = 2
[Indicators]
[indicator]
type = GradientJumpIndicator
variable = zi
[]
[]
[Markers]
[marker]
type = ErrorFractionMarker
indicator = indicator
refine = 0.8
[]
[initial]
type = BoxMarker
bottom_left = '0 1.95 0'
top_right = '2 2 0'
inside = REFINE
outside = DO_NOTHING
[]
[]
[]
[Mesh]
type = GeneratedMesh
dim = 2
ymax = 2
xmax = 2
ny = 40
nx = 40
bias_y = 0.95
[]
[Kernels]
[mass0]
type = PorousFlowMassTimeDerivative
fluid_component = 0
variable = pgas
[]
[flux0]
type = PorousFlowAdvectiveFlux
fluid_component = 0
variable = pgas
[]
[diff0]
type = PorousFlowDispersiveFlux
fluid_component = 0
variable = pgas
disp_long = '0 0'
disp_trans = '0 0'
[]
[mass1]
type = PorousFlowMassTimeDerivative
fluid_component = 1
variable = zi
[]
[flux1]
type = PorousFlowAdvectiveFlux
fluid_component = 1
variable = zi
[]
[diff1]
type = PorousFlowDispersiveFlux
fluid_component = 1
variable = zi
disp_long = '0 0'
disp_trans = '0 0'
[]
[]
[AuxVariables]
[xnacl]
initial_condition = 0.01
[]
[saturation_gas]
order = FIRST
family = MONOMIAL
[]
[xco2l]
order = FIRST
family = MONOMIAL
[]
[density_liquid]
order = FIRST
family = MONOMIAL
[]
[porosity]
order = CONSTANT
family = MONOMIAL
[]
[]
[AuxKernels]
[saturation_gas]
type = PorousFlowPropertyAux
variable = saturation_gas
property = saturation
phase = 1
execute_on = 'timestep_end'
[]
[xco2l]
type = PorousFlowPropertyAux
variable = xco2l
property = mass_fraction
phase = 0
fluid_component = 1
execute_on = 'timestep_end'
[]
[density_liquid]
type = PorousFlowPropertyAux
variable = density_liquid
property = density
phase = 0
execute_on = 'timestep_end'
[]
[]
[Variables]
[pgas]
[]
[zi]
scaling = 1e4
[]
[]
[ICs]
[pressure]
type = FunctionIC
function = 10e6-9.81*1000*y
variable = pgas
[]
[zi]
type = BoundingBoxIC
variable = zi
x1 = 0
x2 = 2
y1 = 1.95
y2 = 2
inside = 0.2
outside = 0
[]
[porosity]
type = RandomIC
variable = porosity
min = 0.25
max = 0.275
seed = 0
[]
[]
[UserObjects]
[dictator]
type = PorousFlowDictator
porous_flow_vars = 'pgas zi'
number_fluid_phases = 2
number_fluid_components = 2
[]
[pc]
type = PorousFlowCapillaryPressureConst
pc = 0
[]
[fs]
type = PorousFlowBrineCO2
brine_fp = brine
co2_fp = co2
capillary_pressure = pc
[]
[]
[Modules]
[FluidProperties]
[co2sw]
type = CO2FluidProperties
[]
[co2]
type = TabulatedFluidProperties
fp = co2sw
[]
[brine]
type = BrineFluidProperties
[]
[]
[]
[Materials]
[temperature]
type = PorousFlowTemperature
temperature = '45'
[]
[brineco2]
type = PorousFlowFluidState
gas_porepressure = 'pgas'
z = 'zi'
temperature_unit = Celsius
xnacl = 'xnacl'
capillary_pressure = pc
fluid_state = fs
[]
[porosity]
type = PorousFlowPorosityConst
porosity = porosity
[]
[permeability]
type = PorousFlowPermeabilityConst
permeability = '1e-11 0 0 0 1e-11 0 0 0 1e-11'
[]
[relperm_water]
type = PorousFlowRelativePermeabilityCorey
phase = 0
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[relperm_gas]
type = PorousFlowRelativePermeabilityCorey
phase = 1
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[diffusivity]
type = PorousFlowDiffusivityConst
diffusion_coeff = '2e-9 2e-9 2e-9 2e-9'
tortuosity = '1 1'
[]
[]
[Preconditioning]
active = basic
[mumps_is_best_for_parallel_jobs]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = ' lu mumps'
[]
[basic]
type = SMP
full = true
petsc_options_iname = '-ksp_type -pc_type -sub_pc_type -sub_pc_factor_shift_type -pc_asm_overlap'
petsc_options_value = 'gmres asm lu NONZERO 2 '
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
end_time = 1e6
nl_max_its = 25
l_max_its = 100
dtmax = 1e4
nl_abs_tol = 1e-6
[TimeStepper]
type = IterationAdaptiveDT
dt = 10
growth_factor = 2
cutback_factor = 0.5
[]
[]
[Functions]
[flux]
type = ParsedFunction
vals = 'delta_xco2 dt'
vars = 'dx dt'
value = 'dx/dt'
[]
[]
[Postprocessors]
[total_co2_in_gas]
type = PorousFlowFluidMass
phase = 1
fluid_component = 1
[]
[total_co2_in_liquid]
type = PorousFlowFluidMass
phase = 0
fluid_component = 1
[]
[numdofs]
type = NumDOFs
[]
[delta_xco2]
type = ChangeOverTimePostprocessor
postprocessor = total_co2_in_liquid
[]
[dt]
type = TimestepSize
[]
[flux]
type = FunctionValuePostprocessor
function = flux
[]
[]
[Outputs]
print_linear_residuals = false
perf_graph = true
exodus = true
csv = true
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)
# Two phase density-driven convection of dissolved CO2 in brine
#
# Initially, the model has a gas phase at the top with a saturation of 0.29
# (which corresponds to an initial value of zi = 0.2).
# Diffusion of the dissolved CO2
# component from the saturated liquid to the unsaturated liquid below reduces the
# amount of CO2 in the gas phase. As the density of the CO2-saturated brine is greater
# than the unsaturated brine, a gravitational instability arises and density-driven
# convection of CO2-rich fingers descend into the unsaturated brine.
#
# The instability is seeded by a random perturbation to the porosity field.
# Mesh adaptivity is used to refine the mesh as the fingers form.
#
# Note: this model is computationally expensive, so should be run with multiple cores,
# preferably on a cluster.
[GlobalParams]
PorousFlowDictator = 'dictator'
gravity = '0 -9.81 0'
[]
[Adaptivity]
max_h_level = 2
marker = marker
initial_marker = initial
initial_steps = 2
[Indicators]
[indicator]
type = GradientJumpIndicator
variable = zi
[]
[]
[Markers]
[marker]
type = ErrorFractionMarker
indicator = indicator
refine = 0.8
[]
[initial]
type = BoxMarker
bottom_left = '0 1.95 0'
top_right = '2 2 0'
inside = REFINE
outside = DO_NOTHING
[]
[]
[]
[Mesh]
type = GeneratedMesh
dim = 2
ymax = 2
xmax = 2
ny = 40
nx = 40
bias_y = 0.95
[]
[Kernels]
[mass0]
type = PorousFlowMassTimeDerivative
fluid_component = 0
variable = pgas
[]
[flux0]
type = PorousFlowAdvectiveFlux
fluid_component = 0
variable = pgas
[]
[diff0]
type = PorousFlowDispersiveFlux
fluid_component = 0
variable = pgas
disp_long = '0 0'
disp_trans = '0 0'
[]
[mass1]
type = PorousFlowMassTimeDerivative
fluid_component = 1
variable = zi
[]
[flux1]
type = PorousFlowAdvectiveFlux
fluid_component = 1
variable = zi
[]
[diff1]
type = PorousFlowDispersiveFlux
fluid_component = 1
variable = zi
disp_long = '0 0'
disp_trans = '0 0'
[]
[]
[AuxVariables]
[xnacl]
initial_condition = 0.01
[]
[saturation_gas]
order = FIRST
family = MONOMIAL
[]
[xco2l]
order = FIRST
family = MONOMIAL
[]
[density_liquid]
order = FIRST
family = MONOMIAL
[]
[porosity]
order = CONSTANT
family = MONOMIAL
[]
[]
[AuxKernels]
[saturation_gas]
type = PorousFlowPropertyAux
variable = saturation_gas
property = saturation
phase = 1
execute_on = 'timestep_end'
[]
[xco2l]
type = PorousFlowPropertyAux
variable = xco2l
property = mass_fraction
phase = 0
fluid_component = 1
execute_on = 'timestep_end'
[]
[density_liquid]
type = PorousFlowPropertyAux
variable = density_liquid
property = density
phase = 0
execute_on = 'timestep_end'
[]
[]
[Variables]
[pgas]
[]
[zi]
scaling = 1e4
[]
[]
[ICs]
[pressure]
type = FunctionIC
function = 10e6-9.81*1000*y
variable = pgas
[]
[zi]
type = BoundingBoxIC
variable = zi
x1 = 0
x2 = 2
y1 = 1.95
y2 = 2
inside = 0.2
outside = 0
[]
[porosity]
type = RandomIC
variable = porosity
min = 0.25
max = 0.275
seed = 0
[]
[]
[UserObjects]
[dictator]
type = PorousFlowDictator
porous_flow_vars = 'pgas zi'
number_fluid_phases = 2
number_fluid_components = 2
[]
[pc]
type = PorousFlowCapillaryPressureConst
pc = 0
[]
[fs]
type = PorousFlowBrineCO2
brine_fp = brine
co2_fp = co2
capillary_pressure = pc
[]
[]
[Modules]
[FluidProperties]
[co2sw]
type = CO2FluidProperties
[]
[co2]
type = TabulatedFluidProperties
fp = co2sw
[]
[brine]
type = BrineFluidProperties
[]
[]
[]
[Materials]
[temperature]
type = PorousFlowTemperature
temperature = '45'
[]
[brineco2]
type = PorousFlowFluidState
gas_porepressure = 'pgas'
z = 'zi'
temperature_unit = Celsius
xnacl = 'xnacl'
capillary_pressure = pc
fluid_state = fs
[]
[porosity]
type = PorousFlowPorosityConst
porosity = porosity
[]
[permeability]
type = PorousFlowPermeabilityConst
permeability = '1e-11 0 0 0 1e-11 0 0 0 1e-11'
[]
[relperm_water]
type = PorousFlowRelativePermeabilityCorey
phase = 0
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[relperm_gas]
type = PorousFlowRelativePermeabilityCorey
phase = 1
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[diffusivity]
type = PorousFlowDiffusivityConst
diffusion_coeff = '2e-9 2e-9 2e-9 2e-9'
tortuosity = '1 1'
[]
[]
[Preconditioning]
active = basic
[mumps_is_best_for_parallel_jobs]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = ' lu mumps'
[]
[basic]
type = SMP
full = true
petsc_options_iname = '-ksp_type -pc_type -sub_pc_type -sub_pc_factor_shift_type -pc_asm_overlap'
petsc_options_value = 'gmres asm lu NONZERO 2 '
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
end_time = 1e6
nl_max_its = 25
l_max_its = 100
dtmax = 1e4
nl_abs_tol = 1e-6
[TimeStepper]
type = IterationAdaptiveDT
dt = 10
growth_factor = 2
cutback_factor = 0.5
[]
[]
[Functions]
[flux]
type = ParsedFunction
vals = 'delta_xco2 dt'
vars = 'dx dt'
value = 'dx/dt'
[]
[]
[Postprocessors]
[total_co2_in_gas]
type = PorousFlowFluidMass
phase = 1
fluid_component = 1
[]
[total_co2_in_liquid]
type = PorousFlowFluidMass
phase = 0
fluid_component = 1
[]
[numdofs]
type = NumDOFs
[]
[delta_xco2]
type = ChangeOverTimePostprocessor
postprocessor = total_co2_in_liquid
[]
[dt]
type = TimestepSize
[]
[flux]
type = FunctionValuePostprocessor
function = flux
[]
[]
[Outputs]
print_linear_residuals = false
perf_graph = true
exodus = true
csv = true
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)
# Two phase density-driven convection of dissolved CO2 in brine
#
# Initially, the model has a gas phase at the top with a saturation of 0.29
# (which corresponds to an initial value of zi = 0.2).
# Diffusion of the dissolved CO2
# component from the saturated liquid to the unsaturated liquid below reduces the
# amount of CO2 in the gas phase. As the density of the CO2-saturated brine is greater
# than the unsaturated brine, a gravitational instability arises and density-driven
# convection of CO2-rich fingers descend into the unsaturated brine.
#
# The instability is seeded by a random perturbation to the porosity field.
# Mesh adaptivity is used to refine the mesh as the fingers form.
#
# Note: this model is computationally expensive, so should be run with multiple cores,
# preferably on a cluster.
[GlobalParams]
PorousFlowDictator = 'dictator'
gravity = '0 -9.81 0'
[]
[Adaptivity]
max_h_level = 2
marker = marker
initial_marker = initial
initial_steps = 2
[Indicators]
[indicator]
type = GradientJumpIndicator
variable = zi
[]
[]
[Markers]
[marker]
type = ErrorFractionMarker
indicator = indicator
refine = 0.8
[]
[initial]
type = BoxMarker
bottom_left = '0 1.95 0'
top_right = '2 2 0'
inside = REFINE
outside = DO_NOTHING
[]
[]
[]
[Mesh]
type = GeneratedMesh
dim = 2
ymax = 2
xmax = 2
ny = 40
nx = 40
bias_y = 0.95
[]
[Kernels]
[mass0]
type = PorousFlowMassTimeDerivative
fluid_component = 0
variable = pgas
[]
[flux0]
type = PorousFlowAdvectiveFlux
fluid_component = 0
variable = pgas
[]
[diff0]
type = PorousFlowDispersiveFlux
fluid_component = 0
variable = pgas
disp_long = '0 0'
disp_trans = '0 0'
[]
[mass1]
type = PorousFlowMassTimeDerivative
fluid_component = 1
variable = zi
[]
[flux1]
type = PorousFlowAdvectiveFlux
fluid_component = 1
variable = zi
[]
[diff1]
type = PorousFlowDispersiveFlux
fluid_component = 1
variable = zi
disp_long = '0 0'
disp_trans = '0 0'
[]
[]
[AuxVariables]
[xnacl]
initial_condition = 0.01
[]
[saturation_gas]
order = FIRST
family = MONOMIAL
[]
[xco2l]
order = FIRST
family = MONOMIAL
[]
[density_liquid]
order = FIRST
family = MONOMIAL
[]
[porosity]
order = CONSTANT
family = MONOMIAL
[]
[]
[AuxKernels]
[saturation_gas]
type = PorousFlowPropertyAux
variable = saturation_gas
property = saturation
phase = 1
execute_on = 'timestep_end'
[]
[xco2l]
type = PorousFlowPropertyAux
variable = xco2l
property = mass_fraction
phase = 0
fluid_component = 1
execute_on = 'timestep_end'
[]
[density_liquid]
type = PorousFlowPropertyAux
variable = density_liquid
property = density
phase = 0
execute_on = 'timestep_end'
[]
[]
[Variables]
[pgas]
[]
[zi]
scaling = 1e4
[]
[]
[ICs]
[pressure]
type = FunctionIC
function = 10e6-9.81*1000*y
variable = pgas
[]
[zi]
type = BoundingBoxIC
variable = zi
x1 = 0
x2 = 2
y1 = 1.95
y2 = 2
inside = 0.2
outside = 0
[]
[porosity]
type = RandomIC
variable = porosity
min = 0.25
max = 0.275
seed = 0
[]
[]
[UserObjects]
[dictator]
type = PorousFlowDictator
porous_flow_vars = 'pgas zi'
number_fluid_phases = 2
number_fluid_components = 2
[]
[pc]
type = PorousFlowCapillaryPressureConst
pc = 0
[]
[fs]
type = PorousFlowBrineCO2
brine_fp = brine
co2_fp = co2
capillary_pressure = pc
[]
[]
[Modules]
[FluidProperties]
[co2sw]
type = CO2FluidProperties
[]
[co2]
type = TabulatedFluidProperties
fp = co2sw
[]
[brine]
type = BrineFluidProperties
[]
[]
[]
[Materials]
[temperature]
type = PorousFlowTemperature
temperature = '45'
[]
[brineco2]
type = PorousFlowFluidState
gas_porepressure = 'pgas'
z = 'zi'
temperature_unit = Celsius
xnacl = 'xnacl'
capillary_pressure = pc
fluid_state = fs
[]
[porosity]
type = PorousFlowPorosityConst
porosity = porosity
[]
[permeability]
type = PorousFlowPermeabilityConst
permeability = '1e-11 0 0 0 1e-11 0 0 0 1e-11'
[]
[relperm_water]
type = PorousFlowRelativePermeabilityCorey
phase = 0
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[relperm_gas]
type = PorousFlowRelativePermeabilityCorey
phase = 1
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[diffusivity]
type = PorousFlowDiffusivityConst
diffusion_coeff = '2e-9 2e-9 2e-9 2e-9'
tortuosity = '1 1'
[]
[]
[Preconditioning]
active = basic
[mumps_is_best_for_parallel_jobs]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = ' lu mumps'
[]
[basic]
type = SMP
full = true
petsc_options_iname = '-ksp_type -pc_type -sub_pc_type -sub_pc_factor_shift_type -pc_asm_overlap'
petsc_options_value = 'gmres asm lu NONZERO 2 '
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
end_time = 1e6
nl_max_its = 25
l_max_its = 100
dtmax = 1e4
nl_abs_tol = 1e-6
[TimeStepper]
type = IterationAdaptiveDT
dt = 10
growth_factor = 2
cutback_factor = 0.5
[]
[]
[Functions]
[flux]
type = ParsedFunction
vals = 'delta_xco2 dt'
vars = 'dx dt'
value = 'dx/dt'
[]
[]
[Postprocessors]
[total_co2_in_gas]
type = PorousFlowFluidMass
phase = 1
fluid_component = 1
[]
[total_co2_in_liquid]
type = PorousFlowFluidMass
phase = 0
fluid_component = 1
[]
[numdofs]
type = NumDOFs
[]
[delta_xco2]
type = ChangeOverTimePostprocessor
postprocessor = total_co2_in_liquid
[]
[dt]
type = TimestepSize
[]
[flux]
type = FunctionValuePostprocessor
function = flux
[]
[]
[Outputs]
print_linear_residuals = false
perf_graph = true
exodus = true
csv = true
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)
# Two phase density-driven convection of dissolved CO2 in brine
#
# Initially, the model has a gas phase at the top with a saturation of 0.29
# (which corresponds to an initial value of zi = 0.2).
# Diffusion of the dissolved CO2
# component from the saturated liquid to the unsaturated liquid below reduces the
# amount of CO2 in the gas phase. As the density of the CO2-saturated brine is greater
# than the unsaturated brine, a gravitational instability arises and density-driven
# convection of CO2-rich fingers descend into the unsaturated brine.
#
# The instability is seeded by a random perturbation to the porosity field.
# Mesh adaptivity is used to refine the mesh as the fingers form.
#
# Note: this model is computationally expensive, so should be run with multiple cores,
# preferably on a cluster.
[GlobalParams]
PorousFlowDictator = 'dictator'
gravity = '0 -9.81 0'
[]
[Adaptivity]
max_h_level = 2
marker = marker
initial_marker = initial
initial_steps = 2
[Indicators]
[indicator]
type = GradientJumpIndicator
variable = zi
[]
[]
[Markers]
[marker]
type = ErrorFractionMarker
indicator = indicator
refine = 0.8
[]
[initial]
type = BoxMarker
bottom_left = '0 1.95 0'
top_right = '2 2 0'
inside = REFINE
outside = DO_NOTHING
[]
[]
[]
[Mesh]
type = GeneratedMesh
dim = 2
ymax = 2
xmax = 2
ny = 40
nx = 40
bias_y = 0.95
[]
[Kernels]
[mass0]
type = PorousFlowMassTimeDerivative
fluid_component = 0
variable = pgas
[]
[flux0]
type = PorousFlowAdvectiveFlux
fluid_component = 0
variable = pgas
[]
[diff0]
type = PorousFlowDispersiveFlux
fluid_component = 0
variable = pgas
disp_long = '0 0'
disp_trans = '0 0'
[]
[mass1]
type = PorousFlowMassTimeDerivative
fluid_component = 1
variable = zi
[]
[flux1]
type = PorousFlowAdvectiveFlux
fluid_component = 1
variable = zi
[]
[diff1]
type = PorousFlowDispersiveFlux
fluid_component = 1
variable = zi
disp_long = '0 0'
disp_trans = '0 0'
[]
[]
[AuxVariables]
[xnacl]
initial_condition = 0.01
[]
[saturation_gas]
order = FIRST
family = MONOMIAL
[]
[xco2l]
order = FIRST
family = MONOMIAL
[]
[density_liquid]
order = FIRST
family = MONOMIAL
[]
[porosity]
order = CONSTANT
family = MONOMIAL
[]
[]
[AuxKernels]
[saturation_gas]
type = PorousFlowPropertyAux
variable = saturation_gas
property = saturation
phase = 1
execute_on = 'timestep_end'
[]
[xco2l]
type = PorousFlowPropertyAux
variable = xco2l
property = mass_fraction
phase = 0
fluid_component = 1
execute_on = 'timestep_end'
[]
[density_liquid]
type = PorousFlowPropertyAux
variable = density_liquid
property = density
phase = 0
execute_on = 'timestep_end'
[]
[]
[Variables]
[pgas]
[]
[zi]
scaling = 1e4
[]
[]
[ICs]
[pressure]
type = FunctionIC
function = 10e6-9.81*1000*y
variable = pgas
[]
[zi]
type = BoundingBoxIC
variable = zi
x1 = 0
x2 = 2
y1 = 1.95
y2 = 2
inside = 0.2
outside = 0
[]
[porosity]
type = RandomIC
variable = porosity
min = 0.25
max = 0.275
seed = 0
[]
[]
[UserObjects]
[dictator]
type = PorousFlowDictator
porous_flow_vars = 'pgas zi'
number_fluid_phases = 2
number_fluid_components = 2
[]
[pc]
type = PorousFlowCapillaryPressureConst
pc = 0
[]
[fs]
type = PorousFlowBrineCO2
brine_fp = brine
co2_fp = co2
capillary_pressure = pc
[]
[]
[Modules]
[FluidProperties]
[co2sw]
type = CO2FluidProperties
[]
[co2]
type = TabulatedFluidProperties
fp = co2sw
[]
[brine]
type = BrineFluidProperties
[]
[]
[]
[Materials]
[temperature]
type = PorousFlowTemperature
temperature = '45'
[]
[brineco2]
type = PorousFlowFluidState
gas_porepressure = 'pgas'
z = 'zi'
temperature_unit = Celsius
xnacl = 'xnacl'
capillary_pressure = pc
fluid_state = fs
[]
[porosity]
type = PorousFlowPorosityConst
porosity = porosity
[]
[permeability]
type = PorousFlowPermeabilityConst
permeability = '1e-11 0 0 0 1e-11 0 0 0 1e-11'
[]
[relperm_water]
type = PorousFlowRelativePermeabilityCorey
phase = 0
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[relperm_gas]
type = PorousFlowRelativePermeabilityCorey
phase = 1
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[diffusivity]
type = PorousFlowDiffusivityConst
diffusion_coeff = '2e-9 2e-9 2e-9 2e-9'
tortuosity = '1 1'
[]
[]
[Preconditioning]
active = basic
[mumps_is_best_for_parallel_jobs]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = ' lu mumps'
[]
[basic]
type = SMP
full = true
petsc_options_iname = '-ksp_type -pc_type -sub_pc_type -sub_pc_factor_shift_type -pc_asm_overlap'
petsc_options_value = 'gmres asm lu NONZERO 2 '
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
end_time = 1e6
nl_max_its = 25
l_max_its = 100
dtmax = 1e4
nl_abs_tol = 1e-6
[TimeStepper]
type = IterationAdaptiveDT
dt = 10
growth_factor = 2
cutback_factor = 0.5
[]
[]
[Functions]
[flux]
type = ParsedFunction
vals = 'delta_xco2 dt'
vars = 'dx dt'
value = 'dx/dt'
[]
[]
[Postprocessors]
[total_co2_in_gas]
type = PorousFlowFluidMass
phase = 1
fluid_component = 1
[]
[total_co2_in_liquid]
type = PorousFlowFluidMass
phase = 0
fluid_component = 1
[]
[numdofs]
type = NumDOFs
[]
[delta_xco2]
type = ChangeOverTimePostprocessor
postprocessor = total_co2_in_liquid
[]
[dt]
type = TimestepSize
[]
[flux]
type = FunctionValuePostprocessor
function = flux
[]
[]
[Outputs]
print_linear_residuals = false
perf_graph = true
exodus = true
csv = true
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)
# Two phase density-driven convection of dissolved CO2 in brine
#
# Initially, the model has a gas phase at the top with a saturation of 0.29
# (which corresponds to an initial value of zi = 0.2).
# Diffusion of the dissolved CO2
# component from the saturated liquid to the unsaturated liquid below reduces the
# amount of CO2 in the gas phase. As the density of the CO2-saturated brine is greater
# than the unsaturated brine, a gravitational instability arises and density-driven
# convection of CO2-rich fingers descend into the unsaturated brine.
#
# The instability is seeded by a random perturbation to the porosity field.
# Mesh adaptivity is used to refine the mesh as the fingers form.
#
# Note: this model is computationally expensive, so should be run with multiple cores,
# preferably on a cluster.
[GlobalParams]
PorousFlowDictator = 'dictator'
gravity = '0 -9.81 0'
[]
[Adaptivity]
max_h_level = 2
marker = marker
initial_marker = initial
initial_steps = 2
[Indicators]
[indicator]
type = GradientJumpIndicator
variable = zi
[]
[]
[Markers]
[marker]
type = ErrorFractionMarker
indicator = indicator
refine = 0.8
[]
[initial]
type = BoxMarker
bottom_left = '0 1.95 0'
top_right = '2 2 0'
inside = REFINE
outside = DO_NOTHING
[]
[]
[]
[Mesh]
type = GeneratedMesh
dim = 2
ymax = 2
xmax = 2
ny = 40
nx = 40
bias_y = 0.95
[]
[Kernels]
[mass0]
type = PorousFlowMassTimeDerivative
fluid_component = 0
variable = pgas
[]
[flux0]
type = PorousFlowAdvectiveFlux
fluid_component = 0
variable = pgas
[]
[diff0]
type = PorousFlowDispersiveFlux
fluid_component = 0
variable = pgas
disp_long = '0 0'
disp_trans = '0 0'
[]
[mass1]
type = PorousFlowMassTimeDerivative
fluid_component = 1
variable = zi
[]
[flux1]
type = PorousFlowAdvectiveFlux
fluid_component = 1
variable = zi
[]
[diff1]
type = PorousFlowDispersiveFlux
fluid_component = 1
variable = zi
disp_long = '0 0'
disp_trans = '0 0'
[]
[]
[AuxVariables]
[xnacl]
initial_condition = 0.01
[]
[saturation_gas]
order = FIRST
family = MONOMIAL
[]
[xco2l]
order = FIRST
family = MONOMIAL
[]
[density_liquid]
order = FIRST
family = MONOMIAL
[]
[porosity]
order = CONSTANT
family = MONOMIAL
[]
[]
[AuxKernels]
[saturation_gas]
type = PorousFlowPropertyAux
variable = saturation_gas
property = saturation
phase = 1
execute_on = 'timestep_end'
[]
[xco2l]
type = PorousFlowPropertyAux
variable = xco2l
property = mass_fraction
phase = 0
fluid_component = 1
execute_on = 'timestep_end'
[]
[density_liquid]
type = PorousFlowPropertyAux
variable = density_liquid
property = density
phase = 0
execute_on = 'timestep_end'
[]
[]
[Variables]
[pgas]
[]
[zi]
scaling = 1e4
[]
[]
[ICs]
[pressure]
type = FunctionIC
function = 10e6-9.81*1000*y
variable = pgas
[]
[zi]
type = BoundingBoxIC
variable = zi
x1 = 0
x2 = 2
y1 = 1.95
y2 = 2
inside = 0.2
outside = 0
[]
[porosity]
type = RandomIC
variable = porosity
min = 0.25
max = 0.275
seed = 0
[]
[]
[UserObjects]
[dictator]
type = PorousFlowDictator
porous_flow_vars = 'pgas zi'
number_fluid_phases = 2
number_fluid_components = 2
[]
[pc]
type = PorousFlowCapillaryPressureConst
pc = 0
[]
[fs]
type = PorousFlowBrineCO2
brine_fp = brine
co2_fp = co2
capillary_pressure = pc
[]
[]
[Modules]
[FluidProperties]
[co2sw]
type = CO2FluidProperties
[]
[co2]
type = TabulatedFluidProperties
fp = co2sw
[]
[brine]
type = BrineFluidProperties
[]
[]
[]
[Materials]
[temperature]
type = PorousFlowTemperature
temperature = '45'
[]
[brineco2]
type = PorousFlowFluidState
gas_porepressure = 'pgas'
z = 'zi'
temperature_unit = Celsius
xnacl = 'xnacl'
capillary_pressure = pc
fluid_state = fs
[]
[porosity]
type = PorousFlowPorosityConst
porosity = porosity
[]
[permeability]
type = PorousFlowPermeabilityConst
permeability = '1e-11 0 0 0 1e-11 0 0 0 1e-11'
[]
[relperm_water]
type = PorousFlowRelativePermeabilityCorey
phase = 0
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[relperm_gas]
type = PorousFlowRelativePermeabilityCorey
phase = 1
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[diffusivity]
type = PorousFlowDiffusivityConst
diffusion_coeff = '2e-9 2e-9 2e-9 2e-9'
tortuosity = '1 1'
[]
[]
[Preconditioning]
active = basic
[mumps_is_best_for_parallel_jobs]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = ' lu mumps'
[]
[basic]
type = SMP
full = true
petsc_options_iname = '-ksp_type -pc_type -sub_pc_type -sub_pc_factor_shift_type -pc_asm_overlap'
petsc_options_value = 'gmres asm lu NONZERO 2 '
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
end_time = 1e6
nl_max_its = 25
l_max_its = 100
dtmax = 1e4
nl_abs_tol = 1e-6
[TimeStepper]
type = IterationAdaptiveDT
dt = 10
growth_factor = 2
cutback_factor = 0.5
[]
[]
[Functions]
[flux]
type = ParsedFunction
vals = 'delta_xco2 dt'
vars = 'dx dt'
value = 'dx/dt'
[]
[]
[Postprocessors]
[total_co2_in_gas]
type = PorousFlowFluidMass
phase = 1
fluid_component = 1
[]
[total_co2_in_liquid]
type = PorousFlowFluidMass
phase = 0
fluid_component = 1
[]
[numdofs]
type = NumDOFs
[]
[delta_xco2]
type = ChangeOverTimePostprocessor
postprocessor = total_co2_in_liquid
[]
[dt]
type = TimestepSize
[]
[flux]
type = FunctionValuePostprocessor
function = flux
[]
[]
[Outputs]
print_linear_residuals = false
perf_graph = true
exodus = true
csv = true
[]
(modules/porous_flow/examples/lava_lamp/2phase_convection.i)
# Two phase density-driven convection of dissolved CO2 in brine
#
# Initially, the model has a gas phase at the top with a saturation of 0.29
# (which corresponds to an initial value of zi = 0.2).
# Diffusion of the dissolved CO2
# component from the saturated liquid to the unsaturated liquid below reduces the
# amount of CO2 in the gas phase. As the density of the CO2-saturated brine is greater
# than the unsaturated brine, a gravitational instability arises and density-driven
# convection of CO2-rich fingers descend into the unsaturated brine.
#
# The instability is seeded by a random perturbation to the porosity field.
# Mesh adaptivity is used to refine the mesh as the fingers form.
#
# Note: this model is computationally expensive, so should be run with multiple cores,
# preferably on a cluster.
[GlobalParams]
PorousFlowDictator = 'dictator'
gravity = '0 -9.81 0'
[]
[Adaptivity]
max_h_level = 2
marker = marker
initial_marker = initial
initial_steps = 2
[Indicators]
[indicator]
type = GradientJumpIndicator
variable = zi
[]
[]
[Markers]
[marker]
type = ErrorFractionMarker
indicator = indicator
refine = 0.8
[]
[initial]
type = BoxMarker
bottom_left = '0 1.95 0'
top_right = '2 2 0'
inside = REFINE
outside = DO_NOTHING
[]
[]
[]
[Mesh]
type = GeneratedMesh
dim = 2
ymax = 2
xmax = 2
ny = 40
nx = 40
bias_y = 0.95
[]
[Kernels]
[mass0]
type = PorousFlowMassTimeDerivative
fluid_component = 0
variable = pgas
[]
[flux0]
type = PorousFlowAdvectiveFlux
fluid_component = 0
variable = pgas
[]
[diff0]
type = PorousFlowDispersiveFlux
fluid_component = 0
variable = pgas
disp_long = '0 0'
disp_trans = '0 0'
[]
[mass1]
type = PorousFlowMassTimeDerivative
fluid_component = 1
variable = zi
[]
[flux1]
type = PorousFlowAdvectiveFlux
fluid_component = 1
variable = zi
[]
[diff1]
type = PorousFlowDispersiveFlux
fluid_component = 1
variable = zi
disp_long = '0 0'
disp_trans = '0 0'
[]
[]
[AuxVariables]
[xnacl]
initial_condition = 0.01
[]
[saturation_gas]
order = FIRST
family = MONOMIAL
[]
[xco2l]
order = FIRST
family = MONOMIAL
[]
[density_liquid]
order = FIRST
family = MONOMIAL
[]
[porosity]
order = CONSTANT
family = MONOMIAL
[]
[]
[AuxKernels]
[saturation_gas]
type = PorousFlowPropertyAux
variable = saturation_gas
property = saturation
phase = 1
execute_on = 'timestep_end'
[]
[xco2l]
type = PorousFlowPropertyAux
variable = xco2l
property = mass_fraction
phase = 0
fluid_component = 1
execute_on = 'timestep_end'
[]
[density_liquid]
type = PorousFlowPropertyAux
variable = density_liquid
property = density
phase = 0
execute_on = 'timestep_end'
[]
[]
[Variables]
[pgas]
[]
[zi]
scaling = 1e4
[]
[]
[ICs]
[pressure]
type = FunctionIC
function = 10e6-9.81*1000*y
variable = pgas
[]
[zi]
type = BoundingBoxIC
variable = zi
x1 = 0
x2 = 2
y1 = 1.95
y2 = 2
inside = 0.2
outside = 0
[]
[porosity]
type = RandomIC
variable = porosity
min = 0.25
max = 0.275
seed = 0
[]
[]
[UserObjects]
[dictator]
type = PorousFlowDictator
porous_flow_vars = 'pgas zi'
number_fluid_phases = 2
number_fluid_components = 2
[]
[pc]
type = PorousFlowCapillaryPressureConst
pc = 0
[]
[fs]
type = PorousFlowBrineCO2
brine_fp = brine
co2_fp = co2
capillary_pressure = pc
[]
[]
[Modules]
[FluidProperties]
[co2sw]
type = CO2FluidProperties
[]
[co2]
type = TabulatedFluidProperties
fp = co2sw
[]
[brine]
type = BrineFluidProperties
[]
[]
[]
[Materials]
[temperature]
type = PorousFlowTemperature
temperature = '45'
[]
[brineco2]
type = PorousFlowFluidState
gas_porepressure = 'pgas'
z = 'zi'
temperature_unit = Celsius
xnacl = 'xnacl'
capillary_pressure = pc
fluid_state = fs
[]
[porosity]
type = PorousFlowPorosityConst
porosity = porosity
[]
[permeability]
type = PorousFlowPermeabilityConst
permeability = '1e-11 0 0 0 1e-11 0 0 0 1e-11'
[]
[relperm_water]
type = PorousFlowRelativePermeabilityCorey
phase = 0
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[relperm_gas]
type = PorousFlowRelativePermeabilityCorey
phase = 1
n = 2
s_res = 0.1
sum_s_res = 0.2
[]
[diffusivity]
type = PorousFlowDiffusivityConst
diffusion_coeff = '2e-9 2e-9 2e-9 2e-9'
tortuosity = '1 1'
[]
[]
[Preconditioning]
active = basic
[mumps_is_best_for_parallel_jobs]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = ' lu mumps'
[]
[basic]
type = SMP
full = true
petsc_options_iname = '-ksp_type -pc_type -sub_pc_type -sub_pc_factor_shift_type -pc_asm_overlap'
petsc_options_value = 'gmres asm lu NONZERO 2 '
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
end_time = 1e6
nl_max_its = 25
l_max_its = 100
dtmax = 1e4
nl_abs_tol = 1e-6
[TimeStepper]
type = IterationAdaptiveDT
dt = 10
growth_factor = 2
cutback_factor = 0.5
[]
[]
[Functions]
[flux]
type = ParsedFunction
vals = 'delta_xco2 dt'
vars = 'dx dt'
value = 'dx/dt'
[]
[]
[Postprocessors]
[total_co2_in_gas]
type = PorousFlowFluidMass
phase = 1
fluid_component = 1
[]
[total_co2_in_liquid]
type = PorousFlowFluidMass
phase = 0
fluid_component = 1
[]
[numdofs]
type = NumDOFs
[]
[delta_xco2]
type = ChangeOverTimePostprocessor
postprocessor = total_co2_in_liquid
[]
[dt]
type = TimestepSize
[]
[flux]
type = FunctionValuePostprocessor
function = flux
[]
[]
[Outputs]
print_linear_residuals = false
perf_graph = true
exodus = true
csv = true
[]