OGS
ProcessLib::AssemblyMixin< Process > Class Template Reference

Detailed Description

template<typename Process>
class ProcessLib::AssemblyMixin< Process >

A mixin providing assembly functionality to a specific Process.

The process must be derived from this class (CRTP).

Definition at line 95 of file AssemblyMixin.h.

#include <AssemblyMixin.h>

Inheritance diagram for ProcessLib::AssemblyMixin< Process >:
[legend]
Collaboration diagram for ProcessLib::AssemblyMixin< Process >:
[legend]

Public Member Functions

void initializeAssemblyOnSubmeshes (const int process_id, std::vector< std::reference_wrapper< MeshLib::Mesh > > const &submeshes, std::vector< std::string > const &residuum_names)
 
void updateActiveElements ()
 
void assemble (const double, double const, std::vector< GlobalVector * > const &, std::vector< GlobalVector * > const &, int const, GlobalMatrix &, GlobalMatrix &, GlobalVector &)
 
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)
 

Private Member Functions

Processderived ()
 
Process const & derived () const
 
template<typename Method , typename... Jac>
void assembleGeneric (Method global_assembler_method, 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, Jac &... jac_or_not_jac)
 
 AssemblyMixinBase (AbstractJacobianAssembler &jacobian_assembler)
 
- Private Member Functions inherited from ProcessLib::AssemblyMixinBase
 AssemblyMixinBase (AbstractJacobianAssembler &jacobian_assembler)
 
void initializeAssemblyOnSubmeshes (const int process_id, MeshLib::Mesh &bulk_mesh, std::vector< std::reference_wrapper< MeshLib::Mesh > > const &submeshes, std::vector< std::string > const &residuum_names, std::vector< std::reference_wrapper< ProcessVariable > > const &pvs)
 
void updateActiveElements (ProcessLib::Process const &process)
 

Private Attributes

friend Process
 
- Private Attributes inherited from ProcessLib::AssemblyMixinBase
std::vector< SubmeshAssemblyDatasubmesh_assembly_data_
 
std::vector< std::reference_wrapper< MeshLib::PropertyVector< double > > > residuum_vectors_bulk_
 
std::size_t b_submesh_id_ = 0
 ID of the b vector on submeshes, cf. NumLib::VectorProvider.
 
Assembly::ParallelVectorMatrixAssembler pvma_
 

Additional Inherited Members

- Static Private Member Functions inherited from ProcessLib::AssemblyMixinBase
static void copyResiduumVectorsToBulkMesh (GlobalVector const &rhs, NumLib::LocalToGlobalIndexMap const &local_to_global_index_map, std::vector< std::reference_wrapper< MeshLib::PropertyVector< double > > > residuum_vectors)
 
static void copyResiduumVectorsToSubmesh (GlobalVector const &rhs, NumLib::LocalToGlobalIndexMap const &local_to_global_index_map, SubmeshAssemblyData const &sad)
 

Member Function Documentation

◆ assemble()

template<typename Process >
void ProcessLib::AssemblyMixin< Process >::assemble ( const double ,
double const ,
std::vector< GlobalVector * > const & ,
std::vector< GlobalVector * > const & ,
int const ,
GlobalMatrix & ,
GlobalMatrix & ,
GlobalVector &  )
inline

Definition at line 140 of file AssemblyMixin.h.

145 {
146 /*
147 DBUG("AssemblyMixin assemble(t={}, dt={}, process_id={}).", t, dt,
148 process_id);
149
150 assembleGeneric(&Assembly::ParallelVectorMatrixAssembler::assemble, t,
151 dt, x, x_prev, process_id, M, K, b);
152 */
153 OGS_FATAL("Not yet implemented.");
154 }
#define OGS_FATAL(...)
Definition Error.h:26

References OGS_FATAL.

◆ assembleGeneric()

template<typename Process >
template<typename Method , typename... Jac>
void ProcessLib::AssemblyMixin< Process >::assembleGeneric ( Method global_assembler_method,
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,
Jac &... jac_or_not_jac )
inlineprivate

Generic assembly routine covering both the case with and without Jacobian assembly.

Definition at line 181 of file AssemblyMixin.h.

186 {
187 // TODO why not getDOFTables(x.size()); ?
188 std::vector<NumLib::LocalToGlobalIndexMap const*> const dof_tables{
190
191 auto const& loc_asms = derived().local_assemblers_;
192
193 if (!submesh_assembly_data_.empty())
194 {
196 b, b_submesh_id_);
197
198 for (auto const& sad : submesh_assembly_data_)
199 {
200 b_submesh.setZero();
201
202 (pvma_.*global_assembler_method)(
203 loc_asms, sad.active_element_ids, dof_tables, t, dt, x,
204 x_prev, process_id, M, K, b_submesh, jac_or_not_jac...);
205
206 MathLib::LinAlg::axpy(b, 1.0, b_submesh);
207
209 b_submesh, *(dof_tables.front()), sad);
210 }
211
213 }
214 else
215 {
216 // convention: process variable 0 governs where assembly takes
217 // place (active element IDs)
219 derived().getProcessVariables(process_id)[0];
220
221 (pvma_.*global_assembler_method)(
222 loc_asms, pv.getActiveElementIDs(), dof_tables, t, dt, x,
223 x_prev, process_id, M, K, b, jac_or_not_jac...);
224 }
225
227 b, *(dof_tables.front()), residuum_vectors_bulk_);
228 }
virtual GlobalVector & getVector(std::size_t &id)=0
Get an uninitialized vector with the given id.
virtual void releaseVector(GlobalVector const &x)=0
static void copyResiduumVectorsToSubmesh(GlobalVector const &rhs, NumLib::LocalToGlobalIndexMap const &local_to_global_index_map, SubmeshAssemblyData const &sad)
std::vector< std::reference_wrapper< MeshLib::PropertyVector< double > > > residuum_vectors_bulk_
Assembly::ParallelVectorMatrixAssembler pvma_
static void copyResiduumVectorsToBulkMesh(GlobalVector const &rhs, NumLib::LocalToGlobalIndexMap const &local_to_global_index_map, std::vector< std::reference_wrapper< MeshLib::PropertyVector< double > > > residuum_vectors)
std::size_t b_submesh_id_
ID of the b vector on submeshes, cf. NumLib::VectorProvider.
std::vector< SubmeshAssemblyData > submesh_assembly_data_
std::vector< std::size_t > const & getActiveElementIDs() const
std::vector< std::reference_wrapper< ProcessVariable > > const & getProcessVariables(const int process_id) const
Definition Process.h:156
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition Process.h:361
void axpy(PETScVector &y, PetscScalar const a, PETScVector const &x)
Definition LinAlg.cpp:57
static NUMLIB_EXPORT VectorProvider & provider

