|
| v0.14.0
|
Go to the documentation of this file.
17 using namespace MoFEM;
19 static char help[] =
"...\n\n";
45 int main(
int argc,
char *argv[]) {
51 enum bases { AINSWORTH, DEMKOWICZ,
LASTOP };
53 const char *list[] = {
"ainsworth",
"demkowicz"};
56 PetscInt choise_value = AINSWORTH;
59 if (flg != PETSC_TRUE) {
63 PetscBool ho_geometry = PETSC_FALSE;
67 DMType dm_name =
"DMMOFEM";
84 switch (choise_value) {
99 if (ho_geometry == PETSC_TRUE)
103 constexpr
int order = 3;
109 if (ho_geometry == PETSC_TRUE)
120 auto material_grad_mat = boost::make_shared<MatrixDouble>();
121 auto material_det_vec = boost::make_shared<VectorDouble>();
122 auto material_inv_grad_mat = boost::make_shared<MatrixDouble>();
123 auto jac_ptr = boost::make_shared<MatrixDouble>();
124 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
125 auto det_ptr = boost::make_shared<VectorDouble>();
127 pipeline_mng->getOpDomainRhsPipeline().push_back(
129 pipeline_mng->getOpDomainRhsPipeline().push_back(
131 pipeline_mng->getOpDomainRhsPipeline().push_back(
133 pipeline_mng->getOpDomainRhsPipeline().push_back(
135 pipeline_mng->getOpDomainRhsPipeline().push_back(
139 pipeline_mng->getOpDomainRhsPipeline().push_back(
143 material_grad_mat, material_det_vec, material_inv_grad_mat));
144 pipeline_mng->getOpDomainRhsPipeline().push_back(
146 pipeline_mng->getOpDomainRhsPipeline().push_back(
148 pipeline_mng->getOpDomainRhsPipeline().push_back(
151 pipeline_mng->getOpDomainRhsPipeline().push_back(
new OpVolCurl(t_curl_vol));
154 pipeline_mng->getOpBoundaryRhsPipeline().push_back(
156 pipeline_mng->getOpBoundaryRhsPipeline().push_back(
158 pipeline_mng->getOpBoundaryRhsPipeline().push_back(
166 if (ho_geometry == PETSC_TRUE) {
172 CHKERR pipeline_mng->loopFiniteElements();
174 std::cout.precision(12);
176 t_curl_vol(
i) -= t_curl_skin(
i);
177 double nrm2 = sqrt(t_curl_vol(
i) * t_curl_vol(
i));
179 constexpr
double eps = 1e-8;
180 if (fabs(nrm2) >
eps)
182 "Curl operator not passed test\n");
196 const unsigned int nb_gauss_pts = data.
getDiffN().size1();
197 const unsigned int nb_dofs = data.
getFieldData().size();
207 for (; gg < nb_gauss_pts; gg++) {
208 double w = getGaussPts()(3, gg) * getVolume();
210 for (
unsigned int dd = 0;
dd != nb_dofs;
dd++) {
212 tCurl(
i) +=
w * t_curl(
i);
227 int nb_gauss_pts = data.
getN().size1();
234 for (
int gg = 0; gg < nb_gauss_pts; gg++) {
235 for (
int dd = 0;
dd < nb_dofs;
dd++) {
236 double w = getGaussPts()(2, gg);
237 const double n0 = getNormalsAtGaussPts(gg)[0];
238 const double n1 = getNormalsAtGaussPts(gg)[1];
239 const double n2 = getNormalsAtGaussPts(gg)[2];
240 if (getFEType() == MBTRI) {
244 tCurl(0) += (n1 * t_curl_base(2) - n2 * t_curl_base(1)) *
w;
245 tCurl(1) += (n2 * t_curl_base(0) - n0 * t_curl_base(2)) *
w;
246 tCurl(2) += (n0 * t_curl_base(1) - n1 * t_curl_base(0)) *
w;
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Data on single entity (This is passed as argument to DataOperator::doWork)
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.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
UBlasMatrix< double > MatrixDouble
const VectorDouble & getFieldData() const
get dofs values
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name, LoadFileFunc loadFunc=defaultLoadFileFunc)
Load mesh file.
PipelineManager interface.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Projection of edge entities with one mid-node on hierarchical basis.
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
MatrixDouble & getDiffN(const FieldApproximationBase base)
get derivatives of base functions
FTensor::Tensor1< double, 3 > & tCurl
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use
Simple interface for fast problem set-up.
Deprecated interface functions.
DeprecatedCoreInterface Interface
MoFEMErrorCode getOptions()
get options
FTensor::Tensor2< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 >, Tensor_Dim0, Tensor_Dim1 > getFTensor2DiffN(FieldApproximationBase base)
Get derivatives of base functions for Hdiv space.
#define CHKERR
Inline error check.
implementation of Data Operators for Forces and Sources
FTensor::Tensor1< double, 3 > & tCurl
default operator for TRI element
int main(int argc, char *argv[])
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.
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Calculate normals at Gauss points of triangle element.
virtual bool check_field(const std::string &name) const =0
check if field is in database
Volume finite element base.
FTensor::Index< 'i', SPACE_DIM > i
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.
OpVolCurl(FTensor::Tensor1< double, 3 > &t_curl)
OpFacesRot(FTensor::Tensor1< double, 3 > &t_curl)
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
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)
#define CATCH_ERRORS
Catch errors.
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
FTensor::Index< 'j', 3 > j
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1N(FieldApproximationBase base)
Get base functions for Hdiv/Hcurl spaces.
Set inverse jacobian to base functions.
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
@ HCURL
field with continuous tangents
PetscErrorCode PetscOptionsGetEList(PetscOptions *, const char pre[], const char name[], const char *const *list, PetscInt next, PetscInt *value, PetscBool *set)
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.
transform local reference derivatives of shape function to global derivatives if higher order geometr...
@ MOFEM_ATOM_TEST_INVALID
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
FTensor::Index< 'k', 3 > k
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)