28static char help[] =
"...\n\n";
41 const int nb_int_pts = getGaussPts().size2();
43 auto t_w = getFTensor0IntegrationWeight();
44 double v = getMeasure();
46 for (
int gg = 0; gg != nb_int_pts; gg++) {
50 CHKERR VecSetValue(
vOl, 0, vol, ADD_VALUES);
80 for (
int gg = 0; gg != nb_int_pts; gg++) {
81 vol += (t_coords(
i) * t_normal(
i)) * t_w;
87 CHKERR VecSetValue(
vOl, 0, vol, ADD_VALUES);
107int main(
int argc,
char *argv[]) {
123 DMType dm_name =
"DMMOFEM";
150 auto dm = simple_interface->
getDM();
153 boost::make_shared<VolumeElementForcesAndSourcesCore>(m_field);
155 boost::make_shared<FaceElementForcesAndSourcesCore>(m_field);
158 domain_fe->getRuleHook =
VolRule();
159 boundary_fe->getRuleHook =
FaceRule();
169 auto material_grad_mat = boost::make_shared<MatrixDouble>();
170 auto material_det_vec = boost::make_shared<VectorDouble>();
172 domain_fe->meshPositionsFieldName =
"none";
173 domain_fe->getOpPtrVector().push_back(
177 material_grad_mat, material_det_vec,
nullptr));
178 domain_fe->getOpPtrVector().push_back(
180 domain_fe->getOpPtrVector().push_back(
182 domain_fe->getOpPtrVector().push_back(
183 new OpVolume(
"MESH_NODE_POSITIONS", vol));
185 boundary_fe->getOpPtrVector().push_back(
187 boundary_fe->getOpPtrVector().push_back(
188 new OpFace(
"MESH_NODE_POSITIONS", surf_vol));
196 auto skeleton_fe = boost::make_shared<FEMethod>();
207 skeleton_fe->x_t = x_t;
208 skeleton_fe->x_tt = x_tt;
210 skeleton_fe->preProcessHook = [&]() {
212 if (
f != skeleton_fe->ts_F)
214 if (
A != skeleton_fe->ts_A)
216 if (
B != skeleton_fe->ts_B)
218 if (x != skeleton_fe->ts_u)
220 if (x_t != skeleton_fe->ts_u_t)
222 if (x_tt != skeleton_fe->ts_u_tt)
227 skeleton_fe->postProcessHook = []() {
return 0; };
228 skeleton_fe->operatorHook = []() {
return 0; };
234 CHKERR VecAssemblyBegin(vol);
235 CHKERR VecAssemblyEnd(vol);
236 CHKERR VecGhostUpdateBegin(vol, ADD_VALUES, SCATTER_REVERSE);
237 CHKERR VecGhostUpdateEnd(vol, ADD_VALUES, SCATTER_REVERSE);
238 CHKERR VecAssemblyBegin(surf_vol);
239 CHKERR VecAssemblyEnd(surf_vol);
240 CHKERR VecGhostUpdateBegin(surf_vol, ADD_VALUES, SCATTER_REVERSE);
241 CHKERR VecGhostUpdateEnd(surf_vol, ADD_VALUES, SCATTER_REVERSE);
244 CHKERR VecGetArray(vol, &a_vol);
246 CHKERR VecGetArray(surf_vol, &a_surf_vol);
247 cout <<
"Volume = " << a_vol[0] << endl;
248 cout <<
"Surf Volume = " << a_surf_vol[0] << endl;
249 if (fabs(a_vol[0] - a_surf_vol[0]) > 1e-12) {
252 CHKERR VecRestoreArray(vol, &a_vol);
253 CHKERR VecRestoreArray(vol, &a_surf_vol);
ForcesAndSourcesCore::UserDataOperator UserDataOperator
#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
@ MOFEM_DATA_INCONSISTENCY
#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 ...
auto smartCreateDMMatrix(DM dm)
Get smart matrix from DM.
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
auto smartCreateDMVector(DM dm)
Get smart vector from DM.
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.
FTensor::Index< 'i', SPACE_DIM > i
double v
phase velocity of light in medium (cm/ns)
const FTensor::Tensor2< T, Dim, Dim > Vec
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
implementation of Data Operators for Forces and Sources
SmartPetscObj< Vec > smartVectorDuplicate(SmartPetscObj< Vec > &vec)
Create duplicate vector of smart vector.
auto createSmartGhostVector(MPI_Comm comm, PetscInt n, PetscInt N, PetscInt nghost, const PetscInt ghosts[])
Create smart ghost vector.
DeprecatedCoreInterface Interface
constexpr auto field_name
int main(int argc, char *argv[])
Set integration rule to boundary elements.
int operator()(int, int, int) const
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.
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
default operator for TRI element
auto getFTensor1NormalsAtGaussPts()
get normal at integration points
auto getFTensor1CoordsAtGaussPts()
Get coordinates at integration points assuming linear geometry.
auto getFTensor0IntegrationWeight()
Get integration weights.
@ OPROW
operator doWork function is executed on FE rows
MatrixDouble & getGaussPts()
matrix of integration (Gauss) points for Volume Element
Calculate HO coordinates at gauss points.
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
Set inverse jacobian to base functions.
Projection of edge entities with one mid-node on hierarchical basis.
Simple interface for fast problem set-up.
const std::string getBoundaryFEName() const
Get the Boundary FE Name.
const std::string getSkeletonFEName() const
Get the Skeleton FE Name.
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.
MoFEMErrorCode getOptions()
get options
MoFEMErrorCode getDM(DM *dm)
Get DM.
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name)
Load mesh file.
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
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.
MoFEMErrorCode addSkeletonField(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 skeleton.
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
const std::string getDomainFEName() const
Get the Domain FE Name.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Volume finite element base.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
Operator for bi-linear form, usually to calculate values on left hand side.
OpFace(const std::string &field_name, Vec vol)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
OpVolume(const std::string &field_name, Vec vol)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
int operator()(int, int, int) const