OGS
ProcessLib::ComponentTransport::ComponentTransportProcess Class Referencefinal

Detailed Description

ComponentTransport process

Governing equations

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.

Process Coupling

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.

Note
At the moment there is not any coupling by source or sink terms, i.e., the coupling is implemented only by density and viscosity changes due to concentration changes as well as by the temporal derivatives of each variable.

Definition at line 94 of file ComponentTransportProcess.h.

#include <ComponentTransportProcess.h>

Inheritance diagram for ProcessLib::ComponentTransport::ComponentTransportProcess:
[legend]
Collaboration diagram for ProcessLib::ComponentTransport::ComponentTransportProcess:
[legend]

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 &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)
 
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 &parameters, 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::MeshgetMesh () 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
 

Additional Inherited Members

- Public Types inherited from ProcessLib::Process
using NonlinearSolver = NumLib::NonlinearSolverBase
 
using TimeDiscretization = NumLib::TimeDiscretization
 
- Public Attributes inherited from ProcessLib::Process
std::string const name
 
- Protected Member Functions inherited from ProcessLib::Process
NumLib::ExtrapolatorgetExtrapolator () const
 
NumLib::LocalToGlobalIndexMap const & getSingleComponentDOFTable () const
 
void initializeProcessBoundaryConditionsAndSourceTerms (const NumLib::LocalToGlobalIndexMap &dof_table, const int process_id)
 
virtual void constructDofTable ()
 
void constructMonolithicProcessDofTable ()
 
void constructDofTableOfSpecifiedProcessStaggeredScheme (const int specified_prosess_id)
 
virtual std::tuple< NumLib::LocalToGlobalIndexMap *, bool > getDOFTableForExtrapolatorData () const
 
- Protected Attributes inherited from ProcessLib::Process
MeshLib::Mesh_mesh
 
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_all_nodes
 
std::unique_ptr< NumLib::LocalToGlobalIndexMap_local_to_global_index_map
 
SecondaryVariableCollection _secondary_variables
 
VectorMatrixAssembler _global_assembler
 
const bool _use_monolithic_scheme
 
CoupledSolutionsForStaggeredScheme_coupled_solutions
 
unsigned const _integration_order
 
std::vector< std::unique_ptr< IntegrationPointWriter > > _integration_point_writer
 
GlobalSparsityPattern _sparsity_pattern
 
std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > _process_variables
 
std::vector< BoundaryConditionCollection_boundary_conditions
 

Constructor & Destructor Documentation

◆ ComponentTransportProcess()

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.

44  : Process(std::move(name), mesh, std::move(jacobian_assembler), parameters,
45  integration_order, std::move(process_variables),
46  std::move(secondary_variables), use_monolithic_scheme),
47  _process_data(std::move(process_data)),
48  _surfaceflux(std::move(surfaceflux)),
49  _chemical_solver_interface(std::move(chemical_solver_interface))
50 {
51 }
std::unique_ptr< ChemistryLib::ChemicalSolverInterface > _chemical_solver_interface
std::unique_ptr< ProcessLib::SurfaceFluxData > _surfaceflux
std::string const name
Definition: Process.h:323
Process(std::string name_, MeshLib::Mesh &mesh, std::unique_ptr< 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, SecondaryVariableCollection &&secondary_variables, const bool use_monolithic_scheme=true)
Definition: Process.cpp:22

Member Function Documentation

◆ assembleConcreteProcess()

void ProcessLib::ComponentTransport::ComponentTransportProcess::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 
)
overrideprivatevirtual

Implements ProcessLib::Process.

Definition at line 143 of file ComponentTransportProcess.cpp.

147 {
148  DBUG("Assemble ComponentTransportProcess.");
149 
150  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
151 
152  std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
153  dof_tables;
155  {
156  dof_tables.push_back(std::ref(*_local_to_global_index_map));
157  }
158  else
159  {
160  std::generate_n(
161  std::back_inserter(dof_tables), _process_variables.size(),
162  [&]() { return std::ref(*_local_to_global_index_map); });
163  }
164  // Call global assembler for each local assembly item.
167  pv.getActiveElementIDs(), dof_tables, t, dt, x, xdot, process_id, M, K,
168  b);
169 }
void DBUG(char const *fmt, Args const &... args)
Definition: Logging.h:27
std::vector< std::unique_ptr< ComponentTransportLocalAssemblerInterface > > _local_assemblers
std::vector< std::size_t > const & getActiveElementIDs() const
std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > _process_variables
Definition: Process.h:360
std::vector< std::reference_wrapper< ProcessVariable > > const & getProcessVariables(const int process_id) const
Definition: Process.h:145
VectorMatrixAssembler _global_assembler
Definition: Process.h:333
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition: Process.h:329
const bool _use_monolithic_scheme
Definition: Process.h:335
void assemble(std::size_t const mesh_item_id, LocalAssemblerInterface &local_assembler, std::vector< std::reference_wrapper< NumLib::LocalToGlobalIndexMap >> const &dof_tables, double const t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &xdot, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b)
static void executeSelectedMemberDereferenced(Object &object, Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)

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().

