7#ifndef __FORMS_INTEGRATORS_HPP__
8#define __FORMS_INTEGRATORS_HPP__
24 map<std::string, std::vector<boost::shared_ptr<EssentialBcStorage>>>;
40 boost::shared_ptr<std::vector<unsigned char>> boundary_marker);
76 const double *ptr, InsertMode iora);
95 const double *ptr, InsertMode iora);
144 boost::function<
double(
const double,
const double,
const double)>;
175using VectorFun = boost::function<FTensor::Tensor1<double, DIM>(
176 const double,
const double,
const double)>;
182 OpBaseImpl(
const std::string row_field_name,
const std::string col_field_name,
183 const OpType type, boost::shared_ptr<Range> ents_ptr =
nullptr)
184 :
EleOp(row_field_name, col_field_name, type, false),
198 EntityType col_type,
EntData &row_data,
285template <AssemblyType A,
typename EleOp>
288 [](ForcesAndSourcesCore::UserDataOperator *op_ptr,
289 const EntitiesFieldData::EntData &row_data,
290 const EntitiesFieldData::EntData &col_data, MatrixDouble &
m) {
292 op_ptr->getKSPB(), row_data, col_data,
m, ADD_VALUES);
295template <
typename OpBase>
349template <AssemblyType A,
typename EleOp>
356 <<
" functions is null, returning 0";
360 auto nb_base_functions = data.
getN().size2();
371 nb_base_functions /= 3;
373 if (data.
getN().size2() % 3) {
375 "Number of base functions is not divisible by 3");
386 return nb_base_functions;
389template <AssemblyType A,
typename EleOp>
398 if (entsPtr->find(this->getFEEntityHandle()) == entsPtr->end())
417 nbIntegrationPts = EleOp::getGaussPts().size2();
419 nbRowBaseFunctions = getNbOfBaseFunctions(row_data);
422 locMat.resize(nbRows, nbCols,
false);
426 CHKERR this->iNtegrate(row_data, col_data);
429 auto check_if_assemble_transpose = [&] {
431 if (row_side != col_side || row_type != col_type)
435 }
else if (assembleTranspose) {
441 CHKERR aSsemble(row_data, col_data, check_if_assemble_transpose());
445template <AssemblyType A,
typename EleOp>
451 if (entsPtr->find(this->getFEEntityHandle()) == entsPtr->end())
463 nbIntegrationPts = EleOp::getGaussPts().size2();
465 nbRowBaseFunctions = getNbOfBaseFunctions(row_data);
470 CHKERR this->iNtegrate(row_data);
472 CHKERR this->aSsemble(row_data);
476template <AssemblyType A,
typename EleOp>
479 const bool transpose) {
482 if (!this->timeScalingFun.empty())
483 this->locMat *= this->timeScalingFun(this->getFEMethod()->ts_t);
484 if (!this->feScalingFun.empty())
485 this->locMat *= this->feScalingFun(this->getFEMethod());
489 this->locMatTranspose.resize(this->locMat.size2(), this->locMat.size1(),
491 noalias(this->locMatTranspose) = trans(this->locMat);
492 CHKERR matSetValuesHook(
this, col_data, row_data, this->locMatTranspose);
495 if (!this->onlyTranspose) {
497 CHKERR matSetValuesHook(
this, row_data, col_data, this->locMat);
503template <AssemblyType A,
typename EleOp>
505 if (!this->timeScalingFun.empty())
506 this->locF *= this->timeScalingFun(this->getFEMethod()->ts_t);
507 if (!this->feScalingFun.empty())
508 this->locF *= this->feScalingFun(this->getFEMethod());
511 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.
UBlasVector< int > VectorInt
implementation of Data Operators for Forces and Sources
FTensor::Tensor2< FTensor::PackPtr< double *, S >, DIM1, DIM2 > getFTensor2FromArray(MatrixDouble &data, const size_t rr, const size_t cc=0)
boost::function< double(double)> TimeFun
Lambda function used to scale with time.
boost::function< double(const FEMethod *fe_ptr)> FEFun
Lambda function used to scale with time.
MoFEMErrorCode MatSetValues(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
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.
MoFEMErrorCode VecSetValues(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
boost::function< double()> ConstantFun
Constant function type.
auto getFTensor1FromArray(VectorDouble &data)
Get FTensor1 from array.
constexpr auto field_name
FTensor::Index< 'm', 3 > m
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)
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 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
Store modifed 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)
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)
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.
boost::function< MoFEMErrorCode( ForcesAndSourcesCore::UserDataOperator *op_ptr, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, MatrixDouble &m)> MatSetValuesHook
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
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.
OpSetBc(std::string field_name, bool yes_set, boost::shared_ptr< std::vector< unsigned char > > boundary_marker)
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.
OpUnSetBc(std::string field_name)