OGS 6.2.1-76-gbb689931b
ProcessLib::VectorMatrixAssembler Class Referencefinal

Detailed Description

Utility class used to assemble global matrices and vectors.

The methods of this class get the global matrices and vectors as input and pass only local data on to the local assemblers.

Definition at line 32 of file VectorMatrixAssembler.h.

#include <VectorMatrixAssembler.h>

Public Member Functions

 VectorMatrixAssembler (std::unique_ptr< AbstractJacobianAssembler > &&jacobian_assembler)
 
void preAssemble (const std::size_t mesh_item_id, LocalAssemblerInterface &local_assembler, const NumLib::LocalToGlobalIndexMap &dof_table, const double t, const GlobalVector &x)
 
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, GlobalVector const &x, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, CoupledSolutionsForStaggeredScheme const *const cpl_xs)
 
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, GlobalVector const &x, GlobalVector const &xdot, const double dxdot_dx, const double dx_dx, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac, CoupledSolutionsForStaggeredScheme const *const cpl_xs)
 

Private Attributes

std::vector< double > _local_M_data
 
std::vector< double > _local_K_data
 
std::vector< double > _local_b_data
 
std::vector< double > _local_Jac_data
 
std::unique_ptr< AbstractJacobianAssembler_jacobian_assembler
 Used to assemble the Jacobian. More...
 

Constructor & Destructor Documentation

◆ VectorMatrixAssembler()

ProcessLib::VectorMatrixAssembler::VectorMatrixAssembler ( std::unique_ptr< AbstractJacobianAssembler > &&  jacobian_assembler)
explicit

Definition at line 24 of file VectorMatrixAssembler.cpp.

26  : _jacobian_assembler(std::move(jacobian_assembler))
27 {
28 }
std::unique_ptr< AbstractJacobianAssembler > _jacobian_assembler
Used to assemble the Jacobian.

Member Function Documentation

◆ assemble()

void ProcessLib::VectorMatrixAssembler::assemble ( std::size_t const  mesh_item_id,
LocalAssemblerInterface local_assembler,
std::vector< std::reference_wrapper< NumLib::LocalToGlobalIndexMap >> const &  dof_tables,
double const  t,
GlobalVector const &  x,
GlobalMatrix &  M,
GlobalMatrix &  K,
GlobalVector &  b,
CoupledSolutionsForStaggeredScheme const *const  cpl_xs 
)

AssemblesM, K, and b.

Remarks
Jacobian is not assembled here, see assembleWithJacobian().

Definition at line 41 of file VectorMatrixAssembler.cpp.

References _local_b_data, _local_K_data, _local_M_data, ProcessLib::LocalAssemblerInterface::assemble(), ProcessLib::LocalAssemblerInterface::assembleForStaggeredScheme(), ProcessLib::CoupledSolutionsForStaggeredScheme::dt, ProcessLib::getCurrentLocalSolutions(), NumLib::getIndices(), ProcessLib::getPreviousLocalSolutions(), ProcessLib::CoupledSolutionsForStaggeredScheme::process_id, and MathLib::toMatrix().

Referenced by ProcessLib::HeatConduction::HeatConductionProcess::assembleConcreteProcess(), ProcessLib::SmallDeformation::SmallDeformationProcess< DisplacementDim >::assembleConcreteProcess(), ProcessLib::RichardsFlow::RichardsFlowProcess::assembleConcreteProcess(), ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::assembleConcreteProcess(), ProcessLib::SmallDeformationNonlocal::SmallDeformationNonlocalProcess< DisplacementDim >::assembleConcreteProcess(), ProcessLib::TwoPhaseFlowWithPrho::TwoPhaseFlowWithPrhoProcess::assembleConcreteProcess(), ProcessLib::TES::TESProcess::assembleConcreteProcess(), ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess::assembleConcreteProcess(), ProcessLib::ThermalTwoPhaseFlowWithPP::ThermalTwoPhaseFlowWithPPProcess::assembleConcreteProcess(), ProcessLib::LIE::SmallDeformation::SmallDeformationProcess< DisplacementDim >::assembleConcreteProcess(), ProcessLib::LIE::HydroMechanics::HydroMechanicsProcess< GlobalDim >::assembleConcreteProcess(), ProcessLib::ThermoMechanics::ThermoMechanicsProcess< DisplacementDim >::assembleConcreteProcess(), ProcessLib::HydroMechanics::HydroMechanicsProcess< DisplacementDim >::assembleConcreteProcess(), ProcessLib::RichardsMechanics::RichardsMechanicsProcess< DisplacementDim >::assembleConcreteProcess(), ProcessLib::ThermoHydroMechanics::ThermoHydroMechanicsProcess< DisplacementDim >::assembleConcreteProcess(), ProcessLib::PhaseField::PhaseFieldProcess< DisplacementDim >::assembleConcreteProcess(), ProcessLib::LiquidFlow::LiquidFlowProcess::assembleConcreteProcess(), ProcessLib::GroundwaterFlow::GroundwaterFlowProcess::assembleConcreteProcess(), ProcessLib::HT::HTProcess::assembleConcreteProcess(), ProcessLib::ThermoMechanicalPhaseField::ThermoMechanicalPhaseFieldProcess< DisplacementDim >::assembleConcreteProcess(), ProcessLib::RichardsComponentTransport::RichardsComponentTransportProcess::assembleConcreteProcess(), and ProcessLib::ComponentTransport::ComponentTransportProcess::assembleConcreteProcess().

