|
| v0.14.0
|
Go to the documentation of this file.
30 using namespace MoFEM;
32 #include <MethodForForceScaling.hpp>
39 static char help[] =
"...\n\n";
41 int main(
int argc,
char *argv[]) {
52 PetscInt nb_sub_steps = 10;
53 PetscBool flg_file = PETSC_FALSE;
54 PetscBool is_partitioned = PETSC_FALSE;
55 PetscBool is_atom_test = PETSC_FALSE;
59 CHKERR PetscOptionsBegin(PETSC_COMM_WORLD,
"",
60 "Minimal surface area config",
"none");
61 CHKERR PetscOptionsString(
"-my_file",
"mesh file name",
"",
"mesh.h5m",
63 CHKERR PetscOptionsInt(
"-my_order",
"default approximation order",
"", 2,
65 CHKERR PetscOptionsInt(
"-my_nb_sub_steps",
"number of substeps",
"", 10,
66 &nb_sub_steps, PETSC_NULL);
67 CHKERR PetscOptionsBool(
"-my_is_partitioned",
68 "set if mesh is partitioned (this result that "
69 "each process keeps only part of the mesh",
70 "", PETSC_FALSE, &is_partitioned, PETSC_NULL);
73 "is used with testing, exit with error when diverged",
"",
74 PETSC_FALSE, &is_atom_test, PETSC_NULL);
76 ierr = PetscOptionsEnd();
80 if (flg_file != PETSC_TRUE) {
81 SETERRQ(PETSC_COMM_SELF, 1,
"*** ERROR -my_file (MESH FILE NEEDED)");
84 if (is_partitioned == PETSC_TRUE) {
87 option =
"PARALLEL=BCAST_DELETE;"
88 "PARALLEL_RESOLVE_SHARED_ENTS;"
89 "PARTITION=PARALLEL_PARTITION;";
96 ParallelComm *pcomm = ParallelComm::get_pcomm(&moab,
MYPCOMM_INDEX);
98 pcomm =
new ParallelComm(&moab, PETSC_COMM_WORLD);
130 "MINIMAL_SURFACE_ELEMENT",
"U");
132 "MINIMAL_SURFACE_ELEMENT",
"U");
134 "MINIMAL_SURFACE_ELEMENT",
"U");
137 root_set, MBTRI,
"MINIMAL_SURFACE_ELEMENT");
144 Range bc_edges, bc_ents;
149 CHKERR moab.get_entities_by_type(root_set, MBTRI, tris,
false);
151 CHKERR skin.find_skin(root_set, tris,
false, bc_edges);
153 CHKERR moab.get_connectivity(bc_edges, bc_nodes);
154 bc_ents.merge(bc_edges);
155 bc_ents.merge(bc_nodes);
181 CHKERR DMCreate(PETSC_COMM_WORLD, &bc_dm);
182 CHKERR DMSetType(bc_dm,
"MoFEM");
187 CHKERR DMSetFromOptions(bc_dm);
194 CHKERR DMCreateGlobalVector(bc_dm, &T);
196 CHKERR DMCreateMatrix(bc_dm, &
A);
199 minimal_surface_element.
feBcEdge.getOpPtrVector().push_back(
201 minimal_surface_element.
feBcEdge.getOpPtrVector().push_back(
207 CHKERR MatAssemblyBegin(
A, MAT_FINAL_ASSEMBLY);
208 CHKERR MatAssemblyEnd(
A, MAT_FINAL_ASSEMBLY);
212 CHKERR KSPCreate(PETSC_COMM_WORLD, &solver);
213 CHKERR KSPSetOperators(solver,
A,
A);
214 CHKERR KSPSetFromOptions(solver);
217 CHKERR KSPDestroy(&solver);
221 CHKERR VecGhostUpdateBegin(T, INSERT_VALUES, SCATTER_FORWARD);
222 CHKERR VecGhostUpdateEnd(T, INSERT_VALUES, SCATTER_FORWARD);
240 CHKERR DMCreate(PETSC_COMM_WORLD, &dm);
241 CHKERR DMSetType(dm,
"MoFEM");
245 CHKERR DMSetFromOptions(dm);
256 CHKERR DMCreateGlobalVector(dm, &T);
266 auto det_ptr = boost::make_shared<VectorDouble>();
267 auto jac_ptr = boost::make_shared<MatrixDouble>();
268 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
269 minimal_surface_element.
feRhs.getOpPtrVector().push_back(
271 minimal_surface_element.
feRhs.getOpPtrVector().push_back(
274 minimal_surface_element.
feRhs.getOpPtrVector().push_back(
277 minimal_surface_element.
feRhs.getOpPtrVector().push_back(
279 minimal_surface_element.
feRhs.getOpPtrVector().push_back(
281 "U", mse_common_data,
true));
282 minimal_surface_element.
feRhs.getOpPtrVector().push_back(
287 minimal_surface_element.
feLhs.getOpPtrVector().push_back(
289 minimal_surface_element.
feLhs.getOpPtrVector().push_back(
292 minimal_surface_element.
feLhs.getOpPtrVector().push_back(
295 minimal_surface_element.
feLhs.getOpPtrVector().push_back(
297 minimal_surface_element.
feLhs.getOpPtrVector().push_back(
299 "U", mse_common_data,
false));
300 minimal_surface_element.
feLhs.getOpPtrVector().push_back(
307 &minimal_surface_element.
feRhs, PETSC_NULL,
316 &minimal_surface_element.
feLhs, PETSC_NULL,
322 SNESConvergedReason snes_reason;
327 CHKERR SNESCreate(PETSC_COMM_WORLD, &snes);
332 CHKERR SNESSetFromOptions(snes);
335 CHKERR VecDuplicate(T, &T0);
338 double step_size = 1. / nb_sub_steps;
339 for (
int ss = 1; ss <= nb_sub_steps; ss++) {
340 CHKERR VecAXPY(T, step_size, T0);
342 CHKERR SNESSolve(snes, PETSC_NULL, T);
343 CHKERR SNESGetConvergedReason(snes, &snes_reason);
345 CHKERR SNESGetIterationNumber(snes, &its);
346 CHKERR PetscPrintf(PETSC_COMM_WORLD,
347 "number of Newton iterations = %D\n\n", its);
348 if (snes_reason < 0) {
351 "atom test diverged");
357 CHKERR VecGhostUpdateBegin(T, INSERT_VALUES, SCATTER_FORWARD);
358 CHKERR VecGhostUpdateEnd(T, INSERT_VALUES, SCATTER_FORWARD);
366 CHKERR SNESDestroy(&snes);
377 "PARALLEL=WRITE_PART");
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
#define MYPCOMM_INDEX
default communicator number PCOMM
moab::Interface & postProcMesh
Implementation of minimal surface element.
virtual MoFEMErrorCode modify_finite_element_add_field_row(const std::string &fe_name, const std::string name_row)=0
set field row which finite element use
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Implementation of minimal area element.
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
virtual MoFEMErrorCode add_ents_to_field_by_type(const Range &ents, const EntityType type, const std::string &name, int verb=DEFAULT_VERBOSITY)=0
Add entities to field meshset.
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
Minimal surface equation.
Deprecated interface functions.
DeprecatedCoreInterface Interface
Integrate vector on rhs,.
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string &name, const bool recursive=true)=0
add entities to finite element
#define CHKERR
Inline error check.
SurfaceElement feLhs
To calculate left hand side.
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
virtual moab::Interface & get_moab()=0
implementation of Data Operators for Forces and Sources
virtual MoFEMErrorCode modify_finite_element_add_field_col(const std::string &fe_name, const std::string name_row)=0
set field col which finite element use
virtual MoFEMErrorCode build_finite_elements(int verb=DEFAULT_VERBOSITY)=0
Build finite elements.
SurfaceElement feRhs
To calculate right hand side.
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx)
get MoFEM::SnesCtx data structure
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
int main(int argc, char *argv[])
PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES Jacobian evaluation function
Post-process fields on refined mesh.
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.
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string name_field)=0
set finite element field data
PetscErrorCode SnesRhs(SNES snes, Vec x, Vec f, void *ctx)
This is MoFEM implementation for the right hand side (residual vector) evaluation in SNES solver.
Evaluate function values and gradients at Gauss Pts.
Integrate vector on lhs,.
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Keep date between operators.
#define CATCH_ERRORS
Catch errors.
Evaluate function values and gradients at Gauss Pts.
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
const FTensor::Tensor2< T, Dim, Dim > Vec
virtual MoFEMErrorCode build_fields(int verb=DEFAULT_VERBOSITY)=0
OpCalculateHOJacForFaceImpl< 2 > OpCalculateHOJacForFace
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
@ MOFEM_ATOM_TEST_INVALID
virtual MoFEMErrorCode build_adjacencies(const Range &ents, int verb=DEFAULT_VERBOSITY)=0
build adjacencies
Interface for nonlinear (SNES) solver.
virtual MoFEMErrorCode set_field_order(const EntityHandle meshset, const EntityType type, const std::string &name, const ApproximationOrder order, int verb=DEFAULT_VERBOSITY)=0
Set order approximation of the entities in the field.
EdgeElement feBcEdge
Used to calculate dofs on boundary.
PetscErrorCode SnesMat(SNES snes, Vec x, Mat A, Mat B, void *ctx)
This is MoFEM implementation for the left hand side (tangent matrix) evaluation in SNES solver.
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
virtual MoFEMErrorCode add_field(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_EXCL, int verb=DEFAULT_VERBOSITY)=0
Add field.
MoFEMErrorCode generateReferenceElementMesh()
PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES residual evaluation function
MoFEMErrorCode addFieldValuesPostProc(const std::string field_name, Vec v=PETSC_NULL)
Add operator to post-process L2, H1, Hdiv, Hcurl field value.