◆ assembleWithJacobianConcreteProcess()

void ProcessLib::ComponentTransport::ComponentTransportProcess::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 
)
overrideprivatevirtual

Implements ProcessLib::Process.

Definition at line 171 of file ComponentTransportProcess.cpp.

176 {
177  DBUG("AssembleWithJacobian ComponentTransportProcess.");
178 
179  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
180  std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
181  dof_table = {std::ref(*_local_to_global_index_map)};
182  // Call global assembler for each local assembly item.
185  _local_assemblers, pv.getActiveElementIDs(), dof_table, t, dt, x, xdot,
186  dxdot_dx, dx_dx, process_id, M, K, b, Jac);
187 }
void assembleWithJacobian(std::size_t const mesh_item_id, LocalAssemblerInterface &local_assembler, std::vector< std::reference_wrapper< NumLib::LocalToGlobalIndexMap >> const &dof_tables, 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)

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().

◆ computeSecondaryVariableConcrete()

void ProcessLib::ComponentTransport::ComponentTransportProcess::computeSecondaryVariableConcrete ( double const  t,
double const  dt,
std::vector< GlobalVector * > const &  x,
GlobalVector const &  x_dot,
int const  process_id 
)
overridevirtual

Reimplemented from ProcessLib::Process.

Definition at line 358 of file ComponentTransportProcess.cpp.

364 {
365  if (process_id != 0)
366  {
367  return;
368  }
369 
370  std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
371  dof_tables.reserve(x.size());
372  std::generate_n(std::back_inserter(dof_tables), x.size(),
373  [&]() { return _local_to_global_index_map.get(); });
374 
375  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
378  _local_assemblers, pv.getActiveElementIDs(), dof_tables, t, dt, x,
379  x_dot, process_id);
380 }
virtual void computeSecondaryVariable(std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_dot, int const process_id)
static void executeSelectedMemberOnDereferenced(Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)

References _local_assemblers, ProcessLib::LocalAssemblerInterface::computeSecondaryVariable(), NumLib::SerialExecutor::executeSelectedMemberOnDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), and ProcessLib::Process::getProcessVariables().

◆ extrapolateIntegrationPointValuesToNodes()

void ProcessLib::ComponentTransport::ComponentTransportProcess::extrapolateIntegrationPointValuesToNodes ( const double  t,
std::vector< GlobalVector * > const &  integration_point_values_vectors,
std::vector< GlobalVector * > &  nodal_values_vectors 
)
overridevirtual

Reimplemented from ProcessLib::Process.

Definition at line 327 of file ComponentTransportProcess.cpp.

331 {
332  auto& extrapolator = getExtrapolator();
333  auto const extrapolatables =
335  &ComponentTransportLocalAssemblerInterface::
336  getInterpolatedLocalSolution);
337 
338  for (unsigned transport_process_id = 0;
339  transport_process_id < integration_point_values_vectors.size();
340  ++transport_process_id)
341  {
342  auto const& pv = _process_variables[transport_process_id + 1][0].get();
343  auto const& int_pt_C =
344  integration_point_values_vectors[transport_process_id];
345 
346  extrapolator.extrapolate(pv.getNumberOfGlobalComponents(),
347  extrapolatables, t, {int_pt_C},
348  {_local_to_global_index_map.get()});
349 
350  auto const& nodal_values = extrapolator.getNodalValues();
351  MathLib::LinAlg::copy(nodal_values,
352  *nodal_values_vectors[transport_process_id + 1]);
354  *nodal_values_vectors[transport_process_id + 1]);
355  }
356 }
NumLib::Extrapolator & getExtrapolator() const
Definition: Process.h:185
void finalizeAssembly(PETScMatrix &A)
Definition: LinAlg.cpp:162
void copy(PETScVector const &x, PETScVector &y)
Definition: LinAlg.cpp:37
ExtrapolatableLocalAssemblerCollection< LocalAssemblerCollection > makeExtrapolatable(LocalAssemblerCollection const &local_assemblers, IntegrationPointValuesMethod integration_point_values_method)

