v0.15.0
Loading...
Searching...
No Matches
OpAssembleP Struct Reference

Assemble P *. More...

#include "tutorials/cor-7/src/ElasticityMixedFormulation.hpp"

Inheritance diagram for OpAssembleP:
[legend]
Collaboration diagram for OpAssembleP:
[legend]

Public Member Functions

 OpAssembleP (DataAtIntegrationPts &common_data, BlockData &data)
 
PetscErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
 

Public Attributes

DataAtIntegrationPtscommonData
 
MatrixDouble locP
 
MatrixDouble translocP
 
BlockDatadAta
 

Detailed Description

Assemble P *.

  • \[ {\bf{P}} = - \int\limits_\Omega {{\bf{N}}_p^T\frac{1}{\lambda }{{\bf{N}}_p}d\Omega } \]

Examples
ElasticityMixedFormulation.hpp, and elasticity_mixed_formulation.cpp.

Definition at line 113 of file ElasticityMixedFormulation.hpp.

Constructor & Destructor Documentation

◆ OpAssembleP()

OpAssembleP::OpAssembleP ( DataAtIntegrationPts & common_data,
BlockData & data )
inline
Examples
ElasticityMixedFormulation.hpp.

Definition at line 121 of file ElasticityMixedFormulation.hpp.

122 : VolumeElementForcesAndSourcesCore::UserDataOperator(
123 "P", "P", UserDataOperator::OPROWCOL),
124 commonData(common_data), dAta(data) {
125 sYmm = true;
126 }
DataAtIntegrationPts & commonData

Member Function Documentation

◆ doWork()

PetscErrorCode OpAssembleP::doWork ( int row_side,
int col_side,
EntityType row_type,
EntityType col_type,
EntitiesFieldData::EntData & row_data,
EntitiesFieldData::EntData & col_data )
inline
Examples
ElasticityMixedFormulation.hpp.

Definition at line 128 of file ElasticityMixedFormulation.hpp.

131 {
132
134 const int row_nb_dofs = row_data.getIndices().size();
135 if (!row_nb_dofs)
137 const int col_nb_dofs = col_data.getIndices().size();
138 if (!col_nb_dofs)
140
141 if (dAta.tEts.find(getNumeredEntFiniteElementPtr()->getEnt()) ==
142 dAta.tEts.end()) {
144 }
146 // Set size can clear local tangent matrix
147 locP.resize(row_nb_dofs, col_nb_dofs, false);
148 locP.clear();
149
150 const int row_nb_gauss_pts = row_data.getN().size1();
151 if (!row_nb_gauss_pts)
153 const int row_nb_base_functions = row_data.getN().size2();
154 auto row_base_functions = row_data.getFTensor0N();
155
156 // get data
157 FTensor::Index<'i', 3> i;
158 FTensor::Index<'j', 3> j;
159 const double lambda = commonData.lAmbda;
160
161 double coefficient = commonData.pOisson == 0.5 ? 0. : 1 / lambda;
162
163 // integration
164 if (coefficient != 0.) {
165 for (int gg = 0; gg != row_nb_gauss_pts; gg++) {
166
167 // Get volume and integration weight
168 double w = getVolume() * getGaussPts()(3, gg);
169
170 // INTEGRATION
171 int row_bb = 0;
172 for (; row_bb != row_nb_dofs; row_bb++) {
173 auto col_base_functions = col_data.getFTensor0N(gg, 0);
174 for (int col_bb = 0; col_bb != col_nb_dofs; col_bb++) {
175
176 locP(row_bb, col_bb) -=
177 w * row_base_functions * col_base_functions * coefficient;
178
179 ++col_base_functions;
180 }
181 ++row_base_functions;
182 }
183 for (; row_bb != row_nb_base_functions; row_bb++) {
184 ++row_base_functions;
185 }
186 }
187 }
188
189
191 getFEMethod()->ksp_B, row_nb_dofs, &*row_data.getIndices().begin(),
192 col_nb_dofs, &*col_data.getIndices().begin(), &locP(0, 0), ADD_VALUES);
193
194 // is symmetric
195 if (row_side != col_side || row_type != col_type) {
196 translocP.resize(col_nb_dofs, row_nb_dofs, false);
197 noalias(translocP) = trans(locP);
198 CHKERR MatSetValues(getFEMethod()->ksp_B, col_nb_dofs,
199 &*col_data.getIndices().begin(), row_nb_dofs,
200 &*row_data.getIndices().begin(), &translocP(0, 0),
201 ADD_VALUES);
202 }
203
205 }
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
FTensor::Index< 'i', SPACE_DIM > i
static double lambda
FTensor::Index< 'j', 3 > j
MoFEMErrorCode MatSetValues(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
MoFEMErrorCode getBlockData(BlockData &data)

Member Data Documentation

◆ commonData

DataAtIntegrationPts& OpAssembleP::commonData

◆ dAta

BlockData& OpAssembleP::dAta

◆ locP

MatrixDouble OpAssembleP::locP

◆ translocP

MatrixDouble OpAssembleP::translocP

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