7#ifndef __LINEAR_FORMS_INTEGRATORS_IMPL_HPP__
8#define __LINEAR_FORMS_INTEGRATORS_IMPL_HPP__
13 template <
typename OpBase>
struct S {
S() =
delete; };
18 template <
typename OpBase>
struct S {
S() =
delete; };
22template <
int BASE_DIM,
int FIELD_DIM, IntegrationType I,
typename OpBase>
30template <
typename OpBase>
42 boost::shared_ptr<Range> ents_ptr =
nullptr)
44 sourceFun(source_fun) {}
56 boost::shared_ptr<Range> ents_ptr =
nullptr)
58 sourceFun(source_fun) {}
65template <
int FIELD_DIM,
typename OpBase>
79 boost::shared_ptr<Range> ents_ptr =
nullptr)
90 boost::shared_ptr<Range> ents_ptr =
nullptr)
92 sourceFun(source_fun) {}
99template <
int FIELD_DIM,
typename OpBase>
105 boost::shared_ptr<Range> ents_ptr =
nullptr)
107 sourceFun(source_fun) {}
110 boost::shared_ptr<Range> ents_ptr =
nullptr)
112 sourceFun(source_fun) {}
119template <
int BASE_DIM,
int S, IntegrationType I,
typename OpBase>
122template <
int S,
typename OpBase>
126 const std::string
field_name, boost::shared_ptr<VectorDouble> vec,
127 ScalarFun beta_coeff = [](
double,
double,
double)
constexpr {
return 1; },
128 boost::shared_ptr<Range> ents_ptr =
nullptr)
130 betaCoeff(beta_coeff) {}
142template <
int FIELD_DIM,
int S,
typename OpBase>
146 const std::string
field_name, boost::shared_ptr<MatrixDouble> vec,
147 ScalarFun beta_coeff = [](
double,
double,
double)
constexpr {
return 1; },
148 boost::shared_ptr<Range> ents_ptr =
nullptr)
150 betaCoeff(beta_coeff) {}
159template <
int FIELD_DIM,
int S,
typename OpBase>
163 const std::string
field_name, boost::shared_ptr<MatrixDouble> vec,
164 ScalarFun beta_coeff = [](
double,
double,
double)
constexpr {
return 1; },
165 boost::shared_ptr<Range> ents_ptr =
nullptr)
167 betaCoeff(beta_coeff) {}
180template <
int SPACE_DIM,
int S,
typename OpBase>
187 const std::string
field_name, boost::shared_ptr<MatrixDouble> mat_vals,
188 ScalarFun beta_coeff = [](
double,
double,
double)
constexpr {
return 1; },
189 boost::shared_ptr<Range> ents_ptr =
nullptr)
191 matVals(mat_vals), betaCoeff(beta_coeff) {}
199template <
int SPACE_DIM,
int S,
typename OpBase>
207 const std::string
field_name, boost::shared_ptr<MatrixDouble> mat_vals,
208 ScalarFun beta_coeff = [](
double,
double,
double)
constexpr {
return 1; },
209 boost::shared_ptr<Range> ents_ptr =
nullptr)
211 matVals(mat_vals), betaCoeff(beta_coeff) {}
223template <
int SPACE_DIM,
int S,
typename OpBase>
228 const std::string
field_name, boost::shared_ptr<MatrixDouble> mat_vals,
229 ScalarFun beta_coeff = [](
double,
double,
double)
constexpr {
return 1; })
231 betaCoeff(beta_coeff) {}
250 const std::string
field_name, boost::shared_ptr<MatrixDouble> mat_vals,
251 ScalarFun beta = [](
double,
double,
double) {
return 1; },
252 boost::shared_ptr<Range> ents_ptr =
nullptr)
254 matVals(mat_vals), betaConst(beta) {}
264template <
int SPACE_DIM,
typename OpBase, CoordinateTypes CoordSys>
268 const std::string
field_name, boost::shared_ptr<VectorDouble> vec_vals,
269 ScalarFun beta = [](
double,
double,
double)
constexpr {
return 1; },
270 boost::shared_ptr<Range> ents_ptr =
nullptr)
272 vecVals(vec_vals), betaConst(beta) {}
281template <
int FIELD_DIM,
typename OpBase, CoordinateTypes CoordSys>
286 const std::string
field_name, boost::shared_ptr<VectorDouble> vec,
287 ScalarFun beta = [](
double,
double,
double)
constexpr {
return 1; },
288 boost::shared_ptr<Range> ents_ptr =
nullptr)
309template <
int SPACE_DIM, IntegrationType I,
typename OpBase>
312template <
int SPACE_DIM,
typename OpBase>
315 boost::shared_ptr<MatrixDouble> mat_vals)
319 boost::shared_ptr<MatrixDouble> mat_vals,
322 betaCoeff(beta_fun) {}
343template <
int SPACE_DIM, IntegrationType I,
typename OpBase>
346template <
int SPACE_DIM,
typename OpBase>
349 boost::shared_ptr<MatrixDouble> mat_vals)
353 boost::shared_ptr<MatrixDouble> mat_vals,
356 betaCoeff(beta_fun) {}
377template <
int SPACE_DIM, IntegrationType I,
typename OpBase>
386template <
int FIELD_DIM, IntegrationType I,
typename OpBase>
391 OpBase>::OpNormalMixVecTimesScalarImpl;
394template <
typename OpBase>
398 ScalarFun source_fun = [](
double,
double,
double)
constexpr {
return 1; },
399 boost::shared_ptr<Range> ents_ptr =
nullptr)
401 sourceFun(source_fun) {}
409template <
typename OpBase>
413 ScalarFun source_fun = [](
double,
double,
double)
constexpr {
return 1; },
414 boost::shared_ptr<Range> ents_ptr =
nullptr)
416 sourceFun(source_fun) {}
434template <
int SPACE_DIM, IntegrationType I,
typename OpBase>
437template <
int SPACE_DIM,
typename OpBase>
441 const std::string
field_name, boost::shared_ptr<MatrixDouble> u_ptr,
442 ScalarFun beta_coeff = [](
double,
double,
double)
constexpr {
return 1; },
443 boost::shared_ptr<Range> ents_ptr =
nullptr)
445 betaCoeff(beta_coeff) {}
448 boost::shared_ptr<MatrixDouble>
uPtr;
457template <
int SPACE_DIM,
typename OpBase>
460 const std::string
field_name, boost::shared_ptr<MatrixDouble> u_ptr,
461 boost::shared_ptr<MatrixDouble> y_grad_ptr,
462 ConstantFun source_fun = []()
constexpr {
return 1; })
464 yGradPtr(y_grad_ptr), alphaConstant(source_fun) {}
467 boost::shared_ptr<MatrixDouble>
uPtr;
473template <
int FIELD_DIM,
int SPACE_DIM,
typename OpBase>
477 const std::string
field_name, boost::shared_ptr<MatrixDouble> u_ptr,
478 boost::shared_ptr<MatrixDouble> y_grad_ptr,
479 ConstantFun source_fun = []()
constexpr {
return 1; })
481 yGradPtr(y_grad_ptr), alphaConstant(source_fun) {}
484 boost::shared_ptr<MatrixDouble>
uPtr;
490template <
typename OpBase>
497 const double vol = OpBase::getMeasure();
499 auto t_w = OpBase::getFTensor0IntegrationWeight();
503 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
508 t_w * vol * sourceFun(t_coords(0), t_coords(1), t_coords(2));
523template <
int FIELD_DIM,
typename OpBase>
531 const double vol = OpBase::getMeasure();
533 auto t_w = OpBase::getFTensor0IntegrationWeight();
537 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
541 auto t_source = sourceFun(t_coords(0), t_coords(1), t_coords(2));
543 const double alpha = t_w * vol;
548 t_nf(
i) += alpha * t_row_base * t_source(
i);
560template <
int FIELD_DIM,
typename OpBase>
567 const size_t nb_base_functions = row_data.
getN().size2() / 3;
569 const double vol = OpBase::getMeasure();
571 auto t_w = OpBase::getFTensor0IntegrationWeight();
575 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
579 auto t_source = sourceFun(t_coords(0), t_coords(1), t_coords(2));
581 const double alpha = t_w * vol;
588 for (; rr < nb_base_functions; ++rr)
596template <
int S,
typename OpBase>
602 const double vol = OpBase::getMeasure();
604 auto t_w = OpBase::getFTensor0IntegrationWeight();
613 "Wrong number of integration points %d != %ld",
619 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
624 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
640template <
int FIELD_DIM,
int S,
typename OpBase>
646 const double vol = OpBase::getMeasure();
648 auto t_w = OpBase::getFTensor0IntegrationWeight();
652 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
659 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
661 auto t_nf = OpBase::template getNf<FIELD_DIM>();
665 t_nf(
i) += alpha * t_row_base * t_vec(
i);
678template <
int FIELD_DIM,
int S,
typename OpBase>
683 const size_t nb_base_functions = row_data.
getN().size2() / 3;
685 const double vol = OpBase::getMeasure();
687 auto t_w = OpBase::getFTensor0IntegrationWeight();
691 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
698 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
705 for (; rr < nb_base_functions; ++rr)
714template <
int SPACE_DIM,
int S,
typename OpBase>
721 const double vol = OpBase::getMeasure();
723 auto t_w = OpBase::getFTensor0IntegrationWeight();
729 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
732 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
734 const double alpha = t_w * beta;
753template <
int SPACE_DIM,
int S,
typename OpBase>
760 const double vol = OpBase::getMeasure();
762 auto t_w = OpBase::getFTensor0IntegrationWeight();
768 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
773 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
775 auto t_nf = OpBase::template getNf<SPACE_DIM>();
780 t_nf(
i) += alpha * (t_row_grad(
j) * t_val_grad(
i,
j));
795template <
int SPACE_DIM,
int S,
typename OpBase>
801 const double vol = OpBase::getMeasure();
803 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
805 auto t_w = OpBase::getFTensor0IntegrationWeight();
814 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
816 auto t_nf = OpBase::template getNf<SPACE_DIM>();
821 t_nf(
j) += alpha * (t_row_grad(
i) * t_val_mat(
i,
j));
842 const size_t nb_base_functions = row_data.
getN().size2() / 3;
843 auto t_w = this->getFTensor0IntegrationWeight();
845 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
852 const double alpha = this->getMeasure() * t_w *
853 betaConst(t_coords(0), t_coords(1), t_coords(2));
854 auto t_nf = OpBase::template getNf<FIELD_DIM>();
857 for (; bb != this->nbRows /
FIELD_DIM; ++bb) {
858 const double t_div_base = t_diff_base(
j,
j);
859 t_nf(
i) += alpha * t_div_base * t_u(
i);
861 t_nf(
i) += t_base(0) * (alpha / t_coords(0)) * t_u(
i);
867 for (; bb < nb_base_functions; ++bb) {
880template <
int SPACE_DIM,
typename OpBase, CoordinateTypes CoordSys>
886 const size_t nb_base_functions = row_data.
getN().size2() / 3;
887 auto t_w = this->getFTensor0IntegrationWeight();
889 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
895 const double alpha = this->getMeasure() * t_w *
896 betaConst(t_coords(0), t_coords(1), t_coords(2));
900 for (; bb != this->nbRows; ++bb) {
901 const double t_div_base = t_diff_base(
j,
j);
905 for (; bb < nb_base_functions; ++bb)
936template <
int FIELD_DIM,
typename OpBase, CoordinateTypes CoordSys>
944 const double vol = OpBase::getMeasure();
946 auto t_w = OpBase::getFTensor0IntegrationWeight();
948 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
957 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
962 t_nf(
i) += alpha * t_row_grad(
i) * t_vec;
975template <
int SPACE_DIM,
typename OpBase>
980 const size_t nb_base_functions = row_data.
getN().size2() / 3;
981 auto t_w = this->getFTensor0IntegrationWeight();
982 auto t_coords = this->getFTensor1CoordsAtGaussPts();
988 const double alpha = this->getMeasure() * t_w;
989 auto t_nf = OpBase::template getNf<SPACE_DIM>();
992 for (; bb != this->nbRows /
SPACE_DIM; ++bb) {
993 t_nf(
i) += alpha * betaCoeff(t_coords(0), t_coords(1), t_coords(2)) *
994 t_base(
j) * t_grad(
i,
j);
998 for (; bb < nb_base_functions; ++bb)
1009template <
int SPACE_DIM,
typename OpBase>
1014 const size_t nb_base_functions = row_data.
getN().size2();
1015 auto t_w = this->getFTensor0IntegrationWeight();
1016 auto t_coords = this->getFTensor1CoordsAtGaussPts();
1021 const double alpha = this->getMeasure() * t_w;
1022 auto t_nf = OpBase::template getNf<SPACE_DIM>();
1025 for (; bb != this->nbRows /
SPACE_DIM; ++bb) {
1026 t_nf(
i) += alpha * t_base *
1027 betaCoeff(t_coords(0), t_coords(1), t_coords(2)) * t_div(
i);
1031 for (; bb < nb_base_functions; ++bb)
1042template <
typename OpBase>
1047 const size_t nb_base_functions = row_data.
getN().size2() / 3;
1050 auto t_w = OpBase::getFTensor0IntegrationWeight();
1054 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1056 auto t_normal = OpBase::getFTensor1NormalsAtGaussPts();
1058 if (this->getNumeredEntFiniteElementPtr()->getEntType() == MBTRI)
1063 const double alpha =
1064 t_w * sourceFun(t_coords(0), t_coords(1), t_coords(2)) /
a;
1071 for (; rr < nb_base_functions; ++rr)
1080template <
typename OpBase>
1085 const size_t nb_base_functions = row_data.
getN().size2() / 3;
1089 auto t_w = OpBase::getFTensor0IntegrationWeight();
1093 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1095 auto t_tangent = OpBase::getFTensor1TangentAtGaussPts();
1099 const double alpha = t_w * sourceFun(t_coords(0), t_coords(1), t_coords(2));
1109 for (; rr < nb_base_functions; ++rr)
1118template <
int SPACE_DIM,
typename OpBase>
1126 const size_t nb_base_functions = row_data.
getN().size2() / 3;
1129 auto t_w = OpBase::getFTensor0IntegrationWeight();
1133 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1135 auto t_normal = OpBase::getFTensor1NormalsAtGaussPts();
1139 auto a = OpBase::getMeasure();
1142 auto l2 = std::sqrt(t_normal(
i) * t_normal(
i));
1143 const double alpha =
1144 t_w * betaCoeff(t_coords(0), t_coords(1), t_coords(2)) * (
a / l2);
1146 auto t_nf = OpBase::template getNf<SPACE_DIM>();
1150 t_nf(
J) += alpha * (t_row_base(
i) * t_normal(
i)) * t_u(
J);
1154 for (; rr < nb_base_functions; ++rr)
1164template <
int SPACE_DIM,
typename OpBase>
1170 auto t_w = this->getFTensor0IntegrationWeight();
1177 const double alpha_constant = alphaConstant();
1181 const double vol = OpBase::getMeasure();
1182 const double c = (t_grad_y(
i) * t_u(
i)) * (t_w * vol * alpha_constant);
1201template <
int FIELD_DIM,
int SPACE_DIM,
typename OpBase>
1207 auto t_w = this->getFTensor0IntegrationWeight();
1215 const double alpha_constant = alphaConstant();
1219 const double vol = OpBase::getMeasure();
1222 t_c(
J) = (t_grad_y(
J,
i) * t_u(
i)) * (t_w * vol * alpha_constant);
1227 t_nf(
J) += t_base * t_c(
J);
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_DATA_INCONSISTENCY
CoordinateTypes
Coodinate system.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
FTensor::Index< 'i', SPACE_DIM > i
const double c
speed of light (cm/ns)
FTensor::Index< 'J', DIM1 > J
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
boost::function< double(double)> TimeFun
Lambda function used to scale with time.
FTensor::Tensor1< FTensor::PackPtr< T *, S >, Tensor_Dim > getFTensor1FromMat(ublas::matrix< T, L, A > &data)
Get tensor rank 1 (vector) form data matrix.
FTensor::Tensor2< FTensor::PackPtr< double *, 1 >, Tensor_Dim1, Tensor_Dim2 > getFTensor2FromMat(MatrixDouble &data)
Get tensor rank 2 (matrix) form data matrix.
static auto getFTensor2SymmetricFromMat(ublas::matrix< T, L, A > &data)
Get symmetric tensor rank 2 (matrix) form data matrix.
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
boost::function< double()> ConstantFun
Constant function type.
auto getFTensor1FromArray(VectorDouble &data)
Get FTensor1 from array.
constexpr IntegrationType I
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
OpBaseTimesScalarImpl(const std::string field_name, boost::shared_ptr< VectorDouble > vec, ScalarFun beta_coeff=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
boost::shared_ptr< VectorDouble > sourceVec
boost::shared_ptr< MatrixDouble > sourceVec
OpBaseTimesVectorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > vec, ScalarFun beta_coeff=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
FTensor::Index< 'i', FIELD_DIM > i
boost::shared_ptr< MatrixDouble > sourceVec
FTensor::Index< 'i', FIELD_DIM > i
OpBaseTimesVectorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > vec, ScalarFun beta_coeff=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
boost::shared_ptr< MatrixDouble > yGradPtr
boost::shared_ptr< MatrixDouble > uPtr
ConstantFun alphaConstant
OpConvectiveTermRhsImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > u_ptr, boost::shared_ptr< MatrixDouble > y_grad_ptr, ConstantFun source_fun=[]() constexpr { return 1;})
boost::shared_ptr< MatrixDouble > uPtr
ConstantFun alphaConstant
OpConvectiveTermRhsImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > u_ptr, boost::shared_ptr< MatrixDouble > y_grad_ptr, ConstantFun source_fun=[]() constexpr { return 1;})
boost::shared_ptr< MatrixDouble > yGradPtr
OpGradTimesSymTensorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals, ScalarFun beta_coeff=[](double, double, double) constexpr { return 1;})
boost::shared_ptr< MatrixDouble > matVals
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'j', SPACE_DIM > j
OpGradTimesTensorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals, ScalarFun beta_coeff=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
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=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
boost::shared_ptr< MatrixDouble > matVals
FTensor::Index< 'j', SPACE_DIM > j
summit Index
FTensor::Index< 'i', SPACE_DIM > i
summit Index
boost::shared_ptr< VectorDouble > sourceVec
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &data)
Class dedicated to integrate operator.
OpMixDivTimesUImpl(const std::string field_name, boost::shared_ptr< VectorDouble > vec, ScalarFun beta=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
boost::shared_ptr< VectorDouble > vecVals
OpMixDivTimesUImpl(const std::string field_name, boost::shared_ptr< VectorDouble > vec_vals, ScalarFun beta=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
FTensor::Index< 'j', SPACE_DIM > j
FTensor::Index< 'i', FIELD_DIM > i
FTensor::Index< 'j', SPACE_DIM > j
boost::shared_ptr< MatrixDouble > matVals
OpMixDivTimesUImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals, ScalarFun beta=[](double, double, double) { return 1;}, 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)
OpMixTensorTimesGradUImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals, ScalarFun beta_fun)
Tensor field time gradient of vector field.
OpMixVecTimesDivLambdaImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals, ScalarFun beta_fun)
OpMixVecTimesDivLambdaImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals)
boost::shared_ptr< MatrixDouble > matVals
FTensor::Index< 'i', SPACE_DIM > i
Vector field time divergence of tensor.
OpNormalMixVecTimesScalarImpl(const std::string field_name, ScalarFun source_fun=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
FTensor::Index< 'i', 3 > i
FTensor::Index< 'i', 3 > i
OpNormalMixVecTimesScalarImpl(const std::string field_name, ScalarFun source_fun=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
Multiply vector times normal on the face times scalar function.
OpNormalMixVecTimesVectorFieldImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > u_ptr, ScalarFun beta_coeff=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
boost::shared_ptr< MatrixDouble > uPtr
Multiply vector times normal on the face times vector field.
OpSourceImpl(const std::string field_name, ScalarFun source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)
Construct a new Op Source Impl object.
OpSourceImpl(const std::string field_name, TimeFun time_fun, ScalarFun source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)
Construct a new Op Source Impl object.
VectorFun< FIELD_DIM > sourceFun
OpSourceImpl(const std::string field_name, VectorFun< FIELD_DIM > source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)
Construct a new Op Source Impl object.
OpSourceImpl(const std::string field_name, TimeFun time_fun, VectorFun< FIELD_DIM > source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)
Construct a new Op Source Impl object.
OpSourceImpl(const std::string field_name, TimeFun time_fun, VectorFun< FIELD_DIM > source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)
OpSourceImpl(const std::string field_name, VectorFun< FIELD_DIM > source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)
VectorFun< FIELD_DIM > sourceFun
SourceBoundaryNormalSpecialization()=delete
SourceFunctionSpecialization()=delete