OGS
ProcessLib::Process Class Referenceabstract

Detailed Description

Definition at line 39 of file Process.h.

#include <Process.h>

Inheritance diagram for ProcessLib::Process:
[legend]
Collaboration diagram for ProcessLib::Process:
[legend]

Public Member Functions

 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.
 
void postTimestep (std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, const double t, const double delta_t, int const process_id)
 Postprocessing after a complete timestep.
 
void postNonLinearSolver (std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, 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_prev, int const process_id)
 compute secondary variables for the coupled equations or for output.
 
NumLib::IterationResult postIteration (GlobalVector const &x) final
 
void initialize (std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media)
 
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 updateDeactivatedSubdomains (double const time, const int process_id)
 
virtual bool isMonolithicSchemeUsed () 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 &x_prev, 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 &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac) final
 
void preOutput (const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id)
 
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< MeshLib::IntegrationPointWriter > > const & getIntegrationPointWriters () 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)
 
- Public Member Functions inherited from ProcessLib::SubmeshAssemblySupport
virtual std::vector< std::string > initializeAssemblyOnSubmeshes (std::vector< std::reference_wrapper< MeshLib::Mesh > > const &meshes)
 
virtual ~SubmeshAssemblySupport ()=default
 

Public Attributes

std::string const name
 

Static Public Attributes

static PROCESSLIB_EXPORT const std::string constant_one_parameter_name = "constant_one"
 

Protected Member Functions

NumLib::ExtrapolatorgetExtrapolator () const
 
NumLib::LocalToGlobalIndexMap const & getSingleComponentDOFTable () const
 
void initializeProcessBoundaryConditionsAndSourceTerms (const NumLib::LocalToGlobalIndexMap &dof_table, const int process_id, std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media)
 
virtual void constructDofTable ()
 
void constructMonolithicProcessDofTable ()
 
void constructDofTableOfSpecifiedProcessStaggeredScheme (const int specified_process_id)
 
virtual std::tuple< NumLib::LocalToGlobalIndexMap *, bool > getDOFTableForExtrapolatorData () const
 
std::vector< GlobalIndexTypegetIndicesOfResiduumWithoutInitialCompensation () const override
 

Protected Attributes

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
 
std::unique_ptr< ProcessLib::AbstractJacobianAssembler_jacobian_assembler
 
VectorMatrixAssembler _global_assembler
 
const bool _use_monolithic_scheme
 
unsigned const _integration_order
 
std::vector< std::unique_ptr< MeshLib::IntegrationPointWriter > > _integration_point_writer
 
GlobalSparsityPattern _sparsity_pattern
 
std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > _process_variables
 
std::vector< BoundaryConditionCollection_boundary_conditions
 

Private Member Functions

virtual void initializeConcreteProcess (NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order)=0
 Process specific initialization called by initialize().
 
virtual void initializeBoundaryConditions (std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media)
 
virtual void setInitialConditionsConcreteProcess (std::vector< GlobalVector * > &, double const, int const)
 
virtual void preAssembleConcreteProcess (const double, double const, GlobalVector const &)
 
virtual void assembleConcreteProcess (const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b)=0
 
virtual void assembleWithJacobianConcreteProcess (const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac)=0
 
virtual void preTimestepConcreteProcess (std::vector< GlobalVector * > const &, const double, const double, const int)
 
virtual void postTimestepConcreteProcess (std::vector< GlobalVector * > const &, std::vector< GlobalVector * > const &, const double, const double, int const)
 
virtual void postNonLinearSolverConcreteProcess (std::vector< GlobalVector * > const &, std::vector< GlobalVector * > const &, const double, double const, int const)
 
virtual void preIterationConcreteProcess (const unsigned, GlobalVector const &)
 
virtual void computeSecondaryVariableConcrete (double const, double const, std::vector< GlobalVector * > const &, GlobalVector const &, int const)
 
virtual NumLib::IterationResult postIterationConcreteProcess (GlobalVector const &)
 
virtual void preOutputConcreteProcess (const double, double const, std::vector< GlobalVector * > const &, std::vector< GlobalVector * > const &, int const)
 
void initializeExtrapolator ()
 
void computeSparsityPattern ()
 

Private Attributes

std::vector< SourceTermCollection_source_term_collections
 
ExtrapolatorData _extrapolator_data
 

Constructor & Destructor Documentation

◆ Process()

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 
)

Definition at line 23 of file Process.cpp.

33 : name(std::move(name_)),
34 _mesh(mesh),
35 _secondary_variables(std::move(secondary_variables)),
36 _jacobian_assembler(std::move(jacobian_assembler)),
38 _use_monolithic_scheme(use_monolithic_scheme),
39 _integration_order(integration_order),
40 _process_variables(std::move(process_variables)),
42 [&](const std::size_t number_of_process_variables)
43 -> std::vector<BoundaryConditionCollection>
44 {
45 std::vector<BoundaryConditionCollection> pcs_BCs;
46 pcs_BCs.reserve(number_of_process_variables);
47 for (std::size_t i = 0; i < number_of_process_variables; i++)
48 {
49 pcs_BCs.emplace_back(BoundaryConditionCollection(parameters));
50 }
51 return pcs_BCs;
52 }(_process_variables.size())),
54 [&](const std::size_t number_of_processes)
55 -> std::vector<SourceTermCollection>
56 {
57 std::vector<SourceTermCollection> pcs_sts;
58 pcs_sts.reserve(number_of_processes);
59 for (std::size_t i = 0; i < number_of_processes; i++)
60 {
61 pcs_sts.emplace_back(SourceTermCollection(parameters));
62 }
63 return pcs_sts;
64 }(_process_variables.size()))
65{
66}
std::string const name
Definition: Process.h:346
std::vector< BoundaryConditionCollection > _boundary_conditions
Definition: Process.h:387
MeshLib::Mesh & _mesh
Definition: Process.h:349
std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > _process_variables
Definition: Process.h:382
SecondaryVariableCollection _secondary_variables
Definition: Process.h:354
VectorMatrixAssembler _global_assembler
Definition: Process.h:359
unsigned const _integration_order
Definition: Process.h:366
std::vector< SourceTermCollection > _source_term_collections
Definition: Process.h:393
std::unique_ptr< ProcessLib::AbstractJacobianAssembler > _jacobian_assembler
Definition: Process.h:358
const bool _use_monolithic_scheme
Definition: Process.h:361

Member Function Documentation

◆ assemble()

void ProcessLib::Process::assemble ( const double  t,
double const  dt,
std::vector< GlobalVector * > const &  x,
std::vector< GlobalVector * > const &  x_prev,
int const  process_id,
GlobalMatrix M,
GlobalMatrix K,
GlobalVector b 
)
final

Definition at line 213 of file Process.cpp.

218{
219 assert(x.size() == x_prev.size());
220 for (std::size_t i = 0; i < x.size(); i++)
221 {
224 }
225
226 assembleConcreteProcess(t, dt, x, x_prev, process_id, M, K, b);
227
228 // the last argument is for the jacobian, nullptr is for a unused jacobian
229 _boundary_conditions[process_id].applyNaturalBC(t, x, process_id, K, b,
230 nullptr);
231
232 // the last argument is for the jacobian, nullptr is for a unused jacobian
233 _source_term_collections[process_id].integrate(t, *x[process_id], b,
234 nullptr);
235}
virtual void assembleConcreteProcess(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b)=0
void setLocalAccessibleVector(PETScVector const &x)
Definition: LinAlg.cpp:27

