 |
| v0.10.0
|
Simple interface for fast problem set-up.
More...
#include <src/interfaces/Simple.hpp>
|
MoFEMErrorCode | query_interface (const MOFEMuuid &uuid, UnknownInterface **iface) const |
|
| Simple (const MoFEM::Core &core) |
|
virtual | ~Simple ()=default |
|
MoFEMErrorCode | getOptions () |
| get options More...
|
|
MoFEMErrorCode | loadFile (const std::string options, const std::string mesh_file_name) |
| Load mesh file. More...
|
|
MoFEMErrorCode | loadFile (const std::string mesh_file_name="") |
| Load mesh file with parallel options if number of cores > 1. More...
|
|
MoFEMErrorCode | addDomainField (const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1) |
| Add field on domain. More...
|
|
MoFEMErrorCode | addBoundaryField (const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1) |
| Add field on boundary. More...
|
|
MoFEMErrorCode | addSkeletonField (const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1) |
| Add field on skeleton. More...
|
|
MoFEMErrorCode | addDataField (const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1) |
| Add field on domain. More...
|
|
MoFEMErrorCode | removeDomainField (const std::string &name) |
| Remove field form domain. More...
|
|
MoFEMErrorCode | removeBoundaryField (const std::string &name) |
| Remove field form boundary. More...
|
|
MoFEMErrorCode | removeSkeletonField (const std::string &name) |
| Remove field form skeleton. More...
|
|
MoFEMErrorCode | defineFiniteElements () |
| Define finite elements. More...
|
|
MoFEMErrorCode | defineProblem (const PetscBool is_partitioned=PETSC_TRUE) |
| define problem More...
|
|
MoFEMErrorCode | setFieldOrder (const std::string field_name, const int order, const Range *ents=NULL) |
| Set field order. More...
|
|
MoFEMErrorCode | buildFields () |
| Build fields. More...
|
|
MoFEMErrorCode | buildFiniteElements () |
| Build finite elements. More...
|
|
MoFEMErrorCode | setSkeletonAdjacency (int dim=-1) |
| Set the skeleton adjacency object. More...
|
|
MoFEMErrorCode | buildProblem () |
| Build problem. More...
|
|
MoFEMErrorCode | setUp (const PetscBool is_partitioned=PETSC_TRUE) |
| Setup problem. More...
|
|
MoFEMErrorCode | reSetUp () |
| Rebuild internal MoFEM data structures. More...
|
|
MoFEMErrorCode | getDM (DM *dm) |
| Get DM. More...
|
|
SmartPetscObj< DM > | getDM () |
| Return smart DM object. More...
|
|
int | getDim () const |
| Get the problem dimensio. More...
|
|
const std::string | getDomainFEName () const |
| Get the Domain FE Name. More...
|
|
const std::string | getBoundaryFEName () const |
| Get the Boundary FE Name. More...
|
|
const std::string | getSkeletonFEName () const |
| Get the Skeleton FE Name. More...
|
|
const std::string | getProblemName () const |
| Get the Problem Name. More...
|
|
std::string & | getDomainFEName () |
| Get the Domain FE Name. More...
|
|
std::string & | getBoundaryFEName () |
| Get the Boundary FE Name. More...
|
|
std::string & | getSkeletonFEName () |
| Get the Skeleton FE Name. More...
|
|
std::string & | getProblemName () |
| Get the Problem Name. More...
|
|
std::vector< std::string > & | getOtherFiniteElements () |
| Get the Other Finite Elements. More...
|
|
template<class IFACE > |
MoFEMErrorCode | registerInterface (const MOFEMuuid &uuid, bool error_if_registration_failed=true) |
| Register interface. More...
|
|
template<class IFACE , bool VERIFY = false> |
MoFEMErrorCode | getInterface (const MOFEMuuid &uuid, IFACE *&iface) const |
| Get interface by uuid and return reference to pointer of interface. More...
|
|
template<class IFACE > |
MoFEMErrorCode | getInterface (IFACE *&iface) const |
| Get interface refernce to pointer of interface. More...
|
|
template<class IFACE > |
MoFEMErrorCode | getInterface (IFACE **const iface) const |
| Get interface pointer to pointer of interface. More...
|
|
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0> |
IFACE | getInterface () const |
| Get interface pointer to pointer of interface. More...
|
|
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0> |
IFACE | getInterface () const |
| Get reference to interface. More...
|
|
template<class IFACE > |
IFACE * | getInterface () const |
| Function returning pointer to interface. More...
|
|
virtual | ~UnknownInterface ()=default |
|
virtual MoFEMErrorCode | getLibVersion (Version &version) const |
| Get library version. More...
|
|
virtual const MoFEMErrorCode | getFileVersion (moab::Interface &moab, Version &version) const |
| Get database major version. More...
|
|
virtual MoFEMErrorCode | getInterfaceVersion (Version &version) const |
| Get database major version. More...
|
|
template<> |
MoFEMErrorCode | getInterface (const MOFEMuuid &uuid, UnknownInterface *&iface) const |
|
Simple interface for fast problem set-up.
- Examples
- boundary_marker.cpp, contact.cpp, continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp, dynamic_elastic.cpp, eigen_elastic.cpp, field_blas_set_vertex_dofs.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, heat_equation.cpp, hello_world.cpp, helmholtz.cpp, loop_entities.cpp, mesh_smoothing.cpp, reaction_diffusion.cpp, scalar_check_approximation_2d.cpp, simple_elasticity.cpp, simple_interface.cpp, test_cache_on_entities.cpp, testing_jacobian_of_hook_element.cpp, testing_jacobian_of_hook_scaled_with_density_element.cpp, and wave_equation.cpp.
Definition at line 36 of file Simple.hpp.
◆ Simple()
Definition at line 185 of file Simple.cpp.
191 PetscLogEventRegister(
"buildFiniteElements", 0,
◆ ~Simple()
virtual MoFEM::Simple::~Simple |
( |
| ) |
|
|
virtualdefault |
◆ addBoundaryField()
Add field on boundary.
- Parameters
-
name | name of the filed |
space | space (L2,H1,Hdiv,Hcurl) |
base | approximation base, see FieldApproximationBase |
nb_of_coefficients | number of field coefficients |
tag_type | type of the tag MB_TAG_DENSE or MB_TAG_SPARSE (DENSE is faster and uses less memory, SPARSE is more flexible if you define field on subdomains) |
bh | if MF_EXCL throws error if field exits, MF_ZERO no error if field exist |
verb | verbosity level |
- Returns
- error code
- Examples
- hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, heat_equation.cpp, mesh_smoothing.cpp, simple_interface.cpp, and wave_equation.cpp.
Definition at line 287 of file Simple.cpp.
294 CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
300 "NOFIELD space for boundary filed not implemented in Simple interface");
◆ addDataField()
Add field on domain.
- Parameters
-
name | name of the filed |
space | space (L2,H1,Hdiv,Hcurl) |
base | approximation base, see FieldApproximationBase |
nb_of_coefficients | number of field coefficients |
tag_type | type of the tag MB_TAG_DENSE or MB_TAG_SPARSE (DENSE is faster and uses less memory, SPARSE is more flexible if you define field on subdomains) |
bh | if MF_EXCL throws error if field exits, MF_ZERO no error if field exist |
verb | verbosity level |
- Returns
- error code
- Examples
- hcurl_curl_operator.cpp, and hdiv_divergence_operator.cpp.
Definition at line 324 of file Simple.cpp.
332 CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
◆ addDomainField()
Add field on domain.
- Parameters
-
name | name of the filed |
space | space (L2,H1,Hdiv,Hcurl) |
base | approximation base, see FieldApproximationBase |
nb_of_coefficients | number of field coefficients |
tag_type | type of the tag MB_TAG_DENSE or MB_TAG_SPARSE (DENSE is faster and uses less memory, SPARSE is more flexible if you define field on subdomains) |
bh | if MF_EXCL throws error if field exits, MF_ZERO no error if field exist |
verb | verbosity level |
- Returns
- error code
- Examples
- continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp, field_blas_set_vertex_dofs.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, heat_equation.cpp, loop_entities.cpp, mesh_smoothing.cpp, reaction_diffusion.cpp, scalar_check_approximation_2d.cpp, simple_elasticity.cpp, simple_interface.cpp, test_cache_on_entities.cpp, testing_jacobian_of_hook_element.cpp, testing_jacobian_of_hook_scaled_with_density_element.cpp, and wave_equation.cpp.
Definition at line 269 of file Simple.cpp.
277 CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
◆ addSkeletonField()
Add field on skeleton.
- Parameters
-
name | name of the filed |
space | space (L2,H1,Hdiv,Hcurl) |
base | approximation base, see FieldApproximationBase |
nb_of_coefficients | number of field coefficients |
tag_type | type of the tag MB_TAG_DENSE or MB_TAG_SPARSE (DENSE is faster and uses less memory, SPARSE is more flexible if you define field on subdomains) |
bh | if MF_EXCL throws error if field exits, MF_ZERO no error if field exist |
verb | verbosity level |
- Returns
- error code
- Examples
- continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, and continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp.
Definition at line 305 of file Simple.cpp.
313 CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
319 "NOFIELD space for boundary filed not implemented in Simple interface");
◆ buildFields()
Build fields.
- Returns
- error code
- Examples
- field_blas_set_vertex_dofs.cpp, mesh_smoothing.cpp, and simple_elasticity.cpp.
Definition at line 490 of file Simple.cpp.
495 auto get_skin = [&](
auto meshset) {
497 CHKERR m_field.get_moab().get_entities_by_dimension(meshset,
dIm,
499 Skinner skin(&m_field.get_moab());
501 CHKERR skin.find_skin(0, domain_ents,
false, domain_skin);
503 ParallelComm *pcomm =
505 Range proc_domain_skin;
506 CHKERR pcomm->filter_pstatus(domain_skin,
507 PSTATUS_SHARED | PSTATUS_MULTISHARED,
508 PSTATUS_NOT, -1, &proc_domain_skin);
509 return proc_domain_skin;
512 auto create_boundary_meshset = [&](
auto &&proc_domain_skin) {
522 CHKERR m_field.get_moab().get_adjacencies(proc_domain_skin,
dd,
false,
523 adj, moab::Interface::UNION);
531 auto comm_interface_ptr = m_field.getInterface<CommInterface>();
532 auto bit_ref_ptr = m_field.getInterface<BitRefManager>();
535 auto add_ents_to_field = [&](
auto meshset,
auto dim,
auto &fields) {
537 for (
auto &field : fields) {
538 CHKERR m_field.add_ents_to_field_by_dim(meshset,
dim, field);
539 CHKERR comm_interface_ptr->synchroniseFieldEntities(field, 0);
544 auto make_no_field_ents = [&](
auto &fields) {
546 for (
auto &field : fields) {
547 std::array<double, 6> coords = {0, 0, 0, 0, 0, 0};
548 CHKERR m_field.create_vertices_and_add_to_field(field, coords.data(), 2);
549 CHKERR comm_interface_ptr->makeFieldEntitiesMultishared(field, 0);
563 auto set_order = [&](
auto meshset,
auto dim,
auto &fields) {
565 for (
auto &field : fields) {
568 "Order for field not set %s", field.c_str());
571 const Field *field_ptr = m_field.get_field_structure(field);
572 switch (field_ptr->getSpace()) {
587 "Glasgow we have a problem");
590 auto set_order = [&](
auto field,
auto &ents) {
594 for (
auto o = range.first;
o != range.second; ++
o) {
595 if (!
o->second.second.empty())
596 ents = intersect(ents,
o->second.second);
597 CHKERR m_field.set_field_order(ents, field,
o->second.first);
603 if (field_ptr->getSpace() ==
H1) {
606 CHKERR m_field.get_field_entities_by_dimension(field, 0, ents);
607 CHKERR set_order(field, ents);
614 CHKERR m_field.get_field_entities_by_dimension(field,
dd, ents);
615 CHKERR set_order(field, ents);
627 CHKERR m_field.build_fields();
◆ buildFiniteElements()
◆ buildProblem()
Build problem.
- Returns
- error code
- Examples
- mesh_smoothing.cpp, and simple_elasticity.cpp.
Definition at line 658 of file Simple.cpp.
665 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
667 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
◆ defineFiniteElements()
Define finite elements.
- Returns
- Error code
- Examples
- mesh_smoothing.cpp, and simple_elasticity.cpp.
Definition at line 387 of file Simple.cpp.
391 auto clear_rows_and_cols = [&](
auto &fe_name) {
393 auto fe_ptr = m_field.get_finite_elements();
395 ->get<FiniteElement_name_mi_tag>();
396 auto it_fe = fe_by_name.find(fe_name);
397 if (it_fe != fe_by_name.end()) {
399 if(!fe_by_name.modify(it_fe, FiniteElement_row_change_bit_reset()))
401 "modification unsuccessful");
403 if(!fe_by_name.modify(it_fe, FiniteElement_col_change_bit_reset()))
405 "modification unsuccessful");
416 auto add_fields = [&](
auto &fe_name,
auto &fields) {
418 for (
auto &field : fields) {
419 CHKERR m_field.modify_finite_element_add_field_row(fe_name, field);
420 CHKERR m_field.modify_finite_element_add_field_col(fe_name, field);
421 CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
426 auto add_data_fields = [&](
auto &fe_name,
auto &fields) {
428 for (
auto &field : fields)
429 CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
◆ defineProblem()
MoFEMErrorCode MoFEM::Simple::defineProblem |
( |
const PetscBool |
is_partitioned = PETSC_TRUE | ) |
|
◆ getBoundaryFEName() [1/2]
std::string& MoFEM::Simple::getBoundaryFEName |
( |
| ) |
|
Get the Boundary FE Name.
- Returns
- std::string&
Definition at line 316 of file Simple.hpp.
◆ getBoundaryFEName() [2/2]
const std::string MoFEM::Simple::getBoundaryFEName |
( |
| ) |
const |
◆ getDim()
int MoFEM::Simple::getDim |
( |
| ) |
const |
Get the problem dimensio.
Problem dimension is determined by highest dimension of entities on the mesh.
- Returns
- int
Definition at line 274 of file Simple.hpp.
◆ getDM() [1/2]
◆ getDM() [2/2]
Get DM.
- Parameters
-
- Returns
- error code
- Examples
- continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, heat_equation.cpp, loop_entities.cpp, mesh_smoothing.cpp, reaction_diffusion.cpp, scalar_check_approximation_2d.cpp, simple_elasticity.cpp, simple_interface.cpp, test_cache_on_entities.cpp, testing_jacobian_of_hook_element.cpp, testing_jacobian_of_hook_scaled_with_density_element.cpp, and wave_equation.cpp.
Definition at line 706 of file Simple.cpp.
◆ getDomainFEName() [1/2]
std::string& MoFEM::Simple::getDomainFEName |
( |
| ) |
|
Get the Domain FE Name.
- Returns
- std::string&
Definition at line 309 of file Simple.hpp.
◆ getDomainFEName() [2/2]
const std::string MoFEM::Simple::getDomainFEName |
( |
| ) |
const |
◆ getOptions()
get options
- Returns
- error code
- Examples
- continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp, field_blas_set_vertex_dofs.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, heat_equation.cpp, loop_entities.cpp, mesh_smoothing.cpp, reaction_diffusion.cpp, scalar_check_approximation_2d.cpp, simple_elasticity.cpp, simple_interface.cpp, test_cache_on_entities.cpp, testing_jacobian_of_hook_element.cpp, testing_jacobian_of_hook_scaled_with_density_element.cpp, and wave_equation.cpp.
Definition at line 198 of file Simple.cpp.
199 PetscBool flg = PETSC_TRUE;
201 ierr = PetscOptionsBegin(PETSC_COMM_WORLD,
"",
"Simple interface options",
204 ierr = PetscOptionsString(
"-file_name",
"file name",
"",
"mesh.h5m",
207 ierr = PetscOptionsEnd();
◆ getOtherFiniteElements()
std::vector<std::string>& MoFEM::Simple::getOtherFiniteElements |
( |
| ) |
|
Get the Other Finite Elements.
User can create finite elements using directly core interface and and them to simple problem by this function
- Returns
- std::vector<std::string>&
- Examples
- mesh_smoothing.cpp.
Definition at line 340 of file Simple.hpp.
◆ getProblemName() [1/2]
std::string& MoFEM::Simple::getProblemName |
( |
| ) |
|
◆ getProblemName() [2/2]
const std::string MoFEM::Simple::getProblemName |
( |
| ) |
const |
◆ getSkeletonFEName() [1/2]
std::string& MoFEM::Simple::getSkeletonFEName |
( |
| ) |
|
Get the Skeleton FE Name.
- Returns
- std::string&
Definition at line 323 of file Simple.hpp.
◆ getSkeletonFEName() [2/2]
const std::string MoFEM::Simple::getSkeletonFEName |
( |
| ) |
const |
◆ loadFile() [1/2]
MoFEMErrorCode MoFEM::Simple::loadFile |
( |
const std::string |
mesh_file_name = "" | ) |
|
Load mesh file with parallel options if number of cores > 1.
- Parameters
-
mesh_file_name | file name if not set default or set by command line is used. |
- Returns
- error code
Definition at line 255 of file Simple.cpp.
258 if (m_field.get_comm_size() == 1)
262 "PARALLEL_RESOLVE_SHARED_ENTS;"
263 "PARTITION=PARALLEL_PARTITION;",
◆ loadFile() [2/2]
MoFEMErrorCode MoFEM::Simple::loadFile |
( |
const std::string |
options, |
|
|
const std::string |
mesh_file_name |
|
) |
| |
Load mesh file.
- Parameters
-
options | file load options |
mesh_file_name | file name if not set default or set by command line is used. |
- Returns
- error code
- Examples
- continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp, field_blas_set_vertex_dofs.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, heat_equation.cpp, loop_entities.cpp, mesh_smoothing.cpp, reaction_diffusion.cpp, scalar_check_approximation_2d.cpp, simple_elasticity.cpp, simple_interface.cpp, test_cache_on_entities.cpp, testing_jacobian_of_hook_element.cpp, testing_jacobian_of_hook_scaled_with_density_element.cpp, and wave_equation.cpp.
Definition at line 212 of file Simple.cpp.
224 CHKERR m_field.rebuild_database();
228 CHKERR m_field.get_moab().get_number_entities_by_dimension(
230 if (nb_ents_3d > 0) {
234 CHKERR m_field.get_moab().get_number_entities_by_dimension(
236 if (nb_ents_2d > 0) {
244 CHKERR m_field.get_moab().get_entities_by_dimension(
meshSet,
dIm, ents,
true);
245 CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevel(ents,
bitLevel,
247 ParallelComm *pcomm =
250 pcomm =
new ParallelComm(&m_field.get_moab(), m_field.get_comm());
◆ query_interface()
◆ removeBoundaryField()
MoFEMErrorCode MoFEM::Simple::removeBoundaryField |
( |
const std::string & |
name | ) |
|
Remove field form boundary.
- Parameters
-
- Returns
- MoFEMErrorCode
Definition at line 357 of file Simple.cpp.
361 auto remove_field_from_list = [&](
auto &vec) {
362 auto it = std::find(vec.begin(), vec.end(), name);
◆ removeDomainField()
MoFEMErrorCode MoFEM::Simple::removeDomainField |
( |
const std::string & |
name | ) |
|
Remove field form domain.
- Parameters
-
- Returns
- MoFEMErrorCode
Definition at line 341 of file Simple.cpp.
345 auto remove_field_from_list = [&](
auto &vec) {
346 auto it = std::find(vec.begin(), vec.end(), name);
◆ removeSkeletonField()
MoFEMErrorCode MoFEM::Simple::removeSkeletonField |
( |
const std::string & |
name | ) |
|
Remove field form skeleton.
- Parameters
-
- Returns
- MoFEMErrorCode
Definition at line 372 of file Simple.cpp.
376 auto remove_field_from_list = [&](
auto &vec) {
377 auto it = std::find(vec.begin(), vec.end(), name);
◆ reSetUp()
Rebuild internal MoFEM data structures.
Call this function after you add field or remove it.
- Note
- If you add field, or remove it, finite element and problem needs to be rebuild. However DM can remain the same.
- Returns
- MoFEMErrorCode
Definition at line 684 of file Simple.cpp.
698 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
700 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
◆ setFieldOrder()
MoFEMErrorCode MoFEM::Simple::setFieldOrder |
( |
const std::string |
field_name, |
|
|
const int |
order, |
|
|
const Range * |
ents = NULL |
|
) |
| |
Set field order.
- Parameters
-
std::field_name | field name |
order | order |
range | of entities to which order is set (If null it sat to all entities) |
- Returns
- error code
- Examples
- continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp, field_blas_set_vertex_dofs.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, heat_equation.cpp, loop_entities.cpp, mesh_smoothing.cpp, reaction_diffusion.cpp, scalar_check_approximation_2d.cpp, simple_elasticity.cpp, simple_interface.cpp, test_cache_on_entities.cpp, testing_jacobian_of_hook_element.cpp, testing_jacobian_of_hook_scaled_with_density_element.cpp, and wave_equation.cpp.
Definition at line 482 of file Simple.cpp.
486 {field_name, {
order, ents == NULL ? Range() : Range(*ents)}});
◆ setSkeletonAdjacency() [1/5]
Definition at line 41 of file Simple.cpp.
45 auto defaultSkeletonEdge = [&](
moab::Interface &moab,
const Field &field,
46 const EntFiniteElement &fe,
56 Range bride_adjacency_edge, bride_adjacency;
57 CHKERR moab.get_adjacencies(&fe_ent, 1, 2,
false, bride_adjacency_edge);
59 switch (field.getSpace()) {
61 CHKERR moab.get_connectivity(bride_adjacency_edge, bride_adjacency,
65 CHKERR moab.get_adjacencies(bride_adjacency_edge, 1,
false,
66 bride_adjacency, moab::Interface::UNION);
68 bride_adjacency.merge(bride_adjacency_edge);
74 "this field is not implemented for TRI finite element");
77 bride_adjacency = subtract(bride_adjacency, adjacency);
79 for (
auto e : bride_adjacency)
81 .insert(boost::shared_ptr<SideNumber>(
new SideNumber(e, -1, 0, 0)));
83 adjacency.merge(bride_adjacency);
◆ setSkeletonAdjacency() [2/5]
Definition at line 96 of file Simple.cpp.
100 auto defaultSkeletonEdge = [&](
moab::Interface &moab,
const Field &field,
101 const EntFiniteElement &fe,
111 Range bride_adjacency_face, bride_adjacency;
112 CHKERR moab.get_adjacencies(&fe_ent, 1, 3,
false, bride_adjacency_face);
114 switch (field.getSpace()) {
116 CHKERR moab.get_connectivity(bride_adjacency_face, bride_adjacency,
119 CHKERR moab.get_adjacencies(bride_adjacency_face, 1,
false,
120 bride_adjacency, moab::Interface::UNION);
122 CHKERR moab.get_adjacencies(bride_adjacency_face, 2,
false,
123 bride_adjacency, moab::Interface::UNION);
125 bride_adjacency.merge(bride_adjacency_face);
131 "this field is not implemented for TRI finite element");
134 bride_adjacency = subtract(bride_adjacency, adjacency);
136 for (
auto e : bride_adjacency)
138 .insert(boost::shared_ptr<SideNumber>(
new SideNumber(e, -1, 0, 0)));
140 adjacency.merge(bride_adjacency);
147 defaultSkeletonEdge);
149 defaultSkeletonEdge);
◆ setSkeletonAdjacency() [3/5]
Definition at line 155 of file Simple.cpp.
161 return setSkeletonAdjacency<2>();
163 return setSkeletonAdjacency<3>();
◆ setSkeletonAdjacency() [4/5]
Definition at line 35 of file Simple.cpp.
36 static_assert(DIM == 2 || DIM == 3,
"not implemented");
◆ setSkeletonAdjacency() [5/5]
Set the skeleton adjacency object.
- Parameters
-
- Returns
- MoFEMErrorCode
Definition at line 170 of file Simple.cpp.
176 return setSkeletonAdjacency<2>();
178 return setSkeletonAdjacency<3>();
◆ setUp()
MoFEMErrorCode MoFEM::Simple::setUp |
( |
const PetscBool |
is_partitioned = PETSC_TRUE | ) |
|
Setup problem.
- Returns
- error code
- Examples
- continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, heat_equation.cpp, loop_entities.cpp, reaction_diffusion.cpp, scalar_check_approximation_2d.cpp, simple_interface.cpp, test_cache_on_entities.cpp, testing_jacobian_of_hook_element.cpp, testing_jacobian_of_hook_scaled_with_density_element.cpp, and wave_equation.cpp.
Definition at line 672 of file Simple.cpp.
◆ bitLevel
BitRefLevel of the probelm.
Definition at line 345 of file Simple.hpp.
◆ boundaryFE
std::string MoFEM::Simple::boundaryFE |
|
private |
boundary finite element
Definition at line 366 of file Simple.hpp.
◆ boundaryFields
std::vector<std::string> MoFEM::Simple::boundaryFields |
|
private |
◆ boundaryMeshset
◆ cOre
◆ dataFields
std::vector<std::string> MoFEM::Simple::dataFields |
|
private |
◆ dIm
◆ dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition at line 375 of file Simple.hpp.
◆ domainFE
std::string MoFEM::Simple::domainFE |
|
private |
◆ domainFields
std::vector<std::string> MoFEM::Simple::domainFields |
|
private |
◆ fieldsOrder
std::multimap<std::string, std::pair<int, Range> > MoFEM::Simple::fieldsOrder |
|
private |
◆ meshFileName
char MoFEM::Simple::meshFileName[255] |
|
private |
◆ meshSet
◆ MOFEM_EVENT_SimpleBuildFields
PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFields |
|
private |
◆ MOFEM_EVENT_SimpleBuildFiniteElements
PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFiniteElements |
|
private |
◆ MOFEM_EVENT_SimpleBuildProblem
PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildProblem |
|
private |
◆ MOFEM_EVENT_SimpleKSPSolve
PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleKSPSolve |
|
private |
◆ MOFEM_EVENT_SimpleLoadMesh
PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleLoadMesh |
|
private |
◆ nameOfProblem
std::string MoFEM::Simple::nameOfProblem |
|
private |
◆ noFieldDataFields
std::vector<std::string> MoFEM::Simple::noFieldDataFields |
|
private |
◆ noFieldFields
std::vector<std::string> MoFEM::Simple::noFieldFields |
|
private |
◆ otherFEs
std::vector<std::string> MoFEM::Simple::otherFEs |
|
private |
Other finite elements.
Definition at line 369 of file Simple.hpp.
◆ skeletonFE
std::string MoFEM::Simple::skeletonFE |
|
private |
skeleton finite element
Definition at line 367 of file Simple.hpp.
◆ skeletonFields
std::vector<std::string> MoFEM::Simple::skeletonFields |
|
private |
fields on the skeleton
Definition at line 357 of file Simple.hpp.
The documentation for this struct was generated from the following files:
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MYPCOMM_INDEX
default communicator number PCOMM
const BitRefLevel bitLevel
BitRefLevel of the probelm.
Simple(const MoFEM::Core &core)
MoFEMErrorCode defineFiniteElements()
Define finite elements.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name)
Load mesh file.
char meshFileName[255]
mesh file name
@ L2
field with C-1 continuity
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
std::vector< std::string > noFieldFields
NOFIELD field name.
static MoFEMErrorCode defaultEdge(Interface &moab, const Field &field, const EntFiniteElement &fe, Range &adjacency)
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
std::vector< std::string > boundaryFields
boundary fields
PetscErrorCode DMSetUp_MoFEM(DM dm)
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
std::vector< std::string > skeletonFields
fields on the skeleton
int dIm
dimension of problem
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
DeprecatedCoreInterface Interface
static const MOFEMuuid IDD_MOFEMSimple
#define CHKERR
Inline error check.
std::string nameOfProblem
problem name
MoFEMErrorCode buildProblem()
Build problem.
auto createSmartDM
Creates smart DM object.
MoFEMErrorCode buildFields()
Build fields.
std::vector< std::string > noFieldDataFields
NOFIELD field name.
EntityHandle meshSet
domain meshset
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
std::vector< std::string > domainFields
domain fields
int getDim() const
Get the problem dimensio.
MoFEMErrorCode defineProblem(const PetscBool is_partitioned=PETSC_TRUE)
define problem
MoFEMErrorCode setSkeletonAdjacency()
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
multi_index_container< boost::shared_ptr< SideNumber >, indexed_by< ordered_unique< member< SideNumber, EntityHandle, &SideNumber::ent > >, ordered_non_unique< composite_key< SideNumber, const_mem_fun< SideNumber, EntityType, &SideNumber::getEntType >, member< SideNumber, char, &SideNumber::side_number > > > > > SideNumber_multiIndex
SideNumber_multiIndex for SideNumber.
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
@ MOFEM_OPERATION_UNSUCCESSFUL
#define THROW_MESSAGE(a)
Throw MoFEM exception.
PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr, const char problem_name[], const MoFEM::BitRefLevel bit_level, const MoFEM::BitRefLevel bit_mask=MoFEM::BitRefLevel().set())
Must be called by user to set MoFEM data structures.
std::vector< std::string > dataFields
Data fields.
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
PetscObject getPetscObject(T obj)
std::string domainFE
domain finite element
MoFEMErrorCode buildFiniteElements()
Build finite elements.
EntityHandle boundaryMeshset
meshset with boundary
static MoFEMErrorCode defaultFace(Interface &moab, const Field &field, const EntFiniteElement &fe, Range &adjacency)
@ HCURL
field with continuous tangents
@ MOFEM_DATA_INCONSISTENCY
std::vector< std::string > otherFEs
Other finite elements.
multi_index_container< boost::shared_ptr< FiniteElement >, indexed_by< hashed_unique< tag< FiniteElement_Meshset_mi_tag >, member< FiniteElement, EntityHandle, &FiniteElement::meshset > >, hashed_unique< tag< BitFEId_mi_tag >, const_mem_fun< FiniteElement, BitFEId, &FiniteElement::getId >, HashBit< BitFEId >, EqBit< BitFEId > >, ordered_unique< tag< FiniteElement_name_mi_tag >, const_mem_fun< FiniteElement, boost::string_ref, &FiniteElement::getNameRef > > > > FiniteElement_multiIndex
MultiIndex for entities for FiniteElement.
std::string boundaryFE
boundary finite element
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
std::string skeletonFE
skeleton finite element
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
@ HDIV
field with continuous normal traction
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
std::multimap< std::string, std::pair< int, Range > > fieldsOrder
fields order
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)
@ NOFIELD
scalar or vector of scalars describe (no true field)
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)