#ifndef __POISSON2DISGALERKIN_HPP__
#define __POISSON2DISGALERKIN_HPP__
std::array<VectorInt, 2>
std::array<VectorInt, 2>
std::fill(&doEntities[MBVERTEX], &doEntities[MBMAXTYPE], false);
for (
auto t = moab::CN::TypeDimensionMap[
SPACE_DIM].first;
t <= moab::CN::TypeDimensionMap[
SPACE_DIM].second; ++
t)
}
if ((CN::Dimension(row_type) ==
SPACE_DIM) &&
const auto nb_in_loop = getFEMethod()->nInTheLoop;
areaMap[nb_in_loop] = getMeasure();
senseMap[nb_in_loop] = getSkeletonSense();
if (!nb_in_loop) {
}
} else {
}
}
};
template <
typename T>
inline auto get_ntensor(
T &base_mat) {
&*base_mat.data().begin());
};
template <
typename T>
inline auto get_ntensor(
T &base_mat,
int gg,
int bb) {
double *ptr = &base_mat(gg, bb);
};
double *ptr = &*base_mat.data().begin();
return getFTensor1FromPtr<2>(ptr);
};
template <typename T>
double *ptr = &base_mat(gg, 2 * bb);
return getFTensor1FromPtr<2>(ptr);
};
public:
OpL2LhsPenalty(boost::shared_ptr<FaceSideEle> side_fe_ptr)
EntitiesFieldData::EntData &data) {
const auto in_the_loop =
#ifndef NDEBUG
const std::array<std::string, 2> ele_type_name = {"BOUNDARY", "SKELETON"};
<< "OpL2LhsPenalty inTheLoop " << ele_type_name[in_the_loop];
#endif
t_normal.normalize();
const size_t nb_integration_pts =
getGaussPts().size2();
const double beta =
static_cast<double>(
nitsche) / (in_the_loop + 1);
if (nb_rows) {
locMat.resize(nb_rows, nb_cols,
false);
for (size_t gg = 0; gg != nb_integration_pts; ++gg) {
auto t_mat =
locMat.data().begin();
size_t rr = 0;
for (; rr != nb_rows; ++rr) {
t_vn_plus(
i) = beta * (
phi * t_diff_row_base(
i) /
p);
t_vn(
i) = t_row_base * t_normal(
i) * sense_row - t_vn_plus(
i);
auto t_diff_col_base =
for (size_t cc = 0; cc != nb_cols; ++cc) {
t_un(
i) = -
p * (t_col_base * t_normal(
i) * sense_col -
beta * t_diff_col_base(
i) /
p);
*t_mat -= alpha * (t_vn(
i) * t_un(
i));
*t_mat -= alpha * (t_vn_plus(
i) * (beta * t_diff_col_base(
i)));
++t_col_base;
++t_diff_col_base;
++t_mat;
}
++t_row_base;
++t_diff_row_base;
}
for (; rr < nb_row_base_functions; ++rr) {
++t_row_base;
++t_diff_row_base;
}
++t_w;
}
&*
locMat.data().begin(), ADD_VALUES);
if (!in_the_loop)
}
}
}
}
private:
boost::shared_ptr<FaceSideEle>
};
public:
OpL2BoundaryRhs(boost::shared_ptr<FaceSideEle> side_fe_ptr, ScalarFun
fun)
const auto in_the_loop =
t_normal.normalize();
if (!nb_rows)
rhsF.resize(nb_rows,
false);
const size_t nb_integration_pts =
getGaussPts().size2();
const double beta =
static_cast<double>(
nitsche) / (in_the_loop + 1);
for (size_t gg = 0; gg != nb_integration_pts; ++gg) {
const double source_val =
-
p *
sourceFun(t_coords(0), t_coords(1), t_coords(2));
auto t_f =
rhsF.data().begin();
size_t rr = 0;
for (; rr != nb_rows; ++rr) {
t_vn_plus(
i) = beta * (
phi * t_diff_row_base(
i) /
p);
t_vn(
i) = t_row_base * t_normal(
i) * sense_row - t_vn_plus(
i);
*t_f -= alpha * t_vn(
i) * (source_val * t_normal(
i));
++t_row_base;
++t_diff_row_base;
++t_f;
}
for (; rr < nb_row_base_functions; ++rr) {
++t_row_base;
++t_diff_row_base;
}
++t_w;
++t_coords;
}
ADD_VALUES);
}
private:
boost::shared_ptr<FaceSideEle>
};
};
#endif
BoundaryEle::UserDataOperator BoundaryEleOp
#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_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 fun
Function to approximate.
ElementsAndOps< SPACE_DIM >::FaceSideOp FaceSideOp
#define MOFEM_LOG(channel, severity)
Log.
auto get_ntensor(T &base_mat)
std::array< MatrixDouble, 2 > rowBaseSideMap
std::array< VectorInt, 2 > indicesColSideMap
indices on columns for left hand-side
FTensor::Index< 'i', SPACE_DIM > i
std::array< MatrixDouble, 2 > colDiffBaseSideMap
std::array< double, 2 > areaMap
std::array< MatrixDouble, 2 > colBaseSideMap
std::array< VectorInt, 2 > indicesRowSideMap
indices on rows for left hand-side
std::array< int, 2 > senseMap
auto get_diff_ntensor(T &base_mat)
std::array< MatrixDouble, 2 > rowDiffBaseSideMap
constexpr double t
plate stiffness
constexpr auto field_name
Data on single entity (This is passed as argument to DataOperator::doWork)
MatrixDouble & getDiffN(const FieldApproximationBase base)
get derivatives of base functions
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
const VectorInt & getIndices() const
Get global indices of dofs on entity.
default operator for TRI element
auto getFTensor1Normal()
get normal as tensor
auto getFTensor1CoordsAtGaussPts()
Get coordinates at integration points assuming linear geometry.
auto getFTensor0IntegrationWeight()
Get integration weights.
double getMeasure() const
get measure of element
@ OPSPACE
operator do Work is execute on space data
MatrixDouble & getGaussPts()
matrix of integration (Gauss) points for Volume Element
Operator tp collect data from elements on the side of Edge/Face.
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
VectorDouble rhsF
vector to store local operator right hand side
ScalarFun sourceFun
pointer to function to evaluate value of function on boundary
boost::shared_ptr< FaceSideEle > sideFEPtr
pointer to element to get data on edge/face sides
MatrixDouble locMat
local operator matrix
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
boost::shared_ptr< FaceSideEle > sideFEPtr
pointer to element to get data on edge/face sides