References _boundary_conditions, _source_term_collections, assembleConcreteProcess(), and MathLib::LinAlg::setLocalAccessibleVector().

◆ assembleConcreteProcess()

virtual void ProcessLib::Process::assembleConcreteProcess ( const double  t,
double const  dt,
std::vector< GlobalVector * > const &  x,
std::vector< GlobalVector * > const &  x_prev,
int const  process_id,
GlobalMatrix M,
GlobalMatrix K,
GlobalVector b 
)
privatepure virtual

Implemented in ProcessLib::ComponentTransport::ComponentTransportProcess, ProcessLib::HeatTransportBHE::HeatTransportBHEProcess, ProcessLib::HT::HTProcess, ProcessLib::LIE::HydroMechanics::HydroMechanicsProcess< GlobalDim >, ProcessLib::LIE::SmallDeformation::SmallDeformationProcess< DisplacementDim >, ProcessLib::LiquidFlow::LiquidFlowProcess, ProcessLib::PhaseField::PhaseFieldProcess< DisplacementDim >, ProcessLib::RichardsComponentTransport::RichardsComponentTransportProcess, ProcessLib::RichardsFlow::RichardsFlowProcess, ProcessLib::RichardsMechanics::RichardsMechanicsProcess< DisplacementDim >, ProcessLib::SmallDeformationNonlocal::SmallDeformationNonlocalProcess< DisplacementDim >, ProcessLib::SteadyStateDiffusion::SteadyStateDiffusion, ProcessLib::StokesFlow::StokesFlowProcess< GlobalDim >, ProcessLib::TES::TESProcess, ProcessLib::TH2M::TH2MProcess< DisplacementDim >, ProcessLib::ThermalTwoPhaseFlowWithPP::ThermalTwoPhaseFlowWithPPProcess, ProcessLib::ThermoHydroMechanics::ThermoHydroMechanicsProcess< DisplacementDim >, ProcessLib::ThermoMechanicalPhaseField::ThermoMechanicalPhaseFieldProcess< DisplacementDim >, ProcessLib::ThermoMechanics::ThermoMechanicsProcess< DisplacementDim >, ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >, ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess, ProcessLib::TwoPhaseFlowWithPrho::TwoPhaseFlowWithPrhoProcess, ProcessLib::HeatConduction::HeatConductionProcess, ProcessLib::HydroMechanics::HydroMechanicsProcess< DisplacementDim >, ProcessLib::LargeDeformation::LargeDeformationProcess< DisplacementDim >, and ProcessLib::SmallDeformation::SmallDeformationProcess< DisplacementDim >.

Referenced by assemble().

◆ assembleWithJacobian()

void ProcessLib::Process::assembleWithJacobian ( const double  t,
double const  dt,
std::vector< GlobalVector * > const &  x,
std::vector< GlobalVector * > const &  x_prev,
int const  process_id,
GlobalMatrix M,
GlobalMatrix K,
GlobalVector b,
GlobalMatrix Jac 
)
final

Definition at line 237 of file Process.cpp.

243{
244 assert(x.size() == x_prev.size());
245 for (std::size_t i = 0; i < x.size(); i++)
246 {
249 }
250
251 assembleWithJacobianConcreteProcess(t, dt, x, x_prev, process_id, M, K, b,
252 Jac);
253
254 // TODO: apply BCs to Jacobian.
255 _boundary_conditions[process_id].applyNaturalBC(t, x, process_id, K, b,
256 &Jac);
257
258 // the last argument is for the jacobian, nullptr is for a unused jacobian
259 _source_term_collections[process_id].integrate(t, *x[process_id], b, &Jac);
260}
virtual void assembleWithJacobianConcreteProcess(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac)=0

References _boundary_conditions, _source_term_collections, assembleWithJacobianConcreteProcess(), and MathLib::LinAlg::setLocalAccessibleVector().

◆ assembleWithJacobianConcreteProcess()

virtual void ProcessLib::Process::assembleWithJacobianConcreteProcess ( const double  t,
double const  dt,
std::vector< GlobalVector * > const &  x,
std::vector< GlobalVector * > const &  x_prev,
int const  process_id,
GlobalMatrix M,
GlobalMatrix K,
GlobalVector b,
GlobalMatrix Jac 
)
privatepure virtual

Implemented in ProcessLib::ComponentTransport::ComponentTransportProcess, ProcessLib::HeatTransportBHE::HeatTransportBHEProcess, ProcessLib::HT::HTProcess, ProcessLib::LargeDeformation::LargeDeformationProcess< DisplacementDim >, ProcessLib::LIE::HydroMechanics::HydroMechanicsProcess< GlobalDim >, ProcessLib::LIE::SmallDeformation::SmallDeformationProcess< DisplacementDim >, ProcessLib::LiquidFlow::LiquidFlowProcess, ProcessLib::PhaseField::PhaseFieldProcess< DisplacementDim >, ProcessLib::RichardsComponentTransport::RichardsComponentTransportProcess, ProcessLib::RichardsFlow::RichardsFlowProcess, ProcessLib::RichardsMechanics::RichardsMechanicsProcess< DisplacementDim >, ProcessLib::SmallDeformation::SmallDeformationProcess< DisplacementDim >, ProcessLib::SmallDeformationNonlocal::SmallDeformationNonlocalProcess< DisplacementDim >, ProcessLib::SteadyStateDiffusion::SteadyStateDiffusion, ProcessLib::StokesFlow::StokesFlowProcess< GlobalDim >, ProcessLib::TES::TESProcess, ProcessLib::TH2M::TH2MProcess< DisplacementDim >, ProcessLib::ThermalTwoPhaseFlowWithPP::ThermalTwoPhaseFlowWithPPProcess, ProcessLib::ThermoHydroMechanics::ThermoHydroMechanicsProcess< DisplacementDim >, ProcessLib::ThermoMechanicalPhaseField::ThermoMechanicalPhaseFieldProcess< DisplacementDim >, ProcessLib::ThermoMechanics::ThermoMechanicsProcess< DisplacementDim >, ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >, ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess, ProcessLib::TwoPhaseFlowWithPrho::TwoPhaseFlowWithPrhoProcess, ProcessLib::HeatConduction::HeatConductionProcess, and ProcessLib::HydroMechanics::HydroMechanicsProcess< DisplacementDim >.

Referenced by assembleWithJacobian().

◆ computeSecondaryVariable()

void ProcessLib::Process::computeSecondaryVariable ( double const  t,
double const  dt,
std::vector< GlobalVector * > const &  x,
GlobalVector const &  x_prev,
int const  process_id 
)

compute secondary variables for the coupled equations or for output.

Definition at line 434 of file Process.cpp.

439{
440 for (auto const* solution : x)
443
444 computeSecondaryVariableConcrete(t, dt, x, x_prev, process_id);
445}
virtual void computeSecondaryVariableConcrete(double const, double const, std::vector< GlobalVector * > const &, GlobalVector const &, int const)
Definition: Process.h:277

References computeSecondaryVariableConcrete(), and MathLib::LinAlg::setLocalAccessibleVector().

◆ computeSecondaryVariableConcrete()

◆ computeSparsityPattern()

void ProcessLib::Process::computeSparsityPattern ( )
private

Computes and stores global matrix' sparsity pattern from given DOF-table.

Definition at line 382 of file Process.cpp.

