v0.9.1
Classes | Functions | Variables
hcurl_curl_operator.cpp File Reference

Testich curl-curl operator by applying Stokes theorem. More...

#include <MoFEM.hpp>

Go to the source code of this file.

Classes

struct  OpTetCurl
 
struct  OpFacesRot
 

Functions

int main (int argc, char *argv[])
 

Variables

static char help [] = "...\n\n"
 

Detailed Description

Testich curl-curl operator by applying Stokes theorem.

Definition in file hcurl_curl_operator.cpp.

Function Documentation

◆ main()

int main ( int  argc,
char *  argv[] 
)
Examples
hcurl_curl_operator.cpp.

Definition at line 57 of file hcurl_curl_operator.cpp.

57  {
58 
59  MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
60 
61  try {
62 
63  enum bases { AINSWORTH, DEMKOWICZ, LASTOP };
64 
65  const char *list[] = {"ainsworth", "demkowicz"};
66 
67  PetscBool flg;
68  PetscInt choise_value = AINSWORTH;
69  CHKERR PetscOptionsGetEList(PETSC_NULL, NULL, "-base", list, LASTOP,
70  &choise_value, &flg);
71  if (flg != PETSC_TRUE) {
72  SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSIBLE_CASE, "base not set");
73  }
74 
75  PetscBool ho_geometry = PETSC_FALSE;
76  CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-ho_geometry", &ho_geometry,
77  PETSC_NULL);
78 
79  DMType dm_name = "DMMOFEM";
80  CHKERR DMRegister_MoFEM(dm_name);
81 
82  // Create MoAB
83  moab::Core mb_instance; ///< database
84  moab::Interface &moab = mb_instance; ///< interface
85 
86  // Create MoFEM
87  MoFEM::Core core(moab);
88  MoFEM::Interface &m_field = core;
89 
90  Simple *simple_interface = m_field.getInterface<Simple>();
91  Basic *basic_interface = m_field.getInterface<Basic>();
92  CHKERR simple_interface->getOptions();
93  CHKERR simple_interface->loadFile("");
94 
95  // fields
96  switch (choise_value) {
97  case AINSWORTH:
98  CHKERR simple_interface->addDomainField("HCURL", HCURL,
100  CHKERR simple_interface->addBoundaryField("HCURL", HCURL,
102  break;
103  case DEMKOWICZ:
104  CHKERR simple_interface->addDomainField("HCURL", HCURL,
106  CHKERR simple_interface->addBoundaryField("HCURL", HCURL,
108  break;
109  }
110 
111  if (ho_geometry == PETSC_TRUE)
112  CHKERR simple_interface->addDataField("MESH_NODE_POSITIONS", H1,
114 
115  constexpr int order = 5;
116  CHKERR simple_interface->setFieldOrder("HCURL", order);
117  if (ho_geometry == PETSC_TRUE)
118  CHKERR simple_interface->setFieldOrder("MESH_NODE_POSITIONS", 2);
119  CHKERR simple_interface->setUp();
120 
121  auto integration_rule = [](int, int, int p_data) { return 2 * p_data; };
122  CHKERR basic_interface->setDomainRhsIntegrationRule(integration_rule);
123  CHKERR basic_interface->setBoundaryRhsIntegrationRule(integration_rule);
124 
125  FTensor::Tensor1<double, 3> t_curl_vol;
126  FTensor::Tensor1<double, 3> t_curl_skin;
127 
128  basic_interface->getOpDomainRhsPipeline().push_back(
129  new OpTetCurl(t_curl_vol));
130  basic_interface->getOpBoundaryRhsPipeline().push_back(
131  new OpFacesRot(t_curl_skin));
132 
134  t_curl_vol(i) = 0;
135  t_curl_skin(i) = 0;
136 
137  // project geometry form 10 node tets on higher order approx. functions
138  if (ho_geometry == PETSC_TRUE) {
139  Projection10NodeCoordsOnField ent_method(m_field, "MESH_NODE_POSITIONS");
140  CHKERR m_field.loop_dofs("MESH_NODE_POSITIONS", ent_method);
141  }
142 
143  // Run pipelines on mesh
144  CHKERR basic_interface->loopFiniteElements();
145 
146  std::cout.precision(12);
147 
148  std::cout << "curl_vol " << t_curl_vol << std::endl;
149  std::cout << "curl_skin " << t_curl_skin << std::endl;
150 
151  t_curl_vol(i) -= t_curl_skin(i);
152  double nrm2 = sqrt(t_curl_vol(i) * t_curl_vol(i));
153 
154  constexpr double eps = 1e-8;
155  if (fabs(nrm2) > eps)
156  SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID,
157  "Curl operator not passed test\n");
158  }
159  CATCH_ERRORS;
160 
162 }
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition: Simple.cpp:660
Deprecated interface functions.
Core (interface) class.
Definition: Core.hpp:50
Simple interface for fast problem set-up.
Definition: Simple.hpp:36
Projection of edge entities with one mid-node on hierarchical basis.
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:51
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.
Definition: Simple.cpp:311
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
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.
Definition: Simple.cpp:274
Basic interface.
Definition: Basic.hpp:36
MoFEMErrorCode loadFile(const std::string options="PARALLEL=READ_PART;" "PARALLEL_RESOLVE_SHARED_ENTS;" "PARTITION=PARALLEL_PARTITION;", const std::string mesh_file_name="")
Load mesh file.
Definition: Simple.cpp:212
static char help[]
PetscErrorCode PetscOptionsGetEList(PetscOptions *, const char pre[], const char name[], const char *const *list, PetscInt next, PetscInt *value, PetscBool *set)
constexpr int order
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:48
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:151
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.
Definition: Simple.cpp:256
virtual MoFEMErrorCode loop_dofs(const Problem *problem_ptr, const std::string &field_name, RowColData rc, DofMethod &method, int lower_rank, int upper_rank, int verb=DEFAULT_VERBOSITY)=0
Make a loop over dofs.
field with continuous tangents
Definition: definitions.h:177
#define CHKERR
Inline error check.
Definition: definitions.h:601
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
MoFEMErrorCode getOptions()
get options
Definition: Simple.cpp:198
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:438
continuous field
Definition: definitions.h:176
static const double eps
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition: Simple.cpp:470
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:61
FTensor::Index< 'i', 2 > i
[Common data]
Definition: ContactOps.hpp:26

Variable Documentation

◆ help

char help[] = "...\n\n"
static
Examples
hcurl_curl_operator.cpp.

Definition at line 31 of file hcurl_curl_operator.cpp.