|
| v0.14.0
|
Go to the documentation of this file.
27 using namespace MoFEM;
29 static char help[] =
"...\n\n";
64 int main(
int argc,
char *argv[]) {
78 enum bases { AINSWORTH, DEMKOWICZ, LASBASETOP };
79 const char *list_bases[] = {
"ainsworth",
"demkowicz"};
81 PetscInt choice_base_value = AINSWORTH;
83 LASBASETOP, &choice_base_value, &flg);
84 if (flg != PETSC_TRUE)
87 if (choice_base_value == AINSWORTH)
89 else if (choice_base_value == DEMKOWICZ)
95 DMType dm_name =
"DMMOFEM";
99 auto core_log = logging::core::get();
127 auto get_skin = [&]() {
133 CHKERR skin.find_skin(0, body_ents,
false, skin_ents);
137 auto filter_true_skin = [&](
auto skin) {
139 ParallelComm *pcomm =
141 CHKERR pcomm->filter_pstatus(skin, PSTATUS_SHARED | PSTATUS_MULTISHARED,
142 PSTATUS_NOT, -1, &boundary_ents);
143 return boundary_ents;
146 auto boundary_ents = filter_true_skin(get_skin());
154 CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
"SYM",
156 CHKERR bc_mng->removeBlockDOFsOnEntities(
159 auto project_ho_geometry = [&]() {
161 return m_field.
loop_dofs(
"GEOMETRY", ent_method);
163 CHKERR project_ho_geometry();
171 pip_mng->getOpDomainLhsPipeline().clear();
172 pip_mng->getOpDomainRhsPipeline().clear();
175 auto rule = [&](
int,
int,
int p) {
return 2 * p + 2; };
176 CHKERR pip_mng->setDomainRhsIntegrationRule(rule);
177 CHKERR pip_mng->setDomainLhsIntegrationRule(rule);
178 CHKERR pip_mng->setBoundaryRhsIntegrationRule(rule);
181 {{
"U", {-1, 1}}, {
"SIGMA", {-1, 1}}});
187 auto jacobian = [&](
double r) {
201 auto post_proc = [&](
auto dm,
auto f_res,
auto out_name) {
204 boost::make_shared<PostProcBrokenMeshInMoab<DomainEle>>(m_field);
208 auto sigma_ptr = boost::make_shared<MatrixDouble>();
209 post_proc_fe->getOpPtrVector().push_back(
212 auto u_ptr = boost::make_shared<MatrixDouble>();
213 post_proc_fe->getOpPtrVector().push_back(
216 post_proc_fe->getOpPtrVector().push_back(
220 post_proc_fe->getPostProcMesh(), post_proc_fe->getMapGaussPts(),
221 {}, {{
"U", u_ptr}}, {{
"SIGMA", sigma_ptr}}, {})
229 post_proc_fe->writeFile(out_name);
233 auto test_consistency_of_domain_and_bdy_integrals = [&]() {
237 GAUSS>::OpMixDivTimesU<3, SPACE_DIM, SPACE_DIM, COORD_TYPE>;
253 auto ops_rhs_interior = [&](
auto &pip) {
257 auto u_ptr = boost::make_shared<MatrixDouble>();
258 auto grad_u_ptr = boost::make_shared<MatrixDouble>();
263 pip.push_back(
new OpMixDivURhs(
"SIGMA", u_ptr, beta_domain));
265 new OpMixLambdaGradURhs(
"SIGMA", grad_u_ptr, beta_domain));
267 auto sigma_ptr = boost::make_shared<MatrixDouble>();
268 auto sigma_div_ptr = boost::make_shared<MatrixDouble>();
271 "SIGMA", sigma_div_ptr));
273 "SIGMA", sigma_ptr));
275 new OpMixUTimesDivLambdaRhs(
"U", sigma_div_ptr, beta_domain));
276 pip.push_back(
new OpMixUTimesLambdaRhs(
"U", sigma_ptr, beta_domain));
281 auto ops_rhs_boundary = [&](
auto &pip) {
285 auto u_ptr = boost::make_shared<MatrixDouble>();
287 auto traction_ptr = boost::make_shared<MatrixDouble>();
289 "SIGMA", traction_ptr));
294 pip.push_back(
new OpMixNormalLambdaURhs(
"SIGMA", u_ptr, beta_bdy));
295 pip.push_back(
new OpUTimeTractionRhs(
"U", traction_ptr, beta_bdy));
300 CHKERR ops_rhs_interior(pip_mng->getOpDomainRhsPipeline());
301 CHKERR ops_rhs_boundary(pip_mng->getOpBoundaryRhsPipeline());
304 pip_mng->getDomainRhsFE()->f =
f;
305 pip_mng->getBoundaryRhsFE()->f =
f;
310 simple->getDomainFEName(),
311 pip_mng->getDomainRhsFE());
313 simple->getBoundaryFEName(),
314 pip_mng->getBoundaryRhsFE());
320 CHKERR VecNorm(
f, NORM_2, &f_nrm2);
322 MOFEM_LOG(
"ATOM", Sev::inform) <<
"f_norm2 = " << f_nrm2;
323 if (std::fabs(f_nrm2) > 1e-10) {
325 "tensor_divergence_operator_res_vec.h5m");
334 auto test_lhs_ops = [&]() {
341 auto op_lhs_domain = [&](
auto &pip) {
346 auto unity = []() {
return 1; };
348 new OpMixDivULhs(
"SIGMA",
"U", unity, beta_domain,
true,
false));
350 new OpLambdaGraULhs(
"SIGMA",
"U", unity, beta_domain,
true,
false));
354 CHKERR op_lhs_domain(pip_mng->getOpDomainLhsPipeline());
356 auto diff_x = opt->setRandomFields(
simple->getDM(),
357 {{
"U", {-1, 1}}, {
"SIGMA", {-1, 1}}});
358 constexpr
double eps = 1e-5;
359 auto diff_res = opt->checkCentralFiniteDifference(
360 simple->getDM(),
simple->getDomainFEName(), pip_mng->getDomainRhsFE(),
367 CHKERR VecNorm(diff_res, NORM_2, &fnorm_res);
368 MOFEM_LOG_C(
"ATOM", Sev::inform,
"Test Lhs OPs %3.4e", fnorm_res);
369 if (std::fabs(fnorm_res) > 1e-8)
375 CHKERR test_consistency_of_domain_and_bdy_integrals();
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Data on single entity (This is passed as argument to DataOperator::doWork)
#define EXECUTABLE_DIMENSION
ElementsAndOps< SPACE_DIM >::BoundaryEle BoundaryEle
#define MYPCOMM_INDEX
default communicator number PCOMM
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.
constexpr IntegrationType I
Get values at integration pts for tensor filed rank 1, i.e. vector field.
PipelineManager interface.
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
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.
OpBaseImpl< PETSC, EdgeEleOp > OpBase
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
CoordinateTypes
Coodinate system.
Simple interface for fast problem set-up.
Deprecated interface functions.
DeprecatedCoreInterface Interface
FieldSpace
approximation spaces
#define CHKERR
Inline error check.
auto createDMVector(DM dm)
Get smart vector from DM.
virtual moab::Interface & get_moab()=0
implementation of Data Operators for Forces and Sources
Simple interface for fast problem set-up.
#define MOFEM_LOG_C(channel, severity, format,...)
int main(int argc, char *argv[])
void simple(double P1[], double P2[], double P3[], double c[], const int N)
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Calculate tenor field using vectorial base, i.e. Hdiv/Hcurl.
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
constexpr CoordinateTypes COORD_TYPE
FormsIntegrators< DomainEleOp >::Assembly< A >::LinearForm< I >::OpGradTimesTensor< 1, FIELD_DIM, SPACE_DIM > OpGradTimesTensor
SmartPetscObj< Vec > setRandomFields(SmartPetscObj< DM > dm, std::vector< RandomFieldData > random_fields, boost::shared_ptr< Range > ents=nullptr)
Generate random fileds.
Calculate trace of vector (Hdiv/Hcurl) space.
Calculate directional derivative of the right hand side and compare it with tangent matrix derivative...
Add operators pushing bases from local to physical configuration.
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
#define MAX_DOFS_ON_ENTITY
Maximal number of DOFs on entity.
Calculate divergence of tonsorial field using vectorial base.
constexpr FieldSpace HDIV_SPACE
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
#define MOFEM_LOG(channel, severity)
Log.
#define CATCH_ERRORS
Catch errors.
ForcesAndSourcesCore::UserDataOperator UserDataOperator
@ 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)
FieldApproximationBase
approximation base
FTensor::Index< 'i', SPACE_DIM > i
EntitiesFieldData::EntData EntData
@ MOFEM_ATOM_TEST_INVALID
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
@ HDIV
field with continuous normal traction
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Post post-proc data at points from hash maps.