v0.14.0
Loading...
Searching...
No Matches
Files | Classes | Functions
Index sets (IS)

Construct index sets for MoFEM problems. More...

Collaboration diagram for Index sets (IS):

Files

file  FieldBlas.hpp
 Field basic algebra.
 
file  ISManager.cpp
 IS creating.
 
file  ISManager.hpp
 Interface managing IS.
 

Classes

struct  MoFEM::ISManager
 Section manager is used to create indexes and sections. More...
 

Functions

MoFEMErrorCode MoFEM::ISManager::sectionCreate (const std::string problem_name, PetscSection *s, const RowColData row_col=COL) const
 Create global selection. More...
 
SmartPetscObj< PetscSection > MoFEM::ISManager::sectionCreate (const std::string problem_name, const RowColData row_col=COL) const
 Create global selection. More...
 
MoFEMErrorCode MoFEM::ISManager::isCreateProblemOrder (const std::string problem_name, RowColData rc, int min_order, int max_order, IS *is) const
 create IS for given order range (collective) More...
 
MoFEMErrorCode MoFEM::ISManager::isCreateProblemFieldAndRank (const std::string problem_name, RowColData rc, const std::string field, int min_coeff_idx, int max_coeff_idx, IS *is, Range *ents=nullptr) const
 create IS for given problem, field and rank range (collective) More...
 
MoFEMErrorCode MoFEM::ISManager::isCreateProblemFieldAndRank (const std::string problem_name, RowColData rc, const std::string field, int min_coeff_idx, int max_coeff_idx, SmartPetscObj< IS > &smart_is, Range *ents=nullptr) const
 IS for given problem, field and rank range (collective) More...
 
MoFEMErrorCode MoFEM::ISManager::isCreateProblemFieldAndRankLocal (const std::string problem_name, RowColData rc, const std::string field, int min_coeff_idx, int max_coeff_idx, IS *is, Range *ents=nullptr) const
 create IS for given problem, field and rank range (collective) More...
 
MoFEMErrorCode MoFEM::ISManager::isCreateProblemFieldAndRankLocal (const std::string problem_name, RowColData rc, const std::string field, int min_coeff_idx, int max_coeff_idx, SmartPetscObj< IS > &smart_is, Range *ents=nullptr) const
 IS for given problem, field and rank range (collective) More...
 
MoFEMErrorCode MoFEM::ISManager::isCreateProblemFieldAndEntityType (const std::string problem_name, RowColData rc, const std::string field, EntityType low_type, EntityType hi_type, int min_coeff_idx, int max_coeff_idx, IS *is, Range *ents=nullptr) const
 create IS for given problem, field and type range (collective) More...
 
MoFEMErrorCode MoFEM::ISManager::isCreateFromProblemFieldToOtherProblemField (const std::string x_problem, const std::string x_field_name, RowColData x_rc, const std::string y_problem, const std::string y_field_name, RowColData y_rc, std::vector< int > &idx, std::vector< int > &idy) const
 create IS for give two problems and field More...
 
MoFEMErrorCode MoFEM::ISManager::isCreateFromProblemFieldToOtherProblemField (const std::string x_problem, const std::string x_field_name, RowColData x_rc, const std::string y_problem, const std::string y_field_name, RowColData y_rc, IS *ix, IS *iy) const
 create IS for give two problems and field More...
 
MoFEMErrorCode MoFEM::ISManager::isCreateFromProblemToOtherProblem (const std::string x_problem, RowColData x_rc, const std::string y_problem, RowColData y_rc, std::vector< int > &idx, std::vector< int > &idy) const
 Create is from one problem to other problem. More...
 
MoFEMErrorCode MoFEM::ISManager::isCreateFromProblemToOtherProblem (const std::string x_problem, RowColData x_rc, const std::string y_problem, RowColData y_rc, IS *ix, IS *iy) const
 Create is from one problem to other problem. More...
 

Detailed Description

Construct index sets for MoFEM problems.

Function Documentation

◆ isCreateFromProblemFieldToOtherProblemField() [1/2]

MoFEMErrorCode MoFEM::ISManager::isCreateFromProblemFieldToOtherProblemField ( const std::string  x_problem,
const std::string  x_field_name,
RowColData  x_rc,
const std::string  y_problem,
const std::string  y_field_name,
RowColData  y_rc,
IS *  ix,
IS *  iy 
) const

create IS for give two problems and field

Indices are sorted by global PETSc index in problem_x.

