8#define FieldCoreFunctionBegin \
10 MOFEM_LOG_CHANNEL("WORLD"); \
11 MOFEM_LOG_CHANNEL("SYNC"); \
12 MOFEM_LOG_FUNCTION(); \
13 MOFEM_LOG_TAG("WORLD", "FieldCore"); \
14 MOFEM_LOG_TAG("SYNC", "FieldCore");
20 auto miit = set.find(name);
21 if (miit == set.end()) {
23 " > not in database (top tip: check spelling)");
25 return (*miit)->getId();
30 auto miit = set.find(name);
31 if (miit == set.end())
33 " > not in database (top tip: check spelling)");
34 return (*miit)->getBitNumber();
39 auto miit = set.find(
id);
40 if (miit == set.end())
41 THROW_MESSAGE(
"field not in database (top tip: check spelling)");
42 return (*miit)->getName();
47 auto miit = set.find(
id);
48 if (miit == set.end())
49 THROW_MESSAGE(
"field not in database (top tip: check spelling)");
50 return (*miit)->meshSet;
54 return get_field_meshset(get_field_id(name));
72 std::string(
"field < " + name +
73 " > not in database (top tip: check spelling)")
87 CHKERR get_moab().get_entities_by_dimension(meshset,
dim, ents,
true);
96 CHKERR get_moab().get_entities_by_type(meshset, type, ents,
true);
104 CHKERR get_moab().get_entities_by_handle(meshset, ents,
true);
111 const TagType tag_type,
const enum MoFEMTypes bh,
123 auto add_meshset_to_partition = [&](
auto meshset) {
125 const void *tag_vals[] = {&rAnk};
126 ParallelComm *pcomm = ParallelComm::get_pcomm(
127 &get_moab(), get_basic_entity_data_ptr()->pcommID);
128 Tag part_tag = pcomm->part_tag();
130 CHKERR get_moab().get_entities_by_type_and_tag(0, MBENTITYSET, &part_tag,
131 tag_vals, 1, tagged_sets,
132 moab::Interface::UNION);
133 for (
auto s : tagged_sets)
134 CHKERR get_moab().add_entities(s, &meshset, 1);
138 auto create_tags = [&](
auto meshset,
auto id) {
141 get_moab().tag_set_data(th_FieldId, &meshset, 1, &
id);
143 CHKERR get_moab().tag_set_data(th_FieldSpace, &meshset, 1, &space);
145 CHKERR get_moab().tag_set_data(th_FieldBase, &meshset, 1, &base);
148 void const *tag_data[] = {name.c_str()};
150 tag_sizes[0] = name.size();
151 CHKERR get_moab().tag_set_by_ptr(th_FieldName, &meshset, 1, tag_data,
156 const std::string tag_rank_name =
"_Field_Rank_" + name;
157 CHKERR get_moab().tag_get_handle(tag_rank_name.c_str(), 1, MB_TYPE_INTEGER,
158 th_rank, MB_TAG_CREAT | MB_TAG_SPARSE,
160 CHKERR get_moab().tag_set_data(th_rank, &meshset, 1, &nb_of_coefficients);
167 const std::string name_data_prefix(
"_App_Data_");
168 void const *tag_prefix_data[] = {name_data_prefix.c_str()};
169 int tag_prefix_sizes[1];
170 tag_prefix_sizes[0] = name_data_prefix.size();
171 CHKERR get_moab().tag_set_by_ptr(th_FieldName_DataNamePrefix, &meshset, 1,
172 tag_prefix_data, tag_prefix_sizes);
173 Tag th_app_order, th_field_data, th_field_data_vert;
177 const std::string tag_approx_order_name =
"_App_Order_" + name;
178 CHKERR get_moab().tag_get_handle(
179 tag_approx_order_name.c_str(), 1, MB_TYPE_INTEGER, th_app_order,
180 MB_TAG_CREAT | MB_TAG_SPARSE, &def_approx_order);
183 std::string tag_data_name = name_data_prefix + name;
184 const int def_len = 0;
185 CHKERR get_moab().tag_get_handle(
186 tag_data_name.c_str(), def_len, MB_TYPE_DOUBLE, th_field_data,
187 MB_TAG_CREAT | MB_TAG_VARLEN | MB_TAG_SPARSE, NULL);
188 std::string tag_data_name_verts = name_data_prefix + name +
"_V";
190 def_vert_data.clear();
191 CHKERR get_moab().tag_get_handle(
192 tag_data_name_verts.c_str(), nb_of_coefficients, MB_TYPE_DOUBLE,
193 th_field_data_vert, MB_TAG_CREAT | tag_type, &*def_vert_data.begin());
196 const std::string name_data_prefix(
"_App_Data");
197 void const *tag_prefix_data[] = {name_data_prefix.c_str()};
198 int tag_prefix_sizes[1];
199 tag_prefix_sizes[0] = name_data_prefix.size();
200 CHKERR get_moab().tag_set_by_ptr(th_FieldName_DataNamePrefix, &meshset, 1,
201 tag_prefix_data, tag_prefix_sizes);
202 Tag th_app_order, th_field_data, th_field_data_vert;
206 const std::string tag_approx_order_name =
"_App_Order_" + name;
207 CHKERR get_moab().tag_get_handle(
208 tag_approx_order_name.c_str(), 1, MB_TYPE_INTEGER, th_app_order,
209 MB_TAG_CREAT | MB_TAG_SPARSE, &def_approx_order);
212 std::string tag_data_name = name_data_prefix + name;
213 const int def_len = 0;
214 CHKERR get_moab().tag_get_handle(
215 tag_data_name.c_str(), def_len, MB_TYPE_DOUBLE, th_field_data,
216 MB_TAG_CREAT | MB_TAG_VARLEN | MB_TAG_SPARSE, NULL);
217 std::string tag_data_name_verts = name_data_prefix + name +
"V";
219 def_vert_data.clear();
220 CHKERR get_moab().tag_get_handle(
221 tag_data_name_verts.c_str(), nb_of_coefficients, MB_TYPE_DOUBLE,
222 th_field_data_vert, MB_TAG_CREAT | tag_type, &*def_vert_data.begin());
235 "field is <%s> in database", name.c_str());
240 CHKERR get_moab().create_meshset(MESHSET_SET, meshset);
241 CHKERR add_meshset_to_partition(meshset);
251 "Maximal number of fields exceeded");
255 CHKERR create_tags(meshset,
id);
257 auto p = fIelds.insert(boost::make_shared<Field>(moab, meshset));
259 MOFEM_LOG(
"WORLD", Sev::inform) <<
"Add field " << **
p.first;
261 <<
"Field " << (*
p.first)->getName() <<
" core value < "
262 << this->
getValue() <<
" > field value ) "
263 << (*
p.first)->getBitNumber() <<
" )";
268 "field not inserted %s (top tip, it could be already "
270 Field(moab, meshset).getName().c_str());
281 const TagType tag_type,
const enum MoFEMTypes bh,
283 return this->addField(name, space, base, nb_of_coefficients, tag_type, bh,
288 const std::string &name,
int verb) {
300 idm = get_field_meshset(name);
302 CHKERR get_moab().tag_get_data(th_FieldSpace, &idm, 1, &space);
303 std::vector<int> nb_ents_on_dim(3, 0);
306 CHKERR get_moab().add_entities(idm, ents);
309 CHKERR get_moab().add_entities(idm, ents);
310 for (
int dd = 0; dd !=
dim; ++dd) {
312 CHKERR get_moab().get_adjacencies(ents, dd,
false, adj_ents,
313 moab::Interface::UNION);
316 CHKERR get_moab().get_connectivity(ents, topo_nodes,
true);
318 CHKERR get_moab().get_connectivity(ents, mid_nodes,
false);
319 mid_nodes = subtract(mid_nodes, topo_nodes);
320 adj_ents = subtract(adj_ents, mid_nodes);
322 CHKERR get_moab().add_entities(idm, adj_ents);
323 nb_ents_on_dim[dd] = adj_ents.size();
327 CHKERR get_moab().add_entities(idm, ents);
328 for (
int dd = 1; dd !=
dim; ++dd) {
330 CHKERR get_moab().get_adjacencies(ents, dd,
false, adj_ents,
331 moab::Interface::UNION);
332 CHKERR get_moab().add_entities(idm, adj_ents);
333 nb_ents_on_dim[dd] = adj_ents.size();
337 CHKERR get_moab().add_entities(idm, ents);
340 CHKERR get_moab().get_adjacencies(ents, 2,
false, adj_ents,
341 moab::Interface::UNION);
342 CHKERR get_moab().add_entities(idm, adj_ents);
343 nb_ents_on_dim[2] = adj_ents.size();
348 "sorry, unknown space added to entity");
351 MOFEM_LOG(
"SYNC", Sev::noisy) <<
"add entities to field " << name;
352 MOFEM_LOG(
"SYNC", Sev::noisy) <<
"\tnb. add ents " << ents.size();
353 MOFEM_LOG(
"SYNC", Sev::noisy) <<
"\tnb. add faces " << nb_ents_on_dim[2];
354 MOFEM_LOG(
"SYNC", Sev::noisy) <<
"\tnb. add edges " << nb_ents_on_dim[1];
355 MOFEM_LOG(
"SYNC", Sev::noisy) <<
"\tnb. add nodes " << nb_ents_on_dim[0];
361 const std::string &name,
364 Range ents_dim = ents.subset_by_dimension(
dim);
365 CHKERR addEntsToFieldByDim(ents_dim,
dim, name, verb);
372 const std::string &name,
375 Range ents_type = ents.subset_by_type(type);
376 if (!ents_type.empty()) {
377 const int dim = get_moab().dimension_from_handle(ents_type[0]);
378 CHKERR addEntsToFieldByDim(ents_type,
dim, name, verb);
386 const std::string &name,
387 const bool recursive,
int verb) {
390 CHKERR get_moab().get_entities_by_dimension(meshset,
dim, ents, recursive);
391 CHKERR addEntsToFieldByDim(ents,
dim, name, verb);
398 const std::string &name,
399 const bool recursive,
int verb) {
402 CHKERR get_moab().get_entities_by_type(meshset, type, ents, recursive);
404 const int dim = get_moab().dimension_from_handle(ents[0]);
405 CHKERR addEntsToFieldByDim(ents,
dim, name, verb);
412 const double coords[],
413 int size,
int verb) {
421 auto create_vertices = [&]() {
424 vector<double *> arrays_coord;
426 ReadUtilIface *iface;
427 CHKERR get_moab().query_interface(iface);
428 CHKERR iface->get_node_coords(3, size, 0, startv, arrays_coord);
429 verts.insert(startv, startv + size - 1);
430 for (
int n = 0;
n != size; ++
n)
431 for (
auto d : {0, 1, 2})
432 arrays_coord[d][
n] = coords[3 *
n + d];
437 auto add_verts_to_field = [&]() {
440 CHKERR get_moab().add_entities(field_meshset, verts);
445 CHKERR add_verts_to_field();
460 const auto field_meshset = field_ptr->getMeshset();
461 const auto bit_number = field_ptr->getBitNumber();
464 Range ents_of_id_meshset;
465 CHKERR get_moab().get_entities_by_handle(field_meshset, ents_of_id_meshset,
467 Range field_ents = intersect(ents, ents_of_id_meshset);
470 "change nb. of ents for order in the field <%s> %d",
471 field_ptr->getName().c_str(), field_ents.size(), ents.size());
481 std::copy(eiit, hi_eiit, std::back_inserter(ents_id_view));
486 "current nb. of ents in the multi index field <%s> %d",
487 field_ptr->getName().c_str(), ents_id_view.size());
490 int nb_ents_set_order_up = 0;
491 int nb_ents_set_order_down = 0;
492 int nb_ents_set_order_new = 0;
497 for (
auto pit = field_ents.const_pair_begin();
498 pit != field_ents.const_pair_end(); pit++) {
505 switch (field_ptr->getSpace()) {
509 if (type == MBVERTEX)
511 "Hcurl space on vertices makes no sense");
515 if (type == MBVERTEX)
517 "Hdiv space on vertices makes no sense");
521 "Hdiv space on edges makes no sense");
529 Range ents_in_database;
530 auto vit = ents_id_view.get<1>().lower_bound(first);
531 auto hi_vit = ents_id_view.get<1>().upper_bound(second);
533 for (; vit != hi_vit; ++vit) {
534 ents_in_database.insert(vit->get()->getEnt());
537 (*vit)->getMaxOrder();
539 if (old_approximation_order !=
order) {
541 FieldEntity_multiIndex::iterator miit =
542 entsFields.get<
Unique_mi_tag>().find((*vit)->getLocalUniqueId());
544 if ((*miit)->getMaxOrder() <
order)
545 nb_ents_set_order_up++;
546 if ((*miit)->getMaxOrder() >
order)
547 nb_ents_set_order_down++;
553 bool can_change_size =
true;
563 can_change_size =
false;
564 for (; dit != hi_dit; dit++) {
565 if ((*dit)->getDofOrder() >
order) {
566 bool success = dofsField.modify(dofsField.project<0>(dit),
570 "modification unsuccessful");
576 entsFields.modify(entsFields.project<0>(miit),
577 can_change_size ? modify_order_size_change
578 : modify_order_no_size_change);
581 "modification unsuccessful");
586 Range new_ents = subtract(
Range(first, second), ents_in_database);
587 for (Range::const_pair_iterator pit = new_ents.const_pair_begin();
588 pit != new_ents.const_pair_end(); ++pit) {
589 const auto first = pit->first;
590 const auto second = pit->second;
591 const auto ent_type = get_moab().type_from_handle(first);
593 if (!field_ptr->getFieldOrderTable()[ent_type])
595 "Number of degrees of freedom for entity %s for %s space on "
598 moab::CN::EntityTypeName(ent_type),
599 field_ptr->getSpaceName().c_str(),
600 field_ptr->getApproxBaseName().c_str());
602 auto get_nb_dofs_on_order = [&](
const int order) {
603 return order >= 0 ? (field_ptr->getFieldOrderTable()[ent_type])(
order)
606 const int nb_dofs_on_order = get_nb_dofs_on_order(
order);
607 if (nb_dofs_on_order ||
order == -1) {
609 const int field_rank = field_ptr->getNbOfCoeffs();
610 const int nb_dofs = nb_dofs_on_order * field_rank;
614 refinedEntities.get<
Ent_mi_tag>().lower_bound(first);
616 auto create_tags_for_max_order = [&](
const Range &ents) {
619 std::vector<ApproximationOrder> o_vec(ents.size(),
order);
620 CHKERR get_moab().tag_set_data(field_ptr->th_AppOrder, ents,
626 auto create_tags_for_data = [&](
const Range &ents) {
631 if (ent_type == MBVERTEX) {
632 std::vector<FieldData> d_vec(nb_dofs * ents.size(), 0);
633 CHKERR get_moab().tag_set_data(field_ptr->th_FieldDataVerts,
634 ents, &*d_vec.begin());
636 std::vector<int> tag_size(ents.size(), nb_dofs);
637 std::vector<FieldData> d_vec(nb_dofs, 0);
638 std::vector<void const *> d_vec_ptr(ents.size(),
640 CHKERR get_moab().tag_set_by_ptr(field_ptr->th_FieldData, ents,
650 auto get_ents_in_ref_ent = [&](
auto miit_ref_ent) {
651 auto hi = refinedEntities.get<
Ent_mi_tag>().upper_bound(second);
653 for (; miit_ref_ent != hi; ++miit_ref_ent)
654 in.insert(miit_ref_ent->get()->getEnt());
658 auto get_ents_max_order = [&](
const Range &ents) {
659 boost::shared_ptr<std::vector<const void *>> vec(
660 new std::vector<const void *>());
661 vec->resize(ents.size());
662 CHKERR get_moab().tag_get_by_ptr(field_ptr->th_AppOrder, ents,
667 auto get_ents_field_data_vector_adaptor =
668 [&](
const Range &ents,
669 boost::shared_ptr<std::vector<const void *>> &ents_max_orders) {
671 boost::shared_ptr<std::vector<double *>> vec(
672 new std::vector<double *>());
673 vec->reserve(ents.size());
675 if (
order >= 0 && nb_dofs == 0) {
677 for (
int i = 0;
i != ents.size(); ++
i)
678 vec->emplace_back(
nullptr);
680 moab::ErrorCode
rval;
681 std::vector<int> tag_size(ents.size());
682 std::vector<const void *> d_vec_ptr(ents.size());
685 if (ent_type == MBVERTEX)
686 rval = get_moab().tag_get_by_ptr(field_ptr->th_FieldDataVerts,
687 ents, &*d_vec_ptr.begin(),
690 rval = get_moab().tag_get_by_ptr(field_ptr->th_FieldData,
691 ents, &*d_vec_ptr.begin(),
694 auto cast = [](
auto p) {
700 if (
rval == MB_SUCCESS) {
702 auto tit = d_vec_ptr.begin();
703 auto oit = ents_max_orders->begin();
704 for (
auto sit = tag_size.begin(); sit != tag_size.end();
706 vec->emplace_back(cast(*tit));
710 for (
int i = 0;
i != ents.size(); ++
i)
711 vec->emplace_back(
nullptr);
716 auto oit = ents_max_orders->begin();
717 auto dit = vec->begin();
718 for (
auto eit = ents.begin(); eit != ents.end();
719 ++eit, ++oit, ++dit) {
721 const int ent_order =
723 const int ent_nb_dofs = get_nb_dofs_on_order(ent_order);
728 if (ent_type == MBVERTEX) {
729 rval = get_moab().tag_get_by_ptr(
730 field_ptr->th_FieldDataVerts, &*eit, 1, &ret_val,
733 rval = get_moab().tag_get_by_ptr(
734 field_ptr->th_FieldData, &*eit, 1, &ret_val,
736 if (
rval != MB_SUCCESS) {
738 const int set_tag_size[] = {ent_nb_dofs * field_rank};
739 std::array<FieldData, MAX_DOFS_ON_ENTITY> set_d_vec;
740 std::fill(set_d_vec.begin(),
741 &set_d_vec[set_tag_size[0]], 0);
742 const void *set_d_vec_ptr[] = {set_d_vec.data()};
743 CHKERR get_moab().tag_set_by_ptr(
744 field_ptr->th_FieldData, &*eit, 1, set_d_vec_ptr,
746 rval = get_moab().tag_get_by_ptr(
747 field_ptr->th_FieldData, &*eit, 1, &ret_val,
750 if (
rval != MB_SUCCESS) {
755 <<
"Error is triggered in MOAB, field tag data "
756 "for same reason can not be for accessed.";
758 <<
"Set order: " <<
order;
760 <<
"Nb. dofs on entity for given order: "
764 << moab::CN::EntityTypeName(ent_type);
766 <<
"Field: " << *field_ptr;
771 const_cast<FieldData *&
>(*dit) = cast(ret_val);
779 auto ents_in_ref_ent = get_ents_in_ref_ent(miit_ref_ent);
781 CHKERR create_tags_for_max_order(ents_in_ref_ent);
782 CHKERR create_tags_for_data(ents_in_ref_ent);
783 auto ents_max_order = get_ents_max_order(ents_in_ref_ent);
784 auto ent_field_data =
785 get_ents_field_data_vector_adaptor(ents_in_ref_ent, ents_max_order);
788 auto ents_array = boost::make_shared<std::vector<FieldEntity>>();
793 ents_array->reserve(second - first + 1);
794 auto vit_max_order = ents_max_order->begin();
795 auto vit_field_data = ent_field_data->begin();
796 for (
int i = 0;
i != ents_in_ref_ent.size(); ++
i) {
798 ents_array->emplace_back(
799 field_ptr, *miit_ref_ent,
800 boost::shared_ptr<double *const>(ent_field_data,
802 boost::shared_ptr<const int>(
803 ents_max_order,
static_cast<const int *
>(*vit_max_order)));
808 if (!ents_array->empty())
809 if ((*ents_array)[0].getFieldRawPtr() != field_ptr.get())
811 "Get field ent poiter and field pointer do not match for "
813 field_ptr->getName().c_str());
814 nb_ents_set_order_new += ents_array->size();
818 if (ents_in_ref_ent.size() < (second - first + 1)) {
819 Range ents_not_in_database =
820 subtract(
Range(first, second), ents_in_ref_ent);
821 std::vector<const void *> vec_bits(ents_not_in_database.size());
822 CHKERR get_moab().tag_get_by_ptr(
823 get_basic_entity_data_ptr()->th_RefBitLevel, ents_not_in_database,
825 auto cast = [](
auto p) {
828 for (
auto v : vec_bits)
831 "Try to add entities which are not seeded or added to "
836 auto hint = entsFields.end();
837 for (
auto &
v : *ents_array)
838 hint = entsFields.emplace_hint(hint, ents_array, &
v);
845 "nb. of entities in field <%s> for which order was "
846 "increased %d (order %d)",
847 field_ptr->getName().c_str(), nb_ents_set_order_up,
order);
849 "nb. of entities in field <%s> for which order was "
850 "reduced %d (order %d)",
851 field_ptr->getName().c_str(), nb_ents_set_order_down,
order);
854 "nb. of entities in field <%s> for which order set %d (order %d)",
855 field_ptr->getName().c_str(), nb_ents_set_order_new,
order);
864 "nb. of ents in the multi index field <%s> %d",
865 field_ptr->getName().c_str(), std::distance(eiit, hi_eiit));
889 <<
"Field " << (*field_it)->getName() <<
" core value < "
890 << this->
getValue() <<
" > field value ( "
891 <<
static_cast<int>((*field_it)->getBitNumber()) <<
" )"
892 <<
" nb. of ents " << ents.size();
907 CHKERR get_moab().get_entities_by_type(meshset, type, ents);
908 CHKERR this->set_field_order(ents,
id,
order, verb);
913 const std::string &name,
919 CHKERR this->set_field_order(meshset, type, get_field_id(name),
order, verb);
928 CHKERR this->set_field_order(ents, get_field_id(name),
order, verb);
941 CHKERR this->set_field_order(ents,
id,
order, verb);
955 CHKERR this->set_field_order(ents, get_field_id(name),
order, verb);
961 std::map<EntityType, int> &dof_counter,
965 const auto bit_number = field_ptr->getBitNumber();
969 CHKERR get_moab().get_entities_by_handle(field_ptr->getMeshset(), ents,
972 MOFEM_LOG_C(
"SYNC", Sev::noisy,
"Ents in field %s meshset %d\n",
973 field_ptr->getName().c_str(), ents.size());
982 std::copy(eiit, hi_eiit, std::back_inserter(ents_id_view));
985 boost::shared_ptr<const int> zero_order(
new const int(0));
987 for (
auto ent : ents) {
989 auto ref_ent_it = refinedEntities.get<
Ent_mi_tag>().find(ent);
990 if (ref_ent_it == refinedEntities.get<
Ent_mi_tag>().end())
992 "Entity is not in MoFEM database, entities in field meshset need "
993 "to be seeded (i.e. bit ref level add to them)");
995 auto add_dofs = [&](
auto field_eit) {
1001 auto p = dofsField.insert(
1002 boost::make_shared<DofEntity>(field_eit, 0, rank, rank));
1004 dof_counter[MBENTITYSET]++;
1007 "Dof expected to be created");
1013 auto field_ent_it = ents_id_view.get<1>().find(ent);
1014 if (field_ent_it == ents_id_view.get<1>().end()) {
1016 auto p = entsFields.insert(
1018 boost::make_shared<FieldEntity>(
1019 field_ptr, *ref_ent_it,
1022 boost::shared_ptr<const int>(zero_order, zero_order.get()))
1026 if ((*
p.first)->getFieldRawPtr() != field_ptr.get())
1028 "Get field ent poiter and field pointer do not match for "
1030 field_ptr->getName().c_str());
1034 "Entity should be created");
1047 CHKERR add_dofs(*field_ent_it);
1056 for (; lo_dof != hi_dof; lo_dof++)
1057 MOFEM_LOG(
"SYNC", Sev::noisy) << **lo_dof;
1066 std::map<EntityType, int> &dof_counter,
int verb) {
1079 <<
"Field " << (*field_it)->getName() <<
" core value < "
1080 << this->
getValue() <<
" > field value () "
1081 << (*field_it)->getBitNumber() <<
" )";
1083 CHKERR this->buildFieldForNoFieldImpl(*
field_it, dof_counter, verb);
1089 const BitFieldId id, std::map<EntityType, int> &dof_counter,
1090 std::map<EntityType, int> &inactive_dof_counter,
int verb) {
1101 const int bit_number =
field_it->get()->getBitNumber();
1102 const int rank =
field_it->get()->getNbOfCoeffs();
1105 Range ents_of_id_meshset;
1106 CHKERR get_moab().get_entities_by_handle((*field_it)->meshSet,
1107 ents_of_id_meshset,
false);
1109 MOFEM_LOG_C(
"SYNC", Sev::noisy,
"Ents in field %s meshset %d",
1110 (*field_it)->getName().c_str(), ents_of_id_meshset.size());
1113 for (
auto p_eit = ents_of_id_meshset.pair_begin();
1114 p_eit != ents_of_id_meshset.pair_end(); ++p_eit) {
1123 auto feit = entsFields.get<
Unique_mi_tag>().lower_bound(lo_uid);
1126 auto hi_feit = entsFields.get<
Unique_mi_tag>().upper_bound(hi_uid);
1130 auto dit = dofsField.get<
Unique_mi_tag>().lower_bound(lo_uid);
1131 auto hi_dit = dofsField.get<
Unique_mi_tag>().upper_bound(hi_uid);
1135 boost::shared_ptr<std::vector<DofEntity>> dofs_array =
1136 boost::make_shared<std::vector<DofEntity>>(std::vector<DofEntity>());
1138 int nb_dofs_on_ents = 0;
1139 for (
auto tmp_feit = feit; tmp_feit != hi_feit; ++tmp_feit) {
1140 nb_dofs_on_ents += rank * tmp_feit->get()->getOrderNbDofs(
1141 tmp_feit->get()->getMaxOrder());
1145 dofs_array->reserve(nb_dofs_on_ents);
1148 for (; feit != hi_feit; ++feit) {
1152 for (
int oo = 0; oo <= feit->get()->getMaxOrder(); ++oo) {
1154 for (
int dd = 0; dd < feit->get()->getOrderNbDofsDiff(oo); ++dd) {
1156 for (
int rr = 0; rr < rank; ++rr, ++DD) {
1157 dofs_array->emplace_back(*feit, oo, rr, DD);
1158 ++dof_counter[feit->get()->getEntType()];
1162 if (DD > feit->get()->getNbDofsOnEnt()) {
1163 std::ostringstream ss;
1164 ss <<
"rank " << rAnk <<
" ";
1165 ss << **feit << std::endl;
1167 "Expected number of DOFs on entity not equal to number added "
1168 "to database (DD = %d != %d = "
1169 "feit->get()->getNbDofsOnEnt())\n"
1171 DD, feit->get()->getNbDofsOnEnt(), ss.str().c_str());
1176 int dofs_field_size0 = dofsField.size();
1177 auto hint = dofsField.end();
1178 for (
auto &
v : *dofs_array)
1179 hint = dofsField.emplace_hint(hint, dofs_array, &
v);
1182 field_it->get()->getDofSequenceContainer().push_back(dofs_array);
1185 if (PetscUnlikely(
static_cast<int>(dofs_array.use_count()) !=
1186 static_cast<int>(dofs_array->size() + 1))) {
1188 "Wrong use count %d != %d", dofs_array.use_count(),
1189 dofs_array->size() + 1);
1191 if (dofs_field_size0 + dofs_array->size() != dofsField.size()) {
1193 "Wrong number of inserted DOFs %d != %d", dofs_array->size(),
1194 dofsField.size() - dofs_field_size0);
1206 MOFEM_LOG(
"SYNC", Sev::verbose) <<
"Build field " << field->getName();
1208 std::map<EntityType, int> dof_counter;
1209 std::map<EntityType, int> inactive_dof_counter;
1213 if (field->getApproxBase() ==
USER_BASE)
1214 CHKERR field->rebuildDofsOrderMap();
1216 switch (field->getSpace()) {
1218 CHKERR this->buildFieldForNoField(field->getId(), dof_counter, verb);
1224 CHKERR this->buildFieldForL2H1HcurlHdiv(field->getId(), dof_counter,
1225 inactive_dof_counter, verb);
1232 int nb_added_dofs = 0;
1233 int nb_inactive_added_dofs = 0;
1234 for (
auto const &it : dof_counter) {
1236 <<
"Nb. of dofs (" << moab::CN::EntityTypeName(it.first) <<
") "
1237 << it.second <<
" (inactive " << inactive_dof_counter[it.first]
1239 nb_added_dofs += it.second;
1240 nb_inactive_added_dofs += inactive_dof_counter[it.first];
1243 <<
"Nb. added dofs " << nb_added_dofs <<
" (number of inactive dofs "
1244 << nb_inactive_added_dofs <<
" )";
1268 CHKERR this->buildField(field, verb);
1270 *buildMoFEM = 1 << 0;
1272 MOFEM_LOG(
"SYNC", Sev::inform) <<
"Number of dofs " << dofsField.size();
1286 MOFEM_LOG(
"SYNC", Sev::inform) <<
"List DOFs:";
1287 for (; dit != hi_dit; dit++)
1296 MOFEM_LOG(
"SYNC", Sev::inform) <<
"List Fields:";
1298 MOFEM_LOG(
"SYNC", Sev::inform) << *miit;
1304FieldEntityByUId::iterator
1309FieldEntityByUId::iterator
1314DofEntityByUId::iterator
1319DofEntityByUId::iterator
1324DofEntityByUId::iterator
1331DofEntityByUId::iterator
1338DofEntityByUId::iterator
1345DofEntityByUId::iterator
1358 "field not found < %s >", name.c_str());
1362 CHKERR get_moab().get_number_entities_by_handle(meshset, num_entities);
1364 auto count_field_ents = [&]() {
1365 auto bit_number = (*it)->getBitNumber();
1370 return std::distance(low_eit, hi_eit);
1373 if (count_field_ents() > (
unsigned int)num_entities) {
1375 "not equal number of entities in meshset and field multiindex "
1384 if (it->getSpace() ==
NOFIELD)
1389 CHKERR get_moab().get_number_entities_by_handle(meshset, num_entities);
1391 auto count_field_ents = [&]() {
1392 auto bit_number = it->getBitNumber();
1397 return std::distance(low_eit, hi_eit);
1400 if (count_field_ents() > (
unsigned int)num_entities) {
1402 "not equal number of entities in meshset and field "
1403 "multiindex < %s >",
1404 it->getName().c_str());
#define FieldCoreFunctionBegin
#define MOFEM_LOG_SEVERITY_SYNC(comm, severity)
Synchronise "SYNC" on curtain severity level.
#define MOFEM_LOG_SYNCHRONISE(comm)
Synchronise "SYNC" channel.
#define MOFEM_LOG_C(channel, severity, format,...)
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
FieldApproximationBase
approximation base
@ USER_BASE
user implemented approximation base
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
FieldSpace
approximation spaces
@ L2
field with C-1 continuity
@ NOFIELD
scalar or vector of scalars describe (no true field)
@ HCURL
field with continuous tangents
@ HDIV
field with continuous normal traction
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_OPERATION_UNSUCCESSFUL
@ MOFEM_DATA_INCONSISTENCY
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define BITFEID_SIZE
max number of finite elements
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
FTensor::Index< 'n', SPACE_DIM > n
MoFEMErrorCode getEntitiesByTypeAndRefLevel(const BitRefLevel bit, const BitRefLevel mask, const EntityType type, const EntityHandle meshset, int verb=0) const
add all ents from ref level given by bit to meshset
#define MOFEM_LOG(channel, severity)
Log.
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
#define MOFEM_LOG_ATTRIBUTES(channel, bit)
Add attributes to channel.
#define MOFEM_LOG_FUNCTION()
Set scope.
FTensor::Index< 'i', SPACE_DIM > i
const double v
phase velocity of light in medium (cm/ns)
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
int ApproximationOrder
Approximation on the entity.
std::bitset< BITFIELDID_SIZE > BitFieldId
Field Id.
int FieldCoefficientsNumber
Number of field coefficients.
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
char FieldBitNumber
Field bit number.
implementation of Data Operators for Forces and Sources
multi_index_container< boost::shared_ptr< FieldEntity >, indexed_by< sequenced<>, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< FieldEntity::interface_type_RefEntity, EntityHandle, &FieldEntity::getEnt > > > > FieldEntity_multiIndex_ent_view
EntityHandle get_id_for_max_type()
EntityHandle get_id_for_min_type()
Field_multiIndex::index< FieldName_mi_tag >::type::iterator field_it
const EntityHandle no_handle
No entity handle is indicated by zero handle, i.e. root meshset.
constexpr auto field_name
std::string get_field_name(const BitFieldId id) const
get field name from id
MoFEMErrorCode set_field_order_by_entity_type_and_bit_ref(const BitRefLevel &bit, const BitRefLevel &mask, const EntityType type, const BitFieldId id, const ApproximationOrder order, int verb=DEFAULT_VERBOSITY)
const Field * get_field_structure(const std::string &name, enum MoFEMTypes bh=MF_EXIST) const
get field structure
MoFEMErrorCode add_ents_to_field_by_dim(const Range &ents, const int dim, const std::string &name, int verb=DEFAULT_VERBOSITY)
Add entities to field meshset.
bool check_field(const std::string &name) const
check if field is in database
DofEntityByUId::iterator get_dofs_by_name_and_type_end(const std::string &field_name, const EntityType ent) const
get begin iterator of filed dofs of given name end ent type(instead you can use IT_GET_DOFS_FIELD_BY_...
MoFEMErrorCode get_field_entities_by_dimension(const std::string name, int dim, Range &ents) const
get entities in the field by dimension
MoFEMErrorCode add_ents_to_field_by_type(const Range &ents, const EntityType type, const std::string &name, int verb=DEFAULT_VERBOSITY)
Add entities to field meshset.
MoFEMErrorCode buildFieldForNoFieldImpl(boost::shared_ptr< Field > field_ptr, std::map< EntityType, int > &dof_counter, int verb)
EntityHandle get_field_meshset(const BitFieldId id) const
MoFEMErrorCode create_vertices_and_add_to_field(const std::string name, const double coords[], int size, int verb=DEFAULT_VERBOSITY)
Create a vertices and add to field object.
MoFEMErrorCode set_field_order(const Range &ents, const BitFieldId id, const ApproximationOrder order, int verb=DEFAULT_VERBOSITY)
DofEntityByUId::iterator get_dofs_by_name_and_type_begin(const std::string &field_name, const EntityType type) const
get begin iterator of filed dofs of given name and ent type (instead you can use IT_GET_DOFS_FIELD_BY...
FieldEntityByUId::iterator get_ent_field_by_name_begin(const std::string &field_name) const
get begin iterator of filed ents of given name (instead you can use IT_GET_ENT_FIELD_BY_NAME_FOR_LOOP...
MoFEMErrorCode addEntsToFieldByDim(const Range &ents, const int dim, const std::string &name, int verb=DEFAULT_VERBOSITY)
FieldBitNumber get_field_bit_number(const std::string name) const
get field bit number
BitFieldId get_field_id(const std::string &name) const
Get field Id.
DofEntityByUId::iterator get_dofs_by_name_and_ent_begin(const std::string &field_name, const EntityHandle ent) const
get begin iterator of filed dofs of given name and ent(instead you can use IT_GET_DOFS_FIELD_BY_NAME_...
DofEntityByUId::iterator get_dofs_by_name_and_ent_end(const std::string &field_name, const EntityHandle ent) const
get begin iterator of filed dofs of given name and ent (instead you can use IT_GET_DOFS_FIELD_BY_NAME...
MoFEMErrorCode check_number_of_ents_in_ents_field() const
check data consistency in entitiesPtr
MoFEMErrorCode build_field(const std::string field_name, int verb=DEFAULT_VERBOSITY)
build field by name
DofEntityByUId::iterator get_dofs_by_name_begin(const std::string &field_name) const
get begin iterator of filed dofs of given name (instead you can use IT_GET_DOFS_FIELD_BY_NAME_FOR_LOO...
virtual MoFEMErrorCode add_field(const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)
Add filed.
MoFEMErrorCode setFieldOrderImpl(boost::shared_ptr< Field > field_ptr, const Range &ents, const ApproximationOrder order, int verb)
MoFEMErrorCode buildField(const boost::shared_ptr< Field > &field, int verb=DEFAULT_VERBOSITY)
MoFEMErrorCode list_dofs_by_field_name(const std::string &name) const
MoFEMErrorCode buildFieldForNoField(const BitFieldId id, std::map< EntityType, int > &dof_counter, int verb=DEFAULT_VERBOSITY)
DofEntityByUId::iterator get_dofs_by_name_end(const std::string &field_name) const
get begin iterator of filed dofs of given name (instead you can use IT_GET_DOFS_FIELD_BY_NAME_FOR_LOO...
FieldEntityByUId::iterator get_ent_field_by_name_end(const std::string &field_name) const
get begin iterator of filed dofs of given name (instead you can use IT_GET_ENT_FIELD_BY_NAME_FOR_LOOP...
MoFEMErrorCode build_fields(int verb=DEFAULT_VERBOSITY)
MoFEMErrorCode get_field_entities_by_handle(const std::string name, Range &ents) const
get entities in the field by handle
MoFEMErrorCode addField(const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_coefficients, const TagType tag_type, const enum MoFEMTypes bh, int verb)
Template for add_field.
MoFEMErrorCode buildFieldForL2H1HcurlHdiv(const BitFieldId id, std::map< EntityType, int > &dof_counter, std::map< EntityType, int > &inactive_dof_counter, int verb=DEFAULT_VERBOSITY)
MoFEMErrorCode list_fields() const
list entities in the field
MoFEMErrorCode get_field_entities_by_type(const std::string name, EntityType type, Range &ents) const
get entities in the field by type
const int getValue() const
structure to change FieldEntity order
static UId getHiBitNumberUId(const FieldBitNumber bit_number)
static boost::shared_ptr< FieldData *const > makeSharedFieldDataAdaptorPtr(const boost::shared_ptr< Field > &field_ptr, const boost::shared_ptr< RefEntity > &ref_ents_ptr)
Return shared pointer to entity field data vector adaptor.
static UId getLoBitNumberUId(const FieldBitNumber bit_number)
static UId getLoLocalEntityBitNumber(const char bit_number, const EntityHandle ent)
static UId getHiLocalEntityBitNumber(const char bit_number, const EntityHandle ent)
Provide data structure for (tensor) field approximation.
MultiIndex Tag for field name.
static MoFEMErrorCode getFileVersion(moab::Interface &moab, Version &version)
Get database major version.