v0.14.0
FormsIntegrators.hpp
Go to the documentation of this file.
1 /** \file FormsIntegrators.hpp
2  * \brief Forms integrators
3  * \ingroup mofem_form
4 
5 */
6 
7 #ifndef __FORMS_INTEGRATORS_HPP__
8 #define __FORMS_INTEGRATORS_HPP__
9 
10 namespace MoFEM {
11 
12 //! [Storage and set boundary conditions]
13 
15  EssentialBcStorage(VectorInt &indices) : entityIndices(indices) {}
17  /**
18  * @brief Store modifed indices by field
19  *
20  * Hash map, key is field name, value is storage.
21  *
22  */
23  using HashVectorStorage =
24  map<std::string, std::vector<boost::shared_ptr<EssentialBcStorage>>>;
26 };
27 
28 /**
29  * @brief Set indices on entities on finite element
30  * @ingroup mofem_forms
31  *
32  * If index is marked, its value is set to -1. DOF with such index is not
33  * assembled into the system.
34  *
35  * Indices are stored on the entity.
36  *
37  */
39  OpSetBc(std::string field_name, bool yes_set,
40  boost::shared_ptr<std::vector<unsigned char>> boundary_marker);
41  MoFEMErrorCode doWork(int side, EntityType type,
43 
44 public:
45  bool yesSet;
46  boost::shared_ptr<std::vector<unsigned char>> boundaryMarker;
47 };
48 
50  OpUnSetBc(std::string field_name);
51  MoFEMErrorCode doWork(int side, EntityType type,
53 };
54 
55 /**
56  * @brief Set values to vector in operator
57  * @ingroup mofem_forms
58  *
59  * MoFEM::FieldEntity provides MoFEM::FieldEntity::getWeakStoragePtr() storage
60  * function which allows to transfer data between FEs or operators processing
61  * the same entities.
62  *
63  * When MoFEM::OpSetBc is pushed in weak storage indices taking in account
64  * indices which are skip to take boundary conditions are stored. Those entities
65  * are used by VecSetValues.
66  *
67  * @param V
68  * @param data
69  * @param ptr
70  * @param iora
71  * @return MoFEMErrorCode
72  */
73 template <>
76  const double *ptr, InsertMode iora);
77 
78 /**
79  * @brief Set values to matrix in operator
80  *
81  * See MoFEM::VecSetValues<EssentialBcStorage> for explanation.
82  *
83  * @param M
84  * @param row_data
85  * @param col_data
86  * @param ptr
87  * @param iora
88  * @return MoFEMErrorCode
89  */
90 template <>
93  const EntitiesFieldData::EntData &row_data,
94  const EntitiesFieldData::EntData &col_data,
95  const double *ptr, InsertMode iora);
96 
97 //! [Storage and set boundary conditions]
98 
99 /**
100  * @brief Form integrator assembly types
101  * @ingroup mofem_forms
102  *
103  */
112 };
113 
114 template <int A> struct AssemblyTypeSelector {};
115 
116 template <>
117 inline MoFEMErrorCode MatSetValues<AssemblyTypeSelector<PETSC>>(
118  Mat M, const EntitiesFieldData::EntData &row_data,
119  const EntitiesFieldData::EntData &col_data, const double *ptr,
120  InsertMode iora) {
121  return MatSetValues<EssentialBcStorage>(M, row_data, col_data, ptr, iora);
122 }
123 
124 template <>
125 inline MoFEMErrorCode VecSetValues<AssemblyTypeSelector<PETSC>>(
126  Vec V, const EntitiesFieldData::EntData &data, const double *ptr,
127  InsertMode iora) {
128  return VecSetValues<EssentialBcStorage>(V, data, ptr, iora);
129 }
130 
131 /**
132  * @brief Form integrator integration types
133  * @ingroup mofem_forms
134  *
135  */
137 
138 /**
139  * @brief Scalar function type
140  * @ingroup mofem_forms
141  *
142  */
143 using ScalarFun =
144  boost::function<double(const double, const double, const double)>;
145 
146 inline double scalar_fun_one(const double, const double, const double) {
147  return 1;
148 }
149 
150 /**
151  * @brief Lambda function used to scale with time
152  *
153  */
154 using TimeFun = boost::function<double(double)>;
155 
156 /**
157  * @brief Lambda function used to scale with time
158  *
159  */
160 using FEFun = boost::function<double(const FEMethod *fe_ptr)>;
161 
162 /**
163  * @brief Constant function type
164  *
165  */
166 using ConstantFun = boost::function<double()>;
167 
168 /**
169  * @brief Vector function type
170  * @ingroup mofem_forms
171  *
172  * @tparam DIM dimension of the return
173  */
174 template <int DIM>
175 using VectorFun = boost::function<FTensor::Tensor1<double, DIM>(
176  const double, const double, const double)>;
177 
178 template <AssemblyType A, typename EleOp> struct OpBaseImpl : public EleOp {
179  using OpType = typename EleOp::OpType;
181 
182  OpBaseImpl(const std::string row_field_name, const std::string col_field_name,
183  const OpType type, boost::shared_ptr<Range> ents_ptr = nullptr)
184  : EleOp(row_field_name, col_field_name, type, false),
185  assembleTranspose(false), onlyTranspose(false), entsPtr(ents_ptr) {}
186 
187  /**
188  * \brief Do calculations for the left hand side
189  * @param row_side row side number (local number) of entity on element
190  * @param col_side column side number (local number) of entity on element
191  * @param row_type type of row entity
192  * @param col_type type of column entity
193  * @param row_data data for row
194  * @param col_data data for column
195  * @return error code
196  */
197  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
198  EntityType col_type, EntData &row_data,
199  EntData &col_data);
200 
201  /**
202  * @brief Do calculations for the right hand side
203  *
204  * @param row_side
205  * @param row_type
206  * @param row_data
207  * @return MoFEMErrorCode
208  */
209  MoFEMErrorCode doWork(int row_side, EntityType row_type, EntData &row_data);
210 
211  TimeFun timeScalingFun; ///< assumes that time variable is set
212  FEFun feScalingFun; ///< assumes that time variable is set
213  boost::shared_ptr<Range> entsPtr; ///< Entities on which element is run
214 
215  using MatSetValuesHook = boost::function<MoFEMErrorCode(
217  const EntitiesFieldData::EntData &row_data,
219 
221 
222 protected:
223  template <int DIM>
225  return getFTensor1FromArray<DIM, DIM>(locF);
226  }
227 
228  template <int DIM>
230  getLocMat(const int rr) {
231  return getFTensor2FromArray<DIM, DIM, DIM>(locMat, rr);
232  }
233 
234  int nbRows; ///< number of dofs on rows
235  int nbCols; ///< number if dof on column
236  int nbIntegrationPts; ///< number of integration points
237  int nbRowBaseFunctions; ///< number or row base functions
238 
239  int rowSide; ///< row side number
240  int colSide; ///< column side number
241  EntityType rowType; ///< row type
242  EntityType colType; ///< column type
243 
246 
247  MatrixDouble locMat; ///< local entity block matrix
248  MatrixDouble locMatTranspose; ///< local entity block matrix
249  VectorDouble locF; ///< local entity vector
250 
251  /**
252  * \brief Integrate grad-grad operator
253  * @param row_data row data (consist base functions on row entity)
254  * @param col_data column data (consist base functions on column entity)
255  * @return error code
256  */
257  virtual MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data) {
258  return MOFEM_NOT_IMPLEMENTED;
259  }
260 
261  virtual MoFEMErrorCode aSsemble(EntData &row_data, EntData &col_data,
262  const bool trans);
263 
264  /**
265  * \brief Class dedicated to integrate operator
266  * @param data entity data on element row
267  * @return error code
268  */
270  return MOFEM_NOT_IMPLEMENTED;
271  }
272 
273  virtual MoFEMErrorCode aSsemble(EntData &data);
274 
275  /** \brief Get number of base functions
276  *
277  * @param data
278  * @return number of base functions
279  */
280  virtual size_t getNbOfBaseFunctions(EntitiesFieldData::EntData &data);
281 };
282 
283 template <AssemblyType A, typename EleOp>
287  const EntitiesFieldData::EntData &row_data,
288  const EntitiesFieldData::EntData &col_data, MatrixDouble &m) {
289  return MatSetValues<AssemblyTypeSelector<A>>(
290  op_ptr->getKSPB(), row_data, col_data, m, ADD_VALUES);
291  };
292 
293 /**
294  * @brief Integrator forms
295  * @ingroup mofem_forms
296  *
297  * @tparam EleOp
298  */
299 template <typename EleOp> struct FormsIntegrators {
300 
302  using OpType = typename EleOp::OpType;
303 
304  /**
305  * @brief Assembly methods
306  * @ingroup mofem_forms
307  *
308  * @tparam A
309  */
310  template <AssemblyType A> struct Assembly {
311 
313 
314  /**
315  * @brief Linear form
316  * @ingroup mofem_forms
317  *
318  * @tparam I
319  */
320  template <IntegrationType I> struct LinearForm;
321 
322  /**
323  * @brief Bi linear form
324  * @ingroup mofem_forms
325  *
326  * @tparam I
327  */
328  template <IntegrationType I> struct BiLinearForm;
329 
330  /**
331  * @brief Tri linear form
332  * @ingroup mofem_forms
333  *
334  * @tparam I
335  */
336  template <IntegrationType I> struct TriLinearForm;
337 
338  }; // Assembly
339 }; // namespace MoFEM
340 
341 template <AssemblyType A, typename EleOp>
342 size_t
344  auto nb_base_functions = data.getN().size2();
345  if (data.getBase() != USER_BASE) {
346  switch (data.getSpace()) {
347  case NOSPACE:
348  break;
349  case NOFIELD:
350  break;
351  case H1:
352  break;
353  case HCURL:
354  case HDIV:
355  nb_base_functions /= 3;
356 #ifndef NDEBUG
357  if (data.getN().size2() % 3) {
358  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
359  "Number of base functions is not divisible by 3");
360  }
361 #endif
362  break;
363  case L2:
364  break;
365  default:
366  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
367  "Space %s not implemented", FieldSpaceNames[data.getSpace()]);
368  }
369  }
370  return nb_base_functions;
371 }
372 
373 template <AssemblyType A, typename EleOp>
375 OpBaseImpl<A, EleOp>::doWork(int row_side, int col_side, EntityType row_type,
376  EntityType col_type,
377  EntitiesFieldData::EntData &row_data,
378  EntitiesFieldData::EntData &col_data) {
380 
381  if (entsPtr) {
382  if (entsPtr->find(this->getFEEntityHandle()) == entsPtr->end())
384  }
385 
386  // get number of dofs on row
387  nbRows = row_data.getIndices().size();
388  // if no dofs on row, exit that work, nothing to do here
389  if (!nbRows)
391  rowSide = row_side;
392  rowType = row_type;
393  // get number of dofs on column
394  nbCols = col_data.getIndices().size();
395  // if no dofs on column, exit nothing to do here
396  if (!nbCols)
398  colSide = col_side;
399  colType = col_type;
400  // get number of integration points
401  nbIntegrationPts = EleOp::getGaussPts().size2();
402  // get row base functions
403  nbRowBaseFunctions = getNbOfBaseFunctions(row_data);
404 
405  // set size of local entity bock
406  locMat.resize(nbRows, nbCols, false);
407  // clear matrix
408  locMat.clear();
409  // integrate local matrix for entity block
410  CHKERR this->iNtegrate(row_data, col_data);
411 
412  // assemble local matrix
413  auto check_if_assemble_transpose = [&] {
414  if (this->sYmm) {
415  if (row_side != col_side || row_type != col_type)
416  return true;
417  else
418  return false;
419  } else if (assembleTranspose) {
420  return true;
421  }
422 
423  return false;
424  };
425  CHKERR aSsemble(row_data, col_data, check_if_assemble_transpose());
427 }
428 
429 template <AssemblyType A, typename EleOp>
430 MoFEMErrorCode OpBaseImpl<A, EleOp>::doWork(int row_side, EntityType row_type,
431  EntData &row_data) {
433 
434  if (entsPtr) {
435  if (entsPtr->find(this->getFEEntityHandle()) == entsPtr->end())
437  }
438 
439  // get number of dofs on row
440  nbRows = row_data.getIndices().size();
441  rowSide = row_side;
442  rowType = row_type;
443 
444  if (!nbRows)
446  // get number of integration points
447  nbIntegrationPts = EleOp::getGaussPts().size2();
448  // get row base functions
449  nbRowBaseFunctions = getNbOfBaseFunctions(row_data);
450  // resize and clear the right hand side vector
451  locF.resize(nbRows);
452  locF.clear();
453  // integrate local vector
454  CHKERR this->iNtegrate(row_data);
455  // assemble local vector
456  CHKERR this->aSsemble(row_data);
458 }
459 
460 template <AssemblyType A, typename EleOp>
462  EntData &col_data,
463  const bool transpose) {
465 
466  if (!this->timeScalingFun.empty())
467  this->locMat *= this->timeScalingFun(this->getFEMethod()->ts_t);
468  if (!this->feScalingFun.empty())
469  this->locMat *= this->feScalingFun(this->getFEMethod());
470 
471  // Assemble transpose
472  if (transpose) {
473  this->locMatTranspose.resize(this->locMat.size2(), this->locMat.size1(),
474  false);
475  noalias(this->locMatTranspose) = trans(this->locMat);
476  CHKERR matSetValuesHook(this, col_data, row_data, this->locMatTranspose);
477  }
478 
479  if (!this->onlyTranspose) {
480  // assemble local matrix
481  CHKERR matSetValuesHook(this, row_data, col_data, this->locMat);
482  }
483 
485 }
486 
487 template <AssemblyType A, typename EleOp>
489  if (!this->timeScalingFun.empty())
490  this->locF *= this->timeScalingFun(this->getFEMethod()->ts_t);
491  if (!this->feScalingFun.empty())
492  this->locF *= this->feScalingFun(this->getFEMethod());
493 
494  return VecSetValues<AssemblyTypeSelector<A>>(this->getKSPf(), data,
495  this->locF, ADD_VALUES);
496 }
497 
498 } // namespace MoFEM
499 
500 /**
501  * \defgroup mofem_forms Forms Integrators
502  *
503  * \brief Classes and functions used to evaluate fields at integration pts,
504  *jacobians, etc..
505  *
506  * \ingroup mofem_forces_and_sources
507  **/
508 
509 #endif //__FORMS_INTEGRATORS_HPP__
MoFEM::OpBaseImpl::feScalingFun
FEFun feScalingFun
assumes that time variable is set
Definition: FormsIntegrators.hpp:212
NOSPACE
@ NOSPACE
Definition: definitions.h:83
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
MoFEM::OpSetBc::OpSetBc
OpSetBc(std::string field_name, bool yes_set, boost::shared_ptr< std::vector< unsigned char >> boundary_marker)
Definition: FormsIntegrators.cpp:16
MoFEM::OpBaseImpl::nbIntegrationPts
int nbIntegrationPts
number of integration points
Definition: FormsIntegrators.hpp:236
MoFEM::EntitiesFieldData::EntData
Data on single entity (This is passed as argument to DataOperator::doWork)
Definition: EntitiesFieldData.hpp:127
MoFEM::OpBaseImpl::iNtegrate
virtual MoFEMErrorCode iNtegrate(EntData &data)
Class dedicated to integrate operator.
Definition: FormsIntegrators.hpp:269
H1
@ H1
continuous field
Definition: definitions.h:85
MoFEM::OpBaseImpl::locMat
MatrixDouble locMat
local entity block matrix
Definition: FormsIntegrators.hpp:247
MoFEM::FormsIntegrators::Assembly::LinearForm
Linear form.
Definition: FormsIntegrators.hpp:320
FTensor::Tensor1
Definition: Tensor1_value.hpp:8
MoFEM::FEMethod
structure for User Loop Methods on finite elements
Definition: LoopMethods.hpp:369
MoFEM::BLOCK_PRECONDITIONER_SCHUR
@ BLOCK_PRECONDITIONER_SCHUR
Definition: FormsIntegrators.hpp:109
L2
@ L2
field with C-1 continuity
Definition: definitions.h:88
MoFEM::BLOCK_MAT
@ BLOCK_MAT
Definition: FormsIntegrators.hpp:107
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::Types::MatrixDouble
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
MoFEM::EssentialBcStorage::HashVectorStorage
map< std::string, std::vector< boost::shared_ptr< EssentialBcStorage > >> HashVectorStorage
Store modifed indices by field.
Definition: FormsIntegrators.hpp:24
MoFEM::PETSC
@ PETSC
Definition: FormsIntegrators.hpp:105
MoFEM::ConstantFun
boost::function< double()> ConstantFun
Constant function type.
Definition: FormsIntegrators.hpp:166
MoFEM::OpBaseImpl::assembleTranspose
bool assembleTranspose
Definition: FormsIntegrators.hpp:244
MoFEM::OpUnSetBc::OpUnSetBc
OpUnSetBc(std::string field_name)
Definition: FormsIntegrators.cpp:54
FTensor::Tensor2
Definition: Tensor2_value.hpp:16
USER_BASE
@ USER_BASE
user implemented approximation base
Definition: definitions.h:68
MoFEM::OpSetBc::boundaryMarker
boost::shared_ptr< std::vector< unsigned char > > boundaryMarker
Definition: FormsIntegrators.hpp:46
MoFEM::scalar_fun_one
double scalar_fun_one(const double, const double, const double)
Definition: FormsIntegrators.hpp:146
EleOp
MoFEM::OpBaseImpl::locF
VectorDouble locF
local entity vector
Definition: FormsIntegrators.hpp:249
PlasticOps::M
FTensor::Index< 'M', 3 > M
Definition: PlasticOps.hpp:117
MoFEM::OpBaseImpl
Definition: FormsIntegrators.hpp:178
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:535
MoFEM::TimeFun
boost::function< double(double)> TimeFun
Lambda function used to scale with time.
Definition: FormsIntegrators.hpp:154
MoFEM::OpBaseImpl::getLocMat
FTensor::Tensor2< FTensor::PackPtr< double *, DIM >, DIM, DIM > getLocMat(const int rr)
Definition: FormsIntegrators.hpp:230
MoFEM::OpBaseImpl::timeScalingFun
TimeFun timeScalingFun
assumes that time variable is set
Definition: FormsIntegrators.hpp:211
MoFEM::OpUnSetBc::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: FormsIntegrators.cpp:58
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
MoFEM::VectorFun
boost::function< FTensor::Tensor1< double, DIM >(const double, const double, const double)> VectorFun
Vector function type.
Definition: FormsIntegrators.hpp:176
MoFEM::EntityStorage
Definition: FieldEntsMultiIndices.hpp:16
MoFEM::MatSetValues< EssentialBcStorage >
MoFEMErrorCode MatSetValues< EssentialBcStorage >(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
Set values to matrix in operator.
Definition: FormsIntegrators.cpp:112
MoFEM::ForcesAndSourcesCore::UserDataOperator
Definition: ForcesAndSourcesCore.hpp:549
MoFEM::OpBaseImpl::getNf
FTensor::Tensor1< FTensor::PackPtr< double *, DIM >, DIM > getNf()
Definition: FormsIntegrators.hpp:224
double
convert.type
type
Definition: convert.py:64
MoFEM::OpBaseImpl::iNtegrate
virtual MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
Integrate grad-grad operator.
Definition: FormsIntegrators.hpp:257
MoFEM::FormsIntegrators::Assembly
Assembly methods.
Definition: FormsIntegrators.hpp:310
MoFEM::OpBaseImpl::doWork
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
Do calculations for the left hand side.
Definition: FormsIntegrators.hpp:375
MoFEM::EntitiesFieldData::EntData::getIndices
const VectorInt & getIndices() const
Get global indices of dofs on entity.
Definition: EntitiesFieldData.hpp:1201
MoFEM::FormsIntegrators::Assembly::BiLinearForm
Bi linear form.
Definition: FormsIntegrators.hpp:328
MoFEM::EntitiesFieldData::EntData::getBase
FieldApproximationBase & getBase()
Get approximation base.
Definition: EntitiesFieldData.hpp:1300
MoFEM::GAUSS
@ GAUSS
Definition: FormsIntegrators.hpp:136
MoFEM::EssentialBcStorage::feStorage
static HashVectorStorage feStorage
[Storage and set boundary conditions]
Definition: FormsIntegrators.hpp:25
MoFEM::OpBaseImpl::colSide
int colSide
column side number
Definition: FormsIntegrators.hpp:240
MoFEM::EssentialBcStorage::EssentialBcStorage
EssentialBcStorage(VectorInt &indices)
Definition: FormsIntegrators.hpp:15
MoFEM::OpBaseImpl::OpType
typename EleOp::OpType OpType
Definition: FormsIntegrators.hpp:179
MoFEM::ScalarFun
boost::function< double(const double, const double, const double)> ScalarFun
Scalar function type.
Definition: FormsIntegrators.hpp:144
MoFEM::USER_INTEGRATION
@ USER_INTEGRATION
Definition: FormsIntegrators.hpp:136
MoFEM::OpBaseImpl::locMatTranspose
MatrixDouble locMatTranspose
local entity block matrix
Definition: FormsIntegrators.hpp:248
MoFEM::OpBaseImpl::onlyTranspose
bool onlyTranspose
Definition: FormsIntegrators.hpp:245
MoFEM::FormsIntegrators::Assembly::TriLinearForm
Tri linear form.
Definition: FormsIntegrators.hpp:336
MoFEM::ForcesAndSourcesCore::UserDataOperator
friend class UserDataOperator
Definition: ForcesAndSourcesCore.hpp:482
MoFEM::OpBaseImpl::MatSetValuesHook
boost::function< MoFEMErrorCode(ForcesAndSourcesCore::UserDataOperator *op_ptr, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, MatrixDouble &m)> MatSetValuesHook
Definition: FormsIntegrators.hpp:218
MoFEM::AssemblyTypeSelector
Definition: FormsIntegrators.hpp:114
EntData
EntitiesFieldData::EntData EntData
Definition: child_and_parent.cpp:37
field_name
constexpr auto field_name
Definition: poisson_2d_homogeneous.cpp:13
MoFEM::OpSetBc::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: FormsIntegrators.cpp:22
MoFEM::IntegrationType
IntegrationType
Form integrator integration types.
Definition: FormsIntegrators.hpp:136
MoFEM::FormsIntegrators::OpType
typename EleOp::OpType OpType
Definition: FormsIntegrators.hpp:302
MoFEM::OpBaseImpl::nbRowBaseFunctions
int nbRowBaseFunctions
number or row base functions
Definition: FormsIntegrators.hpp:237
MoFEM::FEFun
boost::function< double(const FEMethod *fe_ptr)> FEFun
Lambda function used to scale with time.
Definition: FormsIntegrators.hpp:160
MoFEM::OpBaseImpl::rowSide
int rowSide
row side number
Definition: FormsIntegrators.hpp:239
MoFEM::OpUnSetBc
Definition: FormsIntegrators.hpp:49
MoFEM::OpBaseImpl::nbCols
int nbCols
number if dof on column
Definition: FormsIntegrators.hpp:235
MoFEM::EssentialBcStorage
[Storage and set boundary conditions]
Definition: FormsIntegrators.hpp:14
MoFEM::LAST_INTEGRATION
@ LAST_INTEGRATION
Definition: FormsIntegrators.hpp:136
FieldSpaceNames
const static char *const FieldSpaceNames[]
Definition: definitions.h:92
MoFEM::OpBaseImpl::nbRows
int nbRows
number of dofs on rows
Definition: FormsIntegrators.hpp:234
MoFEM::AssemblyType
AssemblyType
[Storage and set boundary conditions]
Definition: FormsIntegrators.hpp:104
MoFEM::OpBaseImpl::entsPtr
boost::shared_ptr< Range > entsPtr
Entities on which element is run.
Definition: FormsIntegrators.hpp:213
MoFEM::EntitiesFieldData::EntData::getN
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
Definition: EntitiesFieldData.hpp:1305
MoFEM::LAST_ASSEMBLE
@ LAST_ASSEMBLE
Definition: FormsIntegrators.hpp:111
MoFEM::OpBaseImpl::matSetValuesHook
static MatSetValuesHook matSetValuesHook
Definition: FormsIntegrators.hpp:220
MoFEM::USER_ASSEMBLE
@ USER_ASSEMBLE
Definition: FormsIntegrators.hpp:110
MoFEM::Types::VectorInt
UBlasVector< int > VectorInt
Definition: Types.hpp:67
MoFEM::OpSetBc::yesSet
bool yesSet
Definition: FormsIntegrators.hpp:45
UBlasVector< double >
MoFEM::OpSetBc
Set indices on entities on finite element.
Definition: FormsIntegrators.hpp:38
MoFEM::EntitiesFieldData::EntData::getSpace
FieldSpace & getSpace()
Get field space.
Definition: EntitiesFieldData.hpp:1302
OpBaseImpl
HCURL
@ HCURL
field with continuous tangents
Definition: definitions.h:86
EigenMatrix::Vec
const FTensor::Tensor2< T, Dim, Dim > Vec
Definition: MatrixFunction.hpp:66
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
MoFEM::OpBaseImpl::aSsemble
virtual MoFEMErrorCode aSsemble(EntData &row_data, EntData &col_data, const bool trans)
Definition: FormsIntegrators.hpp:461
MoFEM::FormsIntegrators
Integrator forms.
Definition: FormsIntegrators.hpp:299
MoFEM::BLOCK_SCHUR
@ BLOCK_SCHUR
Definition: FormsIntegrators.hpp:108
m
FTensor::Index< 'm', 3 > m
Definition: shallow_wave.cpp:80
MoFEM::OpBaseImpl::getNbOfBaseFunctions
virtual size_t getNbOfBaseFunctions(EntitiesFieldData::EntData &data)
Get number of base functions.
Definition: FormsIntegrators.hpp:343
MoFEM::OpBaseImpl::rowType
EntityType rowType
row type
Definition: FormsIntegrators.hpp:241
MoFEM::EssentialBcStorage::entityIndices
VectorInt entityIndices
Definition: FormsIntegrators.hpp:16
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
HDIV
@ HDIV
field with continuous normal traction
Definition: definitions.h:87
MOFEM_NOT_IMPLEMENTED
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
MoFEM::SCHUR
@ SCHUR
Definition: FormsIntegrators.hpp:106
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
MoFEM::OpBaseImpl::colType
EntityType colType
column type
Definition: FormsIntegrators.hpp:242
MoFEM::VecSetValues< EssentialBcStorage >
MoFEMErrorCode VecSetValues< EssentialBcStorage >(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Set values to vector in operator.
Definition: FormsIntegrators.cpp:76
MoFEM::OpBaseImpl::OpBaseImpl
OpBaseImpl(const std::string row_field_name, const std::string col_field_name, const OpType type, boost::shared_ptr< Range > ents_ptr=nullptr)
Definition: FormsIntegrators.hpp:182
NOFIELD
@ NOFIELD
scalar or vector of scalars describe (no true field)
Definition: definitions.h:84