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);
136 boost::function<
double(
const double,
const double,
const double)>;
167using VectorFun = boost::function<FTensor::Tensor1<double, DIM>(
168 const double,
const double,
const double)>;
170template <AssemblyType A,
typename EleOp>
struct OpBaseImpl :
public EleOp {
174 OpBaseImpl(
const std::string row_field_name,
const std::string col_field_name,
175 const OpType type, boost::shared_ptr<Range> ents_ptr =
nullptr)
176 : EleOp(row_field_name, col_field_name, type, false),
230 return this->getKSPA();
232 return this->getKSPB();
234 return this->getKSPB();
240 return getFTensor1FromArray<DIM, DIM>(
locF);
246 return getFTensor2FromArray<DIM, DIM, DIM>(
locMat, rr);
331template <AssemblyType A,
typename EleOp>
340 if (entsPtr->find(this->getFEEntityHandle()) == entsPtr->end())
359 nbIntegrationPts = EleOp::getGaussPts().size2();
361 nbRowBaseFunctions = row_data.
getN().size2();
363 locMat.resize(nbRows, nbCols,
false);
367 CHKERR this->iNtegrate(row_data, col_data);
370 auto check_if_assemble_transpose = [&] {
372 if (row_side != col_side || row_type != col_type)
376 }
else if (assembleTranspose) {
382 CHKERR aSsemble(row_data, col_data, check_if_assemble_transpose());
386template <AssemblyType A,
typename EleOp>
392 if (entsPtr->find(this->getFEEntityHandle()) == entsPtr->end())
404 nbIntegrationPts = EleOp::getGaussPts().size2();
406 nbRowBaseFunctions = row_data.
getN().size2();
411 CHKERR this->iNtegrate(row_data);
413 CHKERR this->aSsemble(row_data);
417template <AssemblyType A,
typename EleOp>
420 const bool transpose) {
423 if (!this->timeScalingFun.empty())
424 this->locMat *= this->timeScalingFun(this->getFEMethod()->ts_t);
425 if (!this->feScalingFun.empty())
426 this->locMat *= this->feScalingFun(this->getFEMethod());
430 this->locMatTranspose.resize(this->locMat.size2(), this->locMat.size1(),
432 noalias(this->locMatTranspose) = trans(this->locMat);
433 CHKERR MatSetValues<AssemblyTypeSelector<A>>(
434 matSelector(), col_data, row_data, this->locMatTranspose, ADD_VALUES);
437 if (!this->onlyTranspose) {
439 CHKERR MatSetValues<AssemblyTypeSelector<A>>(
440 matSelector(), row_data, col_data, this->locMat, ADD_VALUES);
446template <AssemblyType A,
typename EleOp>
448 if (!this->timeScalingFun.empty())
449 this->locF *= this->timeScalingFun(this->getFEMethod()->ts_t);
450 if (!this->feScalingFun.empty())
451 this->locF *= this->feScalingFun(this->getFEMethod());
453 return VecSetValues<AssemblyTypeSelector<A>>(this->getKSPf(), data,
454 this->locF, 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
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
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
void setAssembleTo(AssembleTo a_to)
Where to assemble.
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 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.
auto matSelector()
Select matrix.
MatrixDouble locMat
local entity block matrix
FTensor::Tensor2< FTensor::PackPtr< double *, DIM >, DIM, DIM > getLocMat(const int rr)
FEFun feScalingFun
assumes that time variable is set
enum AssembleTo assembleTo
FTensor::Tensor1< FTensor::PackPtr< double *, DIM >, DIM > getNf()
virtual MoFEMErrorCode aSsemble(EntData &data)
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.
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.