porepy.models.fluid_mass_balance module

Class types:
MassBalanceEquations defines subdomain and interface equations through the

terms entering. Darcy type interface relation is assumed.

Specific ConstitutiveLaws and specific SolutionStrategy for both incompressible and compressible case.

Notes

Apertures and specific volumes are not included.

Refactoring needed for constitutive equations. Modularisation and moving to the library.

Upwind for the mobility of the fluid flux is not complete.

class BoundaryConditionsSinglePhaseFlow[source]

Bases: object

Boundary conditions for single-phase flow.

bc_type_darcy(sd)[source]

Dirichlet conditions on all external boundaries.

Parameters

sd (Grid) – Subdomain grid on which to define boundary conditions.

Returns

Boundary condition object.

Return type

BoundaryCondition

bc_type_mobrho(sd)[source]

Dirichlet conditions on all external boundaries.

Parameters

sd (Grid) – Subdomain grid on which to define boundary conditions.

Returns

Boundary condition object.

Return type

BoundaryCondition

bc_values_darcy(subdomains)[source]

Not sure where this one should reside. Note that we could remove the grid_operator BC and DirBC, probably also ParameterArray/Matrix (unless needed to get rid of pp.ad.Discretization. I don’t see how it would be, though).

Parameters

subdomains (list[porepy.grids.grid.Grid]) – List of subdomains.

Returns

Ad array representing the boundary condition values for the Darcy flux.

Return type

Array

bc_values_mobrho(subdomains)[source]

Boundary condition values for the mobility times density.

Units for Dirichlet: kg * m^-3 * Pa^-1 * s^-1

Parameters:
Value is tricky if ..math:

mobility = rho / mu

with

ho and mu being functions of p (or other variables), since variables

are not defined at the boundary. This may lead to inconsistency between boundary conditions for Darcy flux and mobility. For now, we assume that the mobility is constant. TODO: Better solution. Could involve defining boundary grids.

Parameters:

subdomains: List of subdomains.

Returns:

Array with boundary values for the mobility.

Parameters

subdomains (list[porepy.grids.grid.Grid]) –

Return type

Array

domain_boundary_sides: Callable[[Grid], tuple[numpy.ndarray, numpy.ndarray, numpy.ndarray, numpy.ndarray, numpy.ndarray, numpy.ndarray, numpy.ndarray]]

Boundary sides of the domain. Normally defined in a mixin instance of ModelGeometry.

fluid: FluidConstants

Fluid constant object that takes care of scaling of fluid-related quantities. Normally, this is set by a mixin of instance SolutionStrategy.

class ConstitutiveLawsSinglePhaseFlow[source]

Bases: DarcysLaw, DimensionReduction, AdvectiveFlux, ConstantPorosity, ConstantPermeability, FluidDensityFromPressure, ConstantViscosity, FluidMobility

Constitutive equations for single-phase flow.

The combined laws access the following material constants:

solid:

permeability normal_permeability porosity

fluid:

viscosity density compressibility

basis: Callable[[Sequence[pp.GridLike], int], list[pp.ad.Matrix]]

Basis for the local coordinate system. Normally set by a mixin instance of porepy.models.geometry.ModelGeometry.

bc_values_darcy: Callable[[list[pp.Grid]], pp.ad.Array]

Darcy flux boundary conditions. Normally defined in a mixin instance of BoundaryConditionsSinglePhaseFlow.

darcy_keyword: str

Keyword used to identify the Darcy flux discretization. Normally set by a mixin instance of SolutionStrategySinglePhaseFlow.

fluid: pp.FluidConstants

Fluid constant object that takes care of scaling of fluid-related quantities. Normally, this is set by a mixin of instance SolutionStrategy.

interface_darcy_flux: Callable[[list[pp.MortarGrid]], pp.ad.MixedDimensionalVariable]

Darcy flux variables on interfaces. Normally defined in a mixin instance of VariablesSinglePhaseFlow.

interfaces_to_subdomains: Callable[[list[pp.MortarGrid]], list[pp.Grid]]

