OGS
ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim > Struct Template Reference

Detailed Description

template<typename BcOrStData, typename ShapeFunction, typename LowerOrderShapeFunction, typename IntegrationMethod, int GlobalDim>
struct ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >

Common parts of the implementation of Python boundary conditions and source terms.

The implementation of Python boundary conditions and source terms are almost the same, however, their interfaces differ slightly. This struct contains all common parts of the local assemblers.

Definition at line 34 of file BcAndStLocalAssemblerImpl.h.

#include <BcAndStLocalAssemblerImpl.h>

Collaboration diagram for ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >:
[legend]

Public Types

using Traits = NsAndWeightsTraits< ShapeFunction, LowerOrderShapeFunction, GlobalDim >
 

Public Member Functions

 BcAndStLocalAssemblerImpl (MeshLib::Element const &e, bool is_axially_symmetric, unsigned const integration_order, BcOrStData const &data)
 
void assemble (std::size_t const boundary_element_id, NumLib::LocalToGlobalIndexMap const &dof_table_boundary, double const t, GlobalVector const &x, GlobalVector &b, GlobalMatrix *const Jac) const
 

Public Attributes

BcOrStData const & bc_or_st_data
 
MeshLib::Element const & element
 
IntegrationMethod const integration_method
 
std::vector< typename Traits::NsAndWeight > const nss_and_weights
 

Private Types

using ShapeMatrixPolicy = typename Traits::ShapeMatrixPolicy
 

Private Member Functions

std::array< double, 3 > interpolateCoords (typename Traits::NsAndWeight const &ns_and_weight, NumLib::TemplateIsoparametric< ShapeFunction, ShapeMatrixPolicy > const &fe) const
 
void assembleLocalRhs (double const flux, typename Traits::NsAndWeight const &ns_and_weight, Eigen::VectorXd &local_rhs) const
 
void assembleLocalJacobian (std::vector< double > const &dFlux, typename Traits::NsAndWeight const &ns_and_weight, Eigen::MatrixXd &local_Jac) const
 

Member Typedef Documentation

◆ ShapeMatrixPolicy

template<typename BcOrStData , typename ShapeFunction , typename LowerOrderShapeFunction , typename IntegrationMethod , int GlobalDim>
using ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::ShapeMatrixPolicy = typename Traits::ShapeMatrixPolicy
private

Definition at line 41 of file BcAndStLocalAssemblerImpl.h.

◆ Traits

template<typename BcOrStData , typename ShapeFunction , typename LowerOrderShapeFunction , typename IntegrationMethod , int GlobalDim>
using ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::Traits = NsAndWeightsTraits<ShapeFunction, LowerOrderShapeFunction, GlobalDim>

Definition at line 37 of file BcAndStLocalAssemblerImpl.h.

Constructor & Destructor Documentation

◆ BcAndStLocalAssemblerImpl()

template<typename BcOrStData , typename ShapeFunction , typename LowerOrderShapeFunction , typename IntegrationMethod , int GlobalDim>
ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::BcAndStLocalAssemblerImpl ( MeshLib::Element const &  e,
bool  is_axially_symmetric,
unsigned const  integration_order,
BcOrStData const &  data 
)
inline

Definition at line 44 of file BcAndStLocalAssemblerImpl.h.

