Calculate zero, first and second moments of inertia.
Calculate zero, first and second moments of inertia.Example intended to show how to write user data operator and integrate scalar field.
static char help[] =
"...\n\n";
private:
struct OpZero;
struct OpFirst;
struct OpSecond;
};
: public boost::enable_shared_from_this<Example::CommonData> {
return boost::shared_ptr<VectorDouble>(shared_from_this(),
}
};
};
OpZero(boost::shared_ptr<CommonData> &common_data_ptr)
}
private:
};
OpFirst(boost::shared_ptr<CommonData> &common_data_ptr)
}
private:
};
OpSecond(boost::shared_ptr<CommonData> &common_data_ptr)
}
private:
};
}
}
int local_size;
else
local_size = 0;
}
auto set_density = [&](
VectorAdaptor &&field_data,
double *xcoord,
double *ycoord, double *zcoord) {
field_data[0] = 1;
};
}
}
}
const double *array;
"First moment of inertia [ %6.4e, %6.4e, %6.4e ]",
"Second moment of inertia [ %6.4e, %6.4e, %6.4e; %6.4e %6.4e; "
"%6.4e ]",
}
}
const double *array;
PetscBool test = PETSC_FALSE;
constexpr double eps = 1e-8;
constexpr double expected_volume = 1.;
constexpr double expected_first_moment = 0.;
constexpr double expected_second_moment = 1. / 12.;
"Wrong area %6.4e != %6.4e", expected_volume,
if (std::abs(array[
i] - expected_first_moment) >
eps)
"Wrong first moment %6.4e != %6.4e", expected_first_moment,
}
if (std::abs(array[
i] - expected_second_moment) >
eps)
"Wrong second moment %6.4e != %6.4e", expected_second_moment,
}
}
}
int main(
int argc,
char *argv[]) {
const char param_file[] = "param_file.petsc";
try {
DMType dm_name = "DMMOFEM";
moab::Core mb_instance;
moab::Interface &moab = mb_instance;
}
}
if (type == MBVERTEX) {
const int nb_integration_pts =
double element_local_value = 0;
for (int gg = 0; gg != nb_integration_pts; ++gg) {
element_local_value += t_w * t_rho * volume;
++t_w;
++t_rho;
}
ADD_VALUES);
}
}
const int nb_integration_pts = getGaussPts().size2();
auto t_w = getFTensor0IntegrationWeight();
auto t_x =
getFTensor1CoordsAtGaussPts();
const double volume = getMeasure();
t_s;
for (int gg = 0; gg != nb_integration_pts; ++gg) {
t_s(
i) += t_w * t_rho * volume * t_x(
i);
++t_w;
++t_rho;
++t_x;
}
constexpr std::array<int, 3> indices = {
ADD_VALUES);
}
const int nb_integration_pts = getGaussPts().size2();
auto t_w = getFTensor0IntegrationWeight();
auto t_x = getFTensor1CoordsAtGaussPts();
const double volume = getMeasure();
std::array<double, 6> element_local_value;
std::fill(element_local_value.begin(), element_local_value.end(), 0.0);
&element_local_value[0], &element_local_value[1], &element_local_value[2],
&element_local_value[3], &element_local_value[4],
&element_local_value[5]);
for (int gg = 0; gg != nb_integration_pts; ++gg) {
t_I(
i,
j) += (t_w * t_rho * volume) * (t_x(
i) ^ t_x(
j));
++t_w;
++t_rho;
++t_x;
}
constexpr std::array<int, 6> indices = {
&element_local_value[0], ADD_VALUES);
}
#define MOFEM_LOG_C(channel, severity, format,...)
void simple(double P1[], double P2[], double P3[], double c[], const int N)
#define CATCH_ERRORS
Catch errors.
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
#define MoFEMFunctionReturnHot(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 ...
@ MOFEM_ATOM_TEST_INVALID
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
MoFEMErrorCode loopFiniteElements(SmartPetscObj< DM > dm=nullptr)
Iterate finite elements.
boost::ptr_deque< UserDataOperator > & getOpDomainRhsPipeline()
Get the Op Domain Rhs Pipeline object.
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'j', 3 > j
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
VectorShallowArrayAdaptor< double > VectorAdaptor
UBlasVector< double > VectorDouble
implementation of Data Operators for Forces and Sources
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
auto createVectorMPI(MPI_Comm comm, PetscInt n, PetscInt N)
Create MPI Vector.
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
MoFEMErrorCode VecSetValues(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
boost::shared_ptr< VectorDouble > getRhoAtIntegrationPtsPtr()
VecElements
Vector to indicate indices for storing zero, first and second moments of inertia.
SmartPetscObj< Vec > petscVec
Smart pointer which stores PETSc distributed vector.
VectorDouble rhoAtIntegrationPts
Storing density at integration points.
boost::shared_ptr< CommonData > commonDataPtr
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
[ZeroOp]
OpFirst(boost::shared_ptr< CommonData > &common_data_ptr)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
[FirstOp]
boost::shared_ptr< CommonData > commonDataPtr
OpSecond(boost::shared_ptr< CommonData > &common_data_ptr)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
[main]
boost::shared_ptr< CommonData > commonDataPtr
OpZero(boost::shared_ptr< CommonData > &common_data_ptr)
MoFEMErrorCode pushOperators()
[Set density distribution]
boost::shared_ptr< CommonData > commonDataPtr
MoFEMErrorCode setUp()
[Run all]
MoFEMErrorCode checkResults()
[Postprocess results]
MoFEMErrorCode createCommonData()
[Set up problem]
MoFEMErrorCode runProblem()
[Run problem]
MoFEM::Interface & mField
MoFEMErrorCode postProcess()
[Integrate]
MoFEMErrorCode setFieldValues()
[Create common data]
MoFEMErrorCode integrateElements()
[Push operators to pipeline]
virtual MPI_Comm & get_comm() const =0
virtual int get_comm_rank() const =0
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.
std::array< bool, MBMAXTYPE > doEntities
If true operator is executed for entity.
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
MoFEMErrorCode setVertexDofs(VertexCoordsFunction lambda, const std::string field_name, Range *verts=nullptr)
Set DOFs on vertices using user function.
auto getFTensor0IntegrationWeight()
Get integration weights.
double getMeasure() const
get measure of element
@ OPROW
operator doWork function is executed on FE rows
MatrixDouble & getGaussPts()
matrix of integration (Gauss) points for Volume Element
Specialization for double precision scalar field values calculation.
PipelineManager interface.
MoFEMErrorCode setDomainRhsIntegrationRule(RuleHookFun rule)
Set integration rule for domain right-hand side finite element.
Simple interface for fast problem set-up.
MoFEMErrorCode getOptions()
get options
intrusive_ptr for managing petsc objects
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Volume finite element base.