9#ifndef __LINEAR_FORMS_INTEGRATORS_HPP__
10#define __LINEAR_FORMS_INTEGRATORS_HPP__
14template <
int BASE_DIM,
int FIELD_DIM, IntegrationType I,
typename OpBase>
22template <
typename OpBase>
33 boost::shared_ptr<Range> ents_ptr =
nullptr)
43template <
int FIELD_DIM,
typename OpBase>
46 boost::shared_ptr<Range> ents_ptr =
nullptr)
56template <
int BASE_DIM,
typename OpBase>
59 boost::shared_ptr<Range> ents_ptr =
nullptr)
69template <
int BASE_DIM,
int S, IntegrationType I,
typename OpBase>
72template <
int S,
typename OpBase>
76 boost::shared_ptr<VectorDouble> vec,
78 boost::shared_ptr<Range> ents_ptr =
nullptr)
80 betaCoeff(beta_coeff), entsPtr(ents_ptr) {}
93template <
int FIELD_DIM,
int S,
typename OpBase>
97 boost::shared_ptr<MatrixDouble> vec,
99 boost::shared_ptr<Range> ents_ptr =
nullptr)
101 betaCoeff(beta_coeff), entsPtr(ents_ptr) {}
111template <
int BASE_DIM,
int S,
typename OpBase>
116 boost::shared_ptr<MatrixDouble> vec,
118 boost::shared_ptr<Range> ents_ptr =
nullptr)
120 betaCoeff(beta_coeff), entsPtr(ents_ptr) {}
134template <
int SPACE_DIM,
int S,
typename OpBase>
141 boost::shared_ptr<MatrixDouble> mat_vals,
143 boost::shared_ptr<Range> ents_ptr =
nullptr)
145 betaCoeff(beta_coeff) {}
154template <
int SPACE_DIM,
int S,
typename OpBase>
162 const std::string
field_name, boost::shared_ptr<MatrixDouble> mat_vals,
163 ScalarFun beta_coeff = [](
double,
double,
double) {
return 1; },
164 boost::shared_ptr<Range> ents_ptr =
nullptr)
166 betaCoeff(beta_coeff) {}
179template <
int SPACE_DIM,
int S,
typename OpBase>
184 const std::string
field_name, boost::shared_ptr<MatrixDouble> &mat_vals,
185 ScalarFun beta_coeff = [](
double,
double,
double) {
return 1; })
187 betaCoeff(beta_coeff) {}
201template <
int FIELD_DIM,
int SPACE_DIM,
typename OpBase>
205 boost::shared_ptr<MatrixDouble> mat_vals,
ScalarFun beta,
206 boost::shared_ptr<Range> ents_ptr =
nullptr)
208 betaConst(beta), entsPtr(ents_ptr) {}
219template <
int SPACE_DIM,
typename OpBase>
222 boost::shared_ptr<VectorDouble> vec_vals,
ScalarFun beta,
223 boost::shared_ptr<Range> ents_ptr =
nullptr)
225 betaConst(beta), entsPtr(ents_ptr) {}
235template <
int FIELD_DIM,
typename OpBase>
240 boost::shared_ptr<VectorDouble> vec,
ScalarFun beta,
241 boost::shared_ptr<Range> ents_ptr =
nullptr)
243 betaCoeff(beta), entsPtr(ents_ptr) {}
252template <
int SPACE_DIM, IntegrationType I,
typename OpBase>
255template <
int SPACE_DIM,
typename OpBase>
258 boost::shared_ptr<MatrixDouble> &mat_vals)
267template <
int SPACE_DIM, IntegrationType I,
typename OpBase>
270template <
int SPACE_DIM,
typename OpBase>
273 boost::shared_ptr<MatrixDouble> &mat_vals)
293template <
int SPACE_DIM, IntegrationType I,
typename OpBase>
296template <
typename OpBase>
300 boost::shared_ptr<Range> ents_ptr =
nullptr)
311template <
typename OpBase>
315 boost::shared_ptr<Range> ents_ptr =
nullptr)
330template <
int SPACE_DIM,
typename OpBase>
333 const std::string
field_name, boost::shared_ptr<MatrixDouble> u_ptr,
334 boost::shared_ptr<MatrixDouble> y_grad_ptr,
337 yGradPtr(y_grad_ptr), alphaConstant(source_fun) {}
340 boost::shared_ptr<MatrixDouble>
uPtr;
346template <
int FIELD_DIM,
int SPACE_DIM,
typename OpBase>
350 const std::string
field_name, boost::shared_ptr<MatrixDouble> u_ptr,
351 boost::shared_ptr<MatrixDouble> y_grad_ptr,
354 yGradPtr(y_grad_ptr), alphaConstant(source_fun) {}
357 boost::shared_ptr<MatrixDouble>
uPtr;
371template <
typename EleOp>
372template <AssemblyType A>
373template <IntegrationType I>
385 template <
int BASE_DIM,
int FIELD_DIM>
386 struct OpSource :
public OpSourceImpl<BASE_DIM, FIELD_DIM, I, OpBase> {
396 template <
int BASE_DIM,
int S = 1>
398 :
public OpBaseTimesScalarFieldImpl<BASE_DIM, S, I, OpBase> {
399 using OpBaseTimesScalarFieldImpl<
BASE_DIM, S,
I,
400 OpBase>::OpBaseTimesScalarFieldImpl;
411 template <
int BASE_DIM,
int FIELD_DIM,
int S>
412 struct OpBaseTimesVector
413 :
public OpBaseTimesVectorImpl<BASE_DIM, FIELD_DIM, S, I, OpBase> {
415 OpBase>::OpBaseTimesVectorImpl;
432 template <
int BASE_DIM,
int FIELD_DIM,
int SPACE_DIM,
int S = 1>
433 struct OpGradTimesTensor
437 OpBase>::OpGradTimesTensorImpl;
450 template <
int BASE_DIM,
int FIELD_DIM,
int SPACE_DIM,
int S = 1>
451 struct OpGradTimesSymTensor
455 OpBase>::OpGradTimesSymTensorImpl;
465 template <
int BASE_DIM,
int FIELD_DIM,
int SPACE_DIM>
466 struct OpMixDivTimesU
469 OpBase>::OpMixDivTimesUImpl;
477 template <
int SPACE_DIM>
478 struct OpMixTensorTimesGradU
481 OpBase>::OpMixTensorTimesGradUImpl;
489 template <
int SPACE_DIM>
490 struct OpMixVecTimesDivLambda
491 :
public OpMixVecTimesDivLambdaImpl<SPACE_DIM, I, OpBase> {
492 using OpMixVecTimesDivLambdaImpl<
SPACE_DIM,
I,
493 OpBase>::OpMixVecTimesDivLambdaImpl;
506 template <
int SPACE_DIM>
507 struct OpNormalMixVecTimesScalar
508 :
public OpNormalMixVecTimesScalarImpl<SPACE_DIM, I, OpBase> {
509 using OpNormalMixVecTimesScalarImpl<
SPACE_DIM,
I,
510 OpBase>::OpNormalMixVecTimesScalarImpl;
525 template <
int BASE_DIM,
int FIELD_DIM,
int SPACE_DIM>
526 struct OpConvectiveTermRhs
530 OpBase>::OpConvectiveTermRhsImpl;
534template <
typename OpBase>
540 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
545 const double vol = OpBase::getMeasure();
547 auto t_w = OpBase::getFTensor0IntegrationWeight();
551 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
556 t_w * vol * sourceFun(t_coords(0), t_coords(1), t_coords(2));
571template <
int FIELD_DIM,
typename OpBase>
577 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
582 const double vol = OpBase::getMeasure();
584 auto t_w = OpBase::getFTensor0IntegrationWeight();
588 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
592 auto t_source = sourceFun(t_coords(0), t_coords(1), t_coords(2));
594 const double alpha = t_w * vol;
596 auto t_nf = getFTensor1FromArray<FIELD_DIM, FIELD_DIM>(
OpBase::locF);
599 t_nf(
i) +=
alpha * t_row_base * t_source(
i);
611template <
int BASE_DIM,
typename OpBase>
617 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
621 const size_t nb_base_functions = row_data.
getN().size2() /
BASE_DIM;
623 const double vol = OpBase::getMeasure();
625 auto t_w = OpBase::getFTensor0IntegrationWeight();
629 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
633 auto t_source = sourceFun(t_coords(0), t_coords(1), t_coords(2));
635 const double alpha = t_w * vol;
642 for (; rr < nb_base_functions; ++rr)
650template <
int S,
typename OpBase>
655 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
659 const double vol = OpBase::getMeasure();
661 auto t_w = OpBase::getFTensor0IntegrationWeight();
665 auto t_vec = getFTensor0FromVec<S>(*sourceVec);
670 "Wrong number of integration points %d != %d",
676 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
681 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
697template <
int FIELD_DIM,
int S,
typename OpBase>
702 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
706 const double vol = OpBase::getMeasure();
708 auto t_w = OpBase::getFTensor0IntegrationWeight();
712 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
714 auto t_vec = getFTensor1FromMat<FIELD_DIM, S>(*sourceVec);
719 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
721 auto t_nf = OpBase::template getNf<FIELD_DIM>();
725 t_nf(
i) +=
alpha * t_row_base * t_vec(
i);
738template <
int BASE_DIM,
int S,
typename OpBase>
744 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
747 const size_t nb_base_functions = row_data.
getN().size2() /
BASE_DIM;
749 const double vol = OpBase::getMeasure();
751 auto t_w = OpBase::getFTensor0IntegrationWeight();
755 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
757 auto t_vec = getFTensor1FromMat<BASE_DIM, S>(*sourceVec);
762 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
769 for (; rr < nb_base_functions; ++rr)
778template <
int SPACE_DIM,
int S,
typename OpBase>
784 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
788 const double vol = OpBase::getMeasure();
790 auto t_w = OpBase::getFTensor0IntegrationWeight();
794 auto t_val_grad = getFTensor1FromMat<SPACE_DIM, S>(*(matVals));
796 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
799 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
801 const double alpha = t_w * beta;
820template <
int SPACE_DIM,
int S,
typename OpBase>
826 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
830 const double vol = OpBase::getMeasure();
832 auto t_w = OpBase::getFTensor0IntegrationWeight();
836 auto t_val_grad = getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(*(matVals));
838 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
843 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
845 auto t_nf = OpBase::template getNf<SPACE_DIM>();
850 t_nf(
i) +=
alpha * (t_row_grad(
j) * t_val_grad(
i,
j));
865template <
int SPACE_DIM,
int S,
typename OpBase>
871 const double vol = OpBase::getMeasure();
873 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
875 auto t_w = OpBase::getFTensor0IntegrationWeight();
879 auto t_val_mat = getFTensor2SymmetricFromMat<SPACE_DIM, S>(*(matVals));
884 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
886 auto t_nf = OpBase::template getNf<SPACE_DIM>();
891 t_nf(
j) +=
alpha * (t_row_grad(
i) * t_val_mat(
i,
j));
905template <
int FIELD_DIM,
int SPACE_DIM,
typename OpBase>
912 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
916 const size_t nb_base_functions = row_data.
getN().size2() / 3;
917 auto t_w = this->getFTensor0IntegrationWeight();
919 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
921 auto t_u = getFTensor1FromMat<FIELD_DIM>(*(matVals));
925 const double alpha = this->getMeasure() * t_w *
926 betaConst(t_coords(0), t_coords(1), t_coords(2));
927 auto t_nf = OpBase::template getNf<FIELD_DIM>();
930 for (; bb != this->nbRows /
FIELD_DIM; ++bb) {
931 const double t_div_base = t_diff_base(
j,
j);
932 t_nf(
i) +=
alpha * t_div_base * t_u(
i);
936 for (; bb < nb_base_functions; ++bb)
947template <
int SPACE_DIM,
typename OpBase>
953 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
957 const size_t nb_base_functions = row_data.
getN().size2() / 3;
958 auto t_w = this->getFTensor0IntegrationWeight();
960 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
966 const double alpha = this->getMeasure() * t_w *
967 betaConst(t_coords(0), t_coords(1), t_coords(2));
971 for (; bb != this->nbRows; ++bb) {
972 const double t_div_base = t_diff_base(
j,
j);
976 for (; bb < nb_base_functions; ++bb)
1007template <
int FIELD_DIM,
typename OpBase>
1015 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
1020 const double vol = OpBase::getMeasure();
1022 auto t_w = OpBase::getFTensor0IntegrationWeight();
1024 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1032 const double alpha =
1033 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
1034 auto t_nf = getFTensor1FromArray<FIELD_DIM, FIELD_DIM>(
OpBase::locF);
1038 t_nf(
i) +=
alpha * t_row_grad(
i) * t_vec;
1051template <
int SPACE_DIM,
typename OpBase>
1056 const size_t nb_base_functions = row_data.
getN().size2() / 3;
1057 auto t_w = this->getFTensor0IntegrationWeight();
1059 auto t_grad = getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(*(matVals));
1063 const double alpha = this->getMeasure() * t_w;
1064 auto t_nf = OpBase::template getNf<SPACE_DIM>();
1067 for (; bb != this->nbRows /
SPACE_DIM; ++bb) {
1068 t_nf(
i) +=
alpha * t_base(
j) * t_grad(
i,
j);
1072 for (; bb < nb_base_functions; ++bb)
1082template <
int SPACE_DIM,
typename OpBase>
1087 const size_t nb_base_functions = row_data.
getN().size2();
1088 auto t_w = this->getFTensor0IntegrationWeight();
1090 auto t_div = getFTensor1FromMat<SPACE_DIM>(*(matVals));
1093 const double alpha = this->getMeasure() * t_w;
1094 auto t_nf = OpBase::template getNf<SPACE_DIM>();
1097 for (; bb != this->nbRows /
SPACE_DIM; ++bb) {
1098 t_nf(
i) +=
alpha * t_base * t_div(
i);
1102 for (; bb < nb_base_functions; ++bb)
1112template <
typename OpBase>
1117 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
1120 const size_t nb_base_functions = row_data.
getN().size2() / 3;
1122 const double vol = OpBase::getMeasure();
1124 auto t_w = OpBase::getFTensor0IntegrationWeight();
1128 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1130 auto t_normal = OpBase::getFTensor1NormalsAtGaussPts();
1134 const double alpha = t_w * sourceFun(t_coords(0), t_coords(1), t_coords(2));
1141 for (; rr < nb_base_functions; ++rr)
1150template <
typename OpBase>
1155 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
1158 const size_t nb_base_functions = row_data.
getN().size2() / 3;
1161 const double vol = OpBase::getMeasure();
1163 auto t_w = OpBase::getFTensor0IntegrationWeight();
1167 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1169 auto t_tangent = OpBase::getFTensor1TangentAtGaussPts();
1173 const double alpha = t_w * sourceFun(t_coords(0), t_coords(1), t_coords(2));
1183 for (; rr < nb_base_functions; ++rr)
1192template <
int SPACE_DIM,
typename OpBase>
1198 const size_t nb_base_functions = row_data.
getN().size2();
1199 auto t_w = this->getFTensor0IntegrationWeight();
1202 auto t_u = getFTensor1FromMat<SPACE_DIM>(*uPtr);
1203 auto t_grad_y = getFTensor1FromMat<SPACE_DIM>(*yGradPtr);
1206 const double alpha_constant = alphaConstant();
1210 const double vol = OpBase::getMeasure();
1211 const double c = (t_grad_y(
i) * t_u(
i)) * (t_w * vol * alpha_constant);
1230template <
int FIELD_DIM,
int SPACE_DIM,
typename OpBase>
1236 const size_t nb_base_functions = row_data.
getN().size2();
1237 auto t_w = this->getFTensor0IntegrationWeight();
1240 auto t_u = getFTensor1FromMat<SPACE_DIM>(*uPtr);
1241 auto t_grad_y = getFTensor2FromMat<FIELD_DIM, SPACE_DIM>(*yGradPtr);
1245 const double alpha_constant = alphaConstant();
1249 const double vol = OpBase::getMeasure();
1252 t_c(
J) = (t_grad_y(
J,
i) * t_u(
i)) * (t_w * vol * alpha_constant);
1254 auto t_nf = getFTensor1FromArray<FIELD_DIM, FIELD_DIM>(
OpBase::locF);
1257 t_nf(
J) += t_base * t_c(
J);
#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 ...
@ MOFEM_DATA_INCONSISTENCY
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesScalarField< 1, 1 > OpBaseTimesScalarField
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
implementation of Data Operators for Forces and Sources
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
boost::function< double()> ConstantFun
Constant function type.
constexpr auto field_name
OpBaseImpl< PETSC, EdgeEleOp > OpBase
Data on single entity (This is passed as argument to DataOperator::doWork)
FTensor::Tensor2< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 >, Tensor_Dim0, Tensor_Dim1 > getFTensor2DiffN(FieldApproximationBase base)
Get derivatives of base functions for Hdiv space.
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.
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1N(FieldApproximationBase base)
Get base functions for Hdiv/Hcurl spaces.
VectorDouble locF
local entity vector
int nbRows
number of dofs on rows
int nbIntegrationPts
number of integration points
int nbRowBaseFunctions
number or row base functions
OpBaseTimesScalarFieldImpl(const std::string field_name, boost::shared_ptr< VectorDouble > vec, ScalarFun beta_coeff, boost::shared_ptr< Range > ents_ptr=nullptr)
boost::shared_ptr< VectorDouble > sourceVec
boost::shared_ptr< Range > entsPtr
boost::shared_ptr< MatrixDouble > sourceVec
OpBaseTimesVectorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > vec, ScalarFun beta_coeff, boost::shared_ptr< Range > ents_ptr=nullptr)
FTensor::Index< 'i', FIELD_DIM > i
boost::shared_ptr< Range > entsPtr
boost::shared_ptr< Range > entsPtr
boost::shared_ptr< MatrixDouble > sourceVec
OpBaseTimesVectorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > vec, ScalarFun beta_coeff, boost::shared_ptr< Range > ents_ptr=nullptr)
FTensor::Index< 'i', BASE_DIM > i
boost::shared_ptr< MatrixDouble > yGradPtr
OpConvectiveTermRhsImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > u_ptr, boost::shared_ptr< MatrixDouble > y_grad_ptr, ConstantFun source_fun=[]() { return 1;})
boost::shared_ptr< MatrixDouble > uPtr
ConstantFun alphaConstant
boost::shared_ptr< MatrixDouble > uPtr
OpConvectiveTermRhsImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > u_ptr, boost::shared_ptr< MatrixDouble > y_grad_ptr, ConstantFun source_fun=[]() { return 1;})
ConstantFun alphaConstant
boost::shared_ptr< MatrixDouble > yGradPtr
boost::shared_ptr< MatrixDouble > matVals
OpGradTimesSymTensorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > &mat_vals, ScalarFun beta_coeff=[](double, double, double) { return 1;})
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'j', SPACE_DIM > j
boost::shared_ptr< Range > entsPtr
FTensor::Index< 'i', SPACE_DIM > i
summit Index
boost::shared_ptr< MatrixDouble > matVals
OpGradTimesTensorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals, ScalarFun beta_coeff, boost::shared_ptr< Range > ents_ptr=nullptr)
boost::shared_ptr< MatrixDouble > matVals
FTensor::Index< 'j', SPACE_DIM > j
summit Index
boost::shared_ptr< Range > entsPtr
FTensor::Index< 'i', SPACE_DIM > i
summit Index
OpGradTimesTensorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals, ScalarFun beta_coeff=[](double, double, double) { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
OpMixDivTimesUImpl(const std::string field_name, boost::shared_ptr< VectorDouble > vec, ScalarFun beta, boost::shared_ptr< Range > ents_ptr=nullptr)
boost::shared_ptr< VectorDouble > sourceVec
boost::shared_ptr< Range > entsPtr
OpMixDivTimesUImpl(const std::string field_name, boost::shared_ptr< VectorDouble > vec_vals, ScalarFun beta, boost::shared_ptr< Range > ents_ptr=nullptr)
boost::shared_ptr< Range > entsPtr
FTensor::Index< 'j', SPACE_DIM > j
boost::shared_ptr< VectorDouble > vecVals
FTensor::Index< 'j', SPACE_DIM > j
boost::shared_ptr< MatrixDouble > matVals
FTensor::Index< 'i', FIELD_DIM > i
boost::shared_ptr< Range > entsPtr
OpMixDivTimesUImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals, ScalarFun beta, boost::shared_ptr< Range > ents_ptr=nullptr)
boost::shared_ptr< MatrixDouble > matVals
FTensor::Index< 'j', SPACE_DIM > j
FTensor::Index< 'i', SPACE_DIM > i
OpMixTensorTimesGradUImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > &mat_vals)
boost::shared_ptr< MatrixDouble > matVals
FTensor::Index< 'i', SPACE_DIM > i
OpMixVecTimesDivLambdaImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > &mat_vals)
OpNormalMixVecTimesScalarImpl(const std::string field_name, ScalarFun source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)
FTensor::Index< 'i', 3 > i
boost::shared_ptr< Range > entsPtr
FTensor::Index< 'i', 3 > i
boost::shared_ptr< Range > entsPtr
OpNormalMixVecTimesScalarImpl(const std::string field_name, ScalarFun source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)
Multiply vactor times normal on the face times scalar function.
OpSourceImpl(const std::string field_name, ScalarFun source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)
Construct a new Op Source Impl object.
boost::shared_ptr< Range > entsPtr
VectorFun< FIELD_DIM > sourceFun
OpSourceImpl(const std::string field_name, VectorFun< FIELD_DIM > source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)
boost::shared_ptr< Range > entsPtr
VectorFun< BASE_DIM > sourceFun
boost::shared_ptr< Range > entsPtr
OpSourceImpl(const std::string field_name, VectorFun< BASE_DIM > source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)