9#ifndef __FORMS_INTEGRATORS_HPP__
10#define __FORMS_INTEGRATORS_HPP__
20 map<std::string, std::vector<boost::shared_ptr<EssentialBcStorage>>>;
36 boost::shared_ptr<std::vector<unsigned char>> boundary_marker);
65 const double *ptr, InsertMode iora);
105 boost::function<
double(
const double,
const double,
const double)>;
136using VectorFun = boost::function<FTensor::Tensor1<double, DIM>(
137 const double,
const double,
const double)>;
139template <AssemblyType A,
typename EleOp>
struct OpBaseImpl :
public EleOp {
143 OpBaseImpl(
const std::string row_field_name,
const std::string col_field_name,
144 const OpType type, boost::shared_ptr<Range> ents_ptr =
nullptr)
145 : EleOp(row_field_name, col_field_name, type, false),
179 return getFTensor1FromArray<DIM, DIM>(
locF);
185 return getFTensor2FromArray<DIM, DIM, DIM>(
locMat, rr);
216 const bool trans) = 0;
270template <AssemblyType A,
typename EleOp>
279 if (entsPtr->find(this->getFEEntityHandle()) == entsPtr->end())
298 nbIntegrationPts = EleOp::getGaussPts().size2();
300 nbRowBaseFunctions = row_data.
getN().size2();
302 locMat.resize(nbRows, nbCols,
false);
306 CHKERR this->iNtegrate(row_data, col_data);
309 auto check_if_assemble_transpose = [&] {
311 if (row_side != col_side || row_type != col_type)
315 }
else if (assembleTranspose) {
321 CHKERR aSsemble(row_data, col_data, check_if_assemble_transpose());
325template <AssemblyType A,
typename EleOp>
331 if (entsPtr->find(this->getFEEntityHandle()) == entsPtr->end())
343 nbIntegrationPts = EleOp::getGaussPts().size2();
345 nbRowBaseFunctions = row_data.
getN().size2();
350 CHKERR this->iNtegrate(row_data);
352 CHKERR this->aSsemble(row_data);
356template <
typename EleOp>
363 const bool transpose) {
377 this->getKSPB(), col_data, row_data,
384 this->getKSPB(), row_data, col_data, &*this->
locMat.data().begin(),
399 this->getKSPf(), data, &*this->
locF.data().begin(), ADD_VALUES);
#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 ...
#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.
UBlasVector< int > VectorInt
implementation of Data Operators for Forces and Sources
boost::function< double(const FEMethod *fe_ptr)> FEFun
Lambda function used to scale with time.
boost::function< double(double)> TimeFun
Lambda function used to scale with time.
double scalar_fun_one(const double, const double, const double)
MoFEMErrorCode MatSetValues< EssentialBcStorage >(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
Set values to matrix in operator.
boost::function< double()> ConstantFun
Constant function type.
constexpr auto field_name
Data on single entity (This is passed as argument to DataOperator::doWork)
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.
[Storage and set boundary conditions]
static HashVectorStorage feStorage
[Storage and set boundary conditions]
map< std::string, std::vector< boost::shared_ptr< EssentialBcStorage > > > HashVectorStorage
EssentialBcStorage(VectorInt &indices)
structure for User Loop Methods on finite elements
MoFEMErrorCode aSsemble(EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data, const bool transpose)
MoFEMErrorCode aSsemble(EntitiesFieldData::EntData &data)
OpBaseImpl(const std::string row_field_name, const std::string col_field_name, const OpType type, boost::shared_ptr< Range > ents_ptr=nullptr)
VectorDouble locF
local entity vector
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntData &row_data)
Do calculations for the right hand side.
typename EleOp::OpType OpType
TimeFun timeScalingFun
assumes that time variable is set
MatrixDouble locMatTranspose
local entity block matrix
int rowSide
row side number
int colSide
column side number
boost::shared_ptr< Range > entsPtr
Entities on which element is run.
int nbRows
number of dofs on rows
EntityType colType
column type
int nbIntegrationPts
number of integration points
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
Do calculations for the left hand side.
virtual MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
Integrate grad-grad operator.
virtual MoFEMErrorCode iNtegrate(EntData &data)
Class dedicated to integrate operator.
MatrixDouble locMat
local entity block matrix
virtual MoFEMErrorCode aSsemble(EntData &data)=0
FTensor::Tensor2< FTensor::PackPtr< double *, DIM >, DIM, DIM > getLocMat(const int rr)
FEFun feScalingFun
assumes that time variable is set
FTensor::Tensor1< FTensor::PackPtr< double *, DIM >, DIM > getNf()
int nbCols
number if dof on column
virtual MoFEMErrorCode aSsemble(EntData &row_data, EntData &col_data, const bool trans)=0
int nbRowBaseFunctions
number or row base functions
EntityType rowType
row type
Set indices on entities on finite element.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
boost::shared_ptr< std::vector< unsigned char > > boundaryMarker
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.