v0.13.1
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/* This file is part of MoFEM.
7 * MoFEM is free software: you can redistribute it and/or modify it under
8 * the terms of the GNU Lesser General Public License as published by the
9 * Free Software Foundation, either version 3 of the License, or (at your
10 * option) any later version.
11 *
12 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
13 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
14 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
15 * License for more details.
16 *
17 * You should have received a copy of the GNU Lesser General Public
18 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
19
20namespace MoFEM {
21
22EntitiesFieldData::EntData::EntData(const bool allocate_base_matrices)
23 : sEnse(0), oRder(0), bAse(NOBASE), entDataBitRefLevel(),
24 N(baseFunctionsAndBaseDerivatives[ZeroDerivative]),
25 diffN(baseFunctionsAndBaseDerivatives[FirstDerivative]) {
26 if (allocate_base_matrices) {
27
28 for (auto d = 0; d != LastDerivative; ++d) {
29 for (int b = 0; b != LASTBASE; ++b) {
31 }
32 }
33 }
34}
35
36int EntitiesFieldData::EntData::getSense() const { return sEnse; }
37
38boost::shared_ptr<MatrixDouble> &
40 const BaseDerivatives direvatie) {
41 return baseFunctionsAndBaseDerivatives[direvatie][base];
42}
43
44boost::shared_ptr<MatrixDouble> &
46 return N[base];
47}
48
49boost::shared_ptr<MatrixDouble> &EntitiesFieldData::EntData::getDiffNSharedPtr(
50 const FieldApproximationBase base) {
51 return diffN[base];
52}
53
55
57
58 auto set_default = [&]() {
59 std::array<size_t, MBMAXTYPE> count;
60 std::fill(count.begin(), count.end(), 0);
61 const int dim_type = moab::CN::Dimension(type);
62 data->dataOnEntities[MBVERTEX].resize(1);
63 if (type != MBVERTEX) {
64 for (auto dd = dim_type; dd > 0; --dd) {
65 int nb_ents = moab::CN::NumSubEntities(type, dd);
66 for (int ii = 0; ii != nb_ents; ++ii) {
67 auto sub_ent_type = moab::CN::SubEntityType(type, dd, ii);
68 count[sub_ent_type] = nb_ents;
69 }
70 for (auto tt = moab::CN::TypeDimensionMap[dd].first;
71 tt <= moab::CN::TypeDimensionMap[dd].second; ++tt) {
72 data->dataOnEntities[tt].resize(count[tt]);
73 }
74 }
75 }
76 };
77
78 switch (type) {
79 case MBENTITYSET:
80 break;
81
82 default:
83 set_default();
84 }
85}
86
89}
90
95}
96
97static void
99 const boost::shared_ptr<EntitiesFieldData> &data_ptr) {
100
102 using DerivedEntData = DerivedEntitiesFieldData::DerivedEntData;
103
104 for (int tt = MBVERTEX; tt != MBMAXTYPE; ++tt) {
105 auto &ent_data = data_ptr->dataOnEntities[tt];
106 auto &derived_ent_data = derived_data->dataOnEntities[tt];
107 for (auto c = derived_ent_data.size(); c < ent_data.size(); ++c) {
108 boost::shared_ptr<EntData> ent_data_ptr(data_ptr, &ent_data[c]);
109 derived_ent_data.push_back(new DerivedEntData(ent_data_ptr));
110 }
111 derived_ent_data.resize(ent_data.size());
112 }
113}
114
116 const boost::shared_ptr<EntitiesFieldData> &data_ptr)
117 : EntitiesFieldData(), dataPtr(data_ptr) {
119}
120
125}
126
129 sPace = NOSPACE;
130 bAse = NOBASE;
131 fieldEntities.resize(0, false);
132 iNdices.resize(0, false);
133 localIndices.resize(0, false);
134 dOfs.resize(0, false);
135 fieldData.resize(0, false);
137}
138
141 for (EntityType t = MBVERTEX; t != MBMAXTYPE; t++)
142 for (auto &e : dataOnEntities[t])
143 CHKERR e.resetFieldDependentData();
145}
146
149 const FieldApproximationBase base) {
151 auto make_swap = [](boost::shared_ptr<MatrixDouble> &ptr,
152 boost::shared_ptr<MatrixDouble> &ptrBB,
153 boost::shared_ptr<MatrixDouble> &swap_ptr) {
154 if (swap_ptr) {
155 ptr = swap_ptr;
156 swap_ptr.reset();
157 } else {
158 swap_ptr = ptr;
159 ptr = ptrBB;
160 }
161 };
162 make_swap(getNSharedPtr(base), getBBNSharedPtr(field_name), swapBaseNPtr);
163 make_swap(getDiffNSharedPtr(base), getBBDiffNSharedPtr(field_name),
164 swapBaseDiffNPtr);
166}
167
169 const FieldApproximationBase base) {
171 // Note: Do not swap bases on entities sets
172 for (int tt = MBVERTEX; tt != MBENTITYSET; ++tt) {
173 auto &ent_data = dataOnEntities[tt];
174 for (auto &side_data : ent_data)
175 CHKERR side_data.baseSwap(field_name, base);
176 }
178}
179
181 const std::string &field_name, const FieldApproximationBase base) {
183 auto make_swap = [](boost::shared_ptr<MatrixDouble> &ptr,
184 boost::shared_ptr<MatrixDouble> &ptrBB,
185 boost::shared_ptr<MatrixDouble> &swap_ptr) {
186 if (swap_ptr) {
187 ptr = swap_ptr;
188 swap_ptr.reset();
189 } else {
190 swap_ptr = ptr;
191 ptr = ptrBB;
192 }
193 };
194 make_swap(getDerivedNSharedPtr(base), getBBNSharedPtr(field_name),
195 swapBaseNPtr);
196 make_swap(getDerivedDiffNSharedPtr(base), getBBDiffNSharedPtr(field_name),
197 swapBaseDiffNPtr);
199}
200
202 const boost::shared_ptr<EntitiesFieldData::EntData> &ent_data_ptr)
203 : EntitiesFieldData::EntData(false), entDataPtr(ent_data_ptr) {}
204
206 return entDataPtr->getSense();
207}
208
209boost::shared_ptr<MatrixDouble> &
211 const FieldApproximationBase base, const BaseDerivatives derivative) {
212 if (baseFunctionsAndBaseDerivatives[derivative][base])
213 return baseFunctionsAndBaseDerivatives[derivative][base];
214 else
215 return entDataPtr->getNSharedPtr(base, derivative);
216}
217
218boost::shared_ptr<MatrixDouble> &
220 const FieldApproximationBase base) {
221 if (N[base])
222 return N[base];
223 else
224 return entDataPtr->getNSharedPtr(base);
225}
226
227boost::shared_ptr<MatrixDouble> &
229 const FieldApproximationBase base) {
230 if (diffN[base])
231 return diffN[base];
232 else
233 return entDataPtr->getDiffNSharedPtr(base);
234}
235const boost::shared_ptr<MatrixDouble> &
237 const FieldApproximationBase base) const {
238 if (N[base])
239 return N[base];
240 else
241 return entDataPtr->getNSharedPtr(base);
242}
243const boost::shared_ptr<MatrixDouble> &
245 const FieldApproximationBase base) const {
246 if (diffN[base])
247 return diffN[base];
248 else
249 return entDataPtr->getDiffNSharedPtr(base);
250}
251
252std::ostream &operator<<(std::ostream &os,
254 os << "sEnse: " << e.getSense() << std::endl
255 << "oRder: " << e.getOrder() << std::endl
256 << "global indices: " << e.getIndices() << std::endl
257 << "local indices: " << e.getLocalIndices() << std::endl;
258 // FIXME: precision should not be set here
259 os << "fieldData: " << std::fixed << std::setprecision(2) << e.getFieldData()
260 << std::endl;
261 MatrixDouble base = const_cast<EntitiesFieldData::EntData &>(e).getN();
262 MatrixDouble diff_base =
263 const_cast<EntitiesFieldData::EntData &>(e).getDiffN();
264 const double eps = 1e-6;
265 for (unsigned int ii = 0; ii != base.size1(); ii++) {
266 for (unsigned int jj = 0; jj != base.size2(); jj++) {
267 if (fabs(base(ii, jj)) < eps)
268 base(ii, jj) = 0;
269 }
270 }
271 for (unsigned int ii = 0; ii != diff_base.size1(); ii++) {
272 for (unsigned int jj = 0; jj != diff_base.size2(); jj++) {
273 if (fabs(diff_base(ii, jj)) < eps)
274 diff_base(ii, jj) = 0;
275 }
276 }
277 os << "N: " << std::fixed << base << std::endl
278 << "diffN: " << std::fixed << diff_base;
279 return os;
280}
281
282std::ostream &operator<<(std::ostream &os, const EntitiesFieldData &e) {
283 for (EntityType t = MBVERTEX; t != MBMAXTYPE; ++t) {
284 for (unsigned int nn = 0; nn < e.dataOnEntities[t].size(); nn++) {
285 os << "dataOnEntities[" << moab::CN::EntityTypeName(t) << "][" << nn
286 << "]" << std::endl
287 << e.dataOnEntities[t][nn] << std::endl;
288 }
289 }
290 return os;
291}
292
293/** \name Specializations for H1/L2 */
294
295/**@{*/
296
297template <>
299EntitiesFieldData::EntData::getFTensor1FieldData<3>() {
300 if (dOfs[0]->getNbOfCoeffs() != 3) {
301 std::stringstream s;
302 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
303 s << " but you ask for tensor rank 1 dimension 3";
304 THROW_MESSAGE(s.str());
305 }
306 double *ptr = &*fieldData.data().begin();
308 &ptr[2]);
309}
310
311template <>
313EntitiesFieldData::EntData::getFTensor1FieldData<2>() {
314 if (dOfs[0]->getNbOfCoeffs() != 2) {
315 std::stringstream s;
316 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
317 s << " but you ask for tensor rank 1 dimension 2";
318 THROW_MESSAGE(s.str());
319 }
320 double *ptr = &*fieldData.data().begin();
321 return FTensor::Tensor1<FTensor::PackPtr<double *, 2>, 2>(ptr, &ptr[1]);
322}
323
324template <>
326EntitiesFieldData::EntData::getFTensor1FieldData<1>() {
327 if (dOfs[0]->getNbOfCoeffs() != 1) {
328 std::stringstream s;
329 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
330 s << " but you ask for tensor rank 1 dimension 1";
331 THROW_MESSAGE(s.str());
332 }
333 double *ptr = &*fieldData.data().begin();
335}
336
337template <>
339EntitiesFieldData::EntData::getFTensor2FieldData<3, 3>() {
340 if (dOfs[0]->getNbOfCoeffs() != 9) {
341 std::stringstream s;
342 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
343 s << " but you ask for tensor rank 2 dimensions 3 by 3 so 9 coefficients "
344 "is expected";
345 THROW_MESSAGE(s.str());
346 }
347 double *ptr = &*fieldData.data().begin();
349 ptr, &ptr[1], &ptr[2], &ptr[3], &ptr[4], &ptr[5], &ptr[6], &ptr[7],
350 &ptr[8]);
351}
352
353template <>
355EntitiesFieldData::EntData::getFTensor2SymmetricFieldData<3>() {
356 if (dOfs[0]->getNbOfCoeffs() != 6) {
357 std::stringstream s;
358 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
359 s << " but you ask for symmetric tensor rank 2 dimensions 3 by 3 so 6 "
360 "coefficients "
361 "is expected";
362 THROW_MESSAGE(s.str());
363 }
364 double *ptr = &*fieldData.data().begin();
366 ptr, &ptr[1], &ptr[2], &ptr[3], &ptr[4], &ptr[5]);
367}
368
369template <>
371EntitiesFieldData::EntData::getFTensor2SymmetricFieldData<2>() {
372 if (dOfs[0]->getNbOfCoeffs() != 3) {
373 std::stringstream s;
374 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
375 s << " but you ask for symmetric tensor rank 2 dimensions 2 by 2 so 3 "
376 "coefficients "
377 "is expected";
378 THROW_MESSAGE(s.str());
379 }
380 double *ptr = &*fieldData.data().begin();
382 ptr, &ptr[1], &ptr[2]);
383}
384
387 if (dOfs[0]->getNbOfCoeffs() != 1) {
388 std::stringstream s;
389 s << "Wrong number of coefficients is " << dOfs[0]->getNbOfCoeffs();
390 s << " but expected scalar field, tensor of rank 0";
391 THROW_MESSAGE(s.str());
392 }
394 &*fieldData.data().begin());
395}
396
397template <int Tensor_Dim>
400 const FieldApproximationBase base) {
401 std::stringstream s;
402 s << "Template for tensor dimension " << Tensor_Dim << " not implemented";
403 THROW_MESSAGE(s.str());
405}
406
407template <int Tensor_Dim>
410 return getFTensor1DiffN<Tensor_Dim>(bAse);
411}
412
413/**
414 * \brief Get spatial derivative of base function tensor for dimension 3d
415 */
416template <>
418EntitiesFieldData::EntData::getFTensor1DiffN<3>(
419 const FieldApproximationBase base) {
420 double *ptr = &*getDiffN(base).data().begin();
422 &ptr[2]);
423}
424
425/**
426 * \brief Get spatial derivative of base function tensor for dimension 3d
427 */
428template <>
430EntitiesFieldData::EntData::getFTensor1DiffN<3>() {
431 return getFTensor1DiffN<3>(bAse);
432}
433
434/**
435 * \brief Get spatial derivative of base function tensor for dimension 2d
436 */
437template <>
439EntitiesFieldData::EntData::getFTensor1DiffN<2>(
440 const FieldApproximationBase base) {
441 double *ptr = &*getDiffN(base).data().begin();
442 return FTensor::Tensor1<FTensor::PackPtr<double *, 2>, 2>(ptr, &ptr[1]);
443}
444
445/**
446 * \brief Get spatial derivative of base function tensor for dimension 2d
447 */
448template <>
450EntitiesFieldData::EntData::getFTensor1DiffN<2>() {
451 return getFTensor1DiffN<2>(bAse);
452}
453
454template <int Tensor_Dim>
457 const int gg, const int bb) {
458 std::stringstream s;
459 s << "Template for tensor dimension " << Tensor_Dim << " not implemented";
460 THROW_MESSAGE(s.str());
462}
463
464/**
465 * \brief Get spatial derivative of base function tensor for dimension 3d
466 */
467template <>
469EntitiesFieldData::EntData::getFTensor1DiffN<3>(
470 const FieldApproximationBase base, const int gg, const int bb) {
471 double *ptr = &getDiffN(base)(gg, 3 * bb);
473 &ptr[2]);
474}
475
476/**
477 * \brief Get spatial derivative of base function tensor for dimension 3d
478 */
479template <>
481EntitiesFieldData::EntData::getFTensor1DiffN<3>(const int gg, const int bb) {
482 return getFTensor1DiffN<3>(bAse, gg, bb);
483}
484
485/**
486 * \brief Get spatial derivative of base function tensor for dimension 2d
487 */
488template <>
490EntitiesFieldData::EntData::getFTensor1DiffN<2>(
491 const FieldApproximationBase base, const int gg, const int bb) {
492 double *ptr = &getDiffN(base)(gg, 2 * bb);
493 return FTensor::Tensor1<FTensor::PackPtr<double *, 2>, 2>(ptr, &ptr[1]);
494}
495
496/**
497 * \brief Get spatial derivative of base function tensor for dimension 2d
498 */
499template <>
501EntitiesFieldData::EntData::getFTensor1DiffN<2>(const int gg, const int bb) {
502 return getFTensor1DiffN<2>(bAse, gg, bb);
503}
504
505/**@}*/
506
507/** \name Specializations for HDiv/HCrul */
508
509/**@{*/
510
511template <int Tensor_Dim>
514 std::stringstream s;
515 s << "Template for tensor dimension " << Tensor_Dim << " not implemented";
516 THROW_MESSAGE(s.str());
518}
519
520template <int Tensor_Dim>
523 const int gg, const int bb) {
524 std::stringstream s;
525 s << "Template for tensor dimension " << Tensor_Dim << " not implemented";
526 THROW_MESSAGE(s.str());
528}
529
530template <int Tensor_Dim0, int Tensor_Dim1>
532 Tensor_Dim0, Tensor_Dim1>
534 std::stringstream s;
535 s << "Template for tensor dimension " << Tensor_Dim0 << "x" << Tensor_Dim1
536 << " not implemented";
537 THROW_MESSAGE(s.str());
539}
540
541template <int Tensor_Dim0, int Tensor_Dim1>
543 Tensor_Dim0, Tensor_Dim1>
545 const int gg, const int bb) {
546 std::stringstream s;
547 s << "Template for tensor dimension " << Tensor_Dim0 << "x" << Tensor_Dim1
548 << " not implemented";
549 THROW_MESSAGE(s.str());
551}
552
553template <>
555EntitiesFieldData::EntData::getFTensor1N<3>(FieldApproximationBase base) {
556 double *t_n_ptr = &*getN(base).data().begin();
557 return FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 3>(t_n_ptr, // HVEC0
558 &t_n_ptr[HVEC1],
559 &t_n_ptr[HVEC2]);
560}
561
562template <>
564EntitiesFieldData::EntData::getFTensor1N<3>(FieldApproximationBase base,
565 const int gg, const int bb) {
566 double *t_n_ptr = &getN(base)(gg, 3 * bb);
567 return FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 3>(t_n_ptr, // HVEC0
568 &t_n_ptr[HVEC1],
569 &t_n_ptr[HVEC2]);
570}
571
572template <>
574EntitiesFieldData::EntData::getFTensor2DiffN<3, 3>(
576 double *t_diff_n_ptr = &*getDiffN(base).data().begin();
578 t_diff_n_ptr, &t_diff_n_ptr[HVEC0_1], &t_diff_n_ptr[HVEC0_2],
579 &t_diff_n_ptr[HVEC1_0], &t_diff_n_ptr[HVEC1_1], &t_diff_n_ptr[HVEC1_2],
580 &t_diff_n_ptr[HVEC2_0], &t_diff_n_ptr[HVEC2_1], &t_diff_n_ptr[HVEC2_2]);
581}
582
583template <>
585EntitiesFieldData::EntData::getFTensor2DiffN<3, 3>(FieldApproximationBase base,
586 const int gg, const int bb) {
587 double *t_diff_n_ptr = &getDiffN(base)(gg, 9 * bb);
589 t_diff_n_ptr, &t_diff_n_ptr[HVEC0_1], &t_diff_n_ptr[HVEC0_2],
590 &t_diff_n_ptr[HVEC1_0], &t_diff_n_ptr[HVEC1_1], &t_diff_n_ptr[HVEC1_2],
591 &t_diff_n_ptr[HVEC2_0], &t_diff_n_ptr[HVEC2_1], &t_diff_n_ptr[HVEC2_2]);
592}
593
594template <>
596EntitiesFieldData::EntData::getFTensor2DiffN<3, 2>(
598 double *t_diff_n_ptr = &*getDiffN(base).data().begin();
600 t_diff_n_ptr, &t_diff_n_ptr[HVEC0_1], &t_diff_n_ptr[HVEC1_0],
601 &t_diff_n_ptr[HVEC1_1], &t_diff_n_ptr[HVEC2_0], &t_diff_n_ptr[HVEC2_1]);
602}
603
604template <>
606EntitiesFieldData::EntData::getFTensor2DiffN<3, 2>(FieldApproximationBase base,
607 const int gg, const int bb) {
608 double *t_diff_n_ptr = &getDiffN(base)(gg, 6 * bb);
610 t_diff_n_ptr, &t_diff_n_ptr[HVEC0_1], &t_diff_n_ptr[HVEC1_0],
611 &t_diff_n_ptr[HVEC1_1], &t_diff_n_ptr[HVEC2_0], &t_diff_n_ptr[HVEC2_1]);
612}
613
614template <>
617 double *ptr = &(getN(base, BaseDerivatives::SecondDerivative))(0, 0);
619
620 &ptr[2 * HVEC0_0 + 0], &ptr[2 * HVEC0_0 + 1], &ptr[2 * HVEC0_1 + 0],
621 &ptr[2 * HVEC0_1 + 1],
622
623 &ptr[2 * HVEC1_0 + 0], &ptr[2 * HVEC1_0 + 1], &ptr[2 * HVEC1_1 + 0],
624 &ptr[2 * HVEC1_1 + 1],
625
626 &ptr[2 * HVEC2_0 + 0], &ptr[2 * HVEC2_0 + 1], &ptr[2 * HVEC2_1 + 0],
627 &ptr[2 * HVEC2_1 + 1]
628
629 };
630}
631
632template <int Tensor_Dim0, int Tensor_Dim1>
634 Tensor_Dim0, Tensor_Dim1>
636 std::stringstream s;
637 s << "Template for tensor dimension " << Tensor_Dim0 << ", " << Tensor_Dim1
638 << " not implemented";
639 THROW_MESSAGE(s.str());
641 Tensor_Dim0, Tensor_Dim1>();
642}
643
644template <>
646EntitiesFieldData::EntData::getFTensor2N<3, 3>(FieldApproximationBase base) {
647 double *t_n_ptr = &*(getN(base).data().begin());
649
650 &t_n_ptr[HVEC0], &t_n_ptr[HVEC1], &t_n_ptr[HVEC2],
651
652 &t_n_ptr[3 + HVEC0], &t_n_ptr[3 + HVEC1], &t_n_ptr[3 + HVEC2],
653
654 &t_n_ptr[6 + HVEC0], &t_n_ptr[6 + HVEC1], &t_n_ptr[6 + HVEC2]
655
656 );
657}
658
659template <int Tensor_Dim0, int Tensor_Dim1>
661 Tensor_Dim0, Tensor_Dim1>
663 const int gg, const int bb) {
664 std::stringstream s;
665 s << "Template for tensor dimension " << Tensor_Dim0 << ", " << Tensor_Dim1
666 << " not implemented";
667 THROW_MESSAGE(s.str());
669 Tensor_Dim0, Tensor_Dim1>();
670}
671
672template <>
674EntitiesFieldData::EntData::getFTensor2N<3, 3>(FieldApproximationBase base,
675 const int gg, const int bb) {
676 double *t_n_ptr = &getN(base)(gg, 9 * bb);
678
679 &t_n_ptr[HVEC0], &t_n_ptr[HVEC1], &t_n_ptr[HVEC2],
680
681 &t_n_ptr[3 + HVEC0], &t_n_ptr[3 + HVEC1], &t_n_ptr[3 + HVEC2],
682
683 &t_n_ptr[6 + HVEC0], &t_n_ptr[6 + HVEC1], &t_n_ptr[6 + HVEC2]
684
685 );
686}
687
688/**@}*/
689
690/** \name Bernstein-Bezier base only functions */
691
692/**@{*/
693
694boost::shared_ptr<MatrixInt> &
696 const std::string &field_name) {
697 return bbAlphaIndices[field_name];
698}
699
700boost::shared_ptr<MatrixDouble> &
702 return bbN[field_name];
703}
704
705/**
706 * Get shared pointer to BB base base functions
707 */
708const boost::shared_ptr<MatrixDouble> &
710 const std::string &field_name) const {
711 return bbN.at(field_name);
712}
713
714/**
715 * Get shared pointer to BB derivatives of base base functions
716 */
717boost::shared_ptr<MatrixDouble> &
719 return bbDiffN[field_name];
720}
721
722/**
723 * Get shared pointer to derivatives of BB base base functions
724 */
725const boost::shared_ptr<MatrixDouble> &
727 const std::string &field_name) const {
728 return bbDiffN.at(field_name);
729}
730
731std::map<std::string, boost::shared_ptr<MatrixInt>> &
733 return bbAlphaIndices;
734}
735
736std::map<std::string, boost::shared_ptr<MatrixDouble>> &
738 return bbN;
739}
740
741std::map<std::string, boost::shared_ptr<MatrixDouble>> &
743 return bbDiffN;
744}
745
746boost::shared_ptr<MatrixInt> &
748 return bbAlphaIndicesByOrder[o];
749}
750
751boost::shared_ptr<MatrixDouble> &
753 return bbNByOrder[o];
754}
755
756boost::shared_ptr<MatrixDouble> &
758 return bbDiffNByOrder[o];
759}
760
761std::array<boost::shared_ptr<MatrixInt>,
764 return bbAlphaIndicesByOrder;
765}
766
767std::array<boost::shared_ptr<MatrixDouble>,
770 return bbNByOrder;
771}
772
773std::array<boost::shared_ptr<MatrixDouble>,
776 return bbDiffNByOrder;
777}
778
779boost::shared_ptr<MatrixInt> &
781 const std::string &field_name) {
782 return entDataPtr->getBBAlphaIndicesSharedPtr(field_name);
783}
784
785boost::shared_ptr<MatrixDouble> &
787 const std::string &field_name) {
788 return entDataPtr->getBBNSharedPtr(field_name);
789}
790
791const boost::shared_ptr<MatrixDouble> &
793 const std::string &field_name) const {
794 return entDataPtr->getBBNSharedPtr(field_name);
795}
796
797/**
798 * Get shared pointer to BB derivatives of base base functions
799 */
800boost::shared_ptr<MatrixDouble> &
802 const std::string &field_name) {
803 return entDataPtr->getBBDiffNSharedPtr(field_name);
804}
805
806/**
807 * Get shared pointer to derivatives of BB base base functions
808 */
809const boost::shared_ptr<MatrixDouble> &
811 const std::string &field_name) const {
812 return entDataPtr->getBBDiffNSharedPtr(field_name);
813}
814
815/**@}*/
816
818 return entDataBitRefLevel;
819}
820
822 return entDataPtr->getEntDataBitRefLevel();
823}
824
825
826} // namespace MoFEM
static Index< 'o', 3 > o
static const double eps
EntitiesFieldData::EntData EntData
T data[Tensor_Dim]
T data[Tensor_Dim0][Tensor_Dim1]
FieldApproximationBase
approximation base
Definition: definitions.h:71
@ LASTBASE
Definition: definitions.h:82
@ NOBASE
Definition: definitions.h:72
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
@ NOSPACE
Definition: definitions.h:96
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
@ HVEC0
Definition: definitions.h:199
@ HVEC1
Definition: definitions.h:199
@ HVEC2
Definition: definitions.h:199
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
@ HVEC1_1
Definition: definitions.h:209
@ HVEC0_1
Definition: definitions.h:208
@ HVEC1_0
Definition: definitions.h:206
@ HVEC2_1
Definition: definitions.h:210
@ HVEC1_2
Definition: definitions.h:212
@ HVEC2_2
Definition: definitions.h:213
@ HVEC2_0
Definition: definitions.h:207
@ HVEC0_2
Definition: definitions.h:211
@ HVEC0_0
Definition: definitions.h:205
#define CHKERR
Inline error check.
Definition: definitions.h:548
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
Definition: definitions.h:574
const Tensor1_Expr< const dTensor0< T, Dim, i >, typename promote< T, double >::V, Dim, i > d(const Tensor0< T * > &a, const Index< i, Dim > index, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: dTensor0.hpp:27
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.
Definition: Exceptions.hpp:67
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:88
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
static void constructor_derived_data(DerivedEntitiesFieldData *derived_data, const boost::shared_ptr< EntitiesFieldData > &data_ptr)
static void constructor_data(EntitiesFieldData *data, const EntityType type)
DataForcesAndSourcesCore::EntData EntData
constexpr double t
plate stiffness
Definition: plate.cpp:76
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.
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 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)
virtual BitRefLevel & getEntDataBitRefLevel()
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)