Parameters
x_problemname of problem
x_field_namename of field in problem_x
x_rcthat is ROW or COL
y_problemname of problem
y_field_namename of field in problem_y
y_rcthat is ROW or COL
Return values
ixIS indexes in problem_x (can be PETSC_NULL)
iyIS indexes in problem_y

Definition at line 570 of file ISManager.cpp.

573 {
575 const MoFEM::Interface &m_field = cOre;
576 std::vector<int> idx(0), idy(0);
578 x_problem, x_field_name, x_rc, y_problem, y_field_name, y_rc, idx, idy);
579 if (ix != PETSC_NULL) {
580 CHKERR ISCreateGeneral(m_field.get_comm(), idx.size(), &idx[0],
581 PETSC_COPY_VALUES, ix);
582 }
583 CHKERR ISCreateGeneral(m_field.get_comm(), idy.size(), &idy[0],
584 PETSC_COPY_VALUES, iy);
585 if (dEbug) {
586 ISView(*ix, PETSC_VIEWER_STDOUT_WORLD);
587 ISView(*iy, PETSC_VIEWER_STDOUT_WORLD);
588 }
590}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
MoFEMErrorCode isCreateFromProblemFieldToOtherProblemField(const std::string x_problem, const std::string x_field_name, RowColData x_rc, const std::string y_problem, const std::string y_field_name, RowColData y_rc, std::vector< int > &idx, std::vector< int > &idy) const
create IS for give two problems and field
Definition: ISManager.cpp:467
virtual MPI_Comm & get_comm() const =0
Deprecated interface functions.
const MoFEM::Interface & cOre
Definition: ISManager.hpp:28

◆ isCreateFromProblemFieldToOtherProblemField() [2/2]

MoFEMErrorCode MoFEM::ISManager::isCreateFromProblemFieldToOtherProblemField ( const std::string  x_problem,
const std::string  x_field_name,
RowColData  x_rc,
const std::string  y_problem,
const std::string  y_field_name,
RowColData  y_rc,
std::vector< int > &  idx,
std::vector< int > &  idy 
) const

create IS for give two problems and field

Note that indices are ordered in ascending order of local indices in problem_y

Parameters
x_problemname of problem
x_field_namename of field in problem_x
x_rcthat is ROW or COL
y_problemname of problem
y_field_namename of field in problem_y
y_rcthat is ROW or COL
Return values
idxindexes in problem_x
idyindexes in problem_y

Definition at line 467 of file ISManager.cpp.

