1#ifndef __POISSON2DNONHOMOGENEOUS_HPP__
2#define __POISSON2DNONHOMOGENEOUS_HPP__
19typedef boost::function<
double(
const double,
const double,
const double)>
24 OpDomainLhs(std::string row_field_name, std::string col_field_name)
34 const int nb_row_dofs = row_data.
getIndices().size();
35 const int nb_col_dofs = col_data.
getIndices().size();
37 if (nb_row_dofs && nb_col_dofs) {
39 locLhs.resize(nb_row_dofs, nb_col_dofs,
false);
46 const int nb_integration_points =
getGaussPts().size2();
54 for (
int gg = 0; gg != nb_integration_points; gg++) {
56 const double a = t_w * area;
58 for (
int rr = 0; rr != nb_row_dofs; ++rr) {
62 for (
int cc = 0; cc != nb_col_dofs; cc++) {
64 locLhs(rr, cc) += t_row_diff_base(
i) * t_col_diff_base(
i) *
a;
81 CHKERR MatSetValues<MoFEM::EssentialBcStorage>(
85 if (row_side != col_side || row_type != col_type) {
86 transLocLhs.resize(nb_col_dofs, nb_row_dofs,
false);
88 CHKERR MatSetValues<MoFEM::EssentialBcStorage>(
116 locRhs.resize(nb_dofs,
false);
123 const int nb_integration_points =
getGaussPts().size2();
133 for (
int gg = 0; gg != nb_integration_points; gg++) {
135 const double a = t_w * area;
139 for (
int rr = 0; rr != nb_dofs; rr++) {
141 locRhs[rr] += t_base * body_source *
a;
154 CHKERR VecSetValues<MoFEM::EssentialBcStorage>(
180 const int nb_row_dofs = row_data.
getIndices().size();
181 const int nb_col_dofs = col_data.
getIndices().size();
183 if (nb_row_dofs && nb_col_dofs) {
192 const int nb_integration_points =
getGaussPts().size2();
200 for (
int gg = 0; gg != nb_integration_points; gg++) {
201 const double a = t_w * edge;
203 for (
int rr = 0; rr != nb_row_dofs; ++rr) {
207 for (
int cc = 0; cc != nb_col_dofs; cc++) {
223 CHKERR MatSetValues<MoFEM::EssentialBcStorage>(
225 if (row_side != col_side || row_type != col_type) {
228 CHKERR MatSetValues<MoFEM::EssentialBcStorage>(
263 const int nb_integration_points =
getGaussPts().size2();
273 for (
int gg = 0; gg != nb_integration_points; gg++) {
275 const double a = t_w * edge;
276 double boundary_term =
279 for (
int rr = 0; rr != nb_dofs; rr++) {
294 CHKERR VecSetValues<MoFEM::EssentialBcStorage>(
#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.
boost::function< double(const double, const double, const double)> ScalarFunc
FTensor::Index< 'i', 2 > i
constexpr auto field_name
bool sYmm
If true assume that matrix is symmetric structure.
default operator for EDGE element
Data on single entity (This is passed as argument to DataOperator::doWork)
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1DiffN(const FieldApproximationBase base)
Get derivatives of base functions.
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.
const VectorInt & getIndices() const
Get global indices of dofs on entity.
default operator for TRI element
auto getFTensor1CoordsAtGaussPts()
Get coordinates at integration points assuming linear geometry.
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
MatrixDouble & getGaussPts()
matrix of integration (Gauss) points for Volume Element
MatrixDouble transLocBoundaryLhs
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
Operator for bi-linear form, usually to calculate values on left hand side.
MatrixDouble locBoundaryLhs
OpBoundaryLhs(std::string row_field_name, std::string col_field_name)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
VectorDouble locBoundaryRhs
OpBoundaryRhs(std::string field_name, ScalarFunc boundary_function)
OpDomainLhs(std::string row_field_name, std::string col_field_name)
boost::shared_ptr< std::vector< unsigned char > > boundaryMarker
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
Operator for bi-linear form, usually to calculate values on left hand side.
ScalarFunc sourceTermFunc
OpDomainRhs(std::string field_name, ScalarFunc source_term_function)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.