Map from interfaces to the adjacent subdomains. Normally defined in a mixin instance of ModelGeometry.

internal_boundary_normal_to_outwards: Callable[[list[pp.Grid], int], pp.ad.Matrix]

Switch interface normal vectors to point outwards from the subdomain. Normally set by a mixin instance of porepy.models.geometry.ModelGeometry.

mdg: pp.MixedDimensionalGrid

Mixed dimensional grid for the current model. Normally defined in a mixin instance of ModelGeometry.

nd: int

Ambient dimension of the problem. Normally set by a mixin instance of porepy.models.geometry.ModelGeometry.

outwards_internal_boundary_normals: Callable[[list[pp.MortarGrid], bool], pp.ad.Operator]

Outwards normal vectors on internal boundaries. Normally defined in a mixin instance of ModelGeometry.

pressure: Callable[[list[pp.Grid]], pp.ad.MixedDimensionalVariable]

Pressure variable. Normally defined in a mixin instance of VariablesSinglePhaseFlow.

subdomains_to_interfaces: Callable[[list[pp.Grid], list[int]], list[pp.MortarGrid]]

Map from subdomains to the adjacent interfaces. Normally defined in a mixin instance of ModelGeometry.

wrap_grid_attribute: Callable[[Sequence[pp.GridLike], str, int, bool], pp.ad.Matrix]

Wrap grid attributes as Ad operators. Normally set by a mixin instance of porepy.models.geometry.ModelGeometry.

class MassBalanceEquations[source]

Bases: BalanceEquation

Mixed-dimensional mass balance equation.

Balance equation for all subdomains and Darcy-type flux relation on all interfaces of codimension one.

FIXME: Well equations? Low priority.

fluid_flux(subdomains)[source]

Fluid flux.

Darcy flux times density and mobility.

Parameters

subdomains (list[porepy.grids.grid.Grid]) – List of subdomains.

Returns

Operator representing the fluid flux.

Return type

Operator

fluid_mass(subdomains)[source]

The full measure of cell-wise fluid mass.

The product of fluid density and porosity is assumed constant cell-wise, and integrated over the cell volume.

Note

This implementation assumes constant porosity and must be overridden for variable porosity. This has to do with wrapping of scalars as vectors or matrices and will hopefully be improved in the future. Extension to variable density is straightforward.

Parameters

subdomains (list[porepy.grids.grid.Grid]) – List of subdomains.

Returns

Operator representing the cell-wise fluid mass.

Return type

Operator

fluid_source(subdomains)[source]

Fluid source term.

Includes both external sources and inflow from neighboring subdomains of higher dimension (via interfaces).

Note

When overriding this method to assign internal fluid sources, one is advised to call the base class method and add the new contribution, thus ensuring that the source term includes the contribution from the interface fluxes.

Parameters

subdomains (list[porepy.grids.grid.Grid]) – List of subdomains.

Returns

Operator representing the source term.

Return type

Operator

interface_fluid_flux(interfaces)[source]

Interface fluid flux.

Parameters

interfaces (list[porepy.grids.mortar_grid.MortarGrid]) – List of interface grids.

Returns

Operator representing the interface fluid flux.

Return type

Operator

interface_flux_equation(interfaces)[source]

Interface flux equation.

Parameters

interfaces (list[porepy.grids.mortar_grid.MortarGrid]) – List of interface grids.

Returns

Operator representing the interface flux equation.

Return type

Operator

mass_balance_equation(subdomains)[source]

Mass balance equation for subdomains.

Parameters

subdomains (list[porepy.grids.grid.Grid]) – List of subdomains.

Returns

Operator representing the mass balance equation.

Return type

Operator

set_equations()[source]

Set the equations for the mass balance problem.

A mass balance equation is set for all subdomains and a Darcy-type flux relation is set for all interfaces of codimension one.

advective_flux: Callable[[list[pp.Grid], pp.ad.Operator, pp.ad.UpwindAd, pp.ad.Operator, Callable[[list[pp.MortarGrid]], pp.ad.Operator]], pp.ad.Operator]