471 {
472 const MoFEM::Interface &m_field = cOre;
473 const Problem *px_ptr;
474 const Problem *py_ptr;
476
477 CHKERR m_field.get_problem(x_problem, &px_ptr);
478 CHKERR m_field.get_problem(y_problem, &py_ptr);
479
480 typedef multi_index_container<
481 boost::shared_ptr<NumeredDofEntity>,
482
483 indexed_by<
484
485 sequenced<>,
486
487 ordered_non_unique<
488 tag<Composite_Ent_And_EntDofIdx_mi_tag>,
489 composite_key<
490 NumeredDofEntity,
495
496 >>
497 NumeredDofEntity_view_multiIndex;
498
499 NumeredDofEntity_view_multiIndex dofs_view;
500
501 auto x_bit_number = m_field.get_field_bit_number(x_field_name);
502
503 switch (x_rc) {
504 case ROW:
505 dofs_view.insert(
506 dofs_view.end(),
507 px_ptr->numeredRowDofsPtr->get<Unique_mi_tag>().lower_bound(
508 FieldEntity::getLoBitNumberUId(x_bit_number)),
509 px_ptr->numeredRowDofsPtr->get<Unique_mi_tag>().upper_bound(
510 FieldEntity::getHiBitNumberUId(x_bit_number)));
511 break;
512 case COL:
513 dofs_view.insert(
514 dofs_view.end(),
515 px_ptr->numeredColDofsPtr->get<Unique_mi_tag>().lower_bound(
516 FieldEntity::getLoBitNumberUId(x_bit_number)),
517 px_ptr->numeredColDofsPtr->get<Unique_mi_tag>().upper_bound(
518 FieldEntity::getHiBitNumberUId(x_bit_number)));
519 break;
520 default:
521 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
522 "only makes sense for ROWS and COLS");
523 }
524
525 decltype(py_ptr->numeredRowDofsPtr) dofs_ptr;
526 switch (y_rc) {
527 case ROW:
528 dofs_ptr = py_ptr->numeredRowDofsPtr;
529 break;
530 case COL:
531 dofs_ptr = py_ptr->numeredColDofsPtr;
532 break;
533 default:
534 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
535 "only makes sense for ROWS and COLS");
536 }
537
538 std::map<int, int> global_dofs_map;
539 const auto y_bit_number = m_field.get_field_bit_number(y_field_name);
540 auto dit = dofs_ptr->get<Unique_mi_tag>().lower_bound(
541 FieldEntity::getLoBitNumberUId(y_bit_number));
542 auto hi_dit = dofs_ptr->get<Unique_mi_tag>().upper_bound(
543 FieldEntity::getHiBitNumberUId(y_bit_number));
544 const auto rank = m_field.get_comm_rank();
545 for (; dit != hi_dit; ++dit) {
546 if ((*dit)->getPart() == rank) {
547 auto x_dit = dofs_view.get<Composite_Ent_And_EntDofIdx_mi_tag>().find(
548 boost::make_tuple((*dit)->getEnt(), (*dit)->getEntDofIdx()));
549 if (x_dit != dofs_view.get<Composite_Ent_And_EntDofIdx_mi_tag>().end()) {
550 global_dofs_map[(*x_dit)->getPetscGlobalDofIdx()] =
551 (*dit)->getPetscGlobalDofIdx();
552 }
553 }
554 }
555
556 idx.resize(global_dofs_map.size());
557 idy.resize(global_dofs_map.size());
558 {
559 auto ix = idx.begin();
560 auto iy = idy.begin();
561 for (auto mit = global_dofs_map.begin(); mit != global_dofs_map.end();
562 mit++, ix++, iy++) {
563 *ix = mit->first;
564 *iy = mit->second;
565 }
566 }
568}
@ COL
Definition: definitions.h:123
@ ROW
Definition: definitions.h:123
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
virtual const Problem * get_problem(const std::string problem_name) const =0
Get the problem object.
int DofIdx
Index of DOF.
Definition: Types.hpp:18
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number
virtual int get_comm_rank() const =0
static UId getHiBitNumberUId(const FieldBitNumber bit_number)
static UId getLoBitNumberUId(const FieldBitNumber bit_number)
interface_DofEntity< DofEntity > interface_type_DofEntity

◆ isCreateFromProblemToOtherProblem() [1/2]

MoFEMErrorCode MoFEM::ISManager::isCreateFromProblemToOtherProblem ( const std::string  x_problem,
RowColData  x_rc,
const std::string  y_problem,
RowColData  y_rc,
IS *  ix,
IS *  iy 
) const

Create is from one problem to other problem.

Parameters
x_problem
x_rc
y_problem
y_rc
ix
iy
Returns
MoFEMErrorCode

Definition at line 645 of file ISManager.cpp.

647 {
648 const MoFEM::Interface &m_field = cOre;
650 std::vector<int> idx(0), idy(0);
651 CHKERR isCreateFromProblemToOtherProblem(x_problem, x_rc, y_problem, y_rc,
652 idx, idy);
653 CHKERR ISCreateGeneral(m_field.get_comm(), idx.size(), &idx[0],
654 PETSC_COPY_VALUES, ix);
655 CHKERR ISCreateGeneral(m_field.get_comm(), idy.size(), &idy[0],
656 PETSC_COPY_VALUES, iy);
657 if (dEbug) {
658 ISView(*ix, PETSC_VIEWER_STDOUT_WORLD);
659 ISView(*iy, PETSC_VIEWER_STDOUT_WORLD);
660 }
662}
MoFEMErrorCode isCreateFromProblemToOtherProblem(const std::string x_problem, RowColData x_rc, const std::string y_problem, RowColData y_rc, std::vector< int > &idx, std::vector< int > &idy) const
Create is from one problem to other problem.
Definition: ISManager.cpp:592

◆ isCreateFromProblemToOtherProblem() [2/2]

MoFEMErrorCode MoFEM::ISManager::isCreateFromProblemToOtherProblem ( const std::string  x_problem,
RowColData  x_rc,
const std::string  y_problem,
RowColData  y_rc,
std::vector< int > &  idx,
std::vector< int > &  idy 
) const

Create is from one problem to other problem.

Parameters
x_problem
x_rc
y_problem
y_rc
idx
idy
Returns
MoFEMErrorCode

Definition at line 592 of file ISManager.cpp.

