OGS
ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess Class Referencefinal

Detailed Description

A class to simulate the isothermal two-phase flow process with P-P model in porous media.

The gas and capillary pressure are used as primary variables.

Definition at line 39 of file TwoPhaseFlowWithPPProcess.h.

#include <TwoPhaseFlowWithPPProcess.h>

Inheritance diagram for ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess:
[legend]
Collaboration diagram for ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess:
[legend]

Public Member Functions

 TwoPhaseFlowWithPPProcess (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, TwoPhaseFlowWithPPProcessData &&process_data, SecondaryVariableCollection &&secondary_variables, std::map< std::string, std::unique_ptr< MathLib::PiecewiseLinearInterpolation >> const &curves)
 
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
 
virtual void setCoupledTermForTheStaggeredSchemeToLocalAssemblers (int const)
 
virtual void extrapolateIntegrationPointValuesToNodes (const double, std::vector< GlobalVector * > const &, std::vector< GlobalVector * > &)
 
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
 
virtual Eigen::Vector3d getFlux (std::size_t, MathLib::Point3d const &, double const, std::vector< GlobalVector * > const &) const
 
virtual void solveReactionEquation (std::vector< GlobalVector * > &, std::vector< GlobalVector * > const &, double const, double const, NumLib::EquationSystem &, int 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 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

TwoPhaseFlowWithPPProcessData _process_data
 
std::vector< std::unique_ptr< TwoPhaseFlowWithPPLocalAssemblerInterface > > _local_assemblers
 

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

◆ TwoPhaseFlowWithPPProcess()

ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess::TwoPhaseFlowWithPPProcess ( 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,
TwoPhaseFlowWithPPProcessData &&  process_data,
SecondaryVariableCollection &&  secondary_variables,
std::map< std::string, std::unique_ptr< MathLib::PiecewiseLinearInterpolation >> const &  curves 
)

Definition at line 20 of file TwoPhaseFlowWithPPProcess.cpp.

33  : Process(std::move(name), mesh, std::move(jacobian_assembler), parameters,
34  integration_order, std::move(process_variables),
35  std::move(secondary_variables)),
36  _process_data(std::move(process_data))
37 {
38  DBUG("Create TwoPhaseFlowProcess with PP model.");
39 }
void DBUG(char const *fmt, Args const &... args)
Definition: Logging.h:27
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

References DBUG().

Member Function Documentation

◆ assembleConcreteProcess()

void ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess::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 66 of file TwoPhaseFlowWithPPProcess.cpp.

70 {
71  DBUG("Assemble TwoPhaseFlowWithPPProcess.");
72 
73  std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
74  dof_table = {std::ref(*_local_to_global_index_map)};
75  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
76 
77  // Call global assembler for each local assembly item.
80  pv.getActiveElementIDs(), dof_table, t, dt, x, xdot, process_id, M, K,
81  b);
82 }
std::vector< std::size_t > const & getActiveElementIDs() const
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
std::vector< std::unique_ptr< TwoPhaseFlowWithPPLocalAssemblerInterface > > _local_assemblers
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::VectorMatrixAssembler::assemble(), DBUG(), NumLib::SerialExecutor::executeSelectedMemberDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), and ProcessLib::Process::getProcessVariables().

◆ assembleWithJacobianConcreteProcess()

void ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess::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 84 of file TwoPhaseFlowWithPPProcess.cpp.

89 {
90  DBUG("AssembleWithJacobian TwoPhaseFlowWithPPProcess.");
91 
92  std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
93  dof_table = {std::ref(*_local_to_global_index_map)};
94  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
95 
96  // Call global assembler for each local assembly item.
99  _local_assemblers, pv.getActiveElementIDs(), dof_table, t, dt, x, xdot,
100  dxdot_dx, dx_dx, process_id, M, K, b, Jac);
101 }
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().

◆ initializeConcreteProcess()

void ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess::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 41 of file TwoPhaseFlowWithPPProcess.cpp.

45 {
46  const int process_id = 0;
47  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
48  ProcessLib::createLocalAssemblers<TwoPhaseFlowWithPPLocalAssembler>(
49  mesh.getDimension(), mesh.getElements(), dof_table,
51  mesh.isAxiallySymmetric(), integration_order, _process_data);
52 
54  "saturation",
58 
60  "pressure_wet",
64 }
unsigned getShapeFunctionOrder() const
NumLib::Extrapolator & getExtrapolator() const
Definition: Process.h:185
SecondaryVariableCollection _secondary_variables
Definition: Process.h:331
void addSecondaryVariable(std::string const &internal_name, SecondaryVariableFunctions &&fcts)
virtual std::vector< double > const & getIntPtSaturation(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtWetPressure(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
SecondaryVariableFunctions makeExtrapolator(const unsigned num_components, NumLib::Extrapolator &extrapolator, LocalAssemblerCollection const &local_assemblers, typename NumLib::ExtrapolatableLocalAssemblerCollection< LocalAssemblerCollection >::IntegrationPointValuesMethod integration_point_values_method)

References _local_assemblers, _process_data, ProcessLib::Process::_secondary_variables, ProcessLib::SecondaryVariableCollection::addSecondaryVariable(), MeshLib::Mesh::getDimension(), MeshLib::Mesh::getElements(), ProcessLib::Process::getExtrapolator(), ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPLocalAssemblerInterface::getIntPtSaturation(), ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPLocalAssemblerInterface::getIntPtWetPressure(), ProcessLib::Process::getProcessVariables(), ProcessLib::ProcessVariable::getShapeFunctionOrder(), MeshLib::Mesh::isAxiallySymmetric(), and ProcessLib::makeExtrapolator().

◆ isLinear()

bool ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess::isLinear ( ) const
inlineoverride

Definition at line 57 of file TwoPhaseFlowWithPPProcess.h.

57 { return false; }

Member Data Documentation

◆ _local_assemblers

std::vector<std::unique_ptr<TwoPhaseFlowWithPPLocalAssemblerInterface> > ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess::_local_assemblers
private

◆ _process_data

TwoPhaseFlowWithPPProcessData ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess::_process_data
private

Definition at line 76 of file TwoPhaseFlowWithPPProcess.h.

Referenced by initializeConcreteProcess().


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