References ProcessLib::Process::_local_to_global_index_map, MathLib::LinAlg::axpy(), ProcessLib::AssemblyMixinBase::b_submesh_id_, ProcessLib::AssemblyMixinBase::copyResiduumVectorsToBulkMesh(), ProcessLib::AssemblyMixinBase::copyResiduumVectorsToSubmesh(), ProcessLib::AssemblyMixin< Process >::derived(), ProcessLib::ProcessVariable::getActiveElementIDs(), ProcessLib::Process::getProcessVariables(), NumLib::VectorProvider::getVector(), NumLib::GlobalVectorProvider::provider, ProcessLib::AssemblyMixinBase::pvma_, NumLib::VectorProvider::releaseVector(), ProcessLib::AssemblyMixinBase::residuum_vectors_bulk_, MathLib::EigenVector::setZero(), and ProcessLib::AssemblyMixinBase::submesh_assembly_data_.

Referenced by ProcessLib::AssemblyMixin< Process >::assembleWithJacobian().

◆ assembleWithJacobian()

template<typename Process >
void ProcessLib::AssemblyMixin< 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 )
inline

Definition at line 156 of file AssemblyMixin.h.

162 {
163 DBUG("AssemblyMixin assembleWithJacobian(t={}, dt={}, process_id={}).",
164 t, dt, process_id);
165
168 dt, x, x_prev, process_id, M, K, b, Jac);
169 }
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:30
void assembleGeneric(Method global_assembler_method, 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, Jac &... jac_or_not_jac)
void assembleWithJacobian(BaseLib::PolymorphicRandomAccessContainerView< LocalAssemblerInterface > const &local_assemblers, std::vector< std::size_t > const &active_elements, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, const double t, double const dt, std::vector< GlobalVector * > const &xs, std::vector< GlobalVector * > const &x_prevs, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac)

References ProcessLib::AssemblyMixin< Process >::assembleGeneric(), ProcessLib::Assembly::ParallelVectorMatrixAssembler::assembleWithJacobian(), and DBUG().

◆ AssemblyMixinBase()

template<typename Process >
ProcessLib::AssemblyMixinBase::AssemblyMixinBase ( AbstractJacobianAssembler & jacobian_assembler)
inlineexplicitprivate

Definition at line 49 of file AssemblyMixin.h.

50 : pvma_{jacobian_assembler} {};

◆ derived() [1/2]

◆ derived() [2/2]

template<typename Process >
Process const & ProcessLib::AssemblyMixin< Process >::derived ( ) const
inlineprivate

Definition at line 173 of file AssemblyMixin.h.

174 {
175 return static_cast<Process const&>(*this);
176 }

◆ initializeAssemblyOnSubmeshes()

template<typename Process >
void ProcessLib::AssemblyMixin< Process >::initializeAssemblyOnSubmeshes ( const int process_id,
std::vector< std::reference_wrapper< MeshLib::Mesh > > const & submeshes,
std::vector< std::string > const & residuum_names )
inline

Specifies that the assembly of the process should take place on the passed submeshes.

Note
This method must be called at most once.
Attention
The passed submeshes must be non-overlapping and must cover the entire simulation domain. The caller is responsible to meet this requirement.
The order of the passed residuum_names matters and must match the order of the residuum vectors in the process. Otherwise residuum output on submeshes will be wrong.
Note
As of January 2023 assembly on submeshes has neither been tested with domain decomposition, nor with domain deactivation nor staggered coupling. Use at your own risk in these cases.
Todo
Implement and test with domain decomposition, domain deactivation and staggered coupling.

Definition at line 124 of file AssemblyMixin.h.

128 {
130 process_id, derived().getMesh(), submeshes, residuum_names,
131 derived().getProcessVariables(process_id));
132 }
void initializeAssemblyOnSubmeshes(const int process_id, MeshLib::Mesh &bulk_mesh, std::vector< std::reference_wrapper< MeshLib::Mesh > > const &submeshes, std::vector< std::string > const &residuum_names, std::vector< std::reference_wrapper< ProcessVariable > > const &pvs)
OGSMesh getMesh(std::string const &name)

References ProcessLib::AssemblyMixin< Process >::derived(), getMesh(), and ProcessLib::AssemblyMixinBase::initializeAssemblyOnSubmeshes().

◆ updateActiveElements()

template<typename Process >
void ProcessLib::AssemblyMixin< Process >::updateActiveElements ( )
inline

Member Data Documentation

◆ Process

template<typename Process >
friend ProcessLib::AssemblyMixin< Process >::Process
private

Definition at line 100 of file AssemblyMixin.h.


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