83#ifndef __FORMS_INTEGRATORS_HPP__
84#define __FORMS_INTEGRATORS_HPP__
100 map<std::string, std::vector<boost::shared_ptr<EssentialBcStorage>>>;
116 boost::shared_ptr<std::vector<unsigned char>> boundary_marker);
172 const double *ptr, InsertMode iora);
191 const double *ptr, InsertMode iora);
249 boost::function<
double(
const double,
const double,
const double)>;
289using VectorFun = boost::function<FTensor::Tensor1<double, DIM>(
290 const double,
const double,
const double)>;
303 OpBaseImpl(
const std::string row_field_name,
const std::string col_field_name,
304 const OpType type, boost::shared_ptr<Range> ents_ptr =
nullptr)
305 :
EleOp(row_field_name, col_field_name, type, false),
319 EntityType col_type,
EntData &row_data,
353 return getFTensor1FromArray<DIM, DIM>(
locF);
365 return getFTensor2FromArray<DIM, DIM, DIM>(
locMat, rr);
429template <AssemblyType A,
typename EleOp>
432 [](ForcesAndSourcesCore::UserDataOperator *op_ptr,
433 const EntitiesFieldData::EntData &row_data,
434 const EntitiesFieldData::EntData &col_data,
MatrixDouble &
m) {
435 return MatSetValues<AssemblyTypeSelector<A>>(
436 op_ptr->getKSPB(), row_data, col_data,
m, ADD_VALUES);
439template <
typename OpBase>
440struct OpBrokenBaseImpl;
493template <AssemblyType A,
typename EleOp>
500 <<
" functions is null, returning 0";
504 auto nb_base_functions = data.
getN().size2();
515 nb_base_functions /= 3;
517 if (data.
getN().size2() % 3) {
519 "Number of base functions is not divisible by 3");
530 return nb_base_functions;
533template <AssemblyType A,
typename EleOp>
542 if (entsPtr->find(this->getFEEntityHandle()) == entsPtr->end())
561 nbIntegrationPts = EleOp::getGaussPts().size2();
563 nbRowBaseFunctions = getNbOfBaseFunctions(row_data);
566 locMat.resize(nbRows, nbCols,
false);
570 CHKERR this->iNtegrate(row_data, col_data);
573 auto check_if_assemble_transpose = [&] {
575 if (row_side != col_side || row_type != col_type)
579 }
else if (assembleTranspose) {
585 CHKERR aSsemble(row_data, col_data, check_if_assemble_transpose());
589template <AssemblyType A,
typename EleOp>
595 if (entsPtr->find(this->getFEEntityHandle()) == entsPtr->end())
607 nbIntegrationPts = EleOp::getGaussPts().size2();
609 nbRowBaseFunctions = getNbOfBaseFunctions(row_data);
614 CHKERR this->iNtegrate(row_data);
616 CHKERR this->aSsemble(row_data);
620template <AssemblyType A,
typename EleOp>
623 const bool transpose) {
626 if (!this->timeScalingFun.empty())
627 this->locMat *= this->timeScalingFun(this->getFEMethod()->ts_t);
628 if (!this->feScalingFun.empty())
629 this->locMat *= this->feScalingFun(this->getFEMethod());
633 this->locMatTranspose.resize(this->locMat.size2(), this->locMat.size1(),
635 noalias(this->locMatTranspose) = trans(this->locMat);
636 CHKERR matSetValuesHook(
this, col_data, row_data, this->locMatTranspose);
639 if (!this->onlyTranspose) {
641 CHKERR matSetValuesHook(
this, row_data, col_data, this->locMat);
647template <AssemblyType A,
typename EleOp>
649 if (!this->timeScalingFun.empty())
650 this->locF *= this->timeScalingFun(this->getFEMethod()->ts_t);
651 if (!this->feScalingFun.empty())
652 this->locF *= this->feScalingFun(this->getFEMethod());
654 return VecSetValues<AssemblyTypeSelector<A>>(this->getKSPf(), data,
655 this->locF, ADD_VALUES);
@ USER_BASE
user implemented approximation base
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
@ L2
field with C-1 continuity
@ NOFIELD
scalar or vector of scalars describe (no true field)
@ HCURL
field with continuous tangents
@ HDIV
field with continuous normal traction
static const char *const FieldSpaceNames[]
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_DATA_INCONSISTENCY
static const char *const ApproximationBaseNames[]
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MOFEM_LOG(channel, severity)
Log.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
UBlasMatrix< double > MatrixDouble
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.
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
FTensor::Index< 'm', 3 > m
Template selector for assembly type.
Data on single entity (This is passed as argument to DataOperator::doWork)
FieldApproximationBase & getBase()
Get approximation base.
virtual boost::shared_ptr< MatrixDouble > & getNSharedPtr(const FieldApproximationBase base, const BaseDerivatives derivative)
Get shared pointer to base functions with derivatives.
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
FieldSpace & getSpace()
Get field space.
const VectorInt & getIndices() const
Get global indices of degrees of freedom 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
Store modified indices by field.
EssentialBcStorage(VectorInt &indices)
Structure for user loop methods on finite elements.
OpBaseImpl(const std::string row_field_name, const std::string col_field_name, const OpType type, boost::shared_ptr< Range > ents_ptr=nullptr)
Constructor for base operator implementation.
VectorDouble locF
local entity vector
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.
static MatSetValuesHook matSetValuesHook
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 aSsemble(EntData &row_data, EntData &col_data, const bool trans)
Assemble local matrix into global matrix.
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 size_t getNbOfBaseFunctions(EntitiesFieldData::EntData &data)
Get number of base functions.
FTensor::Tensor2< FTensor::PackPtr< double *, DIM >, DIM, DIM > getLocMat(const int rr)
Get local matrix tensor for assembly.
FEFun feScalingFun
assumes that time variable is set
FTensor::Tensor1< FTensor::PackPtr< double *, DIM >, DIM > getNf()
Get local vector tensor for assembly.
int nbCols
number if dof on column
int nbRowBaseFunctions
number or row base functions
EntityType rowType
row type
boost::function< MoFEMErrorCode(ForcesAndSourcesCore::UserDataOperator *op_ptr, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, MatrixDouble &m)> MatSetValuesHook
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
Unset indices on entities on finite element.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Perform unset operation on entity data.