v0.15.0
Loading...
Searching...
No Matches
FormsBrokenSpaceConstraintImpl.hpp
Go to the documentation of this file.
1/**
2 * @file FormsBrokenSpaceConstraintImpl.hpp
3 * @brief Integrator for broken space constraints
4 * @date 2024-07-01
5 *
6 * @copyright Copyright (c) 2024
7 *
8 */
9
10#ifndef __FORMSBROKENSPACECONSTRAINTIMPL_HPP__
11#define __FORMSBROKENSPACECONSTRAINTIMPL_HPP__
12
13namespace MoFEM {
14
15/**
16 * @brief Operator for broken loop side
17 *
18 * That is used follow pattern,
19 *
20 * First: Iterate over skeleton FEs adjacent to Domain FEs
21 * Note: BoundaryEle, i.e. uses skeleton interation rule
22 *
23 * ans then
24 *
25 * Iterate over domain FEs adjacent to skelton, particularly one
26 * domain element.
27 *
28 * You use this to integrate over skelton, but field which is on adjacent domain
29 * element.
30 *
31 * @tparam E
32 */
33template <typename E> struct OpBrokenLoopSide : public OpLoopSide<E> {
34
37
38 MoFEMErrorCode doWork(int side, EntityType type,
41
42 auto prev_side_fe_ptr = OP::getSidePtrFE();
43 if (OP::sideFEName == prev_side_fe_ptr->getFEName()) {
44 auto prev_fe_uid =
45 prev_side_fe_ptr->numeredEntFiniteElementPtr->getFEUId();
47 CHKERR OP::sideFEPtr->setSideFEPtr(OP::getPtrFE());
48 CHKERR OP::sideFEPtr->copyBasicMethod(*OP::getFEMethod());
49 CHKERR OP::sideFEPtr->copyPetscData(*OP::getFEMethod());
53 OP::sideFEPtr->cacheWeakPtr = prev_side_fe_ptr->cacheWeakPtr;
54 OP::sideFEPtr->loopSize = 1;
55 CHKERR OP::sideFEPtr->preProcess();
56 OP::sideFEPtr->nInTheLoop = 0;
57 OP::sideFEPtr->numeredEntFiniteElementPtr =
58 prev_side_fe_ptr->numeredEntFiniteElementPtr;
60 CHKERR OP::sideFEPtr->postProcess();
61 } else {
62 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
63 "sideFEName is different");
64 }
65
67 };
68};
69
72 inline auto &getSense() { return eleSense; }
73 inline auto &getSide() { return eleSide; }
74 inline auto &getType() { return eleType; }
75 inline auto &getData() { return entData; }
76 inline auto &getFlux() { return fluxMat; }
77
78private:
79 int eleSense = 0;
80 int eleSide = 1;
81 EntityType eleType = MBENTITYSET;
84};
85
86template <typename OpBase> struct OpGetBrokenBaseSideData : public OpBase {
87
88 using OP = OpBase;
89
91 const std::string field_name,
92 boost::shared_ptr<std::vector<BrokenBaseSideData>> broken_base_side_data);
93
94 MoFEMErrorCode doWork(int row_side, EntityType row_type,
96
97private:
98 boost::shared_ptr<std::vector<BrokenBaseSideData>> brokenBaseSideData;
99};
100
101template <typename OpBase>
103 const std::string field_name,
104 boost::shared_ptr<std::vector<BrokenBaseSideData>> broken_base_side_data)
105 : OP(field_name, field_name, OP::OPROW),
106 brokenBaseSideData(broken_base_side_data) {}
107
108template <typename OpBase>
110OpGetBrokenBaseSideData<OpBase>::doWork(int row_side, EntityType row_type,
111 EntitiesFieldData::EntData &row_data) {
113
114 if (row_data.getIndices().size() == 0)
116
117 brokenBaseSideData->resize(OP::getLoopSize());
118
119 const auto n_in_the_loop = OP::getNinTheLoop();
120 const auto face_sense = OP::getSkeletonSense();
121
122#ifndef NDEBUG
123 if (face_sense != -1 && face_sense != 1)
124 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "face sense not set");
125#endif // NDEBUG
126
127 auto set_data = [&](auto &side_data) {
128 side_data.getSide() = row_side;
129 side_data.getType() = row_type;
130 side_data.getSense() = face_sense;
131 side_data.getData().sEnse = row_data.sEnse;
132 side_data.getData().sPace = row_data.sPace;
133 side_data.getData().bAse = row_data.bAse;
134 side_data.getData().iNdices = row_data.iNdices;
135 side_data.getData().localIndices = row_data.localIndices;
136 side_data.getData().dOfs = row_data.dOfs;
137 side_data.getData().fieldEntities = row_data.fieldEntities;
138 side_data.getData().fieldData = row_data.fieldData;
139 };
140
141 auto set_base = [&](auto &side_data) {
142 auto base = side_data.getData().getBase();
143 for (auto dd = 0; dd != BaseDerivatives::LastDerivative; ++dd) {
144 side_data.getData().baseFunctionsAndBaseDerivatives[dd][base].reset();
145 if (auto base_ptr = row_data.baseFunctionsAndBaseDerivatives[dd][base]) {
146 side_data.getData().baseFunctionsAndBaseDerivatives[dd][base] =
147 boost::make_shared<MatrixDouble>(*base_ptr);
148 }
149 }
150 };
151
152 set_data((*brokenBaseSideData)[n_in_the_loop]);
153 set_base((*brokenBaseSideData)[n_in_the_loop]);
154
156}
157
158template <typename OpBase> struct OpSetFlux : public OpBase {
159
160 using OP = OpBase;
161
162 OpSetFlux(
163 boost::shared_ptr<std::vector<BrokenBaseSideData>> broken_base_side_data,
164 boost::shared_ptr<MatrixDouble> flux_ptr);
165
166 MoFEMErrorCode doWork(int row_side, EntityType row_type,
168
169private:
170 boost::shared_ptr<std::vector<BrokenBaseSideData>> brokenBaseSideData;
171 boost::shared_ptr<MatrixDouble> fluxPtr;
172};
173
174template <typename OpBase>
176 boost::shared_ptr<std::vector<BrokenBaseSideData>> broken_base_side_data,
177 boost::shared_ptr<MatrixDouble> flux_ptr)
178 : OP(NOSPACE, OP::OPSPACE), brokenBaseSideData(broken_base_side_data),
179 fluxPtr(flux_ptr) {}
180
181template <typename OpBase>
185 auto swap_flux = [&](auto &side_data) { side_data.getFlux().swap(*fluxPtr); };
186 swap_flux((*brokenBaseSideData)[OP::getNinTheLoop()]);
188}
189
190template <typename OpBase> struct OpBrokenBaseImpl : public OpBase {
191
192 using OP = OpBase;
193
195 boost::shared_ptr<std::vector<BrokenBaseSideData>> broken_base_side_data,
196 boost::shared_ptr<Range> ents_ptr = nullptr)
197 : OP(NOSPACE, OP::OPSPACE), brokenBaseSideData(broken_base_side_data) {
198 OP::entsPtr = ents_ptr;
199 }
200
202 const std::string row_field,
203 boost::shared_ptr<std::vector<BrokenBaseSideData>> broken_base_side_data,
204 const bool assmb_transpose, const bool only_transpose,
205 boost::shared_ptr<Range> ents_ptr = nullptr)
206 : OP(row_field, row_field, OP::OPROW, ents_ptr),
207 brokenBaseSideData(broken_base_side_data) {
208 OP::entsPtr = ents_ptr;
209 OP::assembleTranspose = assmb_transpose;
210 OP::onlyTranspose = only_transpose;
211 OP::sYmm = false;
212 }
213
214 MoFEMErrorCode doWork(int row_side, EntityType row_type,
216
217protected:
218 boost::shared_ptr<std::vector<BrokenBaseSideData>> brokenBaseSideData;
219};
220
221template <typename OpBase>
223OpBrokenBaseImpl<OpBase>::doWork(int row_side, EntityType row_type,
224 EntitiesFieldData::EntData &row_data) {
226
227 if (OP::entsPtr) {
228 if (OP::entsPtr->find(this->getFEEntityHandle()) == OP::entsPtr->end())
230 }
231
232#ifndef NDEBUG
233 if (!brokenBaseSideData) {
234 SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSSIBLE_CASE, "space not set");
235 }
236#endif // NDEBUG
237
238 auto do_work_rhs = [this](int row_side, EntityType row_type,
239 EntitiesFieldData::EntData &row_data, int sense) {
241 // get number of dofs on row
242 OP::nbRows = row_data.getIndices().size();
243 if (!OP::nbRows)
245 // get number of integration points
246 OP::nbIntegrationPts = OP::getGaussPts().size2();
247 // get row base functions
248 OP::nbRowBaseFunctions = OP::getNbOfBaseFunctions(row_data);
249 // resize and clear the right hand side vector
250 OP::locF.resize(OP::nbRows, false);
251 OP::locF.clear();
252 // integrate local vector
253 CHKERR this->iNtegrate(row_data);
254 // assemble local vector
255 OP::locF *= sense;
256 CHKERR this->aSsemble(row_data);
258 };
259
260 auto do_work_lhs = [this](int row_side, int col_side, EntityType row_type,
261 EntityType col_type,
263 EntitiesFieldData::EntData &col_data, int sense) {
265
266 auto check_if_assemble_transpose = [&] {
267 if (this->sYmm) {
268 if (OP::rowSide != OP::colSide || OP::rowType != OP::colType)
269 return true;
270 else
271 return false;
272 } else if (OP::assembleTranspose) {
273 return true;
274 }
275 return false;
276 };
277
278 OP::rowSide = row_side;
279 OP::rowType = row_type;
280 OP::colSide = col_side;
281 OP::colType = col_type;
282 OP::nbCols = col_data.getIndices().size();
283 OP::locMat.resize(OP::nbRows, OP::nbCols, false);
284 OP::locMat.clear();
285 CHKERR this->iNtegrate(row_data, col_data);
286 OP::locMat *= sense;
287 CHKERR this->aSsemble(row_data, col_data, check_if_assemble_transpose());
289 };
290
291 switch (OP::opType) {
292 case OP::OPROW:
293
294 OP::nbRows = row_data.getIndices().size();
295 if (!OP::nbRows)
297 OP::nbIntegrationPts = OP::getGaussPts().size2();
298 OP::nbRowBaseFunctions = OP::getNbOfBaseFunctions(row_data);
299
300 if (!OP::nbRows)
302
303 for (auto &bd : *brokenBaseSideData) {
304
305#ifndef NDEBUG
306 if (bd.getData().getSpace() != HDIV && bd.getData().getSpace() != HCURL) {
307 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "%s",
308 (std::string("Expect Hdiv or Hcurl space but received ") +
309 FieldSpaceNames[bd.getData().getSpace()])
310 .c_str());
311 }
312 if (!bd.getData().getNSharedPtr(bd.getData().getBase())) {
313 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
314 "base functions not set");
315 }
316#endif
317
318 CHKERR do_work_lhs(
319
320 // side
321 row_side, bd.getSide(),
322
323 // type
324 row_type, bd.getType(),
325
326 // row_data
327 row_data, bd.getData(),
328
329 // sense
330 bd.getSense()
331
332 );
333 }
334
335 break;
336 case OP::OPSPACE:
337 for (auto &bd : *brokenBaseSideData) {
338 CHKERR do_work_rhs(bd.getSide(), bd.getType(), bd.getData(),
339 bd.getSense());
340 }
341 break;
342 default:
344 (std::string("wrong op type ") +
346 .c_str());
347 }
348
350}
351
352template <int FIELD_DIM, IntegrationType I, typename OpBrokenBase>
354
355template <int FIELD_DIM, typename OpBrokenBase>
357 : public OpBrokenBase {
358
360
362 const std::string row_field,
363 boost::shared_ptr<std::vector<BrokenBaseSideData>> broken_base_side_data,
364 boost::shared_ptr<double> beta_ptr, const bool assmb_transpose,
365 const bool only_transpose, boost::shared_ptr<Range> ents_ptr = nullptr)
366 : OP(row_field, broken_base_side_data, assmb_transpose, only_transpose,
367 ents_ptr),
368 scalarBetaPtr(beta_ptr) {}
369
371 const std::string row_field,
372 boost::shared_ptr<std::vector<BrokenBaseSideData>> broken_base_side_data,
373 double beta, const bool assmb_transpose, const bool only_transpose,
374 boost::shared_ptr<Range> ents_ptr = nullptr)
375 : OpBrokenSpaceConstrainImpl(row_field, broken_base_side_data,
376 boost::make_shared<double>(beta),
377 assmb_transpose, only_transpose, ents_ptr) {}
378
379protected:
380 boost::shared_ptr<double> scalarBetaPtr;
383};
384
385template <int FIELD_DIM, typename OpBase>
388 EntitiesFieldData::EntData &col_data) {
390
391 auto nb_row_dofs = row_data.getIndices().size();
392 auto nb_col_dofs = col_data.getIndices().size();
393 if (!nb_row_dofs || !nb_col_dofs)
395
397 FTENSOR_INDEX(3, J);
398
399 auto t_w = this->getFTensor0IntegrationWeight();
400 auto t_normal = OpBase::getFTensor1NormalsAtGaussPts();
401 size_t nb_base_functions = col_data.getN().size2() / 3;
402
403#ifndef NDEBUG
404 if (nb_row_dofs % FIELD_DIM != 0 || nb_col_dofs % FIELD_DIM != 0) {
405 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
406 "number of dofs not divisible by field dimension");
407 }
408 if (nb_row_dofs > row_data.getN().size2() * FIELD_DIM ||
409 nb_col_dofs > col_data.getN().size2() * FIELD_DIM) {
410 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
411 "number of dofs exceeds number of base functions");
412 }
413#endif // NDEBUG
414
415 auto t_col_base = col_data.getFTensor1N<3>();
416
417 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
418
419 int cc = 0;
420 for (; cc != nb_col_dofs / FIELD_DIM; cc++) {
421 auto t_row_base = row_data.getFTensor0N(gg, 0);
422 for (auto rr = 0; rr != nb_row_dofs / FIELD_DIM; ++rr) {
423 OP::locMat(FIELD_DIM * rr, FIELD_DIM * cc) +=
424 (0.5 * t_w * t_row_base) * (t_normal(J) * t_col_base(J));
425 ++t_row_base;
426 }
427
428 ++t_col_base;
429 }
430
431 for (; cc < nb_base_functions; ++cc)
432 ++t_col_base;
433
434 ++t_w;
435 ++t_normal;
436 }
437
438 for (auto rr = 0; rr != nb_row_dofs / FIELD_DIM; ++rr) {
439 for (auto cc = 0; cc != nb_col_dofs / FIELD_DIM; ++cc) {
440 for (auto dd = 1; dd < FIELD_DIM; ++dd) {
441 OP::locMat(FIELD_DIM * rr + dd, FIELD_DIM * cc + dd) =
442 OP::locMat(FIELD_DIM * rr, FIELD_DIM * cc);
443 }
444 }
445 }
446
447 if (scalarBetaPtr)
448 OP::locMat *= *scalarBetaPtr;
449
451}
452
453template <int FIELD_DIM, IntegrationType I, typename OpBase>
455
456template <int FIELD_DIM, IntegrationType I, typename OpBase>
458
459template <int FIELD_DIM, typename OpBrokenBase>
461 : public OpBrokenBase {
462
464
466 boost::shared_ptr<std::vector<BrokenBaseSideData>> broken_base_side_data,
467 boost::shared_ptr<MatrixDouble> lagrange_ptr,
468 boost::shared_ptr<double> beta_ptr,
469 boost::shared_ptr<Range> ents_ptr = nullptr)
470 : OP(broken_base_side_data, ents_ptr), scalarBetaPtr(beta_ptr),
471 lagrangePtr(lagrange_ptr) {}
472
474 boost::shared_ptr<std::vector<BrokenBaseSideData>> broken_base_side_data,
475 boost::shared_ptr<MatrixDouble> lagrange_ptr, double beta,
476 boost::shared_ptr<Range> ents_ptr = nullptr)
477 : OpBrokenSpaceConstrainDFluxImpl(broken_base_side_data, lagrange_ptr,
478 boost::make_shared<double>(beta),
479 ents_ptr) {}
480
481private:
483
484 boost::shared_ptr<double> scalarBetaPtr;
485 boost::shared_ptr<MatrixDouble> lagrangePtr;
486};
487
488template <int FIELD_DIM, typename OpBase>
491 EntitiesFieldData::EntData &row_data) {
493
495 FTENSOR_INDEX(3, J);
496
497 auto t_w = this->getFTensor0IntegrationWeight();
498 auto t_normal = OpBase::getFTensor1NormalsAtGaussPts();
499 auto t_lagrange = getFTensor1FromMat<FIELD_DIM>(*lagrangePtr);
500
501 auto t_row_base = row_data.getFTensor1N<3>();
502 auto nb_base_functions = row_data.getN().size2() / 3;
503
504 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
505 auto t_vec = getFTensor1FromPtr<FIELD_DIM>(&*OP::locF.data().begin());
506 size_t rr = 0;
507 for (; rr != row_data.getIndices().size() / FIELD_DIM; ++rr) {
508 t_vec(i) += (0.5 * t_w * (t_row_base(J) * t_normal(J))) * t_lagrange(i);
509 ++t_row_base;
510 ++t_vec;
511 }
512 for (; rr < nb_base_functions; ++rr)
513 ++t_row_base;
514 ++t_w;
515 ++t_normal;
516 ++t_lagrange;
517 }
518
519 if(scalarBetaPtr)
520 OP::locF *= *scalarBetaPtr;
521
523}
524
525template <int FIELD_DIM, typename OpBase>
527 : public OpBase {
528
529 using OP = OpBase;
530
532 const std::string row_field,
533 boost::shared_ptr<std::vector<BrokenBaseSideData>> broken_side_data_ptr,
534 boost::shared_ptr<double> beta_ptr,
535 boost::shared_ptr<Range> ents_ptr = nullptr)
536 : OpBase(row_field, row_field, OpBase::OPROW, ents_ptr),
537 brokenSideDataPtr(broken_side_data_ptr), scalarBetaPtr(beta_ptr) {}
538
540 const std::string row_field,
541 boost::shared_ptr<std::vector<BrokenBaseSideData>> broken_side_data_ptr,
542 double beta, boost::shared_ptr<Range> ents_ptr = nullptr)
543 : OpBrokenSpaceConstrainDHybridImpl(row_field, broken_side_data_ptr,
544 boost::make_shared<double>(beta),
545 ents_ptr) {}
546
547private:
548 boost::shared_ptr<double> scalarBetaPtr;
549 boost::shared_ptr<std::vector<BrokenBaseSideData>> brokenSideDataPtr;
550 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data);
551};
552
553template <int FIELD_DIM, typename OpBase>
556 EntitiesFieldData::EntData &row_data) {
558
560 FTENSOR_INDEX(3, J);
561
562 OP::locF.resize(row_data.getIndices().size(), false);
563 OP::locF.clear();
564
565 for (auto &bd : *brokenSideDataPtr) {
566 auto t_w = this->getFTensor0IntegrationWeight();
567 auto t_normal = OpBase::getFTensor1NormalsAtGaussPts();
568 auto t_row_base = row_data.getFTensor0N();
569 auto t_flux = getFTensor2FromMat<FIELD_DIM, 3>(bd.getFlux());
570 auto nb_base_functions = row_data.getN().size2() / 3;
571 auto sense = bd.getSense();
572 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
573 auto t_vec = getFTensor1FromPtr<FIELD_DIM>(&*OP::locF.data().begin());
574 size_t rr = 0;
575 for (; rr != row_data.getIndices().size() / FIELD_DIM; ++rr) {
576 t_vec(i) +=
577 (sense * 0.5 * t_w) * t_row_base * t_normal(J) * t_flux(i, J);
578 ++t_row_base;
579 ++t_vec;
580 }
581 for (; rr < nb_base_functions; ++rr)
582 ++t_row_base;
583 ++t_w;
584 ++t_normal;
585 ++t_flux;
586 }
587 }
588
589 if(scalarBetaPtr)
590 OP::locF *= *scalarBetaPtr;
591
593}
594
595} // namespace MoFEM
596
597#endif // __FORMSBROKENSPACECONSTRAINTIMPL_HPP__
#define FTENSOR_INDEX(DIM, I)
constexpr int FIELD_DIM
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
@ NOSPACE
Definition definitions.h:83
@ HCURL
field with continuous tangents
Definition definitions.h:86
@ HDIV
field with continuous normal traction
Definition definitions.h:87
static const char *const FieldSpaceNames[]
Definition definitions.h:92
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define CHK_MOAB_THROW(err, msg)
Check error code of MoAB function and throw MoFEM exception.
@ MOFEM_IMPOSSIBLE_CASE
Definition definitions.h:35
@ MOFEM_DATA_INCONSISTENCY
Definition definitions.h:31
#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
FTensor::Index< 'J', DIM1 > J
Definition level_set.cpp:30
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition ddTensor0.hpp:33
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
FTensor::Tensor1< FTensor::PackPtr< T *, S >, Tensor_Dim > getFTensor1FromMat(ublas::matrix< T, L, A > &data)
Get tensor rank 1 (vector) form data matrix.
FTensor::Tensor2< FTensor::PackPtr< double *, 1 >, Tensor_Dim1, Tensor_Dim2 > getFTensor2FromMat(MatrixDouble &data)
Get tensor rank 2 (matrix) form data matrix.
FTensor::Tensor1< FTensor::PackPtr< double *, S >, DIM > getFTensor1FromPtr(double *ptr)
Make Tensor1 from pointer.
constexpr auto field_name
OpBaseImpl< PETSC, EdgeEleOp > OpBase
Definition radiation.cpp:29
Data on single entity (This is passed as argument to DataOperator::doWork)
int sEnse
Entity sense (orientation)
VectorDouble fieldData
Field data on entity.
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.
VectorInt localIndices
Local indices on entity.
VectorInt iNdices
Global indices on entity.
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
virtual int getSense() const
get entity sense, need to calculate base functions with conforming approximation fields
VectorFieldEntities fieldEntities
Field entities.
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1N(FieldApproximationBase base)
Get base functions for Hdiv/Hcurl spaces.
std::array< std::array< boost::shared_ptr< MatrixDouble >, LASTBASE >, LastDerivative > baseFunctionsAndBaseDerivatives
const VectorInt & getIndices() const
Get global indices of dofs on entity.
FieldApproximationBase bAse
Field approximation base.
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
boost::shared_ptr< Range > entsPtr
Entities on which element is run.
int nbIntegrationPts
number of integration points
boost::shared_ptr< std::vector< BrokenBaseSideData > > brokenBaseSideData
OpBrokenBaseImpl(const std::string row_field, boost::shared_ptr< std::vector< BrokenBaseSideData > > broken_base_side_data, const bool assmb_transpose, const bool only_transpose, boost::shared_ptr< Range > ents_ptr=nullptr)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
OpBrokenBaseImpl(boost::shared_ptr< std::vector< BrokenBaseSideData > > broken_base_side_data, boost::shared_ptr< Range > ents_ptr=nullptr)
Operator for broken loop side.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data)
OpBrokenSpaceConstrainDFluxImpl(boost::shared_ptr< std::vector< BrokenBaseSideData > > broken_base_side_data, boost::shared_ptr< MatrixDouble > lagrange_ptr, boost::shared_ptr< double > beta_ptr, boost::shared_ptr< Range > ents_ptr=nullptr)
OpBrokenSpaceConstrainDFluxImpl(boost::shared_ptr< std::vector< BrokenBaseSideData > > broken_base_side_data, boost::shared_ptr< MatrixDouble > lagrange_ptr, double beta, boost::shared_ptr< Range > ents_ptr=nullptr)
OpBrokenSpaceConstrainDHybridImpl(const std::string row_field, boost::shared_ptr< std::vector< BrokenBaseSideData > > broken_side_data_ptr, boost::shared_ptr< double > beta_ptr, boost::shared_ptr< Range > ents_ptr=nullptr)
OpBrokenSpaceConstrainDHybridImpl(const std::string row_field, boost::shared_ptr< std::vector< BrokenBaseSideData > > broken_side_data_ptr, double beta, boost::shared_ptr< Range > ents_ptr=nullptr)
OpBrokenSpaceConstrainImpl(const std::string row_field, boost::shared_ptr< std::vector< BrokenBaseSideData > > broken_base_side_data, boost::shared_ptr< double > beta_ptr, const bool assmb_transpose, const bool only_transpose, boost::shared_ptr< Range > ents_ptr=nullptr)
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
OpBrokenSpaceConstrainImpl(const std::string row_field, boost::shared_ptr< std::vector< BrokenBaseSideData > > broken_base_side_data, double beta, const bool assmb_transpose, const bool only_transpose, boost::shared_ptr< Range > ents_ptr=nullptr)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
boost::shared_ptr< std::vector< BrokenBaseSideData > > brokenBaseSideData
OpGetBrokenBaseSideData(const std::string field_name, boost::shared_ptr< std::vector< BrokenBaseSideData > > broken_base_side_data)
Element used to execute operators on side of the element.
const std::string sideFEName
boost::shared_ptr< E > sideFEPtr
OpLoopSide(MoFEM::Interface &m_field, const std::string fe_name, const int side_dim, boost::shared_ptr< Range > fe_range, const LogManager::SeverityLevel sev=Sev::noisy, boost::shared_ptr< AdjCache > adj_cache=nullptr)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
OpSetFlux(boost::shared_ptr< std::vector< BrokenBaseSideData > > broken_base_side_data, boost::shared_ptr< MatrixDouble > flux_ptr)
boost::shared_ptr< std::vector< BrokenBaseSideData > > brokenBaseSideData
boost::shared_ptr< MatrixDouble > fluxPtr