Ad operator representing the advective flux. Normally provided by a mixin instance of AdvectiveFlux.

bc_values_mobrho: Callable[[list[pp.Grid]], pp.ad.Array]

Mobility times density boundary conditions. Normally defined in a mixin instance of BoundaryConditionsSinglePhaseFlow.

equation_system: pp.ad.EquationSystem

EquationSystem object for the current model. Normally defined in a mixin class defining the solution strategy.

fluid_density: Callable[[list[pp.Grid]], pp.ad.Operator]

Fluid density. Defined in a mixin class with a suitable constitutive relation.

interface_advective_flux: Callable[[list[pp.MortarGrid], pp.ad.Operator, pp.ad.UpwindCouplingAd], pp.ad.Operator]

Ad operator representing the advective flux on internal boundaries. Normally provided by a mixin instance of AdvectiveFlux.

interface_darcy_flux: Callable[[list[pp.MortarGrid]], pp.ad.MixedDimensionalVariable]

Darcy flux variable on interfaces. Normally defined in a mixin instance of VariablesSinglePhaseFlow.

interface_darcy_flux_equation: Callable[[list[pp.MortarGrid]], pp.ad.Operator]

Interface Darcy flux equation. Normally provided by a mixin instance of DarcysLaw.

interface_mobility_discretization: Callable[[list[pp.MortarGrid]], pp.ad.UpwindCouplingAd]

Discretization of the fluid mobility on internal boundaries. Normally provided by a mixin instance of FluidMobility.

interfaces_to_subdomains: Callable[[list[pp.MortarGrid]], list[pp.Grid]]

Map from interfaces to the adjacent subdomains. Normally defined in a mixin instance of ModelGeometry.

mdg: pp.MixedDimensionalGrid

Mixed dimensional grid for the current model. Normally defined in a mixin instance of ModelGeometry.

mobility: Callable[[list[pp.Grid]], pp.ad.Operator]

Fluid mobility. Normally provided by a mixin instance of FluidMobility.

mobility_discretization: Callable[[list[pp.Grid]], pp.ad.UpwindAd]

Discretization of the fluid mobility. Normally provided by a mixin instance of FluidMobility.

porosity: Callable[[list[pp.Grid]], pp.ad.Operator]

Porosity of the rock. Normally provided by a mixin instance of ConstantPorosity or a subclass thereof.

subdomains_to_interfaces: Callable[[list[pp.Grid], list[int]], list[pp.MortarGrid]]

Map from subdomains to the adjacent interfaces. Normally defined in a mixin instance of ModelGeometry.

class SinglePhaseFlow(params=None)[source]

Bases: MassBalanceEquations, VariablesSinglePhaseFlow, ConstitutiveLawsSinglePhaseFlow, BoundaryConditionsSinglePhaseFlow, SolutionStrategySinglePhaseFlow, ModelGeometry, DataSavingMixin

Class for single-phase flow in mixed-dimensional porous media.

Parameters

params (Optional[dict]) –

equation_system: pp.ad.EquationSystem

EquationSystem object for the current model. Normally defined in a mixin class defining the solution strategy.

mdg: pp.MixedDimensionalGrid

Mixed dimensional grid for the current model. Normally defined in a mixin instance of ModelGeometry.

nd: int

Ambient dimension of the problem. Normally set by a mixin instance of porepy.models.geometry.ModelGeometry.

time_manager: pp.TimeManager

Time manager. Normally set by a mixin instance of porepy.models.solution_strategy.SolutionStrategy.

class SolutionStrategySinglePhaseFlow(params=None)[source]

Bases: SolutionStrategy

Setup and numerics-related methods for a single-phase flow problem.

At some point, this will be refined to be a more sophisticated (modularised) solution strategy class. More refactoring may be beneficial.

This is not a full-scale model (in the old sense), but must be mixed with balance equations, constitutive laws etc. See user_examples.

Parameters

params (Optional[dict]) – Parameters for the solution strategy.

