v0.14.0
Loading...
Searching...
No Matches
EntitiesFieldData.cpp
Go to the documentation of this file.
1/** \file EntitiesFieldData.cpp
2\brief Implementation for Data Structures in Forces and Sources
3
4*/
5
6
7
8namespace MoFEM {
9
10EntitiesFieldData::EntData::EntData(const bool allocate_base_matrices)
11 : sEnse(0), oRder(0), bAse(NOBASE), entDataBitRefLevel(),
12 N(baseFunctionsAndBaseDerivatives[ZeroDerivative]),
13 diffN(baseFunctionsAndBaseDerivatives[FirstDerivative]) {
14 if (allocate_base_matrices) {
15
16 for (auto d = 0; d != LastDerivative; ++d) {
17 for (int b = 0; b != LASTBASE; ++b) {
19 }
20 }
21 }
22}
23
24int EntitiesFieldData::EntData::getSense() const { return sEnse; }
25
26boost::shared_ptr<MatrixDouble> &
28 const BaseDerivatives direvatie) {
29 return baseFunctionsAndBaseDerivatives[direvatie][base];
30}
31
32boost::shared_ptr<MatrixDouble> &
34 return N[base];
35}
36
37boost::shared_ptr<MatrixDouble> &EntitiesFieldData::EntData::getDiffNSharedPtr(
38 const FieldApproximationBase base) {
39 return diffN[base];
40}
41
42static void constructor_data(EntitiesFieldData *data, const EntityType type) {
43
45
46 auto set_default = [&]() {
47 std::array<size_t, MBMAXTYPE> count;
48 std::fill(count.begin(), count.end(), 0);
49 const int dim_type = moab::CN::Dimension(type);
50 data->dataOnEntities[MBVERTEX].resize(1);
51 if (type != MBVERTEX) {
52 for (auto dd = dim_type; dd > 0; --dd) {
53 int nb_ents = moab::CN::NumSubEntities(type, dd);
54 for (int ii = 0; ii != nb_ents; ++ii) {
55 auto sub_ent_type = moab::CN::SubEntityType(type, dd, ii);
56 count[sub_ent_type] = nb_ents;
57 }
58 for (auto tt = moab::CN::TypeDimensionMap[dd].first;
59 tt <= moab::CN::TypeDimensionMap[dd].second; ++tt) {
60 data->dataOnEntities[tt].resize(count[tt]);
61 }
62 }
63 }
64 };
65
66 switch (type) {
67 case MBENTITYSET:
68 break;
69
70 default:
71 set_default();
72 }
73}
74
76 constructor_data(this, type);
77}
78
81 constructor_data(this, type);
83}
84
85static void
87 const boost::shared_ptr<EntitiesFieldData> &data_ptr) {
88
90 using DerivedEntData = DerivedEntitiesFieldData::DerivedEntData;
91
92 for (int tt = MBVERTEX; tt != MBMAXTYPE; ++tt) {
93 auto &ent_data = data_ptr->dataOnEntities[tt];
94 auto &derived_ent_data = derived_data->dataOnEntities[tt];
95 for (auto c = derived_ent_data.size(); c < ent_data.size(); ++c) {
96 boost::shared_ptr<EntData> ent_data_ptr(data_ptr, &ent_data[c]);
97 derived_ent_data.push_back(new DerivedEntData(ent_data_ptr));
98 }
99 derived_ent_data.resize(ent_data.size());
100 }
101}
102
104 const boost::shared_ptr<EntitiesFieldData> &data_ptr)
105 : EntitiesFieldData(), dataPtr(data_ptr) {
107}
108
113}
114
117 sPace = NOSPACE;
118 bAse = NOBASE;
119 fieldEntities.resize(0, false);
120 iNdices.resize(0, false);
121 localIndices.resize(0, false);
122 dOfs.resize(0, false);
123 fieldData.resize(0, false);
125}
126
129 for (EntityType t = MBVERTEX; t != MBMAXTYPE; t++)
130 for (auto &e : dataOnEntities[t])
131 CHKERR e.resetFieldDependentData();
133}
134
137 const FieldApproximationBase base) {
139 auto make_swap = [](boost::shared_ptr<MatrixDouble> &ptr,
140 boost::shared_ptr<MatrixDouble> &ptrBB,
141 boost::shared_ptr<MatrixDouble> &swap_ptr) {
142 if (swap_ptr) {
143 ptr = swap_ptr;
144 swap_ptr.reset();
145 } else {
146 swap_ptr = ptr;
147 ptr = ptrBB;
148 }
149 };
150 make_swap(getNSharedPtr(base), getBBNSharedPtr(field_name), swapBaseNPtr);
151 make_swap(getDiffNSharedPtr(base), getBBDiffNSharedPtr(field_name),
152 swapBaseDiffNPtr);
154}
155
157 const FieldApproximationBase base) {
159 // Note: Do not swap bases on entities sets
160 for (int tt = MBVERTEX; tt != MBENTITYSET; ++tt) {
161 auto &ent_data = dataOnEntities[tt];
162 for (auto &side_data : ent_data)
163 CHKERR side_data.baseSwap(field_name, base);
164 }
166}
167
169 const std::string &field_name, const FieldApproximationBase base) {
171 auto make_swap = [](boost::shared_ptr<MatrixDouble> &ptr,
172 boost::shared_ptr<MatrixDouble> &ptrBB,
173 boost::shared_ptr<MatrixDouble> &swap_ptr) {
174 if (swap_ptr) {
175 ptr = swap_ptr;
176 swap_ptr.reset();
177 } else {
178 swap_ptr = ptr;
179 ptr = ptrBB;
180 }
181 };
182 make_swap(getDerivedNSharedPtr(base), getBBNSharedPtr(field_name),
183 swapBaseNPtr);
184 make_swap(getDerivedDiffNSharedPtr(base), getBBDiffNSharedPtr(field_name),
185 swapBaseDiffNPtr);
187}
188
190 const boost::shared_ptr<EntitiesFieldData::EntData> &ent_data_ptr)
191 : EntitiesFieldData::EntData(false), entDataPtr(ent_data_ptr) {}
192
194 return entDataPtr->getSense();
195}
196
197boost::shared_ptr<MatrixDouble> &
199 const FieldApproximationBase base, const BaseDerivatives derivative) {
200 if (baseFunctionsAndBaseDerivatives[derivative][base])
201 return baseFunctionsAndBaseDerivatives[derivative][base];
202 else
203 return entDataPtr->getNSharedPtr(base, derivative);
204}
205
206boost::shared_ptr<MatrixDouble> &
208 const FieldApproximationBase base) {
209 if (N[base])
210 return N[base];
211 else
212 return entDataPtr->getNSharedPtr(base);
213}
214
215boost::shared_ptr<MatrixDouble> &
217 const FieldApproximationBase base) {
218 if (diffN[base])
219 return diffN[base];
220 else
221 return entDataPtr->getDiffNSharedPtr(base);
222}
223const boost::shared_ptr<MatrixDouble> &
225 const FieldApproximationBase base) const {
226 if (N[base])
227 return N[base];
228 else
229 return entDataPtr->getNSharedPtr(base);
230}
231const boost::shared_ptr<MatrixDouble> &
233 const FieldApproximationBase base) const {
234 if (diffN[base])
235 return diffN[base];
236 else
237 return entDataPtr->getDiffNSharedPtr(base);
238}
239
240std::ostream &operator<<(std::ostream &os,
242 os << "sEnse: " << e.getSense() << std::endl
243 << "oRder: " << e.getOrder() << std::endl
244 << "global indices: " << e.getIndices() << std::endl
245 << "local indices: " << e.getLocalIndices() << std::endl;
246 // FIXME: precision should not be set here
247 os << "fieldData: " << std::fixed << std::setprecision(2) << e.getFieldData()
248 << std::endl;
249 MatrixDouble base = const_cast<EntitiesFieldData::EntData &>(e).getN();
250 MatrixDouble diff_base =
251 const_cast<EntitiesFieldData::EntData &>(e).getDiffN();
252 const double eps = 1e-6;
253 for (unsigned int ii = 0; ii != base.size1(); ii++) {
254 for (unsigned int jj = 0; jj != base.size2(); jj++) {
255 if (fabs(base(ii, jj)) < eps)
256 base(ii, jj) = 0;
257 }
258 }
259 for (unsigned int ii = 0; ii != diff_base.size1(); ii++) {
260 for (unsigned int jj = 0; jj != diff_base.size2(); jj++) {
261 if (fabs(diff_base(ii, jj)) < eps)
262 diff_base(ii, jj) = 0;
263 }
264 }
265 os << "N: " << std::fixed << base << std::endl
266 << "diffN: " << std::fixed << diff_base;
267 return os;
268}
269
270std::ostream &operator<<(std::ostream &os, const EntitiesFieldData &e) {
271 for (EntityType t = MBVERTEX; t != MBMAXTYPE; ++t) {
272 for (unsigned int nn = 0; nn < e.dataOnEntities[t].size(); nn++) {
273 os << "dataOnEntities[" << moab::CN::EntityTypeName(t) << "][" << nn
274 << "]" << std::endl
275 << e.dataOnEntities[t][nn] << std::endl;
276 }
277 }
278 return os;
279}
280
281/** \name Specializations for H1/L2 */
282
283/**@{*/
284
285template <>
287EntitiesFieldData::EntData::getFTensor1FieldData<3>() {
288 if (dOfs[0]->getNbOfCoeffs() != 3) {
289 std::stringstream s;
290 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
291 s << " but you ask for tensor rank 1 dimension 3";
292 THROW_MESSAGE(s.str());
293 }
294 double *ptr = &*fieldData.data().begin();
296 &ptr[2]);
297}
298
299template <>
301EntitiesFieldData::EntData::getFTensor1FieldData<2>() {
302 if (dOfs[0]->getNbOfCoeffs() != 2) {
303 std::stringstream s;
304 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
305 s << " but you ask for tensor rank 1 dimension 2";
306 THROW_MESSAGE(s.str());
307 }
308 double *ptr = &*fieldData.data().begin();
309 return FTensor::Tensor1<FTensor::PackPtr<double *, 2>, 2>(ptr, &ptr[1]);
310}
311
312template <>
314EntitiesFieldData::EntData::getFTensor1FieldData<1>() {
315 if (dOfs[0]->getNbOfCoeffs() != 1) {
316 std::stringstream s;
317 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
318 s << " but you ask for tensor rank 1 dimension 1";
319 THROW_MESSAGE(s.str());
320 }
321 double *ptr = &*fieldData.data().begin();
323}
324
325template <>
327EntitiesFieldData::EntData::getFTensor2FieldData<1, 1>() {
328 if (dOfs[0]->getNbOfCoeffs() != 1) {
329 std::stringstream s;
330 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
331 s << " but you ask for tensor rank 2 dimensions 1 by 1 so 1 coefficients "
332 "is expected";
333 THROW_MESSAGE(s.str());
334 }
335 double *ptr = &*fieldData.data().begin();
337}
338
339template <>
341EntitiesFieldData::EntData::getFTensor2FieldData<1, 2>() {
342 if (dOfs[0]->getNbOfCoeffs() != 1) {
343 std::stringstream s;
344 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
345 s << " but you ask for tensor rank 2 dimensions 1 by 2 so 2 coefficients "
346 "is expected";
347 THROW_MESSAGE(s.str());
348 }
349 double *ptr = &*fieldData.data().begin();
350 return FTensor::Tensor2<FTensor::PackPtr<double *, 2>, 1, 2>(ptr, &ptr[1]);
351}
352
353template <>
355EntitiesFieldData::EntData::getFTensor2FieldData<1, 3>() {
356 if (dOfs[0]->getNbOfCoeffs() != 1) {
357 std::stringstream s;
358 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
359 s << " but you ask for tensor rank 2 dimensions 1 by 3 so 3 coefficients "
360 "is expected";
361 THROW_MESSAGE(s.str());
362 }
363 double *ptr = &*fieldData.data().begin();
364 return FTensor::Tensor2<FTensor::PackPtr<double *, 3>, 1, 3>(ptr, &ptr[1],
365 &ptr[2]);
366}
367
368template <>
370EntitiesFieldData::EntData::getFTensor2FieldData<2, 2>() {
371 if (dOfs[0]->getNbOfCoeffs() != 4) {
372 std::stringstream s;
373 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
374 s << " but you ask for tensor rank 2 dimensions 2 by 2 so 4 coefficients "
375 "is expected";
376 THROW_MESSAGE(s.str());
377 }
378 double *ptr = &*fieldData.data().begin();
380 ptr, &ptr[1], &ptr[2], &ptr[3]);
381}
382
383template <>
385EntitiesFieldData::EntData::getFTensor2FieldData<3, 3>() {
386 if (dOfs[0]->getNbOfCoeffs() != 9) {
387 std::stringstream s;
388 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
389 s << " but you ask for tensor rank 2 dimensions 3 by 3 so 9 coefficients "
390 "is expected";
391 THROW_MESSAGE(s.str());
392 }
393 double *ptr = &*fieldData.data().begin();
395 ptr, &ptr[1], &ptr[2], &ptr[3], &ptr[4], &ptr[5], &ptr[6], &ptr[7],
396 &ptr[8]);
397}
398
399template <>
401EntitiesFieldData::EntData::getFTensor2SymmetricFieldData<3>() {
402 if (dOfs[0]->getNbOfCoeffs() != 6) {
403 std::stringstream s;
404 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
405 s << " but you ask for symmetric tensor rank 2 dimensions 3 by 3 so 6 "
406 "coefficients "
407 "is expected";
408 THROW_MESSAGE(s.str());
409 }
410 double *ptr = &*fieldData.data().begin();
412 ptr, &ptr[1], &ptr[2], &ptr[3], &ptr[4], &ptr[5]);
413}
414
415template <>
417EntitiesFieldData::EntData::getFTensor2SymmetricFieldData<2>() {
418 if (dOfs[0]->getNbOfCoeffs() != 3) {
419 std::stringstream s;
420 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
421 s << " but you ask for symmetric tensor rank 2 dimensions 2 by 2 so 3 "
422 "coefficients "
423 "is expected";
424 THROW_MESSAGE(s.str());
425 }
426 double *ptr = &*fieldData.data().begin();
428 ptr, &ptr[1], &ptr[2]);
429}
430
433 if (dOfs[0]->getNbOfCoeffs() != 1) {
434 std::stringstream s;
435 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
436 s << " but expected scalar field, tensor of rank 0";
437 THROW_MESSAGE(s.str());
438 }
440 &*fieldData.data().begin());
441}
442
443template <int Tensor_Dim>
446 const FieldApproximationBase base) {
447 std::stringstream s;
448 s << "Template for tensor dimension " << Tensor_Dim << " not implemented";
449 THROW_MESSAGE(s.str());
451}
452
453template <int Tensor_Dim>
456 return getFTensor1DiffN<Tensor_Dim>(bAse);
457}
458
459/**
460 * \brief Get spatial derivative of base function tensor for dimension 3d
461 */
462template <>
464EntitiesFieldData::EntData::getFTensor1DiffN<3>(
465 const FieldApproximationBase base) {
466 double *ptr = &*getDiffN(base).data().begin();
468 &ptr[2]);
469}
470
471/**
472 * \brief Get spatial derivative of base function tensor for dimension 3d
473 */
474template <>
476EntitiesFieldData::EntData::getFTensor1DiffN<3>() {
477 return getFTensor1DiffN<3>(bAse);
478}
479
480/**
481 * \brief Get spatial derivative of base function tensor for dimension 2d
482 */
483template <>
485EntitiesFieldData::EntData::getFTensor1DiffN<2>(
486 const FieldApproximationBase base) {
487 double *ptr = &*getDiffN(base).data().begin();
488 return FTensor::Tensor1<FTensor::PackPtr<double *, 2>, 2>(ptr, &ptr[1]);
489}
490
491/**
492 * \brief Get spatial derivative of base function tensor for dimension 2d
493 */
494template <>
496EntitiesFieldData::EntData::getFTensor1DiffN<2>() {
497 return getFTensor1DiffN<2>(bAse);
498}
499
500template <int Tensor_Dim>
503 const int gg, const int bb) {
504 std::stringstream s;
505 s << "Template for tensor dimension " << Tensor_Dim << " not implemented";
506 THROW_MESSAGE(s.str());
508}
509
510/**
511 * \brief Get spatial derivative of base function tensor for dimension 3d
512 */
513template <>
515EntitiesFieldData::EntData::getFTensor1DiffN<3>(
516 const FieldApproximationBase base, const int gg, const int bb) {
517 double *ptr = &getDiffN(base)(gg, 3 * bb);
519 &ptr[2]);
520}
521
522/**
523 * \brief Get spatial derivative of base function tensor for dimension 3d
524 */
525template <>
527EntitiesFieldData::EntData::getFTensor1DiffN<3>(const int gg, const int bb) {
528 return getFTensor1DiffN<3>(bAse, gg, bb);
529}
530
531/**
532 * \brief Get spatial derivative of base function tensor for dimension 2d
533 */
534template <>
536EntitiesFieldData::EntData::getFTensor1DiffN<2>(
537 const FieldApproximationBase base, const int gg, const int bb) {
538 double *ptr = &getDiffN(base)(gg, 2 * bb);
539 return FTensor::Tensor1<FTensor::PackPtr<double *, 2>, 2>(ptr, &ptr[1]);
540}
541
542/**
543 * \brief Get spatial derivative of base function tensor for dimension 2d
544 */
545template <>
547EntitiesFieldData::EntData::getFTensor1DiffN<2>(const int gg, const int bb) {
548 return getFTensor1DiffN<2>(bAse, gg, bb);
549}
550
551/**@}*/
552
553/** \name Specializations for HDiv/HCrul */
554
555/**@{*/
556
557template <int Tensor_Dim>
560 std::stringstream s;
561 s << "Template for tensor dimension " << Tensor_Dim << " not implemented";
562 THROW_MESSAGE(s.str());
564}
565
566template <int Tensor_Dim>
569 const int gg, const int bb) {
570 std::stringstream s;
571 s << "Template for tensor dimension " << Tensor_Dim << " not implemented";
572 THROW_MESSAGE(s.str());
574}
575
576template <int Tensor_Dim0, int Tensor_Dim1>
578 Tensor_Dim0, Tensor_Dim1>
580 std::stringstream s;
581 s << "Template for tensor dimension " << Tensor_Dim0 << "x" << Tensor_Dim1
582 << " not implemented";
583 THROW_MESSAGE(s.str());
585}
586
587template <int Tensor_Dim0, int Tensor_Dim1>
589 Tensor_Dim0, Tensor_Dim1>
591 const int gg, const int bb) {
592 std::stringstream s;
593 s << "Template for tensor dimension " << Tensor_Dim0 << "x" << Tensor_Dim1
594 << " not implemented";
595 THROW_MESSAGE(s.str());
597}
598
599template <>
601EntitiesFieldData::EntData::getFTensor1N<3>(FieldApproximationBase base) {
602 double *t_n_ptr = &*getN(base).data().begin();
603 return FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 3>(t_n_ptr, // HVEC0
604 &t_n_ptr[HVEC1],
605 &t_n_ptr[HVEC2]);
606}
607
608template <>
610EntitiesFieldData::EntData::getFTensor1N<3>(FieldApproximationBase base,
611 const int gg, const int bb) {
612 double *t_n_ptr = &getN(base)(gg, 3 * bb);
613 return FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 3>(t_n_ptr, // HVEC0
614 &t_n_ptr[HVEC1],
615 &t_n_ptr[HVEC2]);
616}
617
618template <>
620EntitiesFieldData::EntData::getFTensor2DiffN<3, 3>(
622 double *t_diff_n_ptr = &*getDiffN(base).data().begin();
624 t_diff_n_ptr, &t_diff_n_ptr[HVEC0_1], &t_diff_n_ptr[HVEC0_2],
625 &t_diff_n_ptr[HVEC1_0], &t_diff_n_ptr[HVEC1_1], &t_diff_n_ptr[HVEC1_2],
626 &t_diff_n_ptr[HVEC2_0], &t_diff_n_ptr[HVEC2_1], &t_diff_n_ptr[HVEC2_2]);
627}
628
629template <>
631EntitiesFieldData::EntData::getFTensor2DiffN<3, 3>(FieldApproximationBase base,
632 const int gg, const int bb) {
633 double *t_diff_n_ptr = &getDiffN(base)(gg, 9 * bb);
635 t_diff_n_ptr, &t_diff_n_ptr[HVEC0_1], &t_diff_n_ptr[HVEC0_2],
636 &t_diff_n_ptr[HVEC1_0], &t_diff_n_ptr[HVEC1_1], &t_diff_n_ptr[HVEC1_2],
637 &t_diff_n_ptr[HVEC2_0], &t_diff_n_ptr[HVEC2_1], &t_diff_n_ptr[HVEC2_2]);
638}
639
640template <>
642EntitiesFieldData::EntData::getFTensor2DiffN<3, 2>(
644 double *t_diff_n_ptr = &*getDiffN(base).data().begin();
646 t_diff_n_ptr, &t_diff_n_ptr[HVEC0_1], &t_diff_n_ptr[HVEC1_0],
647 &t_diff_n_ptr[HVEC1_1], &t_diff_n_ptr[HVEC2_0], &t_diff_n_ptr[HVEC2_1]);
648}
649
650template <>
652EntitiesFieldData::EntData::getFTensor2DiffN<3, 2>(FieldApproximationBase base,
653 const int gg, const int bb) {
654 double *t_diff_n_ptr = &getDiffN(base)(gg, 6 * bb);
656 t_diff_n_ptr, &t_diff_n_ptr[HVEC0_1], &t_diff_n_ptr[HVEC1_0],
657 &t_diff_n_ptr[HVEC1_1], &t_diff_n_ptr[HVEC2_0], &t_diff_n_ptr[HVEC2_1]);
658}
659
660template <>
663 double *ptr = &(getN(base, BaseDerivatives::SecondDerivative))(0, 0);
665
666 &ptr[2 * HVEC0_0 + 0], &ptr[2 * HVEC0_0 + 1], &ptr[2 * HVEC0_1 + 0],
667 &ptr[2 * HVEC0_1 + 1],
668
669 &ptr[2 * HVEC1_0 + 0], &ptr[2 * HVEC1_0 + 1], &ptr[2 * HVEC1_1 + 0],
670 &ptr[2 * HVEC1_1 + 1],
671
672 &ptr[2 * HVEC2_0 + 0], &ptr[2 * HVEC2_0 + 1], &ptr[2 * HVEC2_1 + 0],
673 &ptr[2 * HVEC2_1 + 1]
674
675 };
676}
677
678template <int Tensor_Dim0, int Tensor_Dim1>
680 Tensor_Dim0, Tensor_Dim1>
682 std::stringstream s;
683 s << "Template for tensor dimension " << Tensor_Dim0 << ", " << Tensor_Dim1
684 << " not implemented";
685 THROW_MESSAGE(s.str());
687 Tensor_Dim0, Tensor_Dim1>();
688}
689
690template <>
692EntitiesFieldData::EntData::getFTensor2N<3, 3>(FieldApproximationBase base) {
693 double *t_n_ptr = &*(getN(base).data().begin());
695
696 &t_n_ptr[HVEC0], &t_n_ptr[HVEC1], &t_n_ptr[HVEC2],
697
698 &t_n_ptr[3 + HVEC0], &t_n_ptr[3 + HVEC1], &t_n_ptr[3 + HVEC2],
699
700 &t_n_ptr[6 + HVEC0], &t_n_ptr[6 + HVEC1], &t_n_ptr[6 + HVEC2]
701
702 );
703}
704
705template <int Tensor_Dim0, int Tensor_Dim1>
707 Tensor_Dim0, Tensor_Dim1>
709 const int gg, const int bb) {
710 std::stringstream s;
711 s << "Template for tensor dimension " << Tensor_Dim0 << ", " << Tensor_Dim1
712 << " not implemented";
713 THROW_MESSAGE(s.str());
715 Tensor_Dim0, Tensor_Dim1>();
716}
717
718template <>
720EntitiesFieldData::EntData::getFTensor2N<3, 3>(FieldApproximationBase base,
721 const int gg, const int bb) {
722 double *t_n_ptr = &getN(base)(gg, 9 * bb);
724
725 &t_n_ptr[HVEC0], &t_n_ptr[HVEC1], &t_n_ptr[HVEC2],
726
727 &t_n_ptr[3 + HVEC0], &t_n_ptr[3 + HVEC1], &t_n_ptr[3 + HVEC2],
728
729 &t_n_ptr[6 + HVEC0], &t_n_ptr[6 + HVEC1], &t_n_ptr[6 + HVEC2]
730
731 );
732}
733
734/**@}*/
735
736/** \name Bernstein-Bezier base only functions */
737
738/**@{*/
739
740boost::shared_ptr<MatrixInt> &
742 const std::string &field_name) {
743 return bbAlphaIndices[field_name];
744}
745
746boost::shared_ptr<MatrixDouble> &
748 return bbN[field_name];
749}
750
751/**
752 * Get shared pointer to BB base base functions
753 */
754const boost::shared_ptr<MatrixDouble> &
756 const std::string &field_name) const {
757 return bbN.at(field_name);
758}
759
760/**
761 * Get shared pointer to BB derivatives of base base functions
762 */
763boost::shared_ptr<MatrixDouble> &
765 return bbDiffN[field_name];
766}
767
768/**
769 * Get shared pointer to derivatives of BB base base functions
770 */
771const boost::shared_ptr<MatrixDouble> &
773 const std::string &field_name) const {
774 return bbDiffN.at(field_name);
775}
776
777std::map<std::string, boost::shared_ptr<MatrixInt>> &
779 return bbAlphaIndices;
780}
781
782std::map<std::string, boost::shared_ptr<MatrixDouble>> &
784 return bbN;
785}
786
787std::map<std::string, boost::shared_ptr<MatrixDouble>> &
789 return bbDiffN;
790}
791
792boost::shared_ptr<MatrixInt> &
794 return bbAlphaIndicesByOrder[o];
795}
796
797boost::shared_ptr<MatrixDouble> &
799 return bbNByOrder[o];
800}
801
802boost::shared_ptr<MatrixDouble> &
804 return bbDiffNByOrder[o];
805}
806
807std::array<boost::shared_ptr<MatrixInt>,
810 return bbAlphaIndicesByOrder;
811}
812
813std::array<boost::shared_ptr<MatrixDouble>,
816 return bbNByOrder;
817}
818
819std::array<boost::shared_ptr<MatrixDouble>,
822 return bbDiffNByOrder;
823}
824
825boost::shared_ptr<MatrixInt> &
827 const std::string &field_name) {
828 return entDataPtr->getBBAlphaIndicesSharedPtr(field_name);
829}
830
831boost::shared_ptr<MatrixDouble> &
833 const std::string &field_name) {
834 return entDataPtr->getBBNSharedPtr(field_name);
835}
836
837const boost::shared_ptr<MatrixDouble> &
839 const std::string &field_name) const {
840 return entDataPtr->getBBNSharedPtr(field_name);
841}
842
843/**
844 * Get shared pointer to BB derivatives of base base functions
845 */
846boost::shared_ptr<MatrixDouble> &
848 const std::string &field_name) {
849 return entDataPtr->getBBDiffNSharedPtr(field_name);
850}
851
852/**
853 * Get shared pointer to derivatives of BB base base functions
854 */
855const boost::shared_ptr<MatrixDouble> &
857 const std::string &field_name) const {
858 return entDataPtr->getBBDiffNSharedPtr(field_name);
859}
860
861/**@}*/
862
864 return entDataBitRefLevel;
865}
866
867std::vector<BitRefLevel> &
869 return entDataPtr->getEntDataBitRefLevel();
870}
871
872} // namespace MoFEM
static Index< 'o', 3 > o
static const double eps
T data[Tensor_Dim]
T data[Tensor_Dim0][Tensor_Dim1]
FieldApproximationBase
approximation base
Definition: definitions.h:58
@ LASTBASE
Definition: definitions.h:69
@ NOBASE
Definition: definitions.h:59
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
@ NOSPACE
Definition: definitions.h:83
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ HVEC0
Definition: definitions.h:186
@ HVEC1
Definition: definitions.h:186
@ HVEC2
Definition: definitions.h:186
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
@ HVEC1_1
Definition: definitions.h:196
@ HVEC0_1
Definition: definitions.h:195
@ HVEC1_0
Definition: definitions.h:193
@ HVEC2_1
Definition: definitions.h:197
@ HVEC1_2
Definition: definitions.h:199
@ HVEC2_2
Definition: definitions.h:200
@ HVEC2_0
Definition: definitions.h:194
@ HVEC0_2
Definition: definitions.h:198
@ HVEC0_0
Definition: definitions.h:192
#define CHKERR
Inline error check.
Definition: definitions.h:535
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
Definition: definitions.h:561
const double c
speed of light (cm/ns)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
static void constructor_derived_data(DerivedEntitiesFieldData *derived_data, const boost::shared_ptr< EntitiesFieldData > &data_ptr)
static void constructor_data(EntitiesFieldData *data, const EntityType type)
constexpr double t
plate stiffness
Definition: plate.cpp:59
constexpr auto field_name
const int N
Definition: speed_test.cpp:3
Derived ata on single entity (This is passed as argument to DataOperator::doWork)
boost::shared_ptr< MatrixDouble > & getBBNSharedPtr(const std::string &field_name)
MoFEMErrorCode baseSwap(const std::string &field_name, const FieldApproximationBase base)
Used by Bernstein base to keep temporally pointer.
std::vector< BitRefLevel > & getEntDataBitRefLevel()
int getSense() const
get entity sense, need to calculate base functions with conforming approximation fields
boost::shared_ptr< MatrixInt > & getBBAlphaIndicesSharedPtr(const std::string &field_name)
boost::shared_ptr< MatrixDouble > & getBBDiffNSharedPtr(const std::string &field_name)
boost::shared_ptr< MatrixDouble > & getDiffNSharedPtr(const FieldApproximationBase base)
DerivedEntData(const boost::shared_ptr< EntitiesFieldData::EntData > &ent_data_ptr)
boost::shared_ptr< MatrixDouble > & getNSharedPtr(const FieldApproximationBase base, const BaseDerivatives derivative)
this class derive data form other data structure
MoFEMErrorCode setElementType(const EntityType type)
const boost::shared_ptr< EntitiesFieldData > dataPtr
DerivedEntitiesFieldData(const boost::shared_ptr< EntitiesFieldData > &data_ptr)
Data on single entity (This is passed as argument to DataOperator::doWork)
virtual std::array< boost::shared_ptr< MatrixDouble >, MaxBernsteinBezierOrder > & getBBDiffNByOrderArray()
virtual std::array< boost::shared_ptr< MatrixDouble >, MaxBernsteinBezierOrder > & getBBNByOrderArray()
auto getFTensor2N()
Get base functions for Hdiv space.
VectorDouble fieldData
Field data on entity.
ApproximationOrder getOrder() const
get approximation order
VectorInt localIndices
Local indices on entity.
virtual std::array< boost::shared_ptr< MatrixInt >, MaxBernsteinBezierOrder > & getBBAlphaIndicesByOrderArray()
virtual std::map< std::string, boost::shared_ptr< MatrixInt > > & getBBAlphaIndicesMap()
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1DiffN()
Get derivatives of base functions.
static constexpr size_t MaxBernsteinBezierOrder
VectorInt iNdices
Global indices on entity.
virtual boost::shared_ptr< MatrixInt > & getBBAlphaIndicesSharedPtr(const std::string &field_name)
const VectorDouble & getFieldData() const
get dofs values
virtual boost::shared_ptr< MatrixDouble > & getBBDiffNByOrderSharedPtr(const size_t o)
get BB base derivative by order
virtual std::map< std::string, boost::shared_ptr< MatrixDouble > > & getBBNMap()
get hash map of base function for BB base, key is a field name
virtual MoFEMErrorCode baseSwap(const std::string &field_name, const FieldApproximationBase base)
Swap bases functions.
const VectorInt & getLocalIndices() const
get local indices of dofs on entity
virtual std::map< std::string, boost::shared_ptr< MatrixDouble > > & getBBDiffNMap()
get hash map of direvarives base function for BB base, key is a field name
virtual std::vector< BitRefLevel > & getEntDataBitRefLevel()
virtual int getSense() const
get entity sense, need to calculate base functions with conforming approximation fields
VectorFieldEntities fieldEntities
Field entities.
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0FieldData()
Resturn scalar files as a FTensor of rank 0.
virtual boost::shared_ptr< MatrixInt > & getBBAlphaIndicesByOrderSharedPtr(const size_t o)
get ALpha indices for BB base by order
FTensor::Tensor3< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 *Tensor_Dim2 >, Tensor_Dim0, Tensor_Dim1, Tensor_Dim2 > getFTensor3Diff2N()
Get second derivatives of base functions for Hvec space.
auto getFTensor1N()
Get base functions for Hdiv space.
virtual boost::shared_ptr< MatrixDouble > & getBBNSharedPtr(const std::string &field_name)
std::array< std::array< boost::shared_ptr< MatrixDouble >, LASTBASE >, LastDerivative > baseFunctionsAndBaseDerivatives
FTensor::Tensor2< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 >, Tensor_Dim0, Tensor_Dim1 > getFTensor2DiffN()
Get derivatives of base functions for Hdiv space.
virtual boost::shared_ptr< MatrixDouble > & getDiffNSharedPtr(const FieldApproximationBase base)
virtual boost::shared_ptr< MatrixDouble > & getBBNByOrderSharedPtr(const size_t o)
get BB base by order
virtual boost::shared_ptr< MatrixDouble > & getNSharedPtr(const FieldApproximationBase base, const BaseDerivatives direvatie)
const VectorInt & getIndices() const
Get global indices of dofs on entity.
EntData(const bool allocate_base_matrices=true)
FieldApproximationBase bAse
Field approximation base.
virtual boost::shared_ptr< MatrixDouble > & getBBDiffNSharedPtr(const std::string &field_name)
data structure for finite element entity
friend std::ostream & operator<<(std::ostream &os, const EntitiesFieldData &e)
MoFEMErrorCode resetFieldDependentData()
std::bitset< LASTBASE > bAse
bases on element
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
virtual MoFEMErrorCode baseSwap(const std::string &field_name, const FieldApproximationBase base)
Swap approximation base.
virtual MoFEMErrorCode setElementType(const EntityType type)