594 {
595 const MoFEM::Interface &m_field = cOre;
596 const Problem *px_ptr;
597 const Problem *py_ptr;
599 CHKERR m_field.get_problem(x_problem, &px_ptr);
600 CHKERR m_field.get_problem(y_problem, &py_ptr);
601 NumeredDofEntityByLocalIdx::iterator y_dit, hi_y_dit;
602 switch (y_rc) {
603 case ROW:
604 y_dit =
605 py_ptr->numeredRowDofsPtr->get<PetscLocalIdx_mi_tag>().lower_bound(0);
606 hi_y_dit =
607 py_ptr->numeredRowDofsPtr->get<PetscLocalIdx_mi_tag>().lower_bound(
608 py_ptr->getNbLocalDofsRow()); // should be lower
609 break;
610 case COL:
611 y_dit =
612 py_ptr->numeredColDofsPtr->get<PetscLocalIdx_mi_tag>().lower_bound(0);
613 hi_y_dit =
614 py_ptr->numeredColDofsPtr->get<PetscLocalIdx_mi_tag>().lower_bound(
615 py_ptr->getNbLocalDofsCol()); // should be lower
616 break;
617 default:
618 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
619 }
620 const NumeredDofEntityByUId *x_numered_dofs_by_uid;
621 switch (x_rc) {
622 case ROW:
623 x_numered_dofs_by_uid = &(px_ptr->numeredRowDofsPtr->get<Unique_mi_tag>());
624 break;
625 case COL:
626 x_numered_dofs_by_uid = &(px_ptr->numeredColDofsPtr->get<Unique_mi_tag>());
627 break;
628 default:
629 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
630 }
631 for (; y_dit != hi_y_dit; y_dit++) {
632 if ((*y_dit)->getPart() != (unsigned int)m_field.get_comm_rank()) {
633 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
634 }
635 NumeredDofEntityByUId::iterator x_dit;
636 x_dit = x_numered_dofs_by_uid->find((*y_dit)->getLocalUniqueId());
637 if (x_dit == x_numered_dofs_by_uid->end())
638 continue;
639 idx.push_back((*x_dit)->getPetscGlobalDofIdx());
640 idy.push_back((*y_dit)->getPetscGlobalDofIdx());
641 }
643}
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
NumeredDofEntity_multiIndex::index< Unique_mi_tag >::type NumeredDofEntityByUId
Numbered DoF multi-index by UId.

◆ isCreateProblemFieldAndEntityType()

MoFEMErrorCode MoFEM::ISManager::isCreateProblemFieldAndEntityType ( const std::string  problem_name,
RowColData  rc,
const std::string  field,
EntityType  low_type,
EntityType  hi_type,
int  min_coeff_idx,
int  max_coeff_idx,
IS *  is,
Range ents = nullptr 
) const

create IS for given problem, field and type range (collective)

Parameters
problem
rc
field
low_type
hi_type
min_coeff_idx
max_coeff_idx
is
ents
Returns
MoFEMErrorCode

Definition at line 449 of file ISManager.cpp.

452 {
453 const MoFEM::Interface &m_field = cOre;
455 EntityHandle field_meshset = m_field.get_field_meshset(field);
456 Range ents;
457 for (; low_type <= hi_type; ++low_type)
458 CHKERR m_field.get_moab().get_entities_by_type(field_meshset, low_type,
459 ents, true);
460 if (ents_ptr)
461 ents = intersect(ents, *ents_ptr);
462 CHKERR isCreateProblemFieldAndRank(problem_name, rc, field, min_coeff_idx,
463 max_coeff_idx, is, &ents);
465}
MoFEMErrorCode isCreateProblemFieldAndRank(const std::string problem_name, RowColData rc, const std::string field, int min_coeff_idx, int max_coeff_idx, IS *is, Range *ents=nullptr) const
create IS for given problem, field and rank range (collective)
Definition: ISManager.cpp:265
virtual moab::Interface & get_moab()=0
virtual EntityHandle get_field_meshset(const std::string name) const =0
get field meshset

◆ isCreateProblemFieldAndRank() [1/2]

MoFEMErrorCode MoFEM::ISManager::isCreateProblemFieldAndRank ( const std::string  problem_name,
RowColData  rc,
const std::string  field,
int  min_coeff_idx,
int  max_coeff_idx,
IS *  is,
Range ents = nullptr 
) const

create IS for given problem, field and rank range (collective)

