OGS
|
The flow process is described by
\[ \phi \frac{\partial \rho}{\partial p} \frac{\partial p}{\partial t} + \phi \frac{\partial \rho}{\partial C} \frac{\partial C}{\partial t} - \nabla \cdot \left[\frac{\kappa}{\mu(C)} \rho \nabla \left( p + \rho g z \right)\right] + Q_p = 0, \]
where the storage \(S\) has been substituted by \(\phi \frac{\partial \rho}{\partial p}\), \(\phi\) is the porosity, \(C\) is the concentration, \(p\) is the pressure, \(\kappa\) is permeability, \(\mu\) is viscosity of the fluid, \(\rho\) is the density of the fluid, and \(g\) is the gravitational acceleration.
The mass transport process is described by
\[ \phi R C \frac{\partial \rho}{\partial p} \frac{\partial p}{\partial t} + \phi R \left(\rho + C \frac{\partial \rho}{\partial C}\right) \frac{\partial C}{\partial t} - \nabla \cdot \left[\frac{\kappa}{\mu(C)} \rho C \nabla \left( p + \rho g z \right) + \rho D \nabla C\right] + Q_C + R \vartheta \phi \rho C = 0, \]
where \(R\) is the retardation factor, \(\vec{q} = -\frac{\kappa}{\mu(C)} \nabla \left( p + \rho g z \right)\) is the Darcy velocity, \(D\) is the hydrodynamic dispersion tensor, \(\vartheta\) is the decay rate.
For the hydrodynamic dispersion tensor the relation
\[ D = (\phi D_d + \beta_T \|\vec{q}\|) I + (\beta_L - \beta_T) \frac{\vec{q} \vec{q}^T}{\|\vec{q}\|} \]
is implemented, where \(D_d\) is the molecular diffusion coefficient, \(\beta_L\) the longitudinal dispersivity of chemical species, and \(\beta_T\) the transverse dispersivity of chemical species.
The implementation uses a monolithic approach, i.e., both processes are assembled within one global system of equations.
The advective term of the concentration equation is given by the confined groundwater flow process, i.e., the concentration distribution depends on Darcy velocity of the groundwater flow process. On the other hand the concentration dependencies of the viscosity and density in the groundwater flow couples the H process to the C process.
Definition at line 94 of file ComponentTransportProcess.h.
#include <ComponentTransportProcess.h>
Public Member Functions | |
ComponentTransportProcess (std::string name, MeshLib::Mesh &mesh, std::unique_ptr< ProcessLib::AbstractJacobianAssembler > &&jacobian_assembler, std::vector< std::unique_ptr< ParameterLib::ParameterBase >> const ¶meters, unsigned const integration_order, std::vector< std::vector< std::reference_wrapper< ProcessVariable >>> &&process_variables, ComponentTransportProcessData &&process_data, SecondaryVariableCollection &&secondary_variables, bool const use_monolithic_scheme, std::unique_ptr< ProcessLib::SurfaceFluxData > &&surfaceflux, std::unique_ptr< ChemistryLib::ChemicalSolverInterface > &&chemical_solver_interface) | |
Eigen::Vector3d | getFlux (std::size_t const element_id, MathLib::Point3d const &p, double const t, std::vector< GlobalVector * > const &x) const override |
void | setCoupledTermForTheStaggeredSchemeToLocalAssemblers (int const process_id) override |
void | solveReactionEquation (std::vector< GlobalVector * > &x, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, NumLib::EquationSystem &ode_sys, int const process_id) override |
void | extrapolateIntegrationPointValuesToNodes (const double t, std::vector< GlobalVector * > const &integration_point_values_vectors, std::vector< GlobalVector * > &nodal_values_vectors) override |
void | computeSecondaryVariableConcrete (double const, double const, std::vector< GlobalVector * > const &x, GlobalVector const &, int const) override |
void | postTimestepConcreteProcess (std::vector< GlobalVector * > const &x, const double t, const double dt, int const process_id) override |
ODESystem interface | |
bool | isLinear () const override |
Public Member Functions inherited from ProcessLib::Process | |
Process (std::string name_, MeshLib::Mesh &mesh, std::unique_ptr< AbstractJacobianAssembler > &&jacobian_assembler, std::vector< std::unique_ptr< ParameterLib::ParameterBase >> const ¶meters, unsigned const integration_order, std::vector< std::vector< std::reference_wrapper< ProcessVariable >>> &&process_variables, SecondaryVariableCollection &&secondary_variables, const bool use_monolithic_scheme=true) | |
void | preTimestep (std::vector< GlobalVector * > const &x, const double t, const double delta_t, const int process_id) |
Preprocessing before starting assembly for new timestep. More... | |
void | postTimestep (std::vector< GlobalVector * > const &x, const double t, const double delta_t, int const process_id) |
Postprocessing after a complete timestep. More... | |
void | postNonLinearSolver (GlobalVector const &x, GlobalVector const &xdot, const double t, double const dt, int const process_id) |
void | preIteration (const unsigned iter, GlobalVector const &x) final |
void | computeSecondaryVariable (double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_dot, int const process_id) |
compute secondary variables for the coupled equations or for output. More... | |
NumLib::IterationResult | postIteration (GlobalVector const &x) final |
void | initialize () |
void | setInitialConditions (std::vector< GlobalVector * > &process_solutions, std::vector< GlobalVector * > const &process_solutions_prev, double const t, int const process_id) |
MathLib::MatrixSpecifications | getMatrixSpecifications (const int process_id) const override |
void | setCoupledSolutionsForStaggeredScheme (CoupledSolutionsForStaggeredScheme *const coupled_solutions) |
void | updateDeactivatedSubdomains (double const time, const int process_id) |
bool | isMonolithicSchemeUsed () const |
void | preAssemble (const double t, double const dt, GlobalVector const &x) final |
void | assemble (const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &xdot, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b) final |
void | assembleWithJacobian (const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &xdot, const double dxdot_dx, const double dx_dx, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac) final |
std::vector< NumLib::IndexValueVector< GlobalIndexType > > const * | getKnownSolutions (double const t, GlobalVector const &x, int const process_id) const final |
virtual NumLib::LocalToGlobalIndexMap const & | getDOFTable (const int) const |
MeshLib::Mesh & | getMesh () const |
std::vector< std::reference_wrapper< ProcessVariable > > const & | getProcessVariables (const int process_id) const |
SecondaryVariableCollection const & | getSecondaryVariables () const |
std::vector< std::unique_ptr< IntegrationPointWriter > > const * | getIntegrationPointWriter (MeshLib::Mesh const &mesh) const |
Private Member Functions | |
void | initializeConcreteProcess (NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order) override |
Process specific initialization called by initialize(). More... | |
void | setInitialConditionsConcreteProcess (std::vector< GlobalVector * > &x, double const t, int const process_id) override |
void | assembleConcreteProcess (const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &xdot, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b) override |
void | assembleWithJacobianConcreteProcess (const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &xdot, const double dxdot_dx, const double dx_dx, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac) override |
Private Attributes | |
ComponentTransportProcessData | _process_data |
std::vector< std::unique_ptr< ComponentTransportLocalAssemblerInterface > > | _local_assemblers |
std::unique_ptr< ProcessLib::SurfaceFluxData > | _surfaceflux |
std::unique_ptr< ChemistryLib::ChemicalSolverInterface > | _chemical_solver_interface |
ProcessLib::ComponentTransport::ComponentTransportProcess::ComponentTransportProcess | ( | std::string | name, |
MeshLib::Mesh & | mesh, | ||
std::unique_ptr< ProcessLib::AbstractJacobianAssembler > && | jacobian_assembler, | ||
std::vector< std::unique_ptr< ParameterLib::ParameterBase >> const & | parameters, | ||
unsigned const | integration_order, | ||
std::vector< std::vector< std::reference_wrapper< ProcessVariable >>> && | process_variables, | ||
ComponentTransportProcessData && | process_data, | ||
SecondaryVariableCollection && | secondary_variables, | ||
bool const | use_monolithic_scheme, | ||
std::unique_ptr< ProcessLib::SurfaceFluxData > && | surfaceflux, | ||
std::unique_ptr< ChemistryLib::ChemicalSolverInterface > && | chemical_solver_interface | ||
) |
Definition at line 30 of file ComponentTransportProcess.cpp.
|
overrideprivatevirtual |
Implements ProcessLib::Process.
Definition at line 143 of file ComponentTransportProcess.cpp.
References ProcessLib::Process::_global_assembler, _local_assemblers, ProcessLib::Process::_local_to_global_index_map, ProcessLib::Process::_process_variables, ProcessLib::Process::_use_monolithic_scheme, ProcessLib::VectorMatrixAssembler::assemble(), DBUG(), NumLib::SerialExecutor::executeSelectedMemberDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), and ProcessLib::Process::getProcessVariables().
|
overrideprivatevirtual |
Implements ProcessLib::Process.
Definition at line 171 of file ComponentTransportProcess.cpp.
References ProcessLib::Process::_global_assembler, _local_assemblers, ProcessLib::Process::_local_to_global_index_map, ProcessLib::VectorMatrixAssembler::assembleWithJacobian(), DBUG(), NumLib::SerialExecutor::executeSelectedMemberDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), and ProcessLib::Process::getProcessVariables().
|
overridevirtual |
Reimplemented from ProcessLib::Process.
Definition at line 358 of file ComponentTransportProcess.cpp.
References _local_assemblers, ProcessLib::LocalAssemblerInterface::computeSecondaryVariable(), NumLib::SerialExecutor::executeSelectedMemberOnDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), and ProcessLib::Process::getProcessVariables().
|
overridevirtual |
Reimplemented from ProcessLib::Process.
Definition at line 327 of file ComponentTransportProcess.cpp.
References _local_assemblers, ProcessLib::Process::_process_variables, MathLib::LinAlg::copy(), MathLib::LinAlg::finalizeAssembly(), ProcessLib::Process::getExtrapolator(), and NumLib::makeExtrapolatable().
Referenced by setInitialConditionsConcreteProcess().
|
overridevirtual |
Reimplemented from ProcessLib::Process.
Definition at line 189 of file ComponentTransportProcess.cpp.
References _local_assemblers, ProcessLib::Process::_local_to_global_index_map, ProcessLib::getCoupledLocalSolutions(), and NumLib::getRowColumnIndices().
|
overrideprivatevirtual |
Process specific initialization called by initialize().
Implements ProcessLib::Process.
Definition at line 53 of file ComponentTransportProcess.cpp.
References _chemical_solver_interface, _local_assemblers, _process_data, ProcessLib::Process::_process_variables, ProcessLib::Process::_secondary_variables, ProcessLib::Process::_use_monolithic_scheme, ProcessLib::SecondaryVariableCollection::addSecondaryVariable(), MeshLib::Cell, NumLib::SerialExecutor::executeSelectedMemberOnDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), MeshLib::Mesh::getDimension(), MeshLib::Mesh::getElements(), ProcessLib::Process::getExtrapolator(), ProcessLib::ComponentTransport::ComponentTransportLocalAssemblerInterface::getIntPtDarcyVelocity(), ProcessLib::Process::getProcessVariables(), ProcessLib::ProcessVariable::getShapeFunctionOrder(), MeshLib::Mesh::isAxiallySymmetric(), ProcessLib::makeExtrapolator(), ProcessLib::ComponentTransport::ComponentTransportProcessData::mesh_prop_velocity, and ProcessLib::ComponentTransport::ComponentTransportLocalAssemblerInterface::setChemicalSystemID().
|
inlineoverride |
Definition at line 117 of file ComponentTransportProcess.h.
|
overridevirtual |
Reimplemented from ProcessLib::Process.
Definition at line 382 of file ComponentTransportProcess.cpp.
References ProcessLib::Process::_integration_order, _local_assemblers, ProcessLib::Process::_mesh, _surfaceflux, NumLib::SerialExecutor::executeSelectedMemberOnDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), ProcessLib::Process::getProcessVariables(), and ProcessLib::LocalAssemblerInterface::postTimestep().
|
overridevirtual |
Reimplemented from ProcessLib::Process.
Definition at line 207 of file ComponentTransportProcess.cpp.
References ProcessLib::Process::_coupled_solutions, _local_assemblers, DBUG(), NumLib::SerialExecutor::executeSelectedMemberOnDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), and ProcessLib::Process::getProcessVariables().
|
overrideprivatevirtual |
Reimplemented from ProcessLib::Process.
Definition at line 108 of file ComponentTransportProcess.cpp.
References _chemical_solver_interface, _local_assemblers, BaseLib::RunTime::elapsed(), NumLib::SerialExecutor::executeSelectedMemberOnDereferenced(), extrapolateIntegrationPointValuesToNodes(), ProcessLib::ProcessVariable::getActiveElementIDs(), ProcessLib::Process::getProcessVariables(), INFO(), ProcessLib::ComponentTransport::ComponentTransportLocalAssemblerInterface::initializeChemicalSystem(), and BaseLib::RunTime::start().
|
overridevirtual |
Reimplemented from ProcessLib::Process.
Definition at line 219 of file ComponentTransportProcess.cpp.
References _chemical_solver_interface, _local_assemblers, ProcessLib::Process::_mesh, _process_data, ProcessLib::ComponentTransport::ComponentTransportLocalAssemblerInterface::assembleReactionEquation(), MathLib::LinAlg::aypx(), ProcessLib::ComponentTransport::ComponentTransportProcessData::chemical_solver_interface, MathLib::LinAlg::copy(), BaseLib::RunTime::elapsed(), NumLib::SerialExecutor::executeSelectedMemberOnDereferenced(), MathLib::finalizeMatrixAssembly(), MathLib::finalizeVectorAssembly(), ProcessLib::ProcessVariable::getActiveElementIDs(), NumLib::MatrixProvider::getMatrix(), NumLib::EquationSystem::getMatrixSpecifications(), MeshLib::Mesh::getNumberOfNodes(), ProcessLib::Process::getProcessVariables(), NumLib::VectorProvider::getVector(), INFO(), ChemistryLib::ChemicalSolverInterface::linear_solver, ProcessLib::ComponentTransport::ComponentTransportProcessData::lookup_table, MathLib::LinAlg::matMult(), NumLib::GlobalVectorProvider::provider, NumLib::GlobalMatrixProvider::provider, NumLib::MatrixProvider::releaseMatrix(), NumLib::VectorProvider::releaseVector(), ProcessLib::ComponentTransport::ComponentTransportLocalAssemblerInterface::setChemicalSystem(), MathLib::EigenVector::setZero(), MathLib::EigenLisLinearSolver::solve(), and BaseLib::RunTime::start().
|
private |
Definition at line 179 of file ComponentTransportProcess.h.
Referenced by initializeConcreteProcess(), setInitialConditionsConcreteProcess(), and solveReactionEquation().
|
private |
Definition at line 174 of file ComponentTransportProcess.h.
Referenced by assembleConcreteProcess(), assembleWithJacobianConcreteProcess(), computeSecondaryVariableConcrete(), extrapolateIntegrationPointValuesToNodes(), getFlux(), initializeConcreteProcess(), postTimestepConcreteProcess(), setCoupledTermForTheStaggeredSchemeToLocalAssemblers(), setInitialConditionsConcreteProcess(), and solveReactionEquation().
|
private |
Definition at line 171 of file ComponentTransportProcess.h.
Referenced by initializeConcreteProcess(), and solveReactionEquation().
|
private |
Definition at line 176 of file ComponentTransportProcess.h.
Referenced by postTimestepConcreteProcess().