References _local_assemblers, ProcessLib::Process::_process_variables, MathLib::LinAlg::copy(), MathLib::LinAlg::finalizeAssembly(), ProcessLib::Process::getExtrapolator(), and NumLib::makeExtrapolatable().

Referenced by setInitialConditionsConcreteProcess().

◆ getFlux()

Eigen::Vector3d ProcessLib::ComponentTransport::ComponentTransportProcess::getFlux ( std::size_t const  element_id,
MathLib::Point3d const &  p,
double const  t,
std::vector< GlobalVector * > const &  x 
) const
overridevirtual

Reimplemented from ProcessLib::Process.

Definition at line 189 of file ComponentTransportProcess.cpp.

194 {
195  std::vector<GlobalIndexType> indices_cache;
196  auto const r_c_indices = NumLib::getRowColumnIndices(
197  element_id, *_local_to_global_index_map, indices_cache);
198 
199  std::vector<std::vector<GlobalIndexType>> indices_of_all_coupled_processes{
200  x.size(), r_c_indices.rows};
201  auto const local_xs =
202  getCoupledLocalSolutions(x, indices_of_all_coupled_processes);
203 
204  return _local_assemblers[element_id]->getFlux(p, t, local_xs);
205 }
NumLib::LocalToGlobalIndexMap::RowColumnIndices getRowColumnIndices(std::size_t const id, NumLib::LocalToGlobalIndexMap const &dof_table, std::vector< GlobalIndexType > &indices)
std::vector< double > getCoupledLocalSolutions(std::vector< GlobalVector * > const &global_solutions, std::vector< std::vector< GlobalIndexType >> const &indices)

References _local_assemblers, ProcessLib::Process::_local_to_global_index_map, ProcessLib::getCoupledLocalSolutions(), and NumLib::getRowColumnIndices().

◆ initializeConcreteProcess()

void ProcessLib::ComponentTransport::ComponentTransportProcess::initializeConcreteProcess ( NumLib::LocalToGlobalIndexMap const &  dof_table,
MeshLib::Mesh const &  mesh,
unsigned const  integration_order 
)
overrideprivatevirtual

Process specific initialization called by initialize().

Implements ProcessLib::Process.

Definition at line 53 of file ComponentTransportProcess.cpp.

57 {
58  _process_data.mesh_prop_velocity = MeshLib::getOrCreateMeshProperty<double>(
59  const_cast<MeshLib::Mesh&>(mesh), "velocity",
60  MeshLib::MeshItemType::Cell, mesh.getDimension());
61 
62  const int process_id = 0;
63  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
64 
65  std::vector<std::reference_wrapper<ProcessLib::ProcessVariable>>
66  transport_process_variables;
68  {
69  for (auto pv_iter = std::next(_process_variables[process_id].begin());
70  pv_iter != _process_variables[process_id].end();
71  ++pv_iter)
72  {
73  transport_process_variables.push_back(*pv_iter);
74  }
75  }
76  else
77  {
78  for (auto pv_iter = std::next(_process_variables.begin());
79  pv_iter != _process_variables.end();
80  ++pv_iter)
81  {
82  transport_process_variables.push_back((*pv_iter)[0]);
83  }
84  }
85 
86  ProcessLib::createLocalAssemblers<LocalAssemblerData>(
87  mesh.getDimension(), mesh.getElements(), dof_table,
89  mesh.isAxiallySymmetric(), integration_order, _process_data,
90  transport_process_variables);
91 
93  {
97 
98  _chemical_solver_interface->initialize();
99  }
100 
102  "darcy_velocity",
104  mesh.getDimension(), getExtrapolator(), _local_assemblers,
106 }
virtual std::vector< double > const & getIntPtDarcyVelocity(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
unsigned getShapeFunctionOrder() const
SecondaryVariableCollection _secondary_variables
Definition: Process.h:331
void addSecondaryVariable(std::string const &internal_name, SecondaryVariableFunctions &&fcts)
SecondaryVariableFunctions makeExtrapolator(const unsigned num_components, NumLib::Extrapolator &extrapolator, LocalAssemblerCollection const &local_assemblers, typename NumLib::ExtrapolatableLocalAssemblerCollection< LocalAssemblerCollection >::IntegrationPointValuesMethod integration_point_values_method)

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().

◆ isLinear()

bool ProcessLib::ComponentTransport::ComponentTransportProcess::isLinear ( ) const
inlineoverride

Definition at line 117 of file ComponentTransportProcess.h.

117 { return false; }

◆ postTimestepConcreteProcess()

void ProcessLib::ComponentTransport::ComponentTransportProcess::postTimestepConcreteProcess ( std::vector< GlobalVector * > const &  x,
const double  t,
const double  dt,
int const  process_id 
)
overridevirtual

Reimplemented from ProcessLib::Process.

Definition at line 382 of file ComponentTransportProcess.cpp.

387 {
388  if (process_id != 0)
389  {
390  return;
391  }
392 
393  std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
394  dof_tables.reserve(x.size());
395  std::generate_n(std::back_inserter(dof_tables), x.size(),
396  [&]() { return _local_to_global_index_map.get(); });
397 
398  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
401  _local_assemblers, pv.getActiveElementIDs(), dof_tables, x, t, dt);
402 
403  if (!_surfaceflux) // computing the surfaceflux is optional
404  {
405  return;
406  }
407  _surfaceflux->integrate(x, t, *this, process_id, _integration_order, _mesh,
408  pv.getActiveElementIDs());
409 }
virtual void postTimestep(std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, double const t, double const dt)
MeshLib::Mesh & _mesh
Definition: Process.h:326
unsigned const _integration_order
Definition: Process.h:344