Parameters
problemname
rcROW or COL
fieldname
min_coeff_idx
max_coeff_idx
entsif not null get dofs only on given entities
Return values
isout value

Definition at line 265 of file ISManager.cpp.

267 {
268 const MoFEM::Interface &m_field = cOre;
269 const Problem *problem_ptr;
271 CHKERR m_field.get_problem(problem_name, &problem_ptr);
272 const auto bit_number = m_field.get_field_bit_number(field);
273
274 using DofsByUId = NumeredDofEntity_multiIndex::index<Unique_mi_tag>::type;
275 DofsByUId::iterator it, hi_it;
276 int nb_loc_dofs;
277 switch (rc) {
278 case ROW:
279 nb_loc_dofs = problem_ptr->getNbLocalDofsRow();
280 it = problem_ptr->numeredRowDofsPtr->get<Unique_mi_tag>().lower_bound(
282 hi_it = problem_ptr->numeredRowDofsPtr->get<Unique_mi_tag>().upper_bound(
284 break;
285 case COL:
286 nb_loc_dofs = problem_ptr->getNbLocalDofsCol();
287 it = problem_ptr->numeredColDofsPtr->get<Unique_mi_tag>().lower_bound(
289 hi_it = problem_ptr->numeredColDofsPtr->get<Unique_mi_tag>().upper_bound(
291 break;
292 default:
293 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
294 }
295
296 std::vector<int> idx_vec;
297 idx_vec.reserve(std::distance(it, hi_it));
298 for (; it != hi_it; ++it) {
299
300 auto true_if_dof_on_entity = [&]() {
301 if (ents_ptr) {
302 return ents_ptr->find((*it)->getEnt()) != ents_ptr->end();
303 } else {
304 return true;
305 }
306 };
307
308 auto check = [&]() {
309 const auto coeff_idx = (*it)->getDofCoeffIdx();
310 if (
311
312 (*it)->getPetscLocalDofIdx() >= nb_loc_dofs ||
313
314 coeff_idx < min_coeff_idx || coeff_idx > max_coeff_idx
315
316 )
317 return false;
318 else
319 return true;
320 };
321
322 if (check()) {
323 if (true_if_dof_on_entity()) {
324 idx_vec.emplace_back((*it)->getPetscGlobalDofIdx());
325 }
326 }
327 }
328
329 int *id;
330 CHKERR PetscMalloc(idx_vec.size() * sizeof(int), &id);
331 std::copy(idx_vec.begin(), idx_vec.end(), id);
332 CHKERR ISCreateGeneral(m_field.get_comm(), idx_vec.size(), id,
333 PETSC_OWN_POINTER, is);
334
336}

◆ isCreateProblemFieldAndRank() [2/2]

MoFEMErrorCode MoFEM::ISManager::isCreateProblemFieldAndRank ( const std::string  problem_name,
RowColData  rc,
const std::string  field,
int  min_coeff_idx,
int  max_coeff_idx,
SmartPetscObj< IS > &  smart_is,
Range ents = nullptr 
) const

IS for given problem, field and rank range (collective)

Parameters
problemname
rcROW or COL
fieldname
min_coeff_idx
max_coeff_idx
entsif not null get dofs only on given entities
Return values
isout value

Definition at line 338 of file ISManager.cpp.

341 {
343 IS is;
344 CHKERR isCreateProblemFieldAndRank(problem_name, rc, field, min_coeff_idx,
345 max_coeff_idx, &is, ents_ptr);
346 smart_is = SmartPetscObj<IS>(is);
348}

◆ isCreateProblemFieldAndRankLocal() [1/2]

MoFEMErrorCode MoFEM::ISManager::isCreateProblemFieldAndRankLocal ( const std::string  problem_name,
RowColData  rc,
const std::string  field,
int  min_coeff_idx,
int  max_coeff_idx,
IS *  is,
Range ents = nullptr 
) const

create IS for given problem, field and rank range (collective)

Parameters
problemname
rcROW or COL
fieldname
min_coeff_idx
max_coeff_idx
entsif not null get dofs only on given entities
Return values
isout value

Definition at line 350 of file ISManager.cpp.

352 {
353 const MoFEM::Interface &m_field = cOre;
354 const Problem *problem_ptr;
356 CHKERR m_field.get_problem(problem_name, &problem_ptr);
357 const auto bit_number = m_field.get_field_bit_number(field);
358
359 auto get_low_hi_uid = [&]() {
360 return std::make_pair(FieldEntity::getLoBitNumberUId(bit_number),
362 };
363
364 auto get_low_hi_uid_by_entities = [&](auto f, auto s) {
365 return std::make_pair(DofEntity::getLoFieldEntityUId(bit_number, f),
366 DofEntity::getHiFieldEntityUId(bit_number, s));
367 };
368
369 auto get_low_hi = [&](auto lo_uid, auto hi_uid) {
370 using DofsByUId = NumeredDofEntity_multiIndex::index<Unique_mi_tag>::type;
371 DofsByUId::iterator it, hi_it;
372 switch (rc) {
373 case ROW:
374 it = problem_ptr->numeredRowDofsPtr->get<Unique_mi_tag>().lower_bound(
375 lo_uid);
376 hi_it = problem_ptr->numeredRowDofsPtr->get<Unique_mi_tag>().upper_bound(
377 hi_uid);
378 break;
379 case COL:
380 it = problem_ptr->numeredColDofsPtr->get<Unique_mi_tag>().lower_bound(
381 lo_uid);
382 hi_it = problem_ptr->numeredColDofsPtr->get<Unique_mi_tag>().upper_bound(
383 hi_uid);
384 break;
385 default:
386 THROW_MESSAGE("not implemented");
387 }
388 return std::make_pair(it, hi_it);
389 };
390
391 auto check = [&](auto it) {
392 const auto coeff_idx = (*it)->getDofCoeffIdx();
393 if (
394
395 coeff_idx < min_coeff_idx || coeff_idx > max_coeff_idx
396
397 )
398 return false;
399 else
400 return true;
401 };
402
403 auto emplace_indices = [&](auto it, auto hi_it, auto &idx_vec) {
404 for (; it != hi_it; ++it) {
405 if (check(it))
406 idx_vec.emplace_back((*it)->getPetscLocalDofIdx());
407 }
408 };
409
410 auto [lo_uid, hi_uid] = get_low_hi_uid();
411 auto [lo, hi] = get_low_hi(lo_uid, hi_uid);
412 std::vector<int> idx_vec;
413 idx_vec.reserve(std::distance(lo, hi));
414
415 if (ents_ptr) {
416 for (auto pit = ents_ptr->const_pair_begin();
417 pit != ents_ptr->const_pair_end(); ++pit) {
418 auto [lo_uid, hi_uid] =
419 get_low_hi_uid_by_entities(pit->first, pit->second);
420 auto [lo, hi] = get_low_hi(lo_uid, hi_uid);
421 emplace_indices(lo, hi, idx_vec);
422 }
423 } else {
424 emplace_indices(lo, hi, idx_vec);
425 }
426
427 int *id;
428 CHKERR PetscMalloc(idx_vec.size() * sizeof(int), &id);
429 std::copy(idx_vec.begin(), idx_vec.end(), id);
430 CHKERR ISCreateGeneral(PETSC_COMM_SELF, idx_vec.size(), id, PETSC_OWN_POINTER,
431 is);
432
434}
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
Definition: definitions.h:561
static UId getLoFieldEntityUId(const FieldBitNumber bit, const EntityHandle ent)
static UId getHiFieldEntityUId(const FieldBitNumber bit, const EntityHandle ent)

◆ isCreateProblemFieldAndRankLocal() [2/2]

MoFEMErrorCode MoFEM::ISManager::isCreateProblemFieldAndRankLocal ( const std::string  problem_name,
RowColData  rc,
const std::string  field,
int  min_coeff_idx,
int  max_coeff_idx,
SmartPetscObj< IS > &  smart_is,
Range ents = nullptr 
) const

IS for given problem, field and rank range (collective)

Parameters
problemname
rcROW or COL
fieldname
min_coeff_idx
max_coeff_idx
entsif not null get dofs only on given entities
Return values
isout value

Definition at line 436 of file ISManager.cpp.

439 {
441 IS is;
442 CHKERR isCreateProblemFieldAndRankLocal(problem_name, rc, field,
443 min_coeff_idx, max_coeff_idx, &is,
444 ents_ptr);
445 smart_is = SmartPetscObj<IS>(is);
447}
MoFEMErrorCode isCreateProblemFieldAndRankLocal(const std::string problem_name, RowColData rc, const std::string field, int min_coeff_idx, int max_coeff_idx, IS *is, Range *ents=nullptr) const
create IS for given problem, field and rank range (collective)
Definition: ISManager.cpp:350

◆ isCreateProblemOrder()

MoFEMErrorCode MoFEM::ISManager::isCreateProblemOrder ( const std::string  problem_name,
RowColData  rc,
int  min_order,
int  max_order,
IS *  is 
) const

create IS for given order range (collective)

Parameters
problemname
rcROW or COL
min_order
max_order
Return values
isout value

Definition at line 197 of file ISManager.cpp.

199 {
200 const MoFEM::Interface &m_field = cOre;
201 const Problem *problem_ptr;
203 CHKERR m_field.get_problem(problem_name, &problem_ptr);
204
205 typedef multi_index_container<
206 boost::shared_ptr<NumeredDofEntity>,
207
208 indexed_by<
209
210 sequenced<>,
211
212 ordered_non_unique<
213 tag<Order_mi_tag>,
216
217 >>
218 NumeredDofEntity_order_view_multiIndex;
219
220 const int rank = m_field.get_comm_rank();
221
222 NumeredDofEntity_order_view_multiIndex dofs_by_order;
223 auto insert_part_range = [&dofs_by_order, rank](auto &dofs) {
224 dofs_by_order.insert(dofs_by_order.end(), dofs.lower_bound(rank),
225 dofs.upper_bound(rank));
226 };
227
228 switch (rc) {
229 case ROW:
230 insert_part_range(problem_ptr->numeredRowDofsPtr->get<Part_mi_tag>());
231 break;
232 case COL:
233 insert_part_range(problem_ptr->numeredColDofsPtr->get<Part_mi_tag>());
234 break;
235 default:
236 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
237 }
238
239 auto lo = dofs_by_order.get<Order_mi_tag>().lower_bound(min_order);
240 auto hi = dofs_by_order.get<Order_mi_tag>().upper_bound(max_order);
241 const int size = std::distance(lo, hi);
242 int *id;
243 CHKERR PetscMalloc(size * sizeof(int), &id);
244 int *id_it = id;
245 for (; lo != hi; ++lo, ++id_it)
246 *id_it = (*lo)->getPetscGlobalDofIdx();
247 sort(id, &id[size]);
248
249 CHKERR ISCreateGeneral(PETSC_COMM_WORLD, size, id, PETSC_OWN_POINTER, is);
250
252}
int ApproximationOrder
Approximation on the entity.
Definition: Types.hpp:26

◆ sectionCreate() [1/2]

SmartPetscObj< PetscSection > MoFEM::ISManager::sectionCreate ( const std::string  problem_name,
const RowColData  row_col = COL 
) const

Create global selection.

Create section for given problem, such that points are sorted by UId, and number of dofs on point is equal to number of dofs on entity

It takes all fields

Parameters
problem_name
fields_list
s
row_colROE or COL, default is ROW
Returns
error code

Definition at line 140 of file ISManager.cpp.

141 {
142
143 PetscSection s;
144 CHK_THROW_MESSAGE(sectionCreate(problem_name, &s, row_col),
145 "Section not created");
146 return SmartPetscObj<PetscSection>(s, false);
147}
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
Definition: definitions.h:595
MoFEMErrorCode sectionCreate(const std::string problem_name, PetscSection *s, const RowColData row_col=COL) const
Create global selection.
Definition: ISManager.cpp:20

◆ sectionCreate() [2/2]

MoFEMErrorCode MoFEM::ISManager::sectionCreate ( const std::string  problem_name,
PetscSection *  s,
const RowColData  row_col = COL 
) const

Create global selection.

Create section for given problem, such that points are sorted by UId, and number of dofs on point is equal to number of dofs on entity

It takes all fields

Parameters
problem_name
fields_list
s
row_colROE or COL, default is ROW
Returns
error code

Definition at line 20 of file ISManager.cpp.

22 {
23 const MoFEM::Interface &m_field = cOre;
24 const Problem *problem_ptr = m_field.get_problem(problem_name);
25 auto fields_ptr = m_field.get_fields();
26 auto fe_ptr = m_field.get_finite_elements();
28 boost::shared_ptr<NumeredDofEntity_multiIndex> dofs;
29 BitFieldId fields_ids(0);
30 switch (row_col) {
31 case ROW:
32 dofs = problem_ptr->numeredRowDofsPtr;
33 for (auto fit = fe_ptr->begin(); fit != fe_ptr->end(); fit++) {
34 if ((fit->get()->getId() & problem_ptr->getBitFEId()).any()) {
35 fields_ids |= fit->get()->getBitFieldIdRow();
36 }
37 }
38 break;
39 case COL:
40 dofs = problem_ptr->numeredColDofsPtr;
41 for (auto fit = fe_ptr->begin(); fit != fe_ptr->end(); fit++) {
42 if ((fit->get()->getId() & problem_ptr->getBitFEId()).any()) {
43 fields_ids |= fit->get()->getBitFieldIdCol();
44 }
45 }
46 break;
47 default:
48 SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
49 "Has to be ROW or COLUMN");
50 }
51 // get fields names on the problem
52 map<std::string, std::pair<int, int>> fields_map;
53 {
54 int field = 0;
55 for (auto fit = fields_ptr->begin(); fit != fields_ptr->end(); fit++) {
56 if ((fit->get()->getId() & fields_ids).any()) {
57 fields_map[fit->get()->getName()].first = field++;
58 fields_map[fit->get()->getName()].second = fit->get()->getNbOfCoeffs();
59 }
60 }
61 }
62 const int proc = m_field.get_comm_rank();
63 CHKERR PetscSectionCreate(PETSC_COMM_WORLD, s);
64 CHKERR PetscSectionSetNumFields(*s, fields_map.size());
65 for (auto mit = fields_map.begin(); mit != fields_map.end(); mit++) {
66 CHKERR PetscSectionSetFieldName(*s, mit->second.first, mit->first.c_str());
67 CHKERR PetscSectionSetFieldComponents(*s, mit->second.first,
68 mit->second.second);
69 }
70 // determine number of points
71 int nb_charts = 0;
72 {
73 auto dit = dofs->begin();
74 auto hi_dit = dofs->end();
75 for (; dit != hi_dit;) {
76 if (static_cast<int>(dit->get()->getPart()) == proc) {
77 const auto &ent_uid = dit->get()->getEntLocalUniqueId();
78 while (dit != hi_dit && dit->get()->getEntLocalUniqueId() == ent_uid) {
79 ++dit;
80 }
81 ++nb_charts;
82 } else {
83 ++dit;
84 }
85 }
86 }
87 // get layout, i.e. chart
88 PetscLayout layout;
89 CHKERR PetscLayoutCreate(PETSC_COMM_WORLD, &layout);
90 CHKERR PetscLayoutSetBlockSize(layout, 1);
91 CHKERR PetscLayoutSetLocalSize(layout, nb_charts);
92 CHKERR PetscLayoutSetUp(layout);
93 int rstart, rend;
94 CHKERR PetscLayoutGetRange(layout, &rstart, &rend);
95 CHKERR PetscLayoutDestroy(&layout);
96 CHKERR PetscSectionSetChart(*s, rstart, rend);
97
98 // loop of all dofs
99 {
100 auto dit = dofs->begin();
101 auto hi_dit = dofs->end();
102 int point = rstart;
103 for (; dit != hi_dit;) {
104 if (static_cast<int>(dit->get()->getPart()) == proc) {
105
106 const auto &field_name = dit->get()->getName();
107
108 int dd = 0;
109 const auto &ent_uid = dit->get()->getEntLocalUniqueId();
110 while (dit != hi_dit && dit->get()->getEntLocalUniqueId() == ent_uid) {
111 const DofIdx loc_idx = dit->get()->getPetscLocalDofIdx();
112 if (loc_idx >= 0)
113 ++dd;
114 ++dit;
115 }
116
117 if (fields_map.find(field_name) == fields_map.end()) {
118 MOFEM_LOG_C("SELF", Sev::warning, "Warning: Field %s not found",
119 dit->get()->getName().c_str());
120 } else {
121 CHKERR PetscSectionAddDof(*s, point, dd);
122 int field = fields_map.at(field_name).first;
123 CHKERR PetscSectionSetFieldDof(*s, point, field, dd);
124 }
125
126 ++point;
127
128 } else {
129 ++dit;
130 }
131 }
132 }
133 // cerr << "done " << proc << endl;
134 CHKERR PetscSectionSetUp(*s);
135 // cerr << "end " << proc << endl;
137}
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:311
virtual const FiniteElement_multiIndex * get_finite_elements() const =0
Get the finite elements object.
virtual const Field_multiIndex * get_fields() const =0
Get the fields object.
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
std::bitset< BITFIELDID_SIZE > BitFieldId
Field Id.
Definition: Types.hpp:42
constexpr auto field_name