383{
386}
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition: Process.h:352
GlobalSparsityPattern _sparsity_pattern
Definition: Process.h:374
GlobalSparsityPattern computeSparsityPattern(LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh)
Computes a sparsity pattern for the given inputs.

References _local_to_global_index_map, _mesh, _sparsity_pattern, and NumLib::computeSparsityPattern().

Referenced by initialize().

◆ constructDofTable()

void ProcessLib::Process::constructDofTable ( )
protectedvirtual

This function is for general cases, in which all equations of the coupled processes have the same number of unknowns. For the general cases with the staggered scheme, all equations of the coupled processes share one DOF table hold by _local_to_global_index_map. Other cases can be considered by overloading this member function in the derived class.

Reimplemented in ProcessLib::HeatTransportBHE::HeatTransportBHEProcess, ProcessLib::HydroMechanics::HydroMechanicsProcess< DisplacementDim >, ProcessLib::LIE::HydroMechanics::HydroMechanicsProcess< GlobalDim >, ProcessLib::LIE::SmallDeformation::SmallDeformationProcess< DisplacementDim >, ProcessLib::PhaseField::PhaseFieldProcess< DisplacementDim >, ProcessLib::RichardsMechanics::RichardsMechanicsProcess< DisplacementDim >, ProcessLib::StokesFlow::StokesFlowProcess< GlobalDim >, ProcessLib::TH2M::TH2MProcess< DisplacementDim >, ProcessLib::ThermoHydroMechanics::ThermoHydroMechanicsProcess< DisplacementDim >, ProcessLib::ThermoMechanicalPhaseField::ThermoMechanicalPhaseFieldProcess< DisplacementDim >, ProcessLib::ThermoMechanics::ThermoMechanicsProcess< DisplacementDim >, and ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >.

Definition at line 262 of file Process.cpp.

263{
265 {
267
268 return;
269 }
270
271 // For staggered scheme:
272 const int specified_process_id = 0;
274}
void constructDofTableOfSpecifiedProcessStaggeredScheme(const int specified_process_id)
Definition: Process.cpp:309
void constructMonolithicProcessDofTable()
Definition: Process.cpp:276

References _use_monolithic_scheme, constructDofTableOfSpecifiedProcessStaggeredScheme(), and constructMonolithicProcessDofTable().

Referenced by initialize().

◆ constructDofTableOfSpecifiedProcessStaggeredScheme()

void ProcessLib::Process::constructDofTableOfSpecifiedProcessStaggeredScheme ( const int  specified_process_id)
protected

Construct the DOF table for a specified process in the staggered scheme, which is stored in the member of this class, _local_to_global_index_map.

Definition at line 309 of file Process.cpp.