References ProcessLib::Process::_integration_order, _local_assemblers, ProcessLib::Process::_mesh, _surfaceflux, NumLib::SerialExecutor::executeSelectedMemberOnDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), ProcessLib::Process::getProcessVariables(), and ProcessLib::LocalAssemblerInterface::postTimestep().

◆ setCoupledTermForTheStaggeredSchemeToLocalAssemblers()

void ProcessLib::ComponentTransport::ComponentTransportProcess::setCoupledTermForTheStaggeredSchemeToLocalAssemblers ( int const  process_id)
overridevirtual

Reimplemented from ProcessLib::Process.

Definition at line 207 of file ComponentTransportProcess.cpp.

209 {
210  DBUG("Set the coupled term for the staggered scheme to local assembers.");
211 
212  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
214  &ComponentTransportLocalAssemblerInterface::
215  setStaggeredCoupledSolutions,
217 }
CoupledSolutionsForStaggeredScheme * _coupled_solutions
Definition: Process.h:339

References ProcessLib::Process::_coupled_solutions, _local_assemblers, DBUG(), NumLib::SerialExecutor::executeSelectedMemberOnDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), and ProcessLib::Process::getProcessVariables().

◆ setInitialConditionsConcreteProcess()

void ProcessLib::ComponentTransport::ComponentTransportProcess::setInitialConditionsConcreteProcess ( std::vector< GlobalVector * > &  x,
double const  t,
int const  process_id 
)
overrideprivatevirtual

Reimplemented from ProcessLib::Process.

Definition at line 108 of file ComponentTransportProcess.cpp.

110 {
112  {
113  return;
114  }
115 
116  if (process_id != static_cast<int>(x.size() - 1))
117  {
118  return;
119  }
120 
121  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
122 
123  std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
124  dof_tables.reserve(x.size());
125  std::generate_n(std::back_inserter(dof_tables), x.size(),
126  [&]() { return _local_to_global_index_map.get(); });
127 
130  _local_assemblers, pv.getActiveElementIDs(), dof_tables, x, t);
131 
132  BaseLib::RunTime time_phreeqc;
133  time_phreeqc.start();
134 
135  _chemical_solver_interface->executeSpeciationCalculation(0. /*dt*/);
136 
138  t, _chemical_solver_interface->getIntPtProcessSolutions(), x);
139 
140  INFO("[time] Phreeqc took {:g} s.", time_phreeqc.elapsed());
141 }
void INFO(char const *fmt, Args const &... args)
Definition: Logging.h:32
Count the running time.
Definition: RunTime.h:29
double elapsed() const
Get the elapsed time in seconds.
Definition: RunTime.h:42
void start()
Start the timer.
Definition: RunTime.h:32
void initializeChemicalSystem(std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, double const t)
void extrapolateIntegrationPointValuesToNodes(const double t, std::vector< GlobalVector * > const &integration_point_values_vectors, std::vector< GlobalVector * > &nodal_values_vectors) override

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().

◆ solveReactionEquation()

void ProcessLib::ComponentTransport::ComponentTransportProcess::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 
)
overridevirtual

Reimplemented from ProcessLib::Process.

Definition at line 219 of file ComponentTransportProcess.cpp.