48 : bc_or_st_data(data),
49 element(e),
50 integration_method(integration_order),
52 computeNsAndWeights<ShapeFunction, LowerOrderShapeFunction,
53 GlobalDim>(e, is_axially_symmetric,
55 {
56 }
auto computeNsAndWeights(MeshLib::Element const &element, bool const is_axially_symmetric, IntegrationMethod const &integration_method)
Definition: NsAndWeight.h:143

Member Function Documentation

◆ assemble()

template<typename BcOrStData , typename ShapeFunction , typename LowerOrderShapeFunction , typename IntegrationMethod , int GlobalDim>
void ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::assemble ( std::size_t const  boundary_element_id,
NumLib::LocalToGlobalIndexMap const &  dof_table_boundary,
double const  t,
GlobalVector const &  x,
GlobalVector b,
GlobalMatrix *const  Jac 
) const
inline

Assembles the local rhs and (possibly) Jacobian for this BC or ST and adds them to the global rhs vector and Jacobian matrix, respectively.

Definition at line 60 of file BcAndStLocalAssemblerImpl.h.

64 {
65 // Variables for input data
66 Eigen::MatrixXd const primary_variables_mat = ProcessLib::
69 dof_table_boundary, x);
70
71 // Variables for output data
72 auto const& indices_this_component =
73 dof_table_boundary(boundary_element_id,
75 .rows;
76 std::vector<GlobalIndexType> indices_all_components_for_Jac;
77
78 auto const num_dof_this_component = indices_this_component.size();
79
80 Eigen::VectorXd local_rhs =
81 Eigen::VectorXd::Zero(num_dof_this_component);
82 Eigen::MatrixXd local_Jac;
83
84 if (Jac)
85 {
86 indices_all_components_for_Jac =
87 NumLib::getIndices(boundary_element_id, dof_table_boundary);
88
89 auto const num_dof_all_components =
90 indices_all_components_for_Jac.size();
91
92 local_Jac = Eigen::MatrixXd::Zero(num_dof_this_component,
93 num_dof_all_components);
94 }
95
96 // Variables for temporary data
97 auto const num_comp_total =
98 dof_table_boundary.getNumberOfGlobalComponents();
99 std::vector<double> prim_vars_data(num_comp_total);
100 auto prim_vars = MathLib::toVector(prim_vars_data);
101
102 // Variables for the integration loop
103 unsigned const num_integration_points =
104 integration_method.getNumberOfPoints();
106 ShapeFunction, ShapeMatrixPolicy>(element);
107
108 for (unsigned ip = 0; ip < num_integration_points; ip++)
109 {
110 auto const& ns_and_weight = nss_and_weights[ip];
111 auto const coords = interpolateCoords(ns_and_weight, fe);
112
114 primary_variables_mat,
116 nss_and_weights[ip], prim_vars);
117
118 auto const [flag, flux, dFlux] =
119 bc_or_st_data.getFlagAndFluxAndDFlux(t, coords, prim_vars_data);
120
121 if (!flag)
122 {
123 // No flux value for this integration point. Skip assembly of
124 // the entire element.
125 return;
126 }
127
128 assembleLocalRhs(flux, ns_and_weight, local_rhs);
129
130 if (static_cast<int>(dFlux.size()) != num_comp_total)
131 {
132 // This strict check is technically mandatory only if a Jacobian
133 // is assembled. However, it is done as a consistency check also
134 // for cases without Jacobian assembly.
135 OGS_FATAL(
136 "The Python BC or ST must return the derivative of the "
137 "flux w.r.t. each component of each primary variable. "
138 "{:d} components expected. {:d} components returned from "
139 "Python."
140 "The expected number of components depends on multiple "
141 "factors. In the case of vectorial primary variables (such "
142 "as displacement), it will depend on the space dimension. "
143 "In the case of staggered coupling, the number of "
144 "components can be different for each staggered process "
145 "and different from monolithic coupling.",
146 num_comp_total, dFlux.size());
147 }
148
149 if (Jac)
150 {
151 assembleLocalJacobian(dFlux, ns_and_weight, local_Jac);
152 }
153 }
154
155 b.add(indices_this_component, local_rhs);
156
157 if (Jac)
158 {
160 indices_this_component, indices_all_components_for_Jac};
161 Jac->add(rci, local_Jac);
162 }
163 }
#define OGS_FATAL(...)
Definition: Error.h:26
int add(IndexType row, IndexType col, double val)
Definition: EigenMatrix.h:87
void add(IndexType rowId, double v)
add entry
Definition: EigenVector.h:79
std::size_t getID() const
Get id of the mesh.
Definition: Mesh.h:115
Eigen::Map< const Vector > toVector(std::vector< double > const &data, Eigen::VectorXd::Index size)
Creates an Eigen mapped vector from the given data vector.
static const double t
std::vector< GlobalIndexType > getIndices(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table)
NumLib::TemplateIsoparametric< ShapeFunction, ShapeMatricesType > createIsoparametricFiniteElement(MeshLib::Element const &e)
Eigen::MatrixXd collectDofsToMatrix(MeshLib::Element const &element, std::size_t const mesh_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x)
void interpolate(Eigen::MatrixXd const &primary_variables_mat, std::vector< std::reference_wrapper< ProcessVariable > > const &pv_refs, NsAndWeight const &ns_and_weight, Eigen::Ref< Eigen::VectorXd > interpolated_primary_variables)
void assembleLocalJacobian(std::vector< double > const &dFlux, typename Traits::NsAndWeight const &ns_and_weight, Eigen::MatrixXd &local_Jac) const
std::array< double, 3 > interpolateCoords(typename Traits::NsAndWeight const &ns_and_weight, NumLib::TemplateIsoparametric< ShapeFunction, ShapeMatrixPolicy > const &fe) const
void assembleLocalRhs(double const flux, typename Traits::NsAndWeight const &ns_and_weight, Eigen::VectorXd &local_rhs) const
MeshLib::Mesh const & bc_or_st_mesh
The domain, where this BC or ST will be applied.
Definition: BcOrStData.h:46
std::vector< std::reference_wrapper< ProcessVariable > > const & all_process_variables_for_this_process
Definition: BcOrStData.h:51

References MathLib::EigenMatrix::add(), MathLib::EigenVector::add(), ProcessLib::BoundaryConditionAndSourceTerm::Python::BcOrStData< BcOrStPythonSideInterface >::all_process_variables_for_this_process, ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::assembleLocalJacobian(), ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::assembleLocalRhs(), ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::bc_or_st_data, ProcessLib::BoundaryConditionAndSourceTerm::Python::BcOrStData< BcOrStPythonSideInterface >::bc_or_st_mesh, ProcessLib::BoundaryConditionAndSourceTerm::Python::collectDofsToMatrix(), NumLib::createIsoparametricFiniteElement(), ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::element, MeshLib::Mesh::getID(), NumLib::getIndices(), NumLib::LocalToGlobalIndexMap::getNumberOfGlobalComponents(), ProcessLib::BoundaryConditionAndSourceTerm::Python::BcOrStData< BcOrStPythonSideInterface >::global_component_id, ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::integration_method, ProcessLib::BoundaryConditionAndSourceTerm::Python::interpolate(), ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::interpolateCoords(), ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::nss_and_weights, OGS_FATAL, NumLib::LocalToGlobalIndexMap::size(), MathLib::t, and MathLib::toVector().

Referenced by ProcessLib::PythonBoundaryConditionLocalAssembler< ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::assemble(), and ProcessLib::SourceTerms::Python::PythonSourceTermLocalAssembler< ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::assemble().

◆ assembleLocalJacobian()

template<typename BcOrStData , typename ShapeFunction , typename LowerOrderShapeFunction , typename IntegrationMethod , int GlobalDim>
void ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::assembleLocalJacobian ( std::vector< double > const &  dFlux,
typename Traits::NsAndWeight const &  ns_and_weight,
Eigen::MatrixXd &  local_Jac 
) const
inlineprivate

Assembles the term \(N_\mathrm{this}^T \cdot \mathrm{dFlux}_c \cdot \mathrm{weight} \cdot N_c \) for each component \(c\) of each primary variable and adds it to the passed local Jacobian.

The shape functions \(N_\alpha\) of different primary variables \(\alpha\) might differ, e.g., for Taylor-Hood elements.

The index this refers to the primary variable to which this boundary condition or source term belongs. The index \(c\) (called other in the source code in the method body) runs over all components of all primary variables of the ProcessLib::Process to which this boundary condition or source term belongs.

Definition at line 204 of file BcAndStLocalAssemblerImpl.h.

208 {
209 auto const& pv_refs =
211 auto const weight = ns_and_weight.weight();
212
213 auto const this_shp_fct_order = bc_or_st_data.shape_function_order;
214 auto const N_this = ns_and_weight.N(this_shp_fct_order);
215
216 Eigen::Index left = 0;
217
218 for (auto pv_ref : pv_refs)
219 {
220 auto const& pv = pv_ref.get();
221 auto const other_num_comp = pv.getNumberOfGlobalComponents();
222 auto const other_shp_fct_order = pv.getShapeFunctionOrder();
223 auto const N_other = ns_and_weight.N(other_shp_fct_order);
224 auto const width = N_other.size();
225
226 for (decltype(+other_num_comp) other_comp = 0;
227 other_comp < other_num_comp;
228 ++other_comp)
229 {
230 // The assignment -= takes into account the sign convention
231 // of 1st-order in time in ODE systems in OpenGeoSys.
232 local_Jac.middleCols(left, width).noalias() -=
233 N_this.transpose() * (dFlux[other_comp] * weight) * N_other;
234
235 left += width;
236 }
237 }
238
239 assert(left == local_Jac.cols());
240 }

References ProcessLib::BoundaryConditionAndSourceTerm::Python::BcOrStData< BcOrStPythonSideInterface >::all_process_variables_for_this_process, ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::bc_or_st_data, ProcessLib::BoundaryConditionAndSourceTerm::Python::NsAndWeight< ShapeFunction, LowerOrderShapeFunction, ShapeMatrix, LowerOrderShapeMatrix >::N(), ProcessLib::BoundaryConditionAndSourceTerm::Python::BcOrStData< BcOrStPythonSideInterface >::shape_function_order, and ProcessLib::BoundaryConditionAndSourceTerm::Python::NsAndWeight< ShapeFunction, LowerOrderShapeFunction, ShapeMatrix, LowerOrderShapeMatrix >::weight().

Referenced by ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::assemble().

◆ assembleLocalRhs()

template<typename BcOrStData , typename ShapeFunction , typename LowerOrderShapeFunction , typename IntegrationMethod , int GlobalDim>
void ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::assembleLocalRhs ( double const  flux,
typename Traits::NsAndWeight const &  ns_and_weight,
Eigen::VectorXd &  local_rhs 
) const
inlineprivate

◆ interpolateCoords()

template<typename BcOrStData , typename ShapeFunction , typename LowerOrderShapeFunction , typename IntegrationMethod , int GlobalDim>
std::array< double, 3 > ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::interpolateCoords ( typename Traits::NsAndWeight const &  ns_and_weight,
NumLib::TemplateIsoparametric< ShapeFunction, ShapeMatrixPolicy > const &  fe 
) const
inlineprivate

Determines the coordinates that the point associated with the passed shape matrix by interpolating the element's node coordinates

Note
We use the shape function of higher order to interpolate coordinates.

Definition at line 171 of file BcAndStLocalAssemblerImpl.h.

175 {
176 auto const& N_higher = ns_and_weight.NHigherOrder();
177 return fe.interpolateCoordinates(N_higher);
178 }
std::array< double, 3 > interpolateCoordinates(typename ShapeMatrices::ShapeType const &N) const
Interpolates the coordinates of the element with the given shape matrix.

References NumLib::TemplateIsoparametric< ShapeFunctionType_, ShapeMatrixTypes_ >::interpolateCoordinates(), and ProcessLib::BoundaryConditionAndSourceTerm::Python::NsAndWeight< ShapeFunction, LowerOrderShapeFunction, ShapeMatrix, LowerOrderShapeMatrix >::NHigherOrder().

Referenced by ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::assemble().

Member Data Documentation

◆ bc_or_st_data

◆ element

◆ integration_method

template<typename BcOrStData , typename ShapeFunction , typename LowerOrderShapeFunction , typename IntegrationMethod , int GlobalDim>
IntegrationMethod const ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::integration_method

◆ nss_and_weights

template<typename BcOrStData , typename ShapeFunction , typename LowerOrderShapeFunction , typename IntegrationMethod , int GlobalDim>
std::vector<typename Traits::NsAndWeight> const ProcessLib::BoundaryConditionAndSourceTerm::Python::BcAndStLocalAssemblerImpl< BcOrStData, ShapeFunction, LowerOrderShapeFunction, IntegrationMethod, GlobalDim >::nss_and_weights

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