26#include <BasicFiniteElements.hpp>
27#include <PoissonOperators.hpp>
30static char help[] =
"...\n\n";
45 double operator()(
const double x,
const double y,
const double z)
const {
46 return 1 + x * x + y * y + z * z * z;
55 const double z)
const {
77 double operator()(
const double x,
const double y,
const double z)
const {
82int main(
int argc,
char *argv[]) {
95 PetscBool flg_test = PETSC_FALSE;
96 CHKERR PetscOptionsBegin(PETSC_COMM_WORLD,
"",
"Poisson's problem options",
102 CHKERR PetscOptionsBool(
"-test",
"if true is ctest",
"", flg_test,
103 &flg_test, PETSC_NULL);
104 ierr = PetscOptionsEnd();
124 boost::shared_ptr<ForcesAndSourcesCore> domain_lhs_fe;
125 boost::shared_ptr<ForcesAndSourcesCore> boundary_lhs_fe;
126 boost::shared_ptr<ForcesAndSourcesCore> domain_rhs_fe;
127 boost::shared_ptr<ForcesAndSourcesCore> boundary_rhs_fe;
128 boost::shared_ptr<ForcesAndSourcesCore> domain_error;
129 boost::shared_ptr<ForcesAndSourcesCore> post_proc_volume;
130 boost::shared_ptr<ForcesAndSourcesCore> null;
136 boundary_lhs_fe, domain_rhs_fe, boundary_rhs_fe);
141 global_error, domain_error);
149 Simple *simple_interface;
157 CHKERR simple_interface->getOptions();
159 CHKERR simple_interface->loadFile();
185 CHKERR simple_interface->setFieldOrder(
"ERROR",0);
189 CHKERR simple_interface->setUp();
201 CHKERR simple_interface->getDM(&dm);
212 dm, simple_interface->getDomainFEName(), domain_lhs_fe, null, null);
214 dm, simple_interface->getBoundaryFEName(), boundary_lhs_fe, null,
218 domain_rhs_fe, null, null);
220 boundary_rhs_fe, null, null);
228 CHKERR DMCreateGlobalVector(dm,&F);
235 CHKERR KSPCreate(PETSC_COMM_WORLD,&solver);
236 CHKERR KSPSetFromOptions(solver);
237 CHKERR KSPSetDM(solver,dm);
249 CHKERR KSPDestroy(&solver);
262 if (flg_test == PETSC_TRUE) {
273 CHKERR boost::static_pointer_cast<PostProcVolumeOnRefinedMesh>(
275 ->writeFile(
"out_vol.h5m");
282 CHKERR VecDestroy(&global_error);
int main(int argc, char *argv[])
#define CATCH_ERRORS
Catch errors.
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
@ L2
field with C-1 continuity
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
#define CHKERR
Inline error check.
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set KSP right hand side evaluation function
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
Set KSP operators and push mofem finite element methods.
const FTensor::Tensor2< T, Dim, Dim > Vec
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
DeprecatedCoreInterface Interface
const double D
diffusivity
FTensor::Tensor1< double, 3 > operator()(const double x, const double y, const double z) const
double operator()(const double x, const double y, const double z) const
double operator()(const double x, const double y, const double z) const
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Deprecated interface functions.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
MoFEMErrorCode createGhostVec(Vec *ghost_vec) const
MoFEMErrorCode testError(Vec ghost_vec)
Test error.
MoFEMErrorCode assembleGhostVector(Vec ghost_vec) const
Assemble error vector.
MoFEMErrorCode printError(Vec ghost_vec)
Print error.
Create finite elements instances.
MoFEMErrorCode creatFEToPostProcessResults(boost::shared_ptr< ForcesAndSourcesCore > &post_proc_volume) const
Create finite element to post-process results.
MoFEMErrorCode createFEToEvaluateError(boost::function< double(const double, const double, const double)> f_u, boost::function< FTensor::Tensor1< double, 3 >(const double, const double, const double)> g_u, Vec global_error, boost::shared_ptr< ForcesAndSourcesCore > &domain_error) const
Create finite element to calculate error.
MoFEMErrorCode createFEToAssembleMatrixAndVector(boost::function< double(const double, const double, const double)> f_u, boost::function< double(const double, const double, const double)> f_source, boost::shared_ptr< ForcesAndSourcesCore > &domain_lhs_fe, boost::shared_ptr< ForcesAndSourcesCore > &boundary_lhs_fe, boost::shared_ptr< ForcesAndSourcesCore > &domain_rhs_fe, boost::shared_ptr< ForcesAndSourcesCore > &boundary_rhs_fe, bool trans=true) const
Create finite element to calculate matrix and vectors.