223 {
224  // todo (renchao): move chemical calculation to elsewhere.
225  if (_process_data.lookup_table && process_id == 0)
226  {
227  INFO("Update process solutions via the look-up table approach");
228  _process_data.lookup_table->lookup(x, x_prev, _mesh.getNumberOfNodes());
229 
230  return;
231  }
232 
234  {
235  return;
236  }
237 
238  // Sequential non-iterative approach applied here to split the reactive
239  // transport process into the transport stage followed by the reaction
240  // stage.
241  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
242 
243  std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
244  dof_tables.reserve(x.size());
245  std::generate_n(std::back_inserter(dof_tables), x.size(),
246  [&]() { return _local_to_global_index_map.get(); });
247 
248  if (process_id == 0)
249  {
252  _local_assemblers, pv.getActiveElementIDs(), dof_tables, x, t, dt);
253 
254  BaseLib::RunTime time_phreeqc;
255  time_phreeqc.start();
256 
257  _chemical_solver_interface->setAqueousSolutionsPrevFromDumpFile();
258 
259  _chemical_solver_interface->executeSpeciationCalculation(dt);
260 
261  INFO("[time] Phreeqc took {:g} s.", time_phreeqc.elapsed());
262 
264  &ComponentTransportLocalAssemblerInterface::
265  postSpeciationCalculation,
267 
268  return;
269  }
270 
271  auto const matrix_specification =
272  ode_sys.getMatrixSpecifications(process_id);
273 
274  std::size_t matrix_id = 0u;
276  matrix_specification, matrix_id);
278  matrix_specification, matrix_id);
279  auto& b =
281 
282  M.setZero();
283  K.setZero();
284  b.setZero();
285 
288  _local_assemblers, pv.getActiveElementIDs(), dof_tables, x, t, dt, M, K,
289  b, process_id);
290 
291  // todo (renchao): incorporate Neumann boundary condition
295 
297  matrix_specification, matrix_id);
298  auto& rhs =
300 
301  A.setZero();
302  rhs.setZero();
303 
304  // compute A
305  MathLib::LinAlg::copy(M, A);
306  MathLib::LinAlg::aypx(A, 1.0 / dt, K);
307 
308  // compute rhs
309  MathLib::LinAlg::matMult(M, *x[process_id], rhs);
310  MathLib::LinAlg::aypx(rhs, 1.0 / dt, b);
311 
312  using Tag = NumLib::NonlinearSolverTag;
314  auto& equation_system = static_cast<EqSys&>(ode_sys);
315  equation_system.applyKnownSolutionsPicard(A, rhs, *x[process_id]);
316 
317  auto& linear_solver =
319  linear_solver.solve(A, rhs, *x[process_id]);
320 
325 }
bool solve(EigenMatrix &A, EigenVector &b, EigenVector &x)
std::size_t getNumberOfNodes() const
Get the number of nodes.
Definition: Mesh.h:89
virtual MathLib::MatrixSpecifications getMatrixSpecifications(const int process_id) const =0
virtual GlobalMatrix & getMatrix(std::size_t &id)=0
Get an uninitialized matrix with the given id.
virtual void releaseMatrix(GlobalMatrix const &A)=0
virtual GlobalVector & getVector(std::size_t &id)=0
Get an uninitialized vector with the given id.
virtual void releaseVector(GlobalVector const &x)=0
void setChemicalSystem(std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, double const t, double const dt)
void assembleReactionEquation(std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, double const t, double const dt, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, int const process_id)
NonlinearSolverTag
Tag used to specify which nonlinear solver will be used.
Definition: Types.h:20
void matMult(PETScMatrix const &A, PETScVector const &x, PETScVector &y)
Definition: LinAlg.cpp:141
void aypx(PETScVector &y, double const a, PETScVector const &x)
Definition: LinAlg.cpp:50
void finalizeVectorAssembly(VEC_T &)
General function to finalize the vector assembly.
bool finalizeMatrixAssembly(MAT_T &)
static NUMLIB_EXPORT MatrixProvider & provider
static NUMLIB_EXPORT VectorProvider & provider
ChemistryLib::ChemicalSolverInterface *const chemical_solver_interface

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().

Member Data Documentation

◆ _chemical_solver_interface

std::unique_ptr<ChemistryLib::ChemicalSolverInterface> ProcessLib::ComponentTransport::ComponentTransportProcess::_chemical_solver_interface
private

◆ _local_assemblers

◆ _process_data

ComponentTransportProcessData ProcessLib::ComponentTransport::ComponentTransportProcess::_process_data
private

◆ _surfaceflux

std::unique_ptr<ProcessLib::SurfaceFluxData> ProcessLib::ComponentTransport::ComponentTransportProcess::_surfaceflux
private

Definition at line 176 of file ComponentTransportProcess.h.

Referenced by postTimestepConcreteProcess().


The documentation for this class was generated from the following files: