1#ifndef __POISSONOPERTATORS_HPP__
2#define __POISSONOPERTATORS_HPP__
10using OpVolEle = FaceElementForcesAndSourcesCore::UserDataOperator;
15using EntData = DataForcesAndSourcesCore::EntData;
25 OpLhsUU(std::string domain_field1, std::string domain_field2)
35 const int nb_row_dofs = row_data.getIndices().size();
36 const int nb_col_dofs = col_data.getIndices().size();
38 if (nb_row_dofs && nb_col_dofs) {
39 locLhs.resize(nb_row_dofs, nb_col_dofs,
false);
41 const int nb_integration_pts =
getGaussPts().size2();
44 auto t_row_diff_base = row_data.getFTensor1DiffN<2>();
47 for (
int gg = 0; gg != nb_integration_pts; ++gg) {
48 const double a = t_w * vol;
49 for (
int rr = 0; rr != nb_row_dofs; ++rr) {
50 auto t_col_diff_base = col_data.getFTensor1DiffN<2>(gg, 0);
51 for (
int cc = 0; cc != nb_col_dofs; ++cc) {
52 locLhs(rr, cc) +=
a * t_row_diff_base(
i) * t_col_diff_base(
i);
60 &
locLhs(0, 0), ADD_VALUES);
61 if (row_side != col_side || row_type != col_type) {
62 transLocLhs.resize(nb_col_dofs, nb_row_dofs,
false);
82 const int nb_dofs = data.getIndices().size();
84 locRhs.resize(nb_dofs,
false);
86 const int nb_integration_pts =
getGaussPts().size2();
88 auto t_base = data.getFTensor0N();
91 for (
int gg = 0; gg != nb_integration_pts; ++gg) {
92 const double a = vol * t_w;
94 for (
int rr = 0; rr != nb_dofs; ++rr) {
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
MoFEMErrorCode MatSetValues(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
MoFEMErrorCode VecSetValues(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
DataForcesAndSourcesCore::EntData EntData
FTensor::Index< 'i', 2 > i
bool sYmm
If true assume that matrix is symmetric structure.
default operator for EDGE element
default operator for TRI element
auto getFTensor0IntegrationWeight()
Get integration weights.
double getMeasure() const
get measure of element
@ OPROW
operator doWork function is executed on FE rows
@ OPROWCOL
operator doWork is executed on FE rows &columns
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
MatrixDouble & getGaussPts()
matrix of integration (Gauss) points for Volume Element
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
OpLhsUU(std::string domain_field1, std::string domain_field2)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
OpRhsU(std::string domain_field)