311{
312 // Create single component dof in every of the mesh nodes.
314 std::make_unique<MeshLib::MeshSubset>(_mesh, _mesh.getNodes());
315
316 // Vector of mesh subsets.
317 std::vector<MeshLib::MeshSubset> all_mesh_subsets;
318
319 // Vector of the number of variable components
320 std::vector<int> vec_var_n_components;
321 // Collect the mesh subsets in a vector for each variables' components.
322 std::generate_n(std::back_inserter(all_mesh_subsets),
323 _process_variables[specified_process_id][0]
324 .get()
325 .getNumberOfGlobalComponents(),
326 [&]() { return *_mesh_subset_all_nodes; });
327
328 // Create a vector of the number of variable components.
329 vec_var_n_components.push_back(_process_variables[specified_process_id][0]
330 .get()
331 .getNumberOfGlobalComponents());
333 std::make_unique<NumLib::LocalToGlobalIndexMap>(
334 std::move(all_mesh_subsets), vec_var_n_components,
336
338}
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition: Mesh.h:102
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_all_nodes
Definition: Process.h:350
@ BY_LOCATION
Ordering data by spatial location.

References _local_to_global_index_map, _mesh, _mesh_subset_all_nodes, _process_variables, NumLib::BY_LOCATION, and MeshLib::Mesh::getNodes().

Referenced by constructDofTable().

◆ constructMonolithicProcessDofTable()

void ProcessLib::Process::constructMonolithicProcessDofTable ( )
protected

Construct the DOF table for the monolithic scheme, which is stored in the member of this class, _local_to_global_index_map.

Definition at line 276 of file Process.cpp.

277{
278 // Create single component dof in every of the mesh nodes.
280 std::make_unique<MeshLib::MeshSubset>(_mesh, _mesh.getNodes());
281
282 // Vector of mesh subsets.
283 std::vector<MeshLib::MeshSubset> all_mesh_subsets;
284
285 // Collect the mesh subsets in a vector for the components of each
286 // variables.
287 for (ProcessVariable const& pv : _process_variables[0])
288 {
289 std::generate_n(std::back_inserter(all_mesh_subsets),
290 pv.getNumberOfGlobalComponents(),
291 [&]() { return *_mesh_subset_all_nodes; });
292 }
293
294 // Create a vector of the number of variable components
295 std::vector<int> vec_var_n_components;
296 transform(cbegin(_process_variables[0]), cend(_process_variables[0]),
297 back_inserter(vec_var_n_components),
298 [](ProcessVariable const& pv)
299 { return pv.getNumberOfGlobalComponents(); });
300
302 std::make_unique<NumLib::LocalToGlobalIndexMap>(
303 std::move(all_mesh_subsets), vec_var_n_components,
305
307}

References _local_to_global_index_map, _mesh, _mesh_subset_all_nodes, _process_variables, NumLib::BY_LOCATION, MeshLib::Mesh::getNodes(), and ProcessLib::ProcessVariable::getNumberOfGlobalComponents().

Referenced by constructDofTable().

◆ extrapolateIntegrationPointValuesToNodes()

virtual void ProcessLib::Process::extrapolateIntegrationPointValuesToNodes ( const double  ,
std::vector< GlobalVector * > const &  ,
std::vector< GlobalVector * > &   
)
inlinevirtual

Definition at line 106 of file Process.h.

110 {
111 }

◆ getDOFTable()

◆ getDOFTableForExtrapolatorData()

std::tuple< NumLib::LocalToGlobalIndexMap *, bool > ProcessLib::Process::getDOFTableForExtrapolatorData ( ) const
protectedvirtual

Get the address of a LocalToGlobalIndexMap, and the status of its memory. If the LocalToGlobalIndexMap is created as new in this function, the function also returns a true boolean value to let Extrapolator manage the memory by the address returned by this function.

Returns
Address of a LocalToGlobalIndexMap and its memory status.

Reimplemented in ProcessLib::HT::HTProcess, ProcessLib::HydroMechanics::HydroMechanicsProcess< DisplacementDim >, ProcessLib::RichardsMechanics::RichardsMechanicsProcess< DisplacementDim >, ProcessLib::TH2M::TH2MProcess< DisplacementDim >, ProcessLib::ThermoHydroMechanics::ThermoHydroMechanicsProcess< DisplacementDim >, and ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >.

Definition at line 341 of file Process.cpp.

342{
343 if (_local_to_global_index_map->getNumberOfGlobalComponents() == 1)
344 {
345 // For single-variable-single-component processes reuse the existing DOF
346 // table.
347 const bool manage_storage = false;
348 return std::make_tuple(_local_to_global_index_map.get(),
349 manage_storage);
350 }
351
352 // Otherwise construct a new DOF table.
353 std::vector<MeshLib::MeshSubset> all_mesh_subsets_single_component;
354 all_mesh_subsets_single_component.emplace_back(*_mesh_subset_all_nodes);
355
356 const bool manage_storage = true;
357
358 return std::make_tuple(new NumLib::LocalToGlobalIndexMap(
359 std::move(all_mesh_subsets_single_component),
360 // by location order is needed for output
362 manage_storage);
363}

References _local_to_global_index_map, _mesh_subset_all_nodes, and NumLib::BY_LOCATION.

Referenced by initializeExtrapolator().

◆ getExtrapolator()

NumLib::Extrapolator & ProcessLib::Process::getExtrapolator ( ) const
inlineprotected

◆ getFlux()

virtual Eigen::Vector3d ProcessLib::Process::getFlux ( std::size_t  ,
MathLib::Point3d const &  ,
double const  ,
std::vector< GlobalVector * > const &   
) const
inlinevirtual

◆ getIndicesOfResiduumWithoutInitialCompensation()

std::vector< GlobalIndexType > ProcessLib::Process::getIndicesOfResiduumWithoutInitialCompensation ( ) const
overrideprotected
Returns
The global indices for the entries of the global residuum vector that do not need initial non-equilibrium compensation.

Definition at line 471 of file Process.cpp.

472{
473 std::vector<GlobalIndexType> indices;
474
475 for (std::size_t process_id = 0; process_id < _process_variables.size();
476 process_id++)
477 {
478 auto const& dof_table_of_process = getDOFTable(process_id);
479
480 auto const& per_process_variables = _process_variables[process_id];
481 for (std::size_t variable_id = 0;
482 variable_id < per_process_variables.size();
483 variable_id++)
484 {
485 auto const& pv = per_process_variables[variable_id];
486 DBUG("Set the initial condition of variable {:s} of process {:d}.",
487 pv.get().getName().data(), process_id);
488
489 if ((pv.get().compensateNonEquilibriumInitialResiduum()))
490 {
491 continue;
492 }
493
494 auto const num_comp = pv.get().getNumberOfGlobalComponents();
495
496 for (int component_id = 0; component_id < num_comp; ++component_id)
497 {
498 auto const& mesh_subset = dof_table_of_process.getMeshSubset(
499 variable_id, component_id);
500 auto const mesh_id = mesh_subset.getMeshID();
501 for (auto const* node : mesh_subset.getNodes())
502 {
503 MeshLib::Location const l(
504 mesh_id, MeshLib::MeshItemType::Node, node->getID());
505
506 auto global_index =
507 std::abs(dof_table_of_process.getGlobalIndex(
508 l, variable_id, component_id));
509
510 indices.push_back(global_index);
511 }
512 }
513 }
514 }
515
516 return indices;
517}
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Definition: Logging.h:30
virtual NumLib::LocalToGlobalIndexMap const & getDOFTable(const int) const
Definition: Process.h:147

References _process_variables, DBUG(), getDOFTable(), and MeshLib::Node.

◆ getIntegrationPointWriters()

std::vector< std::unique_ptr< MeshLib::IntegrationPointWriter > > const & ProcessLib::Process::getIntegrationPointWriters ( ) const
inline

Definition at line 166 of file Process.h.

167 {
169 }
std::vector< std::unique_ptr< MeshLib::IntegrationPointWriter > > _integration_point_writer
Definition: Process.h:372

References _integration_point_writer.

Referenced by anonymous_namespace{ProcessOutputData.cpp}::getIntegrationPointWriters().

◆ getKnownSolutions()

std::vector< NumLib::IndexValueVector< GlobalIndexType > > const * ProcessLib::Process::getKnownSolutions ( double const  t,
GlobalVector const &  x,
int const  process_id 
) const
inlinefinal

Definition at line 141 of file Process.h.

143 {
144 return _boundary_conditions[process_id].getKnownSolutions(t, x);
145 }

References _boundary_conditions.

◆ getMatrixSpecifications()

MathLib::MatrixSpecifications ProcessLib::Process::getMatrixSpecifications ( const int  process_id) const
override

Definition at line 189 of file Process.cpp.

191{
192 auto const& l = *_local_to_global_index_map;
193 return {l.dofSizeWithoutGhosts(), l.dofSizeWithoutGhosts(),
194 &l.getGhostIndices(), &_sparsity_pattern};
195}

References _local_to_global_index_map, and _sparsity_pattern.

Referenced by ProcessLib::ComponentTransport::ComponentTransportProcess::preOutputConcreteProcess(), and ProcessLib::ComponentTransport::ComponentTransportProcess::solveReactionEquation().

◆ getMesh()

◆ getProcessVariables()

std::vector< std::reference_wrapper< ProcessVariable > > const & ProcessLib::Process::getProcessVariables ( const int  process_id) const
inline

Definition at line 155 of file Process.h.

156 {
157 return _process_variables[process_id];
158 }

References _process_variables.

Referenced by ProcessLib::LIE::HydroMechanics::HydroMechanicsProcess< GlobalDim >::HydroMechanicsProcess(), ProcessLib::ComponentTransport::ComponentTransportProcess::assembleConcreteProcess(), ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::assembleConcreteProcess(), ProcessLib::HT::HTProcess::assembleConcreteProcess(), ProcessLib::LiquidFlow::LiquidFlowProcess::assembleConcreteProcess(), ProcessLib::RichardsComponentTransport::RichardsComponentTransportProcess::assembleConcreteProcess(), ProcessLib::RichardsFlow::RichardsFlowProcess::assembleConcreteProcess(), ProcessLib::SteadyStateDiffusion::SteadyStateDiffusion::assembleConcreteProcess(), ProcessLib::TES::TESProcess::assembleConcreteProcess(), ProcessLib::ThermalTwoPhaseFlowWithPP::ThermalTwoPhaseFlowWithPPProcess::assembleConcreteProcess(), ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::assembleConcreteProcess(), ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess::assembleConcreteProcess(), ProcessLib::TwoPhaseFlowWithPrho::TwoPhaseFlowWithPrhoProcess::assembleConcreteProcess(), ProcessLib::HeatConduction::HeatConductionProcess::assembleConcreteProcess(), ProcessLib::AssemblyMixin< Process >::assembleGeneric(), ProcessLib::ComponentTransport::ComponentTransportProcess::assembleWithJacobianConcreteProcess(), ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::assembleWithJacobianConcreteProcess(), ProcessLib::HT::HTProcess::assembleWithJacobianConcreteProcess(), ProcessLib::LiquidFlow::LiquidFlowProcess::assembleWithJacobianConcreteProcess(), ProcessLib::RichardsComponentTransport::RichardsComponentTransportProcess::assembleWithJacobianConcreteProcess(), ProcessLib::RichardsFlow::RichardsFlowProcess::assembleWithJacobianConcreteProcess(), ProcessLib::SteadyStateDiffusion::SteadyStateDiffusion::assembleWithJacobianConcreteProcess(), ProcessLib::TES::TESProcess::assembleWithJacobianConcreteProcess(), ProcessLib::ThermalTwoPhaseFlowWithPP::ThermalTwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess(), ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::assembleWithJacobianConcreteProcess(), ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess(), ProcessLib::TwoPhaseFlowWithPrho::TwoPhaseFlowWithPrhoProcess::assembleWithJacobianConcreteProcess(), ProcessLib::HeatConduction::HeatConductionProcess::assembleWithJacobianConcreteProcess(), ProcessLib::HeatConduction::HeatConductionProcess::computeSecondaryVariableConcrete(), ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::computeSecondaryVariableConcrete(), ProcessLib::LiquidFlow::LiquidFlowProcess::computeSecondaryVariableConcrete(), ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::computeSecondaryVariableConcrete(), ProcessLib::ComponentTransport::ComponentTransportProcess::computeSecondaryVariableConcrete(), ProcessLib::createConstraintDirichletBoundaryCondition(), anonymous_namespace{ProcessOutputData.cpp}::getProcessVariablesOfAllProcesses(), ProcessLib::RichardsComponentTransport::RichardsComponentTransportProcess::initializeConcreteProcess(), ProcessLib::SteadyStateDiffusion::SteadyStateDiffusion::postTimestepConcreteProcess(), ProcessLib::HT::HTProcess::postTimestepConcreteProcess(), ProcessLib::ComponentTransport::ComponentTransportProcess::postTimestepConcreteProcess(), ProcessLib::LiquidFlow::LiquidFlowProcess::postTimestepConcreteProcess(), ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::postTimestepConcreteProcess(), ProcessLib::ThermalTwoPhaseFlowWithPP::ThermalTwoPhaseFlowWithPPProcess::preTimestepConcreteProcess(), ProcessLib::ComponentTransport::ComponentTransportProcess::solveReactionEquation(), ProcessLib::AssemblyMixin< Process >::updateActiveElements(), and updateDeactivatedSubdomains().

◆ getSecondaryVariables()

SecondaryVariableCollection const & ProcessLib::Process::getSecondaryVariables ( ) const
inline

Definition at line 160 of file Process.h.

161 {
163 }

References _secondary_variables.

Referenced by ProcessLib::createProcessOutputData().

◆ getSingleComponentDOFTable()

NumLib::LocalToGlobalIndexMap const & ProcessLib::Process::getSingleComponentDOFTable ( ) const
inlineprotected

◆ initialize()

void ProcessLib::Process::initialize ( std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &  media)

Definition at line 97 of file Process.cpp.

99{
100 DBUG("Initialize process.");
101
102 DBUG("Construct dof mappings.");
104
105 DBUG("Compute sparsity pattern");
107
108 DBUG("Initialize the extrapolator");
110
113
114 DBUG("Initialize boundary conditions.");
116}
virtual void initializeBoundaryConditions(std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media)
Definition: Process.cpp:83
virtual void initializeConcreteProcess(NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order)=0
Process specific initialization called by initialize().
void initializeExtrapolator()
Definition: Process.cpp:365
virtual void constructDofTable()
Definition: Process.cpp:262
void computeSparsityPattern()
Definition: Process.cpp:382

References _integration_order, _local_to_global_index_map, _mesh, computeSparsityPattern(), constructDofTable(), DBUG(), initializeBoundaryConditions(), initializeConcreteProcess(), and initializeExtrapolator().

◆ initializeBoundaryConditions()

void ProcessLib::Process::initializeBoundaryConditions ( std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &  media)
privatevirtual

Member function to initialize the boundary conditions for all coupled processes. It is called by initialize().

Reimplemented in ProcessLib::HydroMechanics::HydroMechanicsProcess< DisplacementDim >, ProcessLib::PhaseField::PhaseFieldProcess< DisplacementDim >, ProcessLib::RichardsMechanics::RichardsMechanicsProcess< DisplacementDim >, ProcessLib::StokesFlow::StokesFlowProcess< GlobalDim >, ProcessLib::TH2M::TH2MProcess< DisplacementDim >, ProcessLib::ThermoHydroMechanics::ThermoHydroMechanicsProcess< DisplacementDim >, ProcessLib::ThermoMechanicalPhaseField::ThermoMechanicalPhaseFieldProcess< DisplacementDim >, ProcessLib::ThermoMechanics::ThermoMechanicsProcess< DisplacementDim >, and ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >.

Definition at line 83 of file Process.cpp.

85{
86 // The number of processes is identical to the size of _process_variables,
87 // the vector contains variables for different processes. See the
88 // documentation of _process_variables.
89 const std::size_t number_of_processes = _process_variables.size();
90 for (std::size_t pcs_id = 0; pcs_id < number_of_processes; pcs_id++)
91 {
93 *_local_to_global_index_map, pcs_id, media);
94 }
95}
void initializeProcessBoundaryConditionsAndSourceTerms(const NumLib::LocalToGlobalIndexMap &dof_table, const int process_id, std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media)
Definition: Process.cpp:68

References _local_to_global_index_map, _process_variables, and initializeProcessBoundaryConditionsAndSourceTerms().

Referenced by initialize().

◆ initializeConcreteProcess()

virtual void ProcessLib::Process::initializeConcreteProcess ( NumLib::LocalToGlobalIndexMap const &  dof_table,
MeshLib::Mesh const &  mesh,
unsigned const  integration_order 
)
privatepure virtual

Process specific initialization called by initialize().

Implemented in ProcessLib::ComponentTransport::ComponentTransportProcess, ProcessLib::HeatConduction::HeatConductionProcess, ProcessLib::HeatTransportBHE::HeatTransportBHEProcess, ProcessLib::HT::HTProcess, ProcessLib::HydroMechanics::HydroMechanicsProcess< DisplacementDim >, ProcessLib::LargeDeformation::LargeDeformationProcess< DisplacementDim >, ProcessLib::LIE::HydroMechanics::HydroMechanicsProcess< GlobalDim >, ProcessLib::LIE::SmallDeformation::SmallDeformationProcess< DisplacementDim >, ProcessLib::LiquidFlow::LiquidFlowProcess, ProcessLib::PhaseField::PhaseFieldProcess< DisplacementDim >, ProcessLib::RichardsComponentTransport::RichardsComponentTransportProcess, ProcessLib::RichardsFlow::RichardsFlowProcess, ProcessLib::RichardsMechanics::RichardsMechanicsProcess< DisplacementDim >, ProcessLib::SmallDeformation::SmallDeformationProcess< DisplacementDim >, ProcessLib::SmallDeformationNonlocal::SmallDeformationNonlocalProcess< DisplacementDim >, ProcessLib::SteadyStateDiffusion::SteadyStateDiffusion, ProcessLib::StokesFlow::StokesFlowProcess< GlobalDim >, ProcessLib::TES::TESProcess, ProcessLib::TH2M::TH2MProcess< DisplacementDim >, ProcessLib::ThermalTwoPhaseFlowWithPP::ThermalTwoPhaseFlowWithPPProcess, ProcessLib::ThermoHydroMechanics::ThermoHydroMechanicsProcess< DisplacementDim >, ProcessLib::ThermoMechanicalPhaseField::ThermoMechanicalPhaseFieldProcess< DisplacementDim >, ProcessLib::ThermoMechanics::ThermoMechanicsProcess< DisplacementDim >, ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >, ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess, and ProcessLib::TwoPhaseFlowWithPrho::TwoPhaseFlowWithPrhoProcess.

Referenced by initialize().

◆ initializeExtrapolator()

void ProcessLib::Process::initializeExtrapolator ( )
private

Definition at line 365 of file Process.cpp.

366{
367 NumLib::LocalToGlobalIndexMap* dof_table_single_component;
368 bool manage_storage;
369
370 std::tie(dof_table_single_component, manage_storage) =
372
373 std::unique_ptr<NumLib::Extrapolator> extrapolator(
375 *dof_table_single_component));
376
377 // TODO: Later on the DOF table can change during the simulation!
378 _extrapolator_data = ExtrapolatorData(
379 std::move(extrapolator), dof_table_single_component, manage_storage);
380}
virtual std::tuple< NumLib::LocalToGlobalIndexMap *, bool > getDOFTableForExtrapolatorData() const
Definition: Process.cpp:341

References _extrapolator_data, and getDOFTableForExtrapolatorData().

Referenced by initialize().

◆ initializeProcessBoundaryConditionsAndSourceTerms()

void ProcessLib::Process::initializeProcessBoundaryConditionsAndSourceTerms ( const NumLib::LocalToGlobalIndexMap dof_table,
const int  process_id,
std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &  media 
)
protected

Initialize the boundary conditions for a single process or coupled processes modelled by the monolithic scheme. It is called by initializeBoundaryConditions().

Definition at line 68 of file Process.cpp.

71{
72 auto const& per_process_variables = _process_variables[process_id];
73 auto& per_process_BCs = _boundary_conditions[process_id];
74
75 per_process_BCs.addBCsForProcessVariables(per_process_variables, dof_table,
76 _integration_order, *this, media);
77
78 auto& per_process_sts = _source_term_collections[process_id];
79 per_process_sts.addSourceTermsForProcessVariables(
80 per_process_variables, dof_table, _integration_order);
81}

References _boundary_conditions, _integration_order, _process_variables, and _source_term_collections.

Referenced by initializeBoundaryConditions().

◆ isMonolithicSchemeUsed()

virtual bool ProcessLib::Process::isMonolithicSchemeUsed ( ) const
inlinevirtual

◆ postIteration()

NumLib::IterationResult ProcessLib::Process::postIteration ( GlobalVector const &  x)
final

Definition at line 453 of file Process.cpp.

454{
457}
virtual NumLib::IterationResult postIterationConcreteProcess(GlobalVector const &)
Definition: Process.h:286

References postIterationConcreteProcess(), and MathLib::LinAlg::setLocalAccessibleVector().

◆ postIterationConcreteProcess()

virtual NumLib::IterationResult ProcessLib::Process::postIterationConcreteProcess ( GlobalVector const &  )
inlineprivatevirtual

◆ postNonLinearSolver()

void ProcessLib::Process::postNonLinearSolver ( std::vector< GlobalVector * > const &  x,
std::vector< GlobalVector * > const &  x_prev,
const double  t,
double const  dt,
int const  process_id 
)

Calculates secondary variables, e.g. stress and strain for deformation analysis, only after nonlinear solver being successfully conducted.

Definition at line 417 of file Process.cpp.

421{
422 for (auto* const solution : x)
423 {
425 }
426 for (auto* const solution : x_prev)
427 {
429 }
430
431 postNonLinearSolverConcreteProcess(x, x_prev, t, dt, process_id);
432}
virtual void postNonLinearSolverConcreteProcess(std::vector< GlobalVector * > const &, std::vector< GlobalVector * > const &, const double, double const, int const)
Definition: Process.h:265

References postNonLinearSolverConcreteProcess(), and MathLib::LinAlg::setLocalAccessibleVector().

◆ postNonLinearSolverConcreteProcess()

virtual void ProcessLib::Process::postNonLinearSolverConcreteProcess ( std::vector< GlobalVector * > const &  ,
std::vector< GlobalVector * > const &  ,
const double  ,
double const  ,
int const   
)
inlineprivatevirtual

◆ postTimestep()

void ProcessLib::Process::postTimestep ( std::vector< GlobalVector * > const &  x,
std::vector< GlobalVector * > const &  x_prev,
const double  t,
const double  delta_t,
int const  process_id 
)

Postprocessing after a complete timestep.

Definition at line 398 of file Process.cpp.

402{
403 for (auto* const solution : x)
404 {
406 }
407 for (auto* const solution : x_prev)
408 {
410 }
411
412 postTimestepConcreteProcess(x, x_prev, t, delta_t, process_id);
413
414 _boundary_conditions[process_id].postTimestep(t, x, process_id);
415}
virtual void postTimestepConcreteProcess(std::vector< GlobalVector * > const &, std::vector< GlobalVector * > const &, const double, const double, int const)
Definition: Process.h:256

References _boundary_conditions, postTimestepConcreteProcess(), and MathLib::LinAlg::setLocalAccessibleVector().

◆ postTimestepConcreteProcess()

virtual void ProcessLib::Process::postTimestepConcreteProcess ( std::vector< GlobalVector * > const &  ,
std::vector< GlobalVector * > const &  ,
const double  ,
const double  ,
int const   
)
inlineprivatevirtual

◆ preAssemble()

void ProcessLib::Process::preAssemble ( const double  t,
double const  dt,
GlobalVector const &  x 
)
final

Definition at line 207 of file Process.cpp.

209{
211}
virtual void preAssembleConcreteProcess(const double, double const, GlobalVector const &)
Definition: Process.h:231

References preAssembleConcreteProcess().

◆ preAssembleConcreteProcess()

virtual void ProcessLib::Process::preAssembleConcreteProcess ( const double  ,
double const  ,
GlobalVector const &   
)
inlineprivatevirtual

Reimplemented in ProcessLib::SmallDeformationNonlocal::SmallDeformationNonlocalProcess< DisplacementDim >.

Definition at line 231 of file Process.h.

234 {
235 }

Referenced by preAssemble().

◆ preIteration()

void ProcessLib::Process::preIteration ( const unsigned  iter,
GlobalVector const &  x 
)
final

Definition at line 447 of file Process.cpp.

448{
451}
virtual void preIterationConcreteProcess(const unsigned, GlobalVector const &)
Definition: Process.h:272

References preIterationConcreteProcess(), and MathLib::LinAlg::setLocalAccessibleVector().

◆ preIterationConcreteProcess()

virtual void ProcessLib::Process::preIterationConcreteProcess ( const unsigned  ,
GlobalVector const &   
)
inlineprivatevirtual

Reimplemented in ProcessLib::TES::TESProcess.

Definition at line 272 of file Process.h.

274 {
275 }

Referenced by preIteration().

◆ preOutput()

void ProcessLib::Process::preOutput ( const double  t,
double const  dt,
std::vector< GlobalVector * > const &  x,
std::vector< GlobalVector * > const &  x_prev,
int const  process_id 
)

Definition at line 459 of file Process.cpp.

463{
464 for (auto const* solution : x)
466
467 preOutputConcreteProcess(t, dt, x, x_prev, process_id);
468}
virtual void preOutputConcreteProcess(const double, double const, std::vector< GlobalVector * > const &, std::vector< GlobalVector * > const &, int const)
Definition: Process.h:292

References preOutputConcreteProcess(), and MathLib::LinAlg::setLocalAccessibleVector().

◆ preOutputConcreteProcess()

virtual void ProcessLib::Process::preOutputConcreteProcess ( const double  ,
double const  ,
std::vector< GlobalVector * > const &  ,
std::vector< GlobalVector * > const &  ,
int const   
)
inlineprivatevirtual

Reimplemented in ProcessLib::ComponentTransport::ComponentTransportProcess.

Definition at line 292 of file Process.h.

296 {
297 }

Referenced by preOutput().

◆ preTimestep()

void ProcessLib::Process::preTimestep ( std::vector< GlobalVector * > const &  x,
const double  t,
const double  delta_t,
const int  process_id 
)

Preprocessing before starting assembly for new timestep.

Definition at line 388 of file Process.cpp.

390{
391 for (auto* const solution : x)
393 preTimestepConcreteProcess(x, t, delta_t, process_id);
394
395 _boundary_conditions[process_id].preTimestep(t, x, process_id);
396}
virtual void preTimestepConcreteProcess(std::vector< GlobalVector * > const &, const double, const double, const int)
Definition: Process.h:248

References _boundary_conditions, preTimestepConcreteProcess(), and MathLib::LinAlg::setLocalAccessibleVector().

◆ preTimestepConcreteProcess()

◆ setInitialConditions()

void ProcessLib::Process::setInitialConditions ( std::vector< GlobalVector * > &  process_solutions,
std::vector< GlobalVector * > const &  process_solutions_prev,
double const  t,
int const  process_id 
)

Definition at line 118 of file Process.cpp.

123{
124 auto& x = *process_solutions[process_id];
125 auto& x_prev = *process_solutions_prev[process_id];
126
127 // getDOFTableOfProcess can be overloaded by the specific process.
128 auto const& dof_table_of_process = getDOFTable(process_id);
129
130 auto const& per_process_variables = _process_variables[process_id];
131 for (std::size_t variable_id = 0;
132 variable_id < per_process_variables.size();
133 variable_id++)
134 {
137
138 auto const& pv = per_process_variables[variable_id];
139 DBUG("Set the initial condition of variable {:s} of process {:d}.",
140 pv.get().getName().data(), process_id);
141
142 auto const& ic = pv.get().getInitialCondition();
143
144 auto const num_comp = pv.get().getNumberOfGlobalComponents();
145
146 for (int component_id = 0; component_id < num_comp; ++component_id)
147 {
148 auto const& mesh_subset =
149 dof_table_of_process.getMeshSubset(variable_id, component_id);
150 auto const mesh_id = mesh_subset.getMeshID();
151 for (auto const* node : mesh_subset.getNodes())
152 {
154 node->getID());
155
156 pos.setNodeID(node->getID());
157 pos.setCoordinates(*node);
158 auto const& ic_value = ic(t, pos);
159
160 auto global_index =
161 std::abs(dof_table_of_process.getGlobalIndex(l, variable_id,
162 component_id));
163#ifdef USE_PETSC
164 // The global indices of the ghost entries of the global
165 // matrix or the global vectors need to be set as negative
166 // values for equation assembly, however the global indices
167 // start from zero. Therefore, any ghost entry with zero
168 // index is assigned an negative value of the vector size
169 // or the matrix dimension. To assign the initial value for
170 // the ghost entries, the negative indices of the ghost
171 // entries are restored to zero.
172 if (global_index == x.size())
173 global_index = 0;
174#endif
175 x.set(global_index, ic_value[component_id]);
176 }
177 }
178 }
179
181 MathLib::LinAlg::copy(x, x_prev); // pushState
182
185
186 setInitialConditionsConcreteProcess(process_solutions, t, process_id);
187}
void setNodeID(std::size_t node_id)
void setCoordinates(MathLib::Point3d const &coordinates)
virtual void setInitialConditionsConcreteProcess(std::vector< GlobalVector * > &, double const, int const)
Definition: Process.h:224
void finalizeAssembly(PETScMatrix &A)
Definition: LinAlg.cpp:163
void copy(PETScVector const &x, PETScVector &y)
Definition: LinAlg.cpp:37