47 {
48  std::vector<std::vector<GlobalIndexType>> indices_of_processes;
49  indices_of_processes.reserve(dof_tables.size());
50  for (auto dof_table : dof_tables)
51  {
52  indices_of_processes.emplace_back(
53  NumLib::getIndices(mesh_item_id, dof_table.get()));
54  }
55 
56  auto const& indices = (cpl_xs == nullptr)
57  ? indices_of_processes[0]
58  : indices_of_processes[cpl_xs->process_id];
59  _local_M_data.clear();
60  _local_K_data.clear();
61  _local_b_data.clear();
62 
63  if (cpl_xs == nullptr)
64  {
65  auto const local_x = x.get(indices);
66  local_assembler.assemble(t, local_x, _local_M_data, _local_K_data,
68  }
69  else
70  {
71  auto local_coupled_xs0 =
72  getPreviousLocalSolutions(*cpl_xs, indices_of_processes);
73  auto local_coupled_xs =
74  getCurrentLocalSolutions(*cpl_xs, indices_of_processes);
75 
76  ProcessLib::LocalCoupledSolutions local_coupled_solutions(
77  cpl_xs->dt, cpl_xs->process_id, std::move(local_coupled_xs0),
78  std::move(local_coupled_xs));
79 
80  local_assembler.assembleForStaggeredScheme(t, _local_M_data,
82  local_coupled_solutions);
83  }
84 
85  auto const num_r_c = indices.size();
86  auto const r_c_indices =
88 
89  if (!_local_M_data.empty())
90  {
91  auto const local_M = MathLib::toMatrix(_local_M_data, num_r_c, num_r_c);
92  M.add(r_c_indices, local_M);
93  }
94  if (!_local_K_data.empty())
95  {
96  auto const local_K = MathLib::toMatrix(_local_K_data, num_r_c, num_r_c);
97  K.add(r_c_indices, local_K);
98  }
99  if (!_local_b_data.empty())
100  {
101  assert(_local_b_data.size() == num_r_c);
102  b.add(indices, _local_b_data);
103  }
104 }
MathLib::RowColumnIndices< GlobalIndexType > RowColumnIndices
std::vector< GlobalIndexType > getIndices(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table)
std::vector< std::vector< double > > getPreviousLocalSolutions(const CoupledSolutionsForStaggeredScheme &cpl_xs, const std::vector< std::vector< GlobalIndexType >> &indices)
Eigen::Map< const Matrix > toMatrix(std::vector< double > const &data, Eigen::MatrixXd::Index rows, Eigen::MatrixXd::Index cols)
Definition: EigenMapTools.h:71
std::vector< std::vector< double > > getCurrentLocalSolutions(const CoupledSolutionsForStaggeredScheme &cpl_xs, const std::vector< std::vector< GlobalIndexType >> &indices)

◆ assembleWithJacobian()

void ProcessLib::VectorMatrixAssembler::assembleWithJacobian ( std::size_t const  mesh_item_id,
LocalAssemblerInterface local_assembler,
std::vector< std::reference_wrapper< NumLib::LocalToGlobalIndexMap >> const &  dof_tables,
const double  t,
GlobalVector const &  x,
GlobalVector const &  xdot,
const double  dxdot_dx,
const double  dx_dx,
GlobalMatrix &  M,
GlobalMatrix &  K,
GlobalVector &  b,
GlobalMatrix &  Jac,
CoupledSolutionsForStaggeredScheme const *const  cpl_xs 
)

Assembles M, K, b, and the Jacobian Jac of the residual.

Note
The Jacobian must be assembled.

Definition at line 106 of file VectorMatrixAssembler.cpp.

References _jacobian_assembler, _local_b_data, _local_Jac_data, _local_K_data, _local_M_data, ProcessLib::CoupledSolutionsForStaggeredScheme::dt, ProcessLib::getCurrentLocalSolutions(), NumLib::getIndices(), ProcessLib::getPreviousLocalSolutions(), OGS_FATAL, ProcessLib::CoupledSolutionsForStaggeredScheme::process_id, and MathLib::toMatrix().

Referenced by ProcessLib::HeatConduction::HeatConductionProcess::assembleWithJacobianConcreteProcess(), ProcessLib::SmallDeformation::SmallDeformationProcess< DisplacementDim >::assembleWithJacobianConcreteProcess(), ProcessLib::RichardsFlow::RichardsFlowProcess::assembleWithJacobianConcreteProcess(), ProcessLib::TwoPhaseFlowWithPrho::TwoPhaseFlowWithPrhoProcess::assembleWithJacobianConcreteProcess(), ProcessLib::SmallDeformationNonlocal::SmallDeformationNonlocalProcess< DisplacementDim >::assembleWithJacobianConcreteProcess(), ProcessLib::TwoPhaseFlowWithPP::TwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess(), ProcessLib::ThermalTwoPhaseFlowWithPP::ThermalTwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess(), ProcessLib::LIE::SmallDeformation::SmallDeformationProcess< DisplacementDim >::assembleWithJacobianConcreteProcess(), ProcessLib::TES::TESProcess::assembleWithJacobianConcreteProcess(), ProcessLib::LIE::HydroMechanics::HydroMechanicsProcess< GlobalDim >::assembleWithJacobianConcreteProcess(), ProcessLib::ThermoMechanics::ThermoMechanicsProcess< DisplacementDim >::assembleWithJacobianConcreteProcess(), ProcessLib::HydroMechanics::HydroMechanicsProcess< DisplacementDim >::assembleWithJacobianConcreteProcess(), ProcessLib::RichardsMechanics::RichardsMechanicsProcess< DisplacementDim >::assembleWithJacobianConcreteProcess(), ProcessLib::ThermoHydroMechanics::ThermoHydroMechanicsProcess< DisplacementDim >::assembleWithJacobianConcreteProcess(), ProcessLib::PhaseField::PhaseFieldProcess< DisplacementDim >::assembleWithJacobianConcreteProcess(), ProcessLib::LiquidFlow::LiquidFlowProcess::assembleWithJacobianConcreteProcess(), ProcessLib::GroundwaterFlow::GroundwaterFlowProcess::assembleWithJacobianConcreteProcess(), ProcessLib::HT::HTProcess::assembleWithJacobianConcreteProcess(), ProcessLib::ThermoMechanicalPhaseField::ThermoMechanicalPhaseFieldProcess< DisplacementDim >::assembleWithJacobianConcreteProcess(), ProcessLib::RichardsComponentTransport::RichardsComponentTransportProcess::assembleWithJacobianConcreteProcess(), and ProcessLib::ComponentTransport::ComponentTransportProcess::assembleWithJacobianConcreteProcess().

114 {
115  std::vector<std::vector<GlobalIndexType>> indices_of_processes;
116  indices_of_processes.reserve(dof_tables.size());
117  for (auto dof_table : dof_tables)
118  {
119  indices_of_processes.emplace_back(
120  NumLib::getIndices(mesh_item_id, dof_table.get()));
121  }
122 
123  auto const& indices = (cpl_xs == nullptr)
124  ? indices_of_processes[0]
125  : indices_of_processes[cpl_xs->process_id];
126  auto const local_xdot = xdot.get(indices);
127 
128  _local_M_data.clear();
129  _local_K_data.clear();
130  _local_b_data.clear();
131  _local_Jac_data.clear();
132 
133  if (cpl_xs == nullptr)
134  {
135  auto const local_x = x.get(indices);
136  _jacobian_assembler->assembleWithJacobian(
137  local_assembler, t, local_x, local_xdot, dxdot_dx, dx_dx,
139  }
140  else
141  {
142  auto local_coupled_xs0 =
143  getPreviousLocalSolutions(*cpl_xs, indices_of_processes);
144  auto local_coupled_xs =
145  getCurrentLocalSolutions(*cpl_xs, indices_of_processes);
146 
147  ProcessLib::LocalCoupledSolutions local_coupled_solutions(
148  cpl_xs->dt, cpl_xs->process_id, std::move(local_coupled_xs0),
149  std::move(local_coupled_xs));
150 
151  _jacobian_assembler->assembleWithJacobianForStaggeredScheme(
152  local_assembler, t, local_xdot, dxdot_dx, dx_dx, _local_M_data,
154  local_coupled_solutions);
155  }
156 
157  auto const num_r_c = indices.size();
158  auto const r_c_indices =
160 
161  if (!_local_M_data.empty())
162  {
163  auto const local_M = MathLib::toMatrix(_local_M_data, num_r_c, num_r_c);
164  M.add(r_c_indices, local_M);
165  }
166  if (!_local_K_data.empty())
167  {
168  auto const local_K = MathLib::toMatrix(_local_K_data, num_r_c, num_r_c);
169  K.add(r_c_indices, local_K);
170  }
171  if (!_local_b_data.empty())
172  {
173  assert(_local_b_data.size() == num_r_c);
174  b.add(indices, _local_b_data);
175  }
176  if (!_local_Jac_data.empty())
177  {
178  auto const local_Jac =
179  MathLib::toMatrix(_local_Jac_data, num_r_c, num_r_c);
180  Jac.add(r_c_indices, local_Jac);
181  }
182  else
183  {
184  OGS_FATAL(
185  "No Jacobian has been assembled! This might be due to programming "
186  "errors in the local assembler of the current process.");
187  }
188 }
MathLib::RowColumnIndices< GlobalIndexType > RowColumnIndices
std::vector< GlobalIndexType > getIndices(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table)
std::vector< std::vector< double > > getPreviousLocalSolutions(const CoupledSolutionsForStaggeredScheme &cpl_xs, const std::vector< std::vector< GlobalIndexType >> &indices)
Eigen::Map< const Matrix > toMatrix(std::vector< double > const &data, Eigen::MatrixXd::Index rows, Eigen::MatrixXd::Index cols)
Definition: EigenMapTools.h:71
#define OGS_FATAL(fmt,...)
Definition: Error.h:63
std::vector< std::vector< double > > getCurrentLocalSolutions(const CoupledSolutionsForStaggeredScheme &cpl_xs, const std::vector< std::vector< GlobalIndexType >> &indices)
std::unique_ptr< AbstractJacobianAssembler > _jacobian_assembler
Used to assemble the Jacobian.

◆ preAssemble()

void ProcessLib::VectorMatrixAssembler::preAssemble ( const std::size_t  mesh_item_id,
LocalAssemblerInterface local_assembler,
const NumLib::LocalToGlobalIndexMap dof_table,
const double  t,
const GlobalVector &  x 
)

Definition at line 30 of file VectorMatrixAssembler.cpp.

References NumLib::getIndices(), and ProcessLib::LocalAssemblerInterface::preAssemble().

Referenced by ProcessLib::SmallDeformationNonlocal::SmallDeformationNonlocalProcess< DisplacementDim >::preAssembleConcreteProcess().

34 {
35  auto const indices = NumLib::getIndices(mesh_item_id, dof_table);
36  auto const local_x = x.get(indices);
37 
38  local_assembler.preAssemble(t, local_x);
39 }
std::vector< GlobalIndexType > getIndices(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table)

Member Data Documentation

◆ _jacobian_assembler

std::unique_ptr<AbstractJacobianAssembler> ProcessLib::VectorMatrixAssembler::_jacobian_assembler
private

Used to assemble the Jacobian.

Definition at line 75 of file VectorMatrixAssembler.h.

Referenced by assembleWithJacobian().

◆ _local_b_data

std::vector<double> ProcessLib::VectorMatrixAssembler::_local_b_data
private

Definition at line 71 of file VectorMatrixAssembler.h.

Referenced by assemble(), and assembleWithJacobian().

◆ _local_Jac_data

std::vector<double> ProcessLib::VectorMatrixAssembler::_local_Jac_data
private

Definition at line 72 of file VectorMatrixAssembler.h.

Referenced by assembleWithJacobian().

◆ _local_K_data

std::vector<double> ProcessLib::VectorMatrixAssembler::_local_K_data
private

Definition at line 70 of file VectorMatrixAssembler.h.

Referenced by assemble(), and assembleWithJacobian().

◆ _local_M_data

std::vector<double> ProcessLib::VectorMatrixAssembler::_local_M_data
private

Definition at line 69 of file VectorMatrixAssembler.h.

Referenced by assemble(), and assembleWithJacobian().


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