before_nonlinear_iteration()[source]

Evaluate Darcy flux for each subdomain and interface and store in the data dictionary for use in upstream weighting.

initial_condition()[source]

New formulation requires darcy flux (the flux is “advective” with mobilities included).

Return type

None

set_discretization_parameters()[source]

Set default (unitary/zero) parameters for the flow problem.

The parameter fields of the data dictionaries are updated for all subdomains and interfaces. The data to be set is related to:

  • The fluid diffusion, e.g., the permeability and boundary conditions for the pressure. This applies to both subdomains and interfaces.

  • Boundary conditions for the advective flux. This applies to subdomains only.

Return type

None

bc_type_darcy: Callable[[pp.Grid], pp.BoundaryCondition]

Function that returns the boundary condition type for the Darcy flux. Normally provided by a mixin instance of BoundaryConditionsSinglePhaseFlow.

bc_type_mobrho: Callable[[pp.Grid], pp.BoundaryCondition]

Function that returns the boundary condition type for the advective flux. Normally provided by a mixin instance of BoundaryConditionsSinglePhaseFlow.

darcy_keyword: str

Keyword for Darcy flux term.

Used to access discretization parameters and store discretization matrices.

interface_darcy_flux_variable: str

Name of the primary variable representing the Darcy flux on the interface.

mobility_keyword: str

Keyword for mobility factor.

Used to access discretization parameters and store discretization matrices.

nd: int

Ambient dimension of the problem. Normally set by a mixin instance of porepy.models.geometry.ModelGeometry.

permeability: Callable[[list[pp.Grid]], np.ndarray]

Function that returns the permeability of a subdomain. Normally provided by a mixin class with a suitable permeability definition.

pressure_variable: str

Name of the pressure variable.

specific_volume: Callable[[list[pp.Grid]], pp.ad.Operator]

Function that returns the specific volume of a subdomain. Normally provided by a mixin of instance DimensionReduction.

class VariablesSinglePhaseFlow[source]

Bases: VariableMixin

Creates necessary variables (pressure, interface flux) and provides getter methods for these and their reference values. Getters construct mixed-dimensional variables on the fly, and can be called on any subset of the grids where the variable is defined. Setter method (assig_variables), however, must create on all grids where the variable is to be used.

Note

Wrapping in class methods and not calling equation_system directly allows for easier changes of primary variables. As long as all calls to fluid_flux() accept Operators as return values, we can in theory add it as a primary variable and solved mixed form. Similarly for different formulations of the pressure (e.g. pressure head) or enthalpy/ temperature for the energy equation.

create_variables()[source]

Assign primary variables to subdomains and interfaces of the mixed-dimensional grid.

Return type

None

interface_darcy_flux(interfaces)[source]

Interface Darcy flux.

Parameters

interfaces (list[porepy.grids.mortar_grid.MortarGrid]) – List of interface grids.

Returns

Variable representing the interface Darcy flux.

Return type

MixedDimensionalVariable

pressure(subdomains)[source]
Parameters

subdomains (list[porepy.grids.grid.Grid]) –

Return type

MixedDimensionalVariable

reference_pressure(subdomains)[source]

Reference pressure.

Parameters

subdomains (list[porepy.grids.grid.Grid]) – List of subdomains.

Returns

Operator representing the reference pressure.

Return type

Operator

equation_system: EquationSystem

EquationSystem object for the current model. Normally defined in a mixin class defining the solution strategy.

fluid: FluidConstants

Fluid constant object that takes care of scaling of fluid-related quantities. Normally, this is set by a mixin of instance SolutionStrategy.

interface_darcy_flux_variable: str

Name of the primary variable representing the Darcy flux across an interface. Normally defined in a mixin of instance SolutionStrategySinglePhaseFlow.

mdg: MixedDimensionalGrid

Mixed dimensional grid for the current model. Normally defined in a mixin instance of ModelGeometry.

pressure_variable: str

Name of the primary variable representing the pressure. Normally defined in a mixin of instance SolutionStrategySinglePhaseFlow.