References _process_variables, MathLib::LinAlg::copy(), DBUG(), MathLib::LinAlg::finalizeAssembly(), getDOFTable(), MeshLib::Node, ParameterLib::SpatialPosition::setCoordinates(), setInitialConditionsConcreteProcess(), MathLib::LinAlg::setLocalAccessibleVector(), and ParameterLib::SpatialPosition::setNodeID().

◆ setInitialConditionsConcreteProcess()

◆ solveReactionEquation()

virtual void ProcessLib::Process::solveReactionEquation ( std::vector< GlobalVector * > &  ,
std::vector< GlobalVector * > const &  ,
double const  ,
double const  ,
NumLib::EquationSystem ,
int const   
)
inlinevirtual

Reimplemented in ProcessLib::ComponentTransport::ComponentTransportProcess.

Definition at line 182 of file Process.h.

187 {
188 }

◆ updateDeactivatedSubdomains()

void ProcessLib::Process::updateDeactivatedSubdomains ( double const  time,
const int  process_id 
)

Definition at line 197 of file Process.cpp.

199{
200 auto const& variables_per_process = getProcessVariables(process_id);
201 for (auto const& variable : variables_per_process)
202 {
203 variable.get().updateDeactivatedSubdomains(time);
204 }
205}
std::vector< std::reference_wrapper< ProcessVariable > > const & getProcessVariables(const int process_id) const
Definition: Process.h:155

References getProcessVariables().

Member Data Documentation

◆ _boundary_conditions

std::vector<BoundaryConditionCollection> ProcessLib::Process::_boundary_conditions
protected

Vector for boundary conditions. For the monolithic scheme or a single process, the size of the vector is one. For the staggered scheme, the size of vector is the number of the coupled processes.

Definition at line 387 of file Process.h.

Referenced by assemble(), assembleWithJacobian(), ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::createBHEBoundaryConditionTopBottom(), getKnownSolutions(), initializeProcessBoundaryConditionsAndSourceTerms(), postTimestep(), and preTimestep().

◆ _extrapolator_data

ExtrapolatorData ProcessLib::Process::_extrapolator_data
private

Definition at line 395 of file Process.h.

Referenced by getExtrapolator(), getSingleComponentDOFTable(), and initializeExtrapolator().

◆ _global_assembler

VectorMatrixAssembler ProcessLib::Process::_global_assembler
protected

Definition at line 359 of file Process.h.

Referenced by ProcessLib::ComponentTransport::ComponentTransportProcess::assembleConcreteProcess(), ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::assembleConcreteProcess(), ProcessLib::HT::HTProcess::assembleConcreteProcess(), ProcessLib::LiquidFlow::LiquidFlowProcess::assembleConcreteProcess(), ProcessLib::RichardsComponentTransport::RichardsComponentTransportProcess::assembleConcreteProcess(), ProcessLib::RichardsFlow::RichardsFlowProcess::assembleConcreteProcess(), ProcessLib::SteadyStateDiffusion::SteadyStateDiffusion::assembleConcreteProcess(), ProcessLib::TES::TESProcess::assembleConcreteProcess(), ProcessLib::ThermalTwoPhaseFlowWithPP::ThermalTwoPhaseFlowWithPPProcess::assembleConcreteProcess(), ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::assembleConcreteProcess(), ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess::assembleConcreteProcess(), ProcessLib::TwoPhaseFlowWithPrho::TwoPhaseFlowWithPrhoProcess::assembleConcreteProcess(), ProcessLib::HeatConduction::HeatConductionProcess::assembleConcreteProcess(), ProcessLib::ComponentTransport::ComponentTransportProcess::assembleWithJacobianConcreteProcess(), ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::assembleWithJacobianConcreteProcess(), ProcessLib::HT::HTProcess::assembleWithJacobianConcreteProcess(), ProcessLib::LiquidFlow::LiquidFlowProcess::assembleWithJacobianConcreteProcess(), ProcessLib::RichardsComponentTransport::RichardsComponentTransportProcess::assembleWithJacobianConcreteProcess(), ProcessLib::RichardsFlow::RichardsFlowProcess::assembleWithJacobianConcreteProcess(), ProcessLib::SteadyStateDiffusion::SteadyStateDiffusion::assembleWithJacobianConcreteProcess(), ProcessLib::TES::TESProcess::assembleWithJacobianConcreteProcess(), ProcessLib::ThermalTwoPhaseFlowWithPP::ThermalTwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess(), ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess(), ProcessLib::TwoPhaseFlowWithPrho::TwoPhaseFlowWithPrhoProcess::assembleWithJacobianConcreteProcess(), and ProcessLib::HeatConduction::HeatConductionProcess::assembleWithJacobianConcreteProcess().

◆ _integration_order

unsigned const ProcessLib::Process::_integration_order
protected

◆ _integration_point_writer

◆ _jacobian_assembler

◆ _local_to_global_index_map

std::unique_ptr<NumLib::LocalToGlobalIndexMap> ProcessLib::Process::_local_to_global_index_map
protected

Definition at line 352 of file Process.h.

Referenced by ProcessLib::ComponentTransport::ComponentTransportProcess::assembleConcreteProcess(), ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::assembleConcreteProcess(), ProcessLib::HT::HTProcess::assembleConcreteProcess(), ProcessLib::LiquidFlow::LiquidFlowProcess::assembleConcreteProcess(), ProcessLib::RichardsComponentTransport::RichardsComponentTransportProcess::assembleConcreteProcess(), ProcessLib::RichardsFlow::RichardsFlowProcess::assembleConcreteProcess(), ProcessLib::SteadyStateDiffusion::SteadyStateDiffusion::assembleConcreteProcess(), ProcessLib::TES::TESProcess::assembleConcreteProcess(), ProcessLib::ThermalTwoPhaseFlowWithPP::ThermalTwoPhaseFlowWithPPProcess::assembleConcreteProcess(), ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::assembleConcreteProcess(), ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess::assembleConcreteProcess(), ProcessLib::TwoPhaseFlowWithPrho::TwoPhaseFlowWithPrhoProcess::assembleConcreteProcess(), ProcessLib::HeatConduction::HeatConductionProcess::assembleConcreteProcess(), ProcessLib::AssemblyMixin< Process >::assembleGeneric(), ProcessLib::ComponentTransport::ComponentTransportProcess::assembleWithJacobianConcreteProcess(), ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::assembleWithJacobianConcreteProcess(), ProcessLib::HT::HTProcess::assembleWithJacobianConcreteProcess(), ProcessLib::LiquidFlow::LiquidFlowProcess::assembleWithJacobianConcreteProcess(), ProcessLib::RichardsComponentTransport::RichardsComponentTransportProcess::assembleWithJacobianConcreteProcess(), ProcessLib::RichardsFlow::RichardsFlowProcess::assembleWithJacobianConcreteProcess(), ProcessLib::SteadyStateDiffusion::SteadyStateDiffusion::assembleWithJacobianConcreteProcess(), ProcessLib::TES::TESProcess::assembleWithJacobianConcreteProcess(), ProcessLib::ThermalTwoPhaseFlowWithPP::ThermalTwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess(), ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::assembleWithJacobianConcreteProcess(), ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess(), ProcessLib::TwoPhaseFlowWithPrho::TwoPhaseFlowWithPrhoProcess::assembleWithJacobianConcreteProcess(), ProcessLib::HeatConduction::HeatConductionProcess::assembleWithJacobianConcreteProcess(), ProcessLib::TES::TESProcess::computeEquilibriumLoading(), ProcessLib::TES::TESProcess::computeRelativeHumidity(), computeSparsityPattern(), ProcessLib::TES::TESProcess::computeVapourPartialPressure(), ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::constructDofTable(), constructDofTableOfSpecifiedProcessStaggeredScheme(), constructMonolithicProcessDofTable(), ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::createBHEBoundaryConditionTopBottom(), getDOFTable(), getDOFTableForExtrapolatorData(), ProcessLib::HT::HTProcess::getDOFTableForExtrapolatorData(), ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::getDOFTables(), ProcessLib::ComponentTransport::ComponentTransportProcess::getFlux(), ProcessLib::LiquidFlow::LiquidFlowProcess::getFlux(), ProcessLib::HT::HTProcess::getFlux(), ProcessLib::SteadyStateDiffusion::SteadyStateDiffusion::getFlux(), getMatrixSpecifications(), initialize(), initializeBoundaryConditions(), ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::initializeConcreteProcess(), ProcessLib::TES::TESProcess::postIterationConcreteProcess(), ProcessLib::ComponentTransport::ComponentTransportProcess::preOutputConcreteProcess(), ProcessLib::ThermalTwoPhaseFlowWithPP::ThermalTwoPhaseFlowWithPPProcess::preTimestepConcreteProcess(), and ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::setInitialConditionsConcreteProcess().

◆ _mesh

◆ _mesh_subset_all_nodes

◆ _process_variables

◆ _secondary_variables

◆ _source_term_collections

std::vector<SourceTermCollection> ProcessLib::Process::_source_term_collections
private

Vector for nodal source term collections. For the monolithic scheme or a single process, the size of the vector is one. For the staggered scheme, the size of vector is the number of the coupled processes.

Definition at line 393 of file Process.h.

Referenced by assemble(), assembleWithJacobian(), and initializeProcessBoundaryConditionsAndSourceTerms().

◆ _sparsity_pattern

GlobalSparsityPattern ProcessLib::Process::_sparsity_pattern
protected

Definition at line 374 of file Process.h.

Referenced by computeSparsityPattern(), and getMatrixSpecifications().

◆ _use_monolithic_scheme

const bool ProcessLib::Process::_use_monolithic_scheme
protected

Definition at line 361 of file Process.h.

Referenced by ProcessLib::ComponentTransport::ComponentTransportProcess::ComponentTransportProcess(), ProcessLib::HydroMechanics::HydroMechanicsProcess< DisplacementDim >::HydroMechanicsProcess(), ProcessLib::ComponentTransport::ComponentTransportProcess::assembleConcreteProcess(), ProcessLib::HT::HTProcess::assembleConcreteProcess(), ProcessLib::ComponentTransport::ComponentTransportProcess::assembleWithJacobianConcreteProcess(), ProcessLib::HT::HTProcess::assembleWithJacobianConcreteProcess(), constructDofTable(), ProcessLib::HT::HTProcess::getDOFTableForExtrapolatorData(), ProcessLib::RichardsMechanics::RichardsMechanicsProcess< DisplacementDim >::hasMechanicalProcess(), ProcessLib::TH2M::TH2MProcess< DisplacementDim >::hasMechanicalProcess(), ProcessLib::ThermoHydroMechanics::ThermoHydroMechanicsProcess< DisplacementDim >::hasMechanicalProcess(), ProcessLib::ComponentTransport::ComponentTransportProcess::initializeConcreteProcess(), ProcessLib::HT::HTProcess::initializeConcreteProcess(), isMonolithicSchemeUsed(), ProcessLib::HT::HTProcess::postTimestepConcreteProcess(), ProcessLib::ComponentTransport::ComponentTransportProcess::postTimestepConcreteProcess(), ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::postTimestepConcreteProcess(), ProcessLib::ComponentTransport::ComponentTransportProcess::preOutputConcreteProcess(), and ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::setInitialConditionsConcreteProcess().

◆ constant_one_parameter_name

const std::string ProcessLib::Process::constant_one_parameter_name = "constant_one"
static

◆ name

std::string const ProcessLib::Process::name

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