v0.14.0
ContactPrismElementForcesAndSourcesCore.cpp
Go to the documentation of this file.
1 /** \file ContactPrismElementForcesAndSourcesCore.cpp
2 
3 \brief Implementation of contact prism element
4 
5 */
6 
7 namespace MoFEM {
8 
9 ContactPrismElementForcesAndSourcesCore::
10  ContactPrismElementForcesAndSourcesCore(Interface &m_field)
11  : ForcesAndSourcesCore(m_field),
12  dataOnMaster{
13 
14  nullptr,
15  boost::make_shared<EntitiesFieldData>(MBENTITYSET), // NOFIELD
16  boost::make_shared<EntitiesFieldData>(MBENTITYSET), // H1
17  boost::make_shared<EntitiesFieldData>(MBENTITYSET), // HCURL
18  boost::make_shared<EntitiesFieldData>(MBENTITYSET), // HDIV
19  boost::make_shared<EntitiesFieldData>(MBENTITYSET) // L2
20 
21  },
22  dataOnSlave{
23 
24  nullptr,
25  boost::make_shared<EntitiesFieldData>(MBENTITYSET), // NOFIELD
26  boost::make_shared<EntitiesFieldData>(MBENTITYSET), // H1
27  boost::make_shared<EntitiesFieldData>(MBENTITYSET), // HCURL
28  boost::make_shared<EntitiesFieldData>(MBENTITYSET), // HDIV
29  boost::make_shared<EntitiesFieldData>(MBENTITYSET) // L2
30 
31  },
32  derivedDataOnMaster{
33 
34  nullptr,
35  boost::make_shared<DerivedEntitiesFieldData>(dataOnMaster[NOFIELD]),
36  boost::make_shared<DerivedEntitiesFieldData>(dataOnMaster[H1]),
37  boost::make_shared<DerivedEntitiesFieldData>(dataOnMaster[HCURL]),
38  boost::make_shared<DerivedEntitiesFieldData>(dataOnMaster[HDIV]),
39  boost::make_shared<DerivedEntitiesFieldData>(dataOnMaster[L2])
40 
41  },
42  derivedDataOnSlave{
43 
44  nullptr,
45  boost::make_shared<DerivedEntitiesFieldData>(dataOnSlave[NOFIELD]),
46  boost::make_shared<DerivedEntitiesFieldData>(dataOnSlave[H1]),
47  boost::make_shared<DerivedEntitiesFieldData>(dataOnSlave[HCURL]),
48  boost::make_shared<DerivedEntitiesFieldData>(dataOnSlave[HDIV]),
49  boost::make_shared<DerivedEntitiesFieldData>(dataOnSlave[L2])
50 
51  },
52  dataH1Master(*dataOnMaster[H1].get()),
53  dataH1Slave(*dataOnSlave[H1].get()),
54  dataNoFieldSlave(*dataOnSlave[NOFIELD].get()),
55  dataNoFieldMaster(*dataOnMaster[NOFIELD].get()),
56  dataHcurlMaster(*dataOnMaster[HCURL].get()),
57  dataHcurlSlave(*dataOnSlave[HCURL].get()),
58  dataHdivMaster(*dataOnMaster[HDIV].get()),
59  dataL2Master(*dataOnMaster[L2].get()),
60  dataHdivSlave(*dataOnSlave[HDIV].get()),
61  dataL2Slave(*dataOnSlave[L2].get()), opContravariantTransform() {
62 
63  getUserPolynomialBase() =
64  boost::shared_ptr<BaseFunction>(new TriPolynomialBase());
65 
66  // Data on elements for proper spaces
67  dataOnMaster[H1]->setElementType(MBTRI);
68  derivedDataOnMaster[H1]->setElementType(MBTRI);
69  dataOnSlave[H1]->setElementType(MBTRI);
70  derivedDataOnSlave[H1]->setElementType(MBTRI);
71 
72  dataOnMaster[HDIV]->setElementType(MBTRI);
73  derivedDataOnMaster[HDIV]->setElementType(MBTRI);
74  dataOnSlave[HDIV]->setElementType(MBTRI);
75  derivedDataOnSlave[HDIV]->setElementType(MBTRI);
76 
77  dataOnMaster[NOFIELD]->dataOnEntities[MBENTITYSET].push_back(
78  new EntitiesFieldData::EntData());
79  dataOnSlave[NOFIELD]->dataOnEntities[MBENTITYSET].push_back(
80  new EntitiesFieldData::EntData());
81 
82  derivedDataOnMaster[NOFIELD]->dataOnEntities[MBENTITYSET].push_back(
83  new EntitiesFieldData::EntData());
84  derivedDataOnSlave[NOFIELD]->dataOnEntities[MBENTITYSET].push_back(
85  new EntitiesFieldData::EntData());
86 
87  CHK_THROW_MESSAGE(createDataOnElement(MBPRISM),
88  "Problem with creation data on element");
89 }
90 
91 MoFEMErrorCode
92 ContactPrismElementForcesAndSourcesCore::setDefaultGaussPts(const int rule) {
93  MoFEMFunctionBegin;
94 
95  if (rule < QUAD_2D_TABLE_SIZE) {
96  if (QUAD_2D_TABLE[rule]->dim != 2) {
97  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "wrong dimension");
98  }
99  if (QUAD_2D_TABLE[rule]->order < rule) {
100  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
101  "wrong order %d != %d", QUAD_2D_TABLE[rule]->order, rule);
102  }
103  nbGaussPts = QUAD_2D_TABLE[rule]->npoints;
104  // For master and slave
105  gaussPtsMaster.resize(3, nbGaussPts, false);
106  gaussPtsSlave.resize(3, nbGaussPts, false);
107 
108  cblas_dcopy(nbGaussPts, &QUAD_2D_TABLE[rule]->points[1], 3,
109  &gaussPtsMaster(0, 0), 1);
110  cblas_dcopy(nbGaussPts, &QUAD_2D_TABLE[rule]->points[2], 3,
111  &gaussPtsMaster(1, 0), 1);
112  cblas_dcopy(nbGaussPts, QUAD_2D_TABLE[rule]->weights, 1,
113  &gaussPtsMaster(2, 0), 1);
114 
115  gaussPtsSlave = gaussPtsMaster;
116 
117  dataH1Master.dataOnEntities[MBVERTEX][0].getN(NOBASE).resize(nbGaussPts, 3,
118  false);
119 
120  dataH1Slave.dataOnEntities[MBVERTEX][0].getN(NOBASE).resize(nbGaussPts, 3,
121  false);
122 
123  double *shape_ptr_master =
124  &*dataH1Master.dataOnEntities[MBVERTEX][0].getN(NOBASE).data().begin();
125  cblas_dcopy(3 * nbGaussPts, QUAD_2D_TABLE[rule]->points, 1,
126  shape_ptr_master, 1);
127  double *shape_ptr_slave =
128  &*dataH1Slave.dataOnEntities[MBVERTEX][0].getN(NOBASE).data().begin();
129  cblas_dcopy(3 * nbGaussPts, QUAD_2D_TABLE[rule]->points, 1, shape_ptr_slave,
130  1);
131  } else {
132  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
133  "rule > quadrature order %d < %d", rule, QUAD_2D_TABLE_SIZE);
134  nbGaussPts = 0;
135  }
136 
137  MoFEMFunctionReturn(0);
138 }
139 
140 MoFEMErrorCode ContactPrismElementForcesAndSourcesCore::operator()() {
141  MoFEMFunctionBegin;
142 
143  if (numeredEntFiniteElementPtr->getEntType() != MBPRISM)
144  MoFEMFunctionReturnHot(0);
145 
146  EntitiesFieldData &data_div = *dataOnElement[HDIV];
147  EntitiesFieldData &data_curl = *dataOnElement[HCURL];
148  EntitiesFieldData &data_l2 = *dataOnElement[HCURL];
149 
150  EntityHandle ent = numeredEntFiniteElementPtr->getEnt();
151  auto get_coord_and_normal = [&]() {
152  MoFEMFunctionBegin;
153  int num_nodes;
154  const EntityHandle *conn;
155  CHKERR mField.get_moab().get_connectivity(ent, conn, num_nodes, true);
156  coords.resize(num_nodes * 3, false);
157  CHKERR mField.get_moab().get_coords(conn, num_nodes,
158  &*coords.data().begin());
159  normal.resize(6, false);
160 
161  for (auto &v : {&tangentMasterOne, &tangentMasterTwo, &tangentSlaveOne,
162  &tangentSlaveTwo}) {
163  v->resize(3);
164  v->clear();
165  }
166 
167  CHKERR Tools::getTriNormal(&coords[0], &normal[0]);
168  CHKERR Tools::getTriNormal(&coords[9], &normal[3]);
169 
170  auto get_vec_ptr = [](VectorDouble &vec_double, int r = 0) {
171  return FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 3>(
172  &vec_double(r + 0), &vec_double(r + 1), &vec_double(r + 2));
173  };
174 
175  auto t_coords_master = get_vec_ptr(coords);
176  auto t_coords_slave = get_vec_ptr(coords, 9);
177  auto t_normal_master = get_vec_ptr(normal);
178  auto t_normal_slave = get_vec_ptr(normal, 3);
179 
180  auto t_t1_master = get_vec_ptr(tangentMasterOne);
181  auto t_t2_master = get_vec_ptr(tangentMasterTwo);
182  auto t_t1_slave = get_vec_ptr(tangentSlaveOne);
183  auto t_t2_slave = get_vec_ptr(tangentSlaveTwo);
184 
185  const double *diff_ptr = Tools::diffShapeFunMBTRI.data();
186  FTensor::Tensor1<FTensor::PackPtr<const double *, 2>, 2> t_diff(
187  &diff_ptr[0], &diff_ptr[1]);
188 
189  FTensor::Index<'i', 3> i;
190  FTensor::Index<'j', 3> j;
191  FTensor::Index<'k', 3> k;
192 
193  FTensor::Number<0> N0;
194 
195  for (int nn = 0; nn != 3; ++nn) {
196  t_t1_master(i) += t_coords_master(i) * t_diff(N0);
197  t_t1_slave(i) += t_coords_slave(i) * t_diff(N0);
198  ++t_coords_master;
199  ++t_coords_slave;
200  ++t_diff;
201  }
202 
203  aRea[0] = sqrt(t_normal_master(i) * t_normal_master(i)) / 2.;
204  aRea[1] = sqrt(t_normal_slave(i) * t_normal_slave(i)) / 2.;
205 
206  t_t2_master(j) =
207  FTensor::levi_civita(i, j, k) * t_normal_master(k) * t_t1_master(i);
208  t_t2_slave(j) =
209  FTensor::levi_civita(i, j, k) * t_normal_slave(k) * t_t1_slave(i);
210 
211  MoFEMFunctionReturn(0);
212  };
213  CHKERR get_coord_and_normal();
214 
215  CHKERR getSpacesAndBaseOnEntities(dataH1);
216 
217  // H1
218  if ((dataH1.spacesOnEntities[MBVERTEX]).test(H1)) {
219  CHKERR getEntitySense<MBEDGE>(dataH1);
220  CHKERR getEntityDataOrder<MBEDGE>(dataH1, H1);
221  CHKERR getEntitySense<MBTRI>(dataH1);
222  CHKERR getEntityDataOrder<MBTRI>(dataH1, H1);
223  }
224 
225  // Hcurl
226  if ((dataH1.spacesOnEntities[MBEDGE]).test(HCURL)) {
227  CHKERR getEntitySense<MBEDGE>(data_curl);
228  CHKERR getEntityDataOrder<MBEDGE>(data_curl, HCURL);
229  CHKERR getEntitySense<MBTRI>(data_curl);
230  CHKERR getEntityDataOrder<MBTRI>(data_curl, HCURL);
231  }
232 
233  // Hdiv
234  if ((dataH1.spacesOnEntities[MBTRI]).test(HDIV)) {
235  CHKERR getEntitySense<MBTRI>(data_div);
236  CHKERR getEntityDataOrder<MBTRI>(data_div, HDIV);
237  dataHdivSlave.spacesOnEntities[MBTRI].set(HDIV);
238  dataHdivMaster.spacesOnEntities[MBTRI].set(HDIV);
239  data_div.spacesOnEntities[MBTRI].set(HDIV);
240  }
241 
242  // L2
243  if ((dataH1.spacesOnEntities[MBTRI]).test(L2)) {
244  CHKERR getEntitySense<MBTRI>(data_l2);
245  CHKERR getEntityDataOrder<MBTRI>(data_l2, L2);
246  }
247 
248  auto clean_data = [](EntitiesFieldData &data) {
249  MoFEMFunctionBegin;
250  data.bAse.reset();
251  for (EntityType t = MBVERTEX; t != MBMAXTYPE; ++t) {
252  data.spacesOnEntities[t].reset();
253  data.basesOnEntities[t].reset();
254  }
255  for (int s = 0; s != LASTSPACE; ++s)
256  data.basesOnSpaces[s].reset();
257 
258  MoFEMFunctionReturn(0);
259  };
260 
261  auto copy_data = [](EntitiesFieldData &data, EntitiesFieldData &copy_data,
262  const int shift) {
263  MoFEMFunctionBegin;
264 
265  if (shift != 0 && shift != 6) {
266  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
267  "Wrong shift for contact prism element");
268  }
269 
270  data.bAse = copy_data.bAse;
271  for (auto t : {MBVERTEX, MBEDGE, MBTRI}) {
272  data.spacesOnEntities[t] = copy_data.spacesOnEntities[t];
273  data.basesOnEntities[t] = copy_data.basesOnEntities[t];
274  data.basesOnSpaces[t] = copy_data.basesOnSpaces[t];
275  data.brokenBasesOnSpaces[t] = copy_data.brokenBasesOnSpaces[t];
276  }
277 
278  for (int ii = 0; ii != 3; ++ii) {
279  data.dataOnEntities[MBEDGE][ii].getSense() =
280  copy_data.dataOnEntities[MBEDGE][ii + shift].getSense();
281  data.dataOnEntities[MBEDGE][ii].getOrder() =
282  copy_data.dataOnEntities[MBEDGE][ii + shift].getOrder();
283  }
284 
285  if (shift == 0) {
286  data.dataOnEntities[MBTRI][0].getSense() =
287  copy_data.dataOnEntities[MBTRI][3].getSense();
288  data.dataOnEntities[MBTRI][0].getOrder() =
289  copy_data.dataOnEntities[MBTRI][3].getOrder();
290  } else {
291  data.dataOnEntities[MBTRI][0].getSense() =
292  copy_data.dataOnEntities[MBTRI][4].getSense();
293  data.dataOnEntities[MBTRI][0].getOrder() =
294  copy_data.dataOnEntities[MBTRI][4].getOrder();
295  }
296 
297  MoFEMFunctionReturn(0);
298  };
299 
300  auto copy_data_hdiv = [](EntitiesFieldData &data,
301  EntitiesFieldData &copy_data, const int shift) {
302  MoFEMFunctionBegin;
303 
304  if (shift != 3 && shift != 4) {
305  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
306  "Wrong shift for contact prism element");
307  }
308 
309  data.bAse = copy_data.bAse;
310 
311  for (auto t : {MBVERTEX, MBTRI}) {
312  data.spacesOnEntities[t] = copy_data.spacesOnEntities[MBVERTEX];
313  data.basesOnEntities[t] = copy_data.basesOnEntities[MBVERTEX];
314  data.basesOnSpaces[t] = copy_data.basesOnSpaces[MBVERTEX];
315  }
316 
317  auto &cpy_ent_dat = copy_data.dataOnEntities[MBTRI][shift];
318  auto &ent_dat = data.dataOnEntities[MBTRI][0];
319  ent_dat.getBase() = cpy_ent_dat.getBase();
320  ent_dat.getSpace() = cpy_ent_dat.getSpace();
321  ent_dat.getSense() = ent_dat.getSense();
322  ent_dat.getOrder() = cpy_ent_dat.getOrder();
323 
324  MoFEMFunctionReturn(0);
325  };
326 
327  CHKERR clean_data(dataH1Slave);
328  CHKERR copy_data(dataH1Slave, dataH1, 6);
329  CHKERR clean_data(dataH1Master);
330  CHKERR copy_data(dataH1Master, dataH1, 0);
331 
332  int order_data = getMaxDataOrder();
333  int order_row = getMaxRowOrder();
334  int order_col = getMaxColOrder(); // maybe two different rules?
335  int rule = getRule(order_row, order_col, order_data);
336 
337  if (rule >= 0) {
338 
339  CHKERR setDefaultGaussPts(rule);
340 
341  } else {
342 
343  // Master-Slave
344  if (gaussPtsMaster.size2() != gaussPtsSlave.size2())
345  SETERRQ(
346  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
347  "Number of Gauss Points at Master triangle is different than slave");
348 
349  CHKERR setGaussPts(order_row, order_col, order_data);
350  nbGaussPts = gaussPtsMaster.size2();
351  dataH1Master.dataOnEntities[MBVERTEX][0].getN(NOBASE).resize(nbGaussPts, 3,
352  false);
353  dataH1Slave.dataOnEntities[MBVERTEX][0].getN(NOBASE).resize(nbGaussPts, 3,
354  false);
355 
356  if (nbGaussPts) {
357  CHKERR Tools::shapeFunMBTRI<1>(&*dataH1Master.dataOnEntities[MBVERTEX][0]
358  .getN(NOBASE)
359  .data()
360  .begin(),
361  &gaussPtsMaster(0, 0),
362  &gaussPtsMaster(1, 0), nbGaussPts);
363 
364  CHKERR Tools::shapeFunMBTRI<1>(
365  &*dataH1Slave.dataOnEntities[MBVERTEX][0].getN(NOBASE).data().begin(),
366  &gaussPtsSlave(0, 0), &gaussPtsSlave(1, 0), nbGaussPts);
367  }
368  }
369 
370  if (nbGaussPts == 0)
371  MoFEMFunctionReturnHot(0);
372 
373  // Get coordinates on slave and master
374  {
375  coordsAtGaussPtsMaster.resize(nbGaussPts, 3, false);
376  coordsAtGaussPtsSlave.resize(nbGaussPts, 3, false);
377  for (int gg = 0; gg < nbGaussPts; gg++) {
378  for (int dd = 0; dd < 3; dd++) {
379  coordsAtGaussPtsMaster(gg, dd) = cblas_ddot(
380  3, &dataH1Master.dataOnEntities[MBVERTEX][0].getN(NOBASE)(gg, 0), 1,
381  &coords[dd], 3);
382  coordsAtGaussPtsSlave(gg, dd) = cblas_ddot(
383  3, &dataH1Slave.dataOnEntities[MBVERTEX][0].getN(NOBASE)(gg, 0), 1,
384  &coords[9 + dd], 3);
385  }
386  }
387  }
388 
389  for (int space = HCURL; space != LASTSPACE; ++space)
390  if (dataOnElement[space]) {
391 
392  dataH1Master.dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE) =
393  dataOnMaster[H1]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE);
394  dataH1Slave.dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE) =
395  dataOnSlave[H1]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE);
396 
397  if (space == HDIV) {
398 
399  if (dataH1.spacesOnEntities[MBTRI].test(HDIV)) {
400 
401  CHKERR clean_data(dataHdivSlave);
402  CHKERR copy_data_hdiv(dataHdivSlave, data_div, 4);
403  CHKERR clean_data(dataHdivMaster);
404  CHKERR copy_data_hdiv(dataHdivMaster, data_div, 3);
405 
406  dataHdivMaster.dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE) =
407  dataOnMaster[H1]->dataOnEntities[MBVERTEX][0].getNSharedPtr(
408  NOBASE);
409  dataHdivSlave.dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE) =
410  dataOnSlave[H1]->dataOnEntities[MBVERTEX][0].getNSharedPtr(
411  NOBASE);
412  }
413  }
414  }
415 
416  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; b++) {
417  if (dataH1.bAse.test(b)) {
418  switch (static_cast<FieldApproximationBase>(b)) {
419  case AINSWORTH_LEGENDRE_BASE:
420  case AINSWORTH_LOBATTO_BASE:
421  if (dataH1.spacesOnEntities[MBVERTEX].test(H1)) {
422  CHKERR getUserPolynomialBase()
423  -> getValue(
424  gaussPtsMaster,
425  boost::shared_ptr<BaseFunctionCtx>(new EntPolynomialBaseCtx(
426  dataH1Master, H1, CONTINUOUS,
427  static_cast<FieldApproximationBase>(b), NOBASE)));
428 
429  CHKERR getUserPolynomialBase()
430  -> getValue(
431  gaussPtsSlave,
432  boost::shared_ptr<BaseFunctionCtx>(new EntPolynomialBaseCtx(
433  dataH1Slave, H1, CONTINUOUS,
434  static_cast<FieldApproximationBase>(b), NOBASE)));
435  }
436  break;
437  case DEMKOWICZ_JACOBI_BASE:
438  if (dataH1.spacesOnEntities[MBTRI].test(HDIV)) {
439 
440  CHKERR getUserPolynomialBase()
441  -> getValue(
442  gaussPtsMaster,
443  boost::shared_ptr<BaseFunctionCtx>(new EntPolynomialBaseCtx(
444  dataHdivMaster, HDIV, CONTINUOUS,
445  static_cast<FieldApproximationBase>(b), NOBASE)));
446  CHKERR getUserPolynomialBase()
447  -> getValue(
448  gaussPtsSlave,
449  boost::shared_ptr<BaseFunctionCtx>(new EntPolynomialBaseCtx(
450  dataHdivSlave, HDIV, CONTINUOUS,
451  static_cast<FieldApproximationBase>(b), NOBASE)));
452 
453  opContravariantTransform.normalRawPtr = &normal;
454  opContravariantTransform.normalShift = 0;
455  CHKERR opContravariantTransform.opRhs(dataHdivMaster);
456  opContravariantTransform.normalShift = 3;
457  CHKERR opContravariantTransform.opRhs(dataHdivSlave);
458  }
459 
460  if (dataH1.spacesOnEntities[MBEDGE].test(HCURL)) {
461  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
462  "Not yet implemented");
463  }
464  if (dataH1.spacesOnEntities[MBTET].test(L2)) {
465  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
466  "Not yet implemented");
467  }
468  break;
469  default:
470  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
471  "Not yet implemented");
472  }
473  }
474  }
475 
476  // Iterate over operators
477  CHKERR loopOverOperators();
478 
479  MoFEMFunctionReturn(0);
480 }
481 
482 MoFEMErrorCode ContactPrismElementForcesAndSourcesCore::loopOverOperators() {
483  MoFEMFunctionBegin;
484 
485  constexpr std::array<UserDataOperator::OpType, 2> types{
486  UserDataOperator::OPROW, UserDataOperator::OPCOL};
487  std::array<std::string, 2> last_eval_field_name{std::string(), std::string()};
488 
489  auto oit = opPtrVector.begin();
490  auto hi_oit = opPtrVector.end();
491 
492  for (; oit != hi_oit; oit++) {
493 
494  oit->setPtrFE(this);
495 
496  if (oit->opType == UserDataOperator::OPSPACE) {
497 
498  // Set field
499  switch (oit->sPace) {
500  case NOSPACE:
501  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Unknown space");
502  case NOFIELD:
503  case H1:
504  case HCURL:
505  case HDIV:
506  case L2:
507  break;
508  default:
509  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
510  "Not implemented for this space", oit->sPace);
511  }
512 
513  // Reseat all data which all field dependent
514  dataOnMaster[oit->sPace]->resetFieldDependentData();
515  dataOnSlave[oit->sPace]->resetFieldDependentData();
516 
517  // Run operator
518  try {
519  CHKERR oit->opRhs(*dataOnMaster[oit->sPace], false);
520  CHKERR oit->opRhs(*dataOnSlave[oit->sPace], false);
521  }
522  CATCH_OP_ERRORS(*oit);
523 
524  } else {
525  boost::shared_ptr<EntitiesFieldData> op_master_data[2];
526  boost::shared_ptr<EntitiesFieldData> op_slave_data[2];
527 
528  for (int ss = 0; ss != 2; ss++) {
529 
530  const std::string field_name =
531  !ss ? oit->rowFieldName : oit->colFieldName;
532  const Field *field_struture = mField.get_field_structure(field_name);
533  const BitFieldId data_id = field_struture->getId();
534  const FieldSpace space = field_struture->getSpace();
535  op_master_data[ss] =
536  !ss ? dataOnMaster[space] : derivedDataOnMaster[space];
537  op_slave_data[ss] =
538  !ss ? dataOnSlave[space] : derivedDataOnSlave[space];
539 
540  if ((oit->getNumeredEntFiniteElementPtr()->getBitFieldIdData() &
541  data_id)
542  .none())
543  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
544  "no data field < %s > on finite element < %s >",
545  field_name.c_str(), getFEName().c_str());
546 
547  if (oit->getOpType() & types[ss] ||
548  oit->getOpType() & UserDataOperator::OPROWCOL) {
549 
550  switch (space) {
551  case NOSPACE:
552  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown space");
553  break;
554  case NOFIELD:
555  case H1:
556  case HCURL:
557  case HDIV:
558  case L2:
559  break;
560  default:
561  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
562  "Not implemented for this space", space);
563  }
564 
565  if (last_eval_field_name[ss] != field_name) {
566  CHKERR getEntityFieldData(*op_master_data[ss], *op_slave_data[ss],
567  field_name, MBEDGE);
568 
569  if (!ss) {
570  struct Extractor {
571  boost::weak_ptr<EntityCacheNumeredDofs>
572  operator()(boost::shared_ptr<FieldEntity> &e) {
573  return e->entityCacheRowDofs;
574  }
575  };
576 
577  CHKERR getEntityIndices(*op_master_data[ss], *op_slave_data[ss],
578  field_name, getRowFieldEnts(), MBEDGE,
579  MBPRISM, Extractor());
580  } else {
581  struct Extractor {
582  boost::weak_ptr<EntityCacheNumeredDofs>
583  operator()(boost::shared_ptr<FieldEntity> &e) {
584  return e->entityCacheColDofs;
585  }
586  };
587  CHKERR getEntityIndices(*op_master_data[ss], *op_slave_data[ss],
588  field_name, getRowFieldEnts(), MBEDGE,
589  MBPRISM, Extractor());
590  }
591 
592  switch (space) {
593  case NOSPACE:
594  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
595  "unknown space");
596  break;
597  case H1: {
598 
599  auto get_indices = [&](auto &master, auto &slave, auto &ents,
600  auto &&ex) {
601  return getNodesIndices(
602  field_name, ents,
603  master.dataOnEntities[MBVERTEX][0].getIndices(),
604  master.dataOnEntities[MBVERTEX][0].getLocalIndices(),
605  slave.dataOnEntities[MBVERTEX][0].getIndices(),
606  slave.dataOnEntities[MBVERTEX][0].getLocalIndices(), ex);
607  };
608 
609  auto get_data = [&](EntitiesFieldData &master_data,
610  EntitiesFieldData &slave_data) {
611  return getNodesFieldData(
612  field_name,
613  master_data.dataOnEntities[MBVERTEX][0].getFieldData(),
614  slave_data.dataOnEntities[MBVERTEX][0].getFieldData(),
615  master_data.dataOnEntities[MBVERTEX][0].getFieldDofs(),
616  slave_data.dataOnEntities[MBVERTEX][0].getFieldDofs(),
617  master_data.dataOnEntities[MBVERTEX][0].getFieldEntities(),
618  slave_data.dataOnEntities[MBVERTEX][0].getFieldEntities(),
619  master_data.dataOnEntities[MBVERTEX][0].getSpace(),
620  slave_data.dataOnEntities[MBVERTEX][0].getSpace(),
621  master_data.dataOnEntities[MBVERTEX][0].getBase(),
622  slave_data.dataOnEntities[MBVERTEX][0].getBase());
623  };
624 
625  if (!ss) {
626 
627  struct Extractor {
628  boost::weak_ptr<EntityCacheNumeredDofs>
629  operator()(boost::shared_ptr<FieldEntity> &e) {
630  return e->entityCacheRowDofs;
631  }
632  };
633 
634  CHKERR get_indices(*op_master_data[ss], *op_slave_data[ss],
635  getRowFieldEnts(), Extractor());
636  } else {
637  struct Extractor {
638  boost::weak_ptr<EntityCacheNumeredDofs>
639  operator()(boost::shared_ptr<FieldEntity> &e) {
640  return e->entityCacheColDofs;
641  }
642  };
643 
644  CHKERR get_indices(*op_master_data[ss], *op_slave_data[ss],
645  getColFieldEnts(), Extractor());
646  }
647 
648  CHKERR get_data(*op_master_data[ss], *op_slave_data[ss]);
649 
650  } break;
651  case HCURL:
652  case HDIV:
653  case L2:
654  break;
655  case NOFIELD:
656  if (!getNinTheLoop()) {
657  // NOFIELD data are the same for each element, can be
658  // retrieved only once
659  const auto bit_number = mField.get_field_bit_number(field_name);
660  if (!ss) {
661  CHKERR getNoFieldRowIndices(*op_master_data[ss], bit_number);
662  CHKERR getNoFieldRowIndices(*op_slave_data[ss], bit_number);
663  } else {
664  CHKERR getNoFieldColIndices(*op_master_data[ss], bit_number);
665  CHKERR getNoFieldColIndices(*op_slave_data[ss], bit_number);
666  }
667  CHKERR getNoFieldFieldData(*op_master_data[ss], bit_number);
668  CHKERR getNoFieldFieldData(*op_slave_data[ss], bit_number);
669  }
670  break;
671  default:
672  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
673  "Not implemented for this space", space);
674  }
675  last_eval_field_name[ss] = field_name;
676  }
677  }
678  }
679 
680  int type;
681 
682  if (UserDataOperator *cast_oit =
683  dynamic_cast<UserDataOperator *>(&*oit)) {
684  type = cast_oit->getFaceType();
685  if (((oit->getOpType() & UserDataOperator::OPROW) ||
686  (oit->getOpType() & UserDataOperator::OPCOL)) &&
687  ((type & UserDataOperator::FACEMASTERMASTER) ||
688  (type & UserDataOperator::FACEMASTERSLAVE) ||
689  (type & UserDataOperator::FACESLAVEMASTER) ||
690  (type & UserDataOperator::FACESLAVESLAVE))) {
691  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
692  "Wrong combination of FaceType and OpType, OPROW or OPCOL "
693  "combined with face-face OpType");
694  }
695 
696  if ((oit->getOpType() & UserDataOperator::OPROWCOL) &&
697  ((type & UserDataOperator::FACEMASTER) ||
698  (type & UserDataOperator::FACESLAVE))) {
699  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
700  "Wrong combination of FaceType and OpType, OPROWCOL "
701  "combined with face-face OpType");
702  }
703 
704  if (!type) {
705  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
706  "Face type is not set");
707  }
708  } else {
709  type = UserDataOperator::FACEMASTER | UserDataOperator::FACESLAVE |
710  UserDataOperator::FACEMASTERMASTER |
711  UserDataOperator::FACEMASTERSLAVE |
712  UserDataOperator::FACESLAVEMASTER |
713  UserDataOperator::FACESLAVESLAVE;
714  }
715 
716  if (oit->getOpType() & UserDataOperator::OPROW &&
717  (type & UserDataOperator::FACEMASTER)) {
718  try {
719  CHKERR oit->opRhs(*op_master_data[0], false);
720  }
721  CATCH_OP_ERRORS(*oit);
722  }
723 
724  if (oit->getOpType() & UserDataOperator::OPROW &&
725  (type & UserDataOperator::FACESLAVE)) {
726  try {
727  CHKERR oit->opRhs(*op_slave_data[0], false);
728  }
729  CATCH_OP_ERRORS(*oit);
730  }
731 
732  if (oit->getOpType() & UserDataOperator::OPCOL &&
733  (type & UserDataOperator::FACEMASTER)) {
734  try {
735  CHKERR oit->opRhs(*op_master_data[1], false);
736  }
737  CATCH_OP_ERRORS(*oit);
738  }
739 
740  if (oit->getOpType() & UserDataOperator::OPCOL &&
741  (type & UserDataOperator::FACESLAVE)) {
742  try {
743  CHKERR oit->opRhs(*op_slave_data[1], false);
744  }
745  CATCH_OP_ERRORS(*oit);
746  }
747 
748  if (oit->getOpType() & UserDataOperator::OPROWCOL &&
749  (type & UserDataOperator::FACEMASTERMASTER)) {
750  try {
751  CHKERR oit->opLhs(*op_master_data[0], *op_master_data[1]);
752  }
753  CATCH_OP_ERRORS(*oit);
754  }
755 
756  if (oit->getOpType() & UserDataOperator::OPROWCOL &&
757  (type & UserDataOperator::FACEMASTERSLAVE)) {
758  try {
759  CHKERR oit->opLhs(*op_master_data[0], *op_slave_data[1]);
760  }
761  CATCH_OP_ERRORS(*oit);
762  }
763 
764  if (oit->getOpType() & UserDataOperator::OPROWCOL &&
765  (type & UserDataOperator::FACESLAVEMASTER)) {
766  try {
767  CHKERR oit->opLhs(*op_slave_data[0], *op_master_data[1]);
768  }
769  CATCH_OP_ERRORS(*oit);
770  }
771 
772  if (oit->getOpType() & UserDataOperator::OPROWCOL &&
773  (type & UserDataOperator::FACESLAVESLAVE)) {
774  try {
775  CHKERR oit->opLhs(*op_slave_data[0], *op_slave_data[1]);
776  }
777  CATCH_OP_ERRORS(*oit);
778  }
779  }
780  }
781  MoFEMFunctionReturn(0);
782 }
783 
784 MoFEMErrorCode ContactPrismElementForcesAndSourcesCore::getEntityFieldData(
785  EntitiesFieldData &master_data, EntitiesFieldData &slave_data,
786  const std::string &field_name, const EntityType type_lo,
787  const EntityType type_hi) const {
788  MoFEMFunctionBegin;
789 
790  auto reset_data = [type_lo, type_hi](auto &data) {
791  for (EntityType t = type_lo; t != type_hi; ++t) {
792  for (auto &dat : data.dataOnEntities[t]) {
793  dat.getOrder() = 0;
794  dat.getBase() = NOBASE;
795  dat.getSpace() = NOSPACE;
796  dat.getFieldData().resize(0, false);
797  dat.getFieldDofs().resize(0, false);
798  dat.getFieldEntities().resize(0, false);
799  }
800  }
801  };
802  reset_data(master_data);
803  reset_data(slave_data);
804 
805  auto &field_ents = getDataFieldEnts();
806  auto bit_number = mField.get_field_bit_number(field_name);
807  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
808  bit_number, get_id_for_min_type(type_lo));
809  auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
810  cmp_uid_lo);
811  if (lo != field_ents.end()) {
812  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
813  bit_number, get_id_for_max_type(type_hi));
814  auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
815  if (lo != hi) {
816  for (auto it = lo; it != hi; ++it)
817  if (auto e = it->lock()) {
818 
819  auto get_data = [&](auto &data, auto type, auto side) {
820  auto &dat = data.dataOnEntities[type][side];
821  auto &ent_field_dofs = dat.getFieldDofs();
822  auto &ent_field_data = dat.getFieldData();
823  dat.getFieldEntities().resize(1, false);
824  dat.getFieldEntities()[0] = e.get();
825  dat.getBase() = e->getApproxBase();
826  dat.getSpace() = e->getSpace();
827  const int ent_order = e->getMaxOrder();
828  dat.getOrder() =
829  dat.getOrder() > ent_order ? dat.getOrder() : ent_order;
830  const auto dof_ent_field_data = e->getEntFieldData();
831  const int nb_dofs_on_ent = e->getNbDofsOnEnt();
832  ent_field_data.resize(nb_dofs_on_ent, false);
833  noalias(ent_field_data) = e->getEntFieldData();
834  ent_field_dofs.resize(nb_dofs_on_ent, false);
835  std::fill(ent_field_dofs.begin(), ent_field_dofs.end(), nullptr);
836  if (auto cache = e->entityCacheDataDofs.lock()) {
837  for (auto dit = cache->loHi[0]; dit != cache->loHi[1]; ++dit) {
838  ent_field_dofs[(*dit)->getEntDofIdx()] =
839  reinterpret_cast<FEDofEntity *>((*dit).get());
840  }
841  }
842  };
843 
844  const EntityType type = e->getEntType();
845  const int side = e->getSideNumberPtr()->side_number;
846 
847  switch (type) {
848  case MBEDGE:
849 
850  if (side < 3)
851  get_data(master_data, type, side);
852  else if (side > 5)
853  get_data(slave_data, type, side - 6);
854 
855  break;
856  case MBTRI:
857 
858  if (side == 3)
859  get_data(master_data, type, 0);
860  if (side == 4)
861  get_data(slave_data, type, 0);
862 
863  break;
864  default:
865  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
866  "Entity type not implemented (FIXME)");
867  };
868 
869  const int brother_side = e->getSideNumberPtr()->brother_side_number;
870  if (brother_side != -1)
871  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
872  "Case with brother side not implemented (FIXME)");
873  }
874  }
875  }
876 
877  MoFEMFunctionReturn(0);
878 }
879 
880 MoFEMErrorCode ContactPrismElementForcesAndSourcesCore::getNodesFieldData(
881  const std::string field_name, VectorDouble &master_nodes_data,
882  VectorDouble &slave_nodes_data, VectorDofs &master_nodes_dofs,
883  VectorDofs &slave_nodes_dofs,
884 
885  VectorFieldEntities &master_field_entities,
886  VectorFieldEntities &slave_field_entities,
887 
888  FieldSpace &master_space, FieldSpace &slave_space,
889  FieldApproximationBase &master_base,
890  FieldApproximationBase &slave_base) const {
891  MoFEMFunctionBegin;
892 
893  auto set_zero = [](auto &nodes_data, auto &nodes_dofs, auto &field_entities) {
894  nodes_data.resize(0, false);
895  nodes_dofs.resize(0, false);
896  field_entities.resize(0, false);
897  };
898  set_zero(master_nodes_data, master_nodes_dofs, master_field_entities);
899  set_zero(slave_nodes_data, slave_nodes_dofs, slave_field_entities);
900 
901  auto field_it = fieldsPtr->get<FieldName_mi_tag>().find(field_name);
902  if (field_it != fieldsPtr->get<FieldName_mi_tag>().end()) {
903 
904  auto bit_number = (*field_it)->getBitNumber();
905  const int nb_dofs_on_vert = (*field_it)->getNbOfCoeffs();
906  master_space = slave_space = (*field_it)->getSpace();
907  master_base = slave_base = (*field_it)->getApproxBase();
908 
909  auto &field_ents = getDataFieldEnts();
910  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
911  bit_number, get_id_for_min_type<MBVERTEX>());
912  auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
913  cmp_uid_lo);
914  if (lo != field_ents.end()) {
915  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
916  bit_number, get_id_for_max_type<MBVERTEX>());
917  auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
918  if (lo != hi) {
919 
920  int nb_dofs = 0;
921  for (auto it = lo; it != hi; ++it) {
922  if (auto e = it->lock()) {
923  nb_dofs += e->getEntFieldData().size();
924  }
925  }
926 
927  if (nb_dofs) {
928 
929  auto init_set = [&](auto &nodes_data, auto &nodes_dofs,
930  auto &field_entities) {
931  constexpr int num_nodes = 3;
932  const int max_nb_dofs = nb_dofs_on_vert * num_nodes;
933  nodes_data.resize(max_nb_dofs, false);
934  nodes_dofs.resize(max_nb_dofs, false);
935  field_entities.resize(num_nodes, false);
936  nodes_data.clear();
937  fill(nodes_dofs.begin(), nodes_dofs.end(), nullptr);
938  fill(field_entities.begin(), field_entities.end(), nullptr);
939  };
940 
941  init_set(master_nodes_data, master_nodes_dofs, master_field_entities);
942  init_set(slave_nodes_data, slave_nodes_dofs, slave_field_entities);
943 
944  for (auto it = lo; it != hi; ++it) {
945  if (auto e = it->lock()) {
946 
947  const auto &sn = e->getSideNumberPtr();
948  int side = sn->side_number;
949 
950  auto set_data = [&](auto &nodes_data, auto &nodes_dofs,
951  auto &field_entities, int side, int pos) {
952  field_entities[side] = e.get();
953  if (auto cache = e->entityCacheDataDofs.lock()) {
954  for (auto dit = cache->loHi[0]; dit != cache->loHi[1];
955  ++dit) {
956  const auto dof_idx = (*dit)->getEntDofIdx();
957  nodes_data[pos + dof_idx] = (*dit)->getFieldData();
958  nodes_dofs[pos + dof_idx] =
959  reinterpret_cast<FEDofEntity *>((*dit).get());
960  }
961  }
962  };
963 
964  if (side < 3)
965  set_data(master_nodes_data, master_nodes_dofs,
966  master_field_entities, side, side * nb_dofs_on_vert);
967  else
968  set_data(slave_nodes_data, slave_nodes_dofs,
969  slave_field_entities, (side - 3),
970  (side - 3) * nb_dofs_on_vert);
971 
972  const int brother_side = sn->brother_side_number;
973  if (brother_side != -1)
974  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
975  "Not implemented (FIXME please)");
976  }
977  }
978  }
979  }
980  }
981  }
982 
983  MoFEMFunctionReturn(0);
984 }
985 
986 template <typename EXTRACTOR>
987 MoFEMErrorCode ContactPrismElementForcesAndSourcesCore::getEntityIndices(
988  EntitiesFieldData &master_data, EntitiesFieldData &slave_data,
989  const std::string &field_name, FieldEntity_vector_view &ents_field,
990  const EntityType type_lo, const EntityType type_hi,
991  EXTRACTOR &&extractor) const {
992  MoFEMFunctionBegin;
993 
994  auto clear_data = [type_lo, type_hi](auto &data) {
995  for (EntityType t = type_lo; t != type_hi; ++t) {
996  for (auto &d : data.dataOnEntities[t]) {
997  d.getIndices().resize(0, false);
998  d.getLocalIndices().resize(0, false);
999  }
1000  }
1001  };
1002 
1003  clear_data(master_data);
1004  clear_data(slave_data);
1005 
1006  auto bit_number = mField.get_field_bit_number(field_name);
1007  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
1008  bit_number, get_id_for_min_type(type_lo));
1009  auto lo = std::lower_bound(ents_field.begin(), ents_field.end(), lo_uid,
1010  cmp_uid_lo);
1011  if (lo != ents_field.end()) {
1012  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
1013  bit_number, get_id_for_max_type(type_hi));
1014  auto hi = std::upper_bound(lo, ents_field.end(), hi_uid, cmp_uid_hi);
1015  if (lo != hi) {
1016 
1017  std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
1018 
1019  for (auto it = lo; it != hi; ++it)
1020  if (auto e = it->lock()) {
1021 
1022  const EntityType type = e->getEntType();
1023  const int side = e->getSideNumberPtr()->side_number;
1024  const int brother_side = e->getSideNumberPtr()->brother_side_number;
1025  if (brother_side != -1)
1026  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1027  "Not implemented case");
1028 
1029  auto get_indices = [&](auto &data, const auto type, const auto side) {
1030  if (auto cache = extractor(e).lock()) {
1031  for (auto dit = cache->loHi[0]; dit != cache->loHi[1]; ++dit) {
1032  auto &dof = (**dit);
1033  auto &dat = data.dataOnEntities[type][side];
1034  auto &ent_field_indices = dat.getIndices();
1035  auto &ent_field_local_indices = dat.getLocalIndices();
1036  if (ent_field_indices.empty()) {
1037  const int nb_dofs_on_ent = dof.getNbDofsOnEnt();
1038  ent_field_indices.resize(nb_dofs_on_ent, false);
1039  ent_field_local_indices.resize(nb_dofs_on_ent, false);
1040  std::fill(ent_field_indices.data().begin(),
1041  ent_field_indices.data().end(), -1);
1042  std::fill(ent_field_local_indices.data().begin(),
1043  ent_field_local_indices.data().end(), -1);
1044  }
1045  const int idx = dof.getEntDofIdx();
1046  ent_field_indices[idx] = dof.getPetscGlobalDofIdx();
1047  ent_field_local_indices[idx] = dof.getPetscLocalDofIdx();
1048  }
1049  }
1050  };
1051 
1052  switch (type) {
1053  case MBEDGE:
1054 
1055  if (side < 3)
1056  get_indices(master_data, type, side);
1057  else if (side > 5)
1058  get_indices(slave_data, type, side - 6);
1059 
1060  break;
1061  case MBTRI:
1062 
1063  if (side == 3)
1064  get_indices(master_data, type, 0);
1065  if (side == 4)
1066  get_indices(slave_data, type, 0);
1067 
1068  break;
1069  default:
1070  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1071  "Entity type not implemented");
1072  }
1073  }
1074  }
1075  }
1076 
1077  MoFEMFunctionReturn(0);
1078 }
1079 
1080 template <typename EXTRACTOR>
1081 MoFEMErrorCode ContactPrismElementForcesAndSourcesCore::getNodesIndices(
1082  const std::string field_name, FieldEntity_vector_view &ents_field,
1083  VectorInt &master_nodes_indices, VectorInt &master_local_nodes_indices,
1084  VectorInt &slave_nodes_indices, VectorInt &slave_local_nodes_indices,
1085  EXTRACTOR &&extractor) const {
1086  MoFEMFunctionBegin;
1087 
1088  master_nodes_indices.resize(0, false);
1089  master_local_nodes_indices.resize(0, false);
1090  slave_nodes_indices.resize(0, false);
1091  slave_local_nodes_indices.resize(0, false);
1092 
1093  auto field_it = fieldsPtr->get<FieldName_mi_tag>().find(field_name);
1094  if (field_it != fieldsPtr->get<FieldName_mi_tag>().end()) {
1095 
1096  auto bit_number = (*field_it)->getBitNumber();
1097  const int nb_dofs_on_vert = (*field_it)->getNbOfCoeffs();
1098 
1099  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
1100  bit_number, get_id_for_min_type<MBVERTEX>());
1101  auto lo = std::lower_bound(ents_field.begin(), ents_field.end(), lo_uid,
1102  cmp_uid_lo);
1103  if (lo != ents_field.end()) {
1104  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
1105  bit_number, get_id_for_max_type<MBVERTEX>());
1106  auto hi = std::upper_bound(lo, ents_field.end(), hi_uid, cmp_uid_hi);
1107  if (lo != hi) {
1108 
1109  int nb_dofs = 0;
1110  for (auto it = lo; it != hi; ++it) {
1111  if (auto e = it->lock()) {
1112  if (auto cache = extractor(e).lock()) {
1113  nb_dofs += std::distance(cache->loHi[0], cache->loHi[1]);
1114  }
1115  }
1116  }
1117 
1118  if (nb_dofs) {
1119 
1120  constexpr int num_nodes = 3;
1121  const int max_nb_dofs = nb_dofs_on_vert * num_nodes;
1122 
1123  auto set_vec_size = [&](auto &nodes_indices,
1124  auto &local_nodes_indices) {
1125  nodes_indices.resize(max_nb_dofs, false);
1126  local_nodes_indices.resize(max_nb_dofs, false);
1127  std::fill(nodes_indices.begin(), nodes_indices.end(), -1);
1128  std::fill(local_nodes_indices.begin(), local_nodes_indices.end(),
1129  -1);
1130  };
1131 
1132  set_vec_size(master_nodes_indices, master_local_nodes_indices);
1133  set_vec_size(slave_nodes_indices, slave_local_nodes_indices);
1134 
1135  for (auto it = lo; it != hi; ++it) {
1136  if (auto e = it->lock()) {
1137 
1138  const int side = e->getSideNumberPtr()->side_number;
1139  const int brother_side =
1140  e->getSideNumberPtr()->brother_side_number;
1141  if (brother_side != -1)
1142  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1143  "Not implemented case");
1144 
1145  auto get_indices = [&](auto &nodes_indices,
1146  auto &local_nodes_indices,
1147  const auto side) {
1148  if (auto cache = extractor(e).lock()) {
1149  for (auto dit = cache->loHi[0]; dit != cache->loHi[1];
1150  ++dit) {
1151  const int idx = (*dit)->getPetscGlobalDofIdx();
1152  const int local_idx = (*dit)->getPetscLocalDofIdx();
1153  const int pos =
1154  side * nb_dofs_on_vert + (*dit)->getDofCoeffIdx();
1155  nodes_indices[pos] = idx;
1156  local_nodes_indices[pos] = local_idx;
1157  }
1158  }
1159  };
1160 
1161  if (side < 3)
1162  get_indices(master_nodes_indices, master_local_nodes_indices,
1163  side);
1164  else if (side > 2)
1165  get_indices(slave_nodes_indices, slave_local_nodes_indices,
1166  side - 3);
1167  else
1168  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1169  "Impossible case");
1170  }
1171  }
1172  }
1173  }
1174  }
1175  }
1176 
1177  MoFEMFunctionReturn(0);
1178 }
1179 
1180 MoFEMErrorCode
1181 ContactPrismElementForcesAndSourcesCore::UserDataOperator::loopSideVolumes(
1182  const string fe_name,
1183  VolumeElementForcesAndSourcesCoreOnContactPrismSide &fe_method,
1184  const int side_type, const EntityHandle ent_for_side) {
1185  return loopSide(fe_name, &fe_method, side_type, ent_for_side);
1186 }
1187 
1188 } // namespace MoFEM
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:460
MoFEM::FEMethod::getRowFieldEnts
auto & getRowFieldEnts() const
Definition: LoopMethods.hpp:423
MoFEM::BasicMethod::getNinTheLoop
int getNinTheLoop() const
get number of evaluated element in the loop
Definition: LoopMethods.hpp:207
MoFEM::ContactPrismElementForcesAndSourcesCore::tangentSlaveOne
VectorDouble tangentSlaveOne
Definition: ContactPrismElementForcesAndSourcesCore.hpp:68
MoFEM::ContactPrismElementForcesAndSourcesCore::tangentMasterOne
VectorDouble tangentMasterOne
Definition: ContactPrismElementForcesAndSourcesCore.hpp:69
MoFEM::ContactPrismElementForcesAndSourcesCore::derivedDataOnSlave
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > derivedDataOnSlave
Definition: ContactPrismElementForcesAndSourcesCore.hpp:92
MoFEM::ContactPrismElementForcesAndSourcesCore::UserDataOperator::FACEMASTERSLAVE
@ FACEMASTERSLAVE
Definition: ContactPrismElementForcesAndSourcesCore.hpp:255
CATCH_OP_ERRORS
#define CATCH_OP_ERRORS(OP)
Definition: ForcesAndSourcesCore.cpp:1355
H1
@ H1
continuous field
Definition: definitions.h:85
MoFEM::ContactPrismElementForcesAndSourcesCore::UserDataOperator::FACESLAVESLAVE
@ FACESLAVESLAVE
Definition: ContactPrismElementForcesAndSourcesCore.hpp:257
MoFEM::ForcesAndSourcesCore::getNoFieldFieldData
MoFEMErrorCode getNoFieldFieldData(const int bit_number, VectorDouble &ent_field_data, VectorDofs &ent_field_dofs, VectorFieldEntities &ent_field) const
Get field data on nodes.
Definition: ForcesAndSourcesCore.cpp:855
MoFEM::DataOperator::opRhs
virtual MoFEMErrorCode opRhs(EntitiesFieldData &data, const bool error_if_no_base=false)
Definition: DataOperators.cpp:140
LASTBASE
@ LASTBASE
Definition: definitions.h:69
MoFEM::ForcesAndSourcesCore::getMaxRowOrder
int getMaxRowOrder() const
Get max order of approximation for field in rows.
Definition: ForcesAndSourcesCore.cpp:131
FTensor::Tensor1
Definition: Tensor1_value.hpp:8
MoFEM::ContactPrismElementForcesAndSourcesCore::opContravariantTransform
OpSetContravariantPiolaTransformOnFace opContravariantTransform
Definition: ContactPrismElementForcesAndSourcesCore.hpp:70
EntityHandle
MoFEM::ContactPrismElementForcesAndSourcesCore::derivedDataOnMaster
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > derivedDataOnMaster
Entity data on element entity columns fields.
Definition: ContactPrismElementForcesAndSourcesCore.hpp:90
MoFEM::ForcesAndSourcesCore::UserDataOperator::OPSPACE
@ OPSPACE
operator do Work is execute on space data
Definition: ForcesAndSourcesCore.hpp:570
MoFEM::EntitiesFieldData::basesOnEntities
std::array< std::bitset< LASTBASE >, MBMAXTYPE > basesOnEntities
bases on entity types
Definition: EntitiesFieldData.hpp:51
MoFEM::EntPolynomialBaseCtx
Class used to pass element data to calculate base functions on tet,triangle,edge.
Definition: EntPolynomialBaseCtx.hpp:22
MoFEM::ContactPrismElementForcesAndSourcesCore::ContactPrismElementForcesAndSourcesCore
ContactPrismElementForcesAndSourcesCore(Interface &m_field)
Definition: ContactPrismElementForcesAndSourcesCore.cpp:10
MoFEM::cmp_uid_hi
static auto cmp_uid_hi(const UId &b, const boost::weak_ptr< FieldEntity > &a)
Definition: ForcesAndSourcesCore.cpp:29
NOBASE
@ NOBASE
Definition: definitions.h:59
CHK_THROW_MESSAGE
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
Definition: definitions.h:609
L2
@ L2
field with C-1 continuity
Definition: definitions.h:88
QUAD_2D_TABLE
static QUAD *const QUAD_2D_TABLE[]
Definition: quad.h:175
MoFEM::CoreInterface::get_field_structure
virtual const Field * get_field_structure(const std::string &name, enum MoFEMTypes bh=MF_EXIST) const =0
get field structure
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::OpSetContravariantPiolaTransformOnFace::normalRawPtr
const VectorDouble * normalRawPtr
Definition: DataOperators.hpp:419
MoFEM::Field::getSpace
FieldSpace getSpace() const
Get field approximation space.
Definition: FieldMultiIndices.hpp:151
MoFEM::FieldEntity_vector_view
std::vector< boost::weak_ptr< FieldEntity > > FieldEntity_vector_view
Definition: FieldEntsMultiIndices.hpp:478
MoFEM::ForcesAndSourcesCore::setGaussPts
virtual MoFEMErrorCode setGaussPts(int order_row, int order_col, int order_data)
set user specific integration rule
Definition: ForcesAndSourcesCore.cpp:1961
MoFEM::ContactPrismElementForcesAndSourcesCore::tangentMasterTwo
VectorDouble tangentMasterTwo
Definition: ContactPrismElementForcesAndSourcesCore.hpp:69
MoFEM::CoreInterface::get_field_bit_number
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number
MoFEM::FEMethod::getDataFieldEnts
const FieldEntity_vector_view & getDataFieldEnts() const
Definition: LoopMethods.hpp:412
MoFEM::ForcesAndSourcesCore::getMaxDataOrder
int getMaxDataOrder() const
Get max order of approximation for data fields.
Definition: ForcesAndSourcesCore.cpp:120
MoFEM::Types::BitFieldId
std::bitset< BITFIELDID_SIZE > BitFieldId
Field Id.
Definition: Types.hpp:42
MoFEM::ContactPrismElementForcesAndSourcesCore::dataHdivMaster
EntitiesFieldData & dataHdivMaster
Definition: ContactPrismElementForcesAndSourcesCore.hpp:101
MoFEM::ContactPrismElementForcesAndSourcesCore::getEntityFieldData
MoFEMErrorCode getEntityFieldData(EntitiesFieldData &master_data, EntitiesFieldData &slave_data, const std::string &field_name, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
function that gets entity field data.
Definition: ContactPrismElementForcesAndSourcesCore.cpp:784
MoFEM::FieldName_mi_tag
MultiIndex Tag for field name.
Definition: TagMultiIndices.hpp:67
QUAD_2D_TABLE_SIZE
#define QUAD_2D_TABLE_SIZE
Definition: quad.h:174
MoFEM::ForcesAndSourcesCore::getMaxColOrder
int getMaxColOrder() const
Get max order of approximation for field in columns.
Definition: ForcesAndSourcesCore.cpp:135
FTensor::levi_civita
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use
Definition: Levi_Civita.hpp:617
MoFEM::ForcesAndSourcesCore::getNoFieldRowIndices
MoFEMErrorCode getNoFieldRowIndices(EntitiesFieldData &data, const int bit_number) const
get col NoField indices
Definition: ForcesAndSourcesCore.cpp:459
MoFEM::ForcesAndSourcesCore::UserDataOperator::OPROWCOL
@ OPROWCOL
operator doWork is executed on FE rows &columns
Definition: ForcesAndSourcesCore.hpp:569
MoFEM::EntitiesFieldData::spacesOnEntities
std::array< std::bitset< LASTSPACE >, MBMAXTYPE > spacesOnEntities
spaces on entity types
Definition: EntitiesFieldData.hpp:49
sdf.r
int r
Definition: sdf.py:8
MoFEM::Field
Provide data structure for (tensor) field approximation.
Definition: FieldMultiIndices.hpp:51
order
constexpr int order
Definition: dg_projection.cpp:18
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
MoFEM::EntitiesFieldData::basesOnSpaces
std::array< std::bitset< LASTBASE >, LASTSPACE > basesOnSpaces
base on spaces
Definition: EntitiesFieldData.hpp:53
MoFEM::BasicMethod::fieldsPtr
const Field_multiIndex * fieldsPtr
raw pointer to fields container
Definition: LoopMethods.hpp:252
MoFEM::ContactPrismElementForcesAndSourcesCore::dataOnSlave
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > dataOnSlave
Definition: ContactPrismElementForcesAndSourcesCore.hpp:81
FieldSpace
FieldSpace
approximation spaces
Definition: definitions.h:82
MoFEM::Field::getId
const BitFieldId & getId() const
Get unique field id.
Definition: FieldMultiIndices.hpp:129
MoFEM::ContactPrismElementForcesAndSourcesCore::UserDataOperator::loopSideVolumes
MoFEMErrorCode loopSideVolumes(const string fe_name, VolumeElementForcesAndSourcesCoreOnContactPrismSide &fe_method, const int side_type, const EntityHandle ent_for_side)
Definition: ContactPrismElementForcesAndSourcesCore.cpp:1181
FTensor::Number< 0 >
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::ContactPrismElementForcesAndSourcesCore::getNodesIndices
MoFEMErrorCode getNodesIndices(const std::string field_name, FieldEntity_vector_view &ents_field, VectorInt &master_nodes_indices, VectorInt &master_local_nodes_indices, VectorInt &slave_nodes_indices, VectorInt &slave_local_nodes_indices, EXTRACTOR &&extractor) const
function that gets nodes indices.
Definition: ContactPrismElementForcesAndSourcesCore.cpp:1081
MoFEM::ForcesAndSourcesCore::getNoFieldColIndices
MoFEMErrorCode getNoFieldColIndices(EntitiesFieldData &data, const int bit_number) const
get col NoField indices
Definition: ForcesAndSourcesCore.cpp:474
MoFEM::ContactPrismElementForcesAndSourcesCore::UserDataOperator::FACESLAVEMASTER
@ FACESLAVEMASTER
Definition: ContactPrismElementForcesAndSourcesCore.hpp:256
MoFEM::CoreInterface::get_moab
virtual moab::Interface & get_moab()=0
MoFEM::ContactPrismElementForcesAndSourcesCore::setDefaultGaussPts
MoFEMErrorCode setDefaultGaussPts(const int rule)
Definition: ContactPrismElementForcesAndSourcesCore.cpp:92
MoFEM::ForcesAndSourcesCore::getUserPolynomialBase
auto & getUserPolynomialBase()
Get the User Polynomial Base object.
Definition: ForcesAndSourcesCore.hpp:97
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
MoFEM::VolumeElementForcesAndSourcesCoreOnContactPrismSide
Base volume element used to integrate on contact surface (could be extended to other volume elements)
Definition: VolumeElementForcesAndSourcesCoreOnContactPrismSide.hpp:23
MoFEM::ContactPrismElementForcesAndSourcesCore::getNodesFieldData
MoFEMErrorCode getNodesFieldData(const std::string field_name, VectorDouble &master_nodes_data, VectorDouble &slave_nodes_data, VectorDofs &master_nodes_dofs, VectorDofs &slave_nodes_dofs, VectorFieldEntities &master_field_entities, VectorFieldEntities &slave_field_entities, FieldSpace &master_space, FieldSpace &slave_space, FieldApproximationBase &master_base, FieldApproximationBase &slave_base) const
function that gets nodes field data.
Definition: ContactPrismElementForcesAndSourcesCore.cpp:880
MoFEM::ContactPrismElementForcesAndSourcesCore::normal
VectorDouble normal
vector storing vector normal to master or slave element
Definition: ContactPrismElementForcesAndSourcesCore.hpp:56
convert.type
type
Definition: convert.py:64
MoFEM::FieldEntity::getLocalUniqueIdCalculate
UId getLocalUniqueIdCalculate()
Get the Local Unique Id Calculate object.
Definition: FieldEntsMultiIndices.hpp:136
MoFEM::ForcesAndSourcesCore::UserDataOperator::OPCOL
@ OPCOL
operator doWork function is executed on FE columns
Definition: ForcesAndSourcesCore.hpp:568
MoFEM::ForcesAndSourcesCore::opPtrVector
boost::ptr_deque< UserDataOperator > opPtrVector
Vector of finite element users data operators.
Definition: ForcesAndSourcesCore.hpp:480
MoFEM::ForcesAndSourcesCore::dataOnElement
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > dataOnElement
Entity data on element entity rows fields.
Definition: ForcesAndSourcesCore.hpp:461
MoFEM::ContactPrismElementForcesAndSourcesCore::tangentSlaveTwo
VectorDouble tangentSlaveTwo
Definition: ContactPrismElementForcesAndSourcesCore.hpp:68
MoFEM::cmp_uid_lo
static auto cmp_uid_lo(const boost::weak_ptr< FieldEntity > &a, const UId &b)
Definition: ForcesAndSourcesCore.cpp:18
MoFEM::VectorDofs
ublas::vector< FEDofEntity *, DofsAllocator > VectorDofs
Definition: EntitiesFieldData.hpp:23
MoFEM::ContactPrismElementForcesAndSourcesCore::gaussPtsMaster
MatrixDouble gaussPtsMaster
Definition: ContactPrismElementForcesAndSourcesCore.hpp:63
MoFEM::ContactPrismElementForcesAndSourcesCore::UserDataOperator
default operator for Contact Prism element
Definition: ContactPrismElementForcesAndSourcesCore.hpp:208
QUAD_::npoints
int npoints
Definition: quad.h:29
MoFEM::ContactPrismElementForcesAndSourcesCore::coordsAtGaussPtsSlave
MatrixDouble coordsAtGaussPtsSlave
Definition: ContactPrismElementForcesAndSourcesCore.hpp:60
LASTSPACE
@ LASTSPACE
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:89
MoFEM::ForcesAndSourcesCore::getRule
virtual int getRule(int order_row, int order_col, int order_data)
another variant of getRule
Definition: ForcesAndSourcesCore.cpp:1955
MoFEM::get_id_for_max_type
EntityHandle get_id_for_max_type()
Definition: RefEntsMultiIndices.hpp:13
AINSWORTH_LOBATTO_BASE
@ AINSWORTH_LOBATTO_BASE
Definition: definitions.h:62
MoFEM::Tools::diffShapeFunMBTRI
static constexpr std::array< double, 6 > diffShapeFunMBTRI
Definition: Tools.hpp:104
MoFEM::FEMethod::getFEName
auto getFEName() const
get finite element name
Definition: LoopMethods.hpp:391
MoFEM::ContactPrismElementForcesAndSourcesCore::dataH1Slave
EntitiesFieldData & dataH1Slave
Definition: ContactPrismElementForcesAndSourcesCore.hpp:95
MoFEM::ContactPrismElementForcesAndSourcesCore::getEntityIndices
MoFEMErrorCode getEntityIndices(EntitiesFieldData &master_data, EntitiesFieldData &slave_data, const std::string &field_name, FieldEntity_vector_view &ents_field, const EntityType type_lo, const EntityType type_hi, EXTRACTOR &&extractor) const
function that gets entity indices.
Definition: ContactPrismElementForcesAndSourcesCore.cpp:987
MoFEM::EntitiesFieldData::brokenBasesOnSpaces
std::array< std::bitset< LASTBASE >, LASTSPACE > brokenBasesOnSpaces
base on spaces
Definition: EntitiesFieldData.hpp:55
MoFEM::ContactPrismElementForcesAndSourcesCore::coords
VectorDouble coords
Definition: ContactPrismElementForcesAndSourcesCore.hpp:57
t
constexpr double t
plate stiffness
Definition: plate.cpp:58
MoFEM::EntitiesFieldData::bAse
std::bitset< LASTBASE > bAse
bases on element
Definition: EntitiesFieldData.hpp:42
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
MoFEM::ForcesAndSourcesCore::getSpacesAndBaseOnEntities
MoFEMErrorCode getSpacesAndBaseOnEntities(EntitiesFieldData &data) const
Get field approximation space and base on entities.
Definition: ForcesAndSourcesCore.cpp:990
EntData
EntitiesFieldData::EntData EntData
Definition: child_and_parent.cpp:37
field_name
constexpr auto field_name
Definition: poisson_2d_homogeneous.cpp:13
FTensor::Index< 'i', 3 >
MoFEM::ForcesAndSourcesCore
structure to get information form mofem into EntitiesFieldData
Definition: ForcesAndSourcesCore.hpp:22
MoFEM::ForcesAndSourcesCore::dataH1
EntitiesFieldData & dataH1
Definition: ForcesAndSourcesCore.hpp:471
v
const double v
phase velocity of light in medium (cm/ns)
Definition: initial_diffusion.cpp:40
FTensor::dd
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
DEMKOWICZ_JACOBI_BASE
@ DEMKOWICZ_JACOBI_BASE
Definition: definitions.h:66
MoFEM::ContactPrismElementForcesAndSourcesCore::nbGaussPts
int nbGaussPts
Definition: ContactPrismElementForcesAndSourcesCore.hpp:202
MoFEM::ForcesAndSourcesCore::UserDataOperator::loopSide
MoFEMErrorCode loopSide(const string &fe_name, ForcesAndSourcesCore *side_fe, const size_t dim, const EntityHandle ent_for_side=0, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy, AdjCache *adj_cache=nullptr)
User calls this function to loop over elements on the side of face. This function calls finite elemen...
Definition: ForcesAndSourcesCore.cpp:1700
MoFEM::Types::VectorInt
UBlasVector< int > VectorInt
Definition: Types.hpp:67
j
FTensor::Index< 'j', 3 > j
Definition: matrix_function.cpp:19
MoFEM::ContactPrismElementForcesAndSourcesCore::dataH1Master
EntitiesFieldData & dataH1Master
Definition: ContactPrismElementForcesAndSourcesCore.hpp:94
MoFEM::ContactPrismElementForcesAndSourcesCore::operator()
MoFEMErrorCode operator()()
function is run for every finite element
Definition: ContactPrismElementForcesAndSourcesCore.cpp:140
MoFEM::ContactPrismElementForcesAndSourcesCore::coordsAtGaussPtsMaster
MatrixDouble coordsAtGaussPtsMaster
Definition: ContactPrismElementForcesAndSourcesCore.hpp:58
MoFEM::ContactPrismElementForcesAndSourcesCore::dataOnMaster
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > dataOnMaster
Entity data on element entity rows fields.
Definition: ContactPrismElementForcesAndSourcesCore.hpp:80
MoFEM::OpSetContravariantPiolaTransformOnFace::normalShift
int normalShift
Shift in vector for linear geometry.
Definition: DataOperators.hpp:431
UBlasVector< double >
AINSWORTH_LEGENDRE_BASE
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
HCURL
@ HCURL
field with continuous tangents
Definition: definitions.h:86
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
MoFEM::FEMethod::getColFieldEnts
auto & getColFieldEnts() const
Definition: LoopMethods.hpp:431
sdf_hertz_2d_axisymm_plane.d
float d
Definition: sdf_hertz_2d_axisymm_plane.py:4
FieldApproximationBase
FieldApproximationBase
approximation base
Definition: definitions.h:58
MoFEM::ContactPrismElementForcesAndSourcesCore::dataHdivSlave
EntitiesFieldData & dataHdivSlave
Definition: ContactPrismElementForcesAndSourcesCore.hpp:102
MoFEM::EntitiesFieldData::dataOnEntities
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
Definition: EntitiesFieldData.hpp:57
MoFEM::VectorFieldEntities
ublas::vector< FieldEntity *, FieldEntAllocator > VectorFieldEntities
Definition: EntitiesFieldData.hpp:31
MoFEM::FEMethod::numeredEntFiniteElementPtr
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
Definition: LoopMethods.hpp:383
MoFEM::ContactPrismElementForcesAndSourcesCore::UserDataOperator::FACEMASTERMASTER
@ FACEMASTERMASTER
Definition: ContactPrismElementForcesAndSourcesCore.hpp:254
MoFEM::field_it
Field_multiIndex::index< FieldName_mi_tag >::type::iterator field_it
Definition: Projection10NodeCoordsOnField.cpp:123
MoFEM::FEDofEntity
keeps information about indexed dofs for the finite element
Definition: DofsMultiIndices.hpp:267
MoFEM::EntitiesFieldData
data structure for finite element entity
Definition: EntitiesFieldData.hpp:40
MoFEM::ContactPrismElementForcesAndSourcesCore::loopOverOperators
MoFEMErrorCode loopOverOperators()
Iterate user data operators.
Definition: ContactPrismElementForcesAndSourcesCore.cpp:482
CONTINUOUS
@ CONTINUOUS
Regular field.
Definition: definitions.h:100
k
FTensor::Index< 'k', 3 > k
Definition: matrix_function.cpp:20
MoFEM::get_id_for_min_type
EntityHandle get_id_for_min_type()
Definition: RefEntsMultiIndices.hpp:18
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
HDIV
@ HDIV
field with continuous normal traction
Definition: definitions.h:87
MoFEM::ContactPrismElementForcesAndSourcesCore::aRea
std::array< double, 2 > aRea
Array storing master and slave faces areas.
Definition: ContactPrismElementForcesAndSourcesCore.hpp:53
MoFEM::Tools::getTriNormal
static MoFEMErrorCode getTriNormal(const double *coords, double *normal, double *d_normal=nullptr)
Get the Tri Normal objectGet triangle normal.
Definition: Tools.cpp:353
MOFEM_NOT_IMPLEMENTED
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
MoFEM::ForcesAndSourcesCore::mField
Interface & mField
Definition: ForcesAndSourcesCore.hpp:24
MoFEM::ContactPrismElementForcesAndSourcesCore::UserDataOperator::FACESLAVE
@ FACESLAVE
Definition: ContactPrismElementForcesAndSourcesCore.hpp:253
NOFIELD
@ NOFIELD
scalar or vector of scalars describe (no true field)
Definition: definitions.h:84
MoFEM::ForcesAndSourcesCore::UserDataOperator::OPROW
@ OPROW
operator doWork function is executed on FE rows
Definition: ForcesAndSourcesCore.hpp:567
MoFEM::ContactPrismElementForcesAndSourcesCore::gaussPtsSlave
MatrixDouble gaussPtsSlave
Definition: ContactPrismElementForcesAndSourcesCore.hpp:65
MoFEM::ContactPrismElementForcesAndSourcesCore::UserDataOperator::FACEMASTER
@ FACEMASTER
Definition: ContactPrismElementForcesAndSourcesCore.hpp:252
Generated by Doxygen 1.8.17 and hosted at University of Glasgow