13#ifndef __BILINEAR_FORMS_INTEGRATORS_HPP__
14#define __BILINEAR_FORMS_INTEGRATORS_HPP__
22template <
int SPACE_DIM,
typename OpBase>
26 const std::string col_field_name,
28 boost::shared_ptr<Range> ents_ptr =
nullptr)
29 :
OpBase(row_field_name, col_field_name,
OpBase::OPROWCOL, ents_ptr),
31 if (row_field_name == col_field_name)
41template <
int FIELD_DIM,
int SPACE_DIM,
typename OpBase>
45 const std::string col_field_name,
47 boost::shared_ptr<Range> ents_ptr =
nullptr)
48 :
OpBase(row_field_name, col_field_name,
OpBase::OPROWCOL, ents_ptr),
50 if (row_field_name == col_field_name)
60template <
int BASE_DIM,
int FIELD_DIM, IntegrationType I,
typename OpBase>
63template <
typename OpBase>
67 const std::string row_field_name,
const std::string col_field_name,
68 ScalarFun beta = [](
double,
double,
double)
constexpr {
return 1; },
69 boost::shared_ptr<Range> ents_ptr =
nullptr)
70 :
OpBase(row_field_name, col_field_name, OpBase::OPROWCOL, ents_ptr),
72 if (row_field_name == col_field_name)
82template <
int FIELD_DIM,
typename OpBase>
93template <
int FIELD_DIM,
typename OpBase>
96 OpMassImpl(
const std::string row_field_name,
const std::string col_field_name,
98 boost::shared_ptr<Range> ents_ptr =
nullptr)
99 :
OpBase(row_field_name, col_field_name,
OpBase::OPROWCOL, ents_ptr),
101 if (row_field_name == col_field_name)
111template <
typename OpBase>
113 OpMassImpl(
const std::string row_field_name,
const std::string col_field_name,
115 boost::shared_ptr<Range> ents_ptr =
nullptr)
116 :
OpBase(row_field_name, col_field_name,
OpBase::OPROWCOL, ents_ptr),
118 if (row_field_name == col_field_name)
128template <
typename OpBase>
130 OpMassImpl(
const std::string row_field_name,
const std::string col_field_name,
132 boost::shared_ptr<Range> ents_ptr =
nullptr)
133 :
OpBase(row_field_name, col_field_name,
OpBase::OPROWCOL, ents_ptr),
135 if (row_field_name == col_field_name)
149template <
int SPACE_DIM,
int S,
typename OpBase>
154 const std::string col_field_name,
155 boost::shared_ptr<MatrixDouble> mat_D,
156 boost::shared_ptr<Range> ents_ptr =
nullptr,
158 :
OpBase(row_field_name, col_field_name,
OpBase::OPROWCOL, ents_ptr),
159 matD(mat_D), betaCoeff(beta) {
160 if (row_field_name == col_field_name)
165 boost::shared_ptr<MatrixDouble>
matD;
175template <
int SPACE_DIM,
int S,
typename OpBase>
180 const std::string col_field_name,
181 boost::shared_ptr<MatrixDouble> mat_D,
183 :
OpBase(row_field_name, col_field_name,
OpBase::OPROWCOL), matD(mat_D),
187 boost::shared_ptr<MatrixDouble>
matD;
197template <
int SPACE_DIM,
int S,
typename OpBase>
202 const std::string col_field_name,
203 boost::shared_ptr<MatrixDouble> mat_D,
204 boost::shared_ptr<Range> ents_ptr =
nullptr)
205 :
OpBase(row_field_name, col_field_name,
OpBase::OPROWCOL, ents_ptr),
207 if (row_field_name == col_field_name)
212 boost::shared_ptr<MatrixDouble>
matD;
217template <
int SPACE_DIM, IntegrationType I,
typename OpBase>
220template <
int SPACE_DIM,
typename OpBase>
223 const std::string col_field_name,
225 const bool assemble_transpose =
false,
226 const bool only_transpose =
false)
227 :
OpBase(row_field_name, col_field_name,
OpBase::OPROWCOL),
228 alphaConstant(alpha_fun) {
229 this->assembleTranspose = assemble_transpose;
230 this->onlyTranspose = only_transpose;
244template <
int SPACE_DIM,
typename OpBase, CoordinateTypes CoordSys>
249 const std::string col_field_name,
ConstantFun alpha_fun,
250 const bool assemble_transpose,
251 const bool only_transpose =
false)
252 :
OpBase(row_field_name, col_field_name,
OpBase::OPROWCOL),
253 alphaConstant(alpha_fun) {
254 this->assembleTranspose = assemble_transpose;
255 this->onlyTranspose = only_transpose;
259 const std::string col_field_name,
ConstantFun alpha_fun,
260 ScalarFun beta_fun,
const bool assemble_transpose,
261 const bool only_transpose =
false)
262 :
OpBase(row_field_name, col_field_name,
OpBase::OPROWCOL),
263 alphaConstant(alpha_fun), betaCoeff(beta_fun) {
264 this->assembleTranspose = assemble_transpose;
265 this->onlyTranspose = only_transpose;
282template <
int SPACE_DIM,
typename OpBase, CoordinateTypes COORDINATE_SYSTEM>
286 const std::string row_field_name,
const std::string col_field_name,
287 ScalarFun alpha_fun = [](
double,
double,
double)
constexpr {
return 1; },
288 const bool assemble_transpose =
false,
const bool only_transpose =
false)
289 :
OpBase(row_field_name, col_field_name, OpBase::OPROWCOL),
290 alphaConstant(alpha_fun) {
291 this->assembleTranspose = assemble_transpose;
292 this->onlyTranspose = only_transpose;
307template <
int SPACE_DIM,
typename OpBase>
311 const std::string col_field_name,
313 const bool assemble_transpose =
false,
314 const bool only_transpose =
false)
315 :
OpBase(row_field_name, col_field_name,
OpBase::OPROWCOL),
316 alphaConstant(alpha_fun) {
317 this->assembleTranspose = assemble_transpose;
318 this->onlyTranspose = only_transpose;
328template <
int SPACE_DIM,
typename OpBase>
332 const std::string col_field_name,
334 const bool assemble_transpose =
false,
335 const bool only_transpose =
false)
336 :
OpBase(row_field_name, col_field_name,
OpBase::OPROWCOL),
337 alphaConstant(alpha_fun) {
338 this->assembleTranspose = assemble_transpose;
339 this->onlyTranspose = only_transpose;
349template <
int SPACE_DIM, IntegrationType I,
typename OpBase>
352template <
int SPACE_DIM,
typename OpBase>
356 const std::string col_field_name,
357 ConstantFun alpha_fun,
const bool assemble_transpose,
358 const bool only_transpose =
false)
359 :
OpBase(row_field_name, col_field_name,
OpBase::OPROWCOL),
360 alphaConstant(alpha_fun) {
361 this->assembleTranspose = assemble_transpose;
362 this->onlyTranspose = only_transpose;
366 const std::string col_field_name,
368 const bool assemble_transpose,
369 const bool only_transpose =
false)
370 :
OpBase(row_field_name, col_field_name,
OpBase::OPROWCOL),
371 alphaConstant(alpha_fun), betaCoeff(beta_coeff) {
372 this->assembleTranspose = assemble_transpose;
373 this->onlyTranspose = only_transpose;
393template <
typename EleOp>
394template <AssemblyType A>
395template <IntegrationType I>
404 template <
int BASE_DIM,
int FIELD_DIM,
int SPACE_DIM>
420 template <
int BASE_DIM,
int FIELD_DIM>
433 template <
int BASE_DIM,
int FIELD_DIM,
int SPACE_DIM,
int S = 0>
447 template <
int BASE_DIM,
int FIELD_DIM,
int SPACE_DIM,
int S = 0>
462 template <
int BASE_DIM,
int FIELD_DIM,
int SPACE_DIM,
int S = 0>
471 template <
int SPACE_DIM>
479 template <
int SPACE_DIM, CoordinateTypes CoordSys = CARTESIAN>
487 template <
int SPACE_DIM, CoordinateTypes COORDINATE_SYSTEM = CARTESIAN>
496 template <
int BASE_DIM,
int FIELD_DIM,
int SPACE_DIM>
505 template <
int SPACE_DIM>
509template <
int SPACE_DIM,
typename OpBase>
516 const double vol = OpBase::getMeasure();
518 auto t_w = OpBase::getFTensor0IntegrationWeight();
522 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
526 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
528 const double alpha = t_w * beta;
553template <
int FIELD_DIM,
int SPACE_DIM,
typename OpBase>
561 const double vol = OpBase::getMeasure();
563 auto t_w = OpBase::getFTensor0IntegrationWeight();
567 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
572 auto get_t_vec = [&](
const int rr) {
573 std::array<double *, FIELD_DIM> ptrs;
582 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
584 const double alpha = t_w * beta;
595 t_vec(
j) += alpha * (t_row_grad(
i) * t_col_grad(
i));
612template <
typename OpBase>
619 auto log_error = [&]() {
629 "Wrong number of base functions on rows %d < %d",
635 "Wrong number of base functions on cols %d < %d",
641 "Wrong number of integration points on rows %d != %d",
647 "Wrong number of integration points on cols %d != %d",
653 const double vol = OpBase::getMeasure();
655 auto t_w = OpBase::getFTensor0IntegrationWeight();
659 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
662 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
664 const double alpha = t_w * beta;
674 *a_mat_ptr += alpha * (t_row_base * t_col_base);
688template <
int FIELD_DIM,
typename OpBase>
694 const double vol = OpBase::getMeasure();
696 auto t_w = OpBase::getFTensor0IntegrationWeight();
700 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
703 auto get_t_vec = [&](
const int rr) {
704 std::array<double *, FIELD_DIM> ptrs;
713 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
715 const double alpha = t_w * beta;
726 t_vec(
i) += alpha * (t_row_base * t_col_base);
740template <
int FIELD_DIM,
typename OpBase>
746 size_t nb_base_functions = row_data.
getN().size2() / 3;
748 const double vol = OpBase::getMeasure();
750 auto t_w = OpBase::getFTensor0IntegrationWeight();
754 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
757 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
759 const double alpha = t_w * beta;
769 (*a_mat_ptr) += alpha * (t_row_base(
i) * t_col_base(
i));
775 for (; rr < nb_base_functions; ++rr)
783template <
typename OpBase>
790 auto get_t_vec = [&](
const int rr) {
794 size_t nb_base_functions = row_data.
getN().size2() / 3;
796 const double vol = OpBase::getMeasure();
798 auto t_w = OpBase::getFTensor0IntegrationWeight();
802 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
805 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
807 const double alpha = t_w * beta;
813 auto t_vec = get_t_vec(2 * rr);
817 t_vec(
I) += alpha * (t_row_base(
k) * t_col_base(
k));
823 for (; rr < nb_base_functions; ++rr)
831template <
typename OpBase>
838 auto get_t_vec = [&](
const int rr) {
843 size_t nb_base_functions = row_data.
getN().size2() / 3;
845 const double vol = OpBase::getMeasure();
847 auto t_w = OpBase::getFTensor0IntegrationWeight();
851 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
854 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
856 const double alpha = t_w * beta;
862 auto t_vec = get_t_vec(3 * rr);
866 t_vec(
i) += alpha * (t_row_base(
k) * t_col_base(
k));
872 for (; rr < nb_base_functions; ++rr)
880template <
int SPACE_DIM,
int S,
typename OpBase>
887 const size_t nb_row_dofs = row_data.
getIndices().size();
888 const size_t nb_col_dofs = col_data.
getIndices().size();
890 if (nb_row_dofs && nb_col_dofs) {
901 double vol = OpBase::getMeasure();
904 auto t_w = OpBase::getFTensor0IntegrationWeight();
911 auto t_D = getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, S>(*matD);
914 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
920 double a = t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
927 auto t_m = OpBase::template getLocMat<SPACE_DIM>(
SPACE_DIM * rr);
933 t_rowD(
l,
j,
k) = t_D(
i,
j,
k,
l) * (
a * t_row_diff_base(
i));
940 for (
int cc = 0; cc <= nb_cols; ++cc) {
943 t_m(
i,
j) += t_rowD(
i,
j,
k) * t_col_diff_base(
k);
968 auto t_m_rr = getFTensor2FromArray<SPACE_DIM, SPACE_DIM, SPACE_DIM>(
970 auto t_m_cc = getFTensor2FromArray<SPACE_DIM, SPACE_DIM>(
974 t_m_rr(
i,
j) = t_m_cc(
j,
i);
985template <
int SPACE_DIM,
int S,
typename OpBase>
999 auto &row_hessian = row_data.
getN(BaseDerivatives::SecondDerivative);
1000 auto &col_hessian = col_data.
getN(BaseDerivatives::SecondDerivative);
1005 "Wrong number of integration pts (%d != %d)",
1008 if (row_hessian.size2() !=
1011 "Wrong number of base functions (%d != %d)",
1017 "Wrong number of base functions (%d < %d)", row_hessian.size2(),
1022 "Wrong number of integration pts (%d != %d)",
1027 "Wrong number of base functions (%d < %d)", col_hessian.size2(),
1033 double vol = OpBase::getMeasure();
1036 auto t_w = OpBase::getFTensor0IntegrationWeight();
1038 auto t_row_diff2 = getFTensor2SymmetricLowerFromPtr<SPACE_DIM>(
1039 &*row_hessian.data().begin());
1043 auto t_D = getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, S>(*matD);
1049 double a = t_w * vol;
1059 t_rowD(
k,
l) = t_D(
i,
j,
k,
l) * (
a * t_row_diff2(
i,
j));
1063 getFTensor2SymmetricLowerFromPtr<SPACE_DIM>(&col_hessian(gg, 0));
1069 *m_ptr += t_rowD(
i,
j) * t_col_diff2(
i,
j);
1094template <
int SPACE_DIM,
int S,
typename OpBase>
1101 const size_t nb_row_dofs = row_data.
getIndices().size();
1102 const size_t nb_col_dofs = col_data.
getIndices().size();
1104 if (nb_row_dofs && nb_col_dofs) {
1112 double vol = OpBase::getMeasure();
1115 auto t_w = OpBase::getFTensor0IntegrationWeight();
1122 getFTensor4FromMat<SPACE_DIM, SPACE_DIM, SPACE_DIM, SPACE_DIM, S>(
1126 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1132 double a = t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
1139 auto t_m = OpBase::template getLocMat<SPACE_DIM>(
SPACE_DIM * rr);
1143 t_row(
i,
k,
l) = t_D(
i,
j,
k,
l) * (
a * t_row_diff_base(
j));
1152 t_m(
i,
k) += t_row(
i,
k,
l) * t_col_diff_base(
l);
1178template <
int SPACE_DIM,
typename OpBase>
1184 auto t_w = this->getFTensor0IntegrationWeight();
1186 size_t nb_base_functions = row_data.
getN().size2() / 3;
1189 const double alpha_constant = alphaConstant();
1192 const double alpha = alpha_constant * this->getMeasure() * t_w;
1196 const double t_row_div_base = t_row_diff_base(
i,
i);
1199 this->locMat(rr, cc) += alpha * t_row_div_base * t_col_base;
1204 for (; rr < nb_base_functions; ++rr)
1213template <
int SPACE_DIM,
typename OpBase, CoordinateTypes CoordSys>
1220 auto t_w = this->getFTensor0IntegrationWeight();
1221 auto t_coords = this->getFTensor1CoordsAtGaussPts();
1223 size_t nb_base_functions = row_data.
getN().size2() / 3;
1226 const double alpha_constant = alphaConstant() * this->getMeasure();
1229 const double alpha =
1230 alpha_constant * t_w * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
1234 auto t_mat_diag = getFTensor1FromArrayDiag<SPACE_DIM, SPACE_DIM>(
1236 const double t_row_div_base = t_row_diff_base(
i,
i);
1240 t_mat_diag(
i) += alpha * t_row_div_base * t_col_base;
1242 t_mat_diag(
i) += alpha * (t_row_base(0) / t_coords(0)) * t_col_base;
1250 for (; rr < nb_base_functions; ++rr) {
1262template <
int SPACE_DIM,
typename OpBase, CoordinateTypes COORDINATE_SYSTEM>
1272 "Number of rows in matrix should be multiple of space dimensions");
1276 if constexpr (COORDINATE_SYSTEM ==
POLAR || COORDINATE_SYSTEM ==
SPHERICAL)
1278 "%s coordiante not implemented",
1281 auto t_w = this->getFTensor0IntegrationWeight();
1282 auto t_coords = this->getFTensor1CoordsAtGaussPts();
1283 size_t nb_base_functions_row = row_data.
getN().size2();
1285 const double vol = this->getMeasure();
1288 const double alpha =
1289 alphaConstant(t_coords(0), t_coords(1), t_coords(2)) * t_w * vol;
1292 auto t_m = getFTensor1FromPtr<SPACE_DIM>(
OpBase::locMat.data().data());
1295 if constexpr (COORDINATE_SYSTEM ==
CARTESIAN) {
1297 const double r_val = alpha * t_row_base;
1300 t_m(
i) += r_val * t_col_diff_base(
i);
1311 const double r_val = alpha * t_row_base;
1315 t_m(
i) += r_val * t_col_diff_base(
i);
1316 t_m(0) += (r_val / t_coords(0)) * t_col_base;
1325 for (; rr < nb_base_functions_row; ++rr)
1335template <
int SPACE_DIM,
typename OpBase>
1342 auto t_w = this->getFTensor0IntegrationWeight();
1344 size_t nb_base_functions = row_data.
getN().size2() / 3;
1346 auto &mat = this->locMat;
1347 const double alpha_constant = alphaConstant();
1350 const double alpha = alpha_constant * this->getMeasure() * t_w;
1356 mat(rr, cc) += alpha * t_row_base(
i) * t_col_diff_base(
i);
1361 for (; rr < nb_base_functions; ++rr)
1370template <
int SPACE_DIM,
typename OpBase>
1377 auto t_w = this->getFTensor0IntegrationWeight();
1379 size_t nb_base_functions = row_data.
getN().size2();
1382 auto get_t_vec = [&](
const int rr) {
1383 std::array<double *, SPACE_DIM> ptrs;
1389 const double alpha_constant = alphaConstant();
1392 const double alpha = alpha_constant * this->getMeasure() * t_w;
1399 t_vec(
i) += alpha * t_row_base * t_col_diff_base(
i);
1405 for (; rr < nb_base_functions; ++rr)
1414template <
int SPACE_DIM,
typename OpBase>
1420 auto t_w = this->getFTensor0IntegrationWeight();
1421 auto t_coords = this->getFTensor1CoordsAtGaussPts();
1423 size_t nb_base_functions = row_data.
getN().size2() / 3;
1425 const double alpha_constant = alphaConstant() * this->getMeasure();
1428 const double alpha =
1429 alpha_constant * betaCoeff(t_coords(0), t_coords(1), t_coords(2)) * t_w;
1433 auto t_mat_diag = getFTensor1FromArrayDiag<SPACE_DIM, SPACE_DIM>(
1438 t_mat_diag(
i) += alpha * t_row_base(
j) * t_col_diff_base(
j);
1445 for (; rr < nb_base_functions; ++rr)
static const char *const CoordinateTypesNames[]
Coordinate system names.
#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()
#define MOFEM_LOG(channel, severity)
Log.
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
implementation of Data Operators for Forces and Sources
double scalar_fun_one(const double, const double, const double)
boost::function< double()> ConstantFun
Constant function type.
constexpr IntegrationType I
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.
const VectorFieldEntities & getFieldEntities() const
get field entities
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.
const VectorInt & getIndices() const
Get global indices of dofs on entity.
int rowSide
row side number
int colSide
column side number
int nbRows
number of dofs on rows
EntityType colType
column type
int nbIntegrationPts
number of integration points
MatrixDouble locMat
local entity block matrix
int nbCols
number if dof on column
int nbRowBaseFunctions
number or row base functions
EntityType rowType
row type
OpGradGradImpl(const std::string row_field_name, const std::string col_field_name, ScalarFun beta=scalar_fun_one, boost::shared_ptr< Range > ents_ptr=nullptr)
FTensor::Index< 'i', SPACE_DIM > i
summit Index
OpGradGradImpl(const std::string row_field_name, const std::string col_field_name, ScalarFun beta=scalar_fun_one, boost::shared_ptr< Range > ents_ptr=nullptr)
FTensor::Index< 'i', SPACE_DIM > i
summit Index
boost::shared_ptr< MatrixDouble > matD
OpGradGradSymTensorGradGradImpl(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< MatrixDouble > mat_D, boost::shared_ptr< Range > ents_ptr=nullptr)
FTensor::Index< 'i', SPACE_DIM > i
summit Index
FTensor::Index< 'i', SPACE_DIM > i
summit Index
boost::shared_ptr< MatrixDouble > matD
OpGradSymTensorGradImpl(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< MatrixDouble > mat_D, boost::shared_ptr< Range > ents_ptr=nullptr, ScalarFun beta=scalar_fun_one)
boost::shared_ptr< MatrixDouble > matD
FTensor::Index< 'i', SPACE_DIM > i
summit Index
OpGradTensorGradImpl(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< MatrixDouble > mat_D, ScalarFun beta=scalar_fun_one)
OpMassImpl(const std::string row_field_name, const std::string col_field_name, ScalarFun beta=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
OpMassImpl(const std::string row_field_name, const std::string col_field_name, ScalarFun beta=scalar_fun_one, boost::shared_ptr< Range > ents_ptr=nullptr)
OpMassImpl(const std::string row_field_name, const std::string col_field_name, ScalarFun beta=scalar_fun_one, boost::shared_ptr< Range > ents_ptr=nullptr)
OpMassImpl(const std::string row_field_name, const std::string col_field_name, ScalarFun beta=scalar_fun_one, boost::shared_ptr< Range > ents_ptr=nullptr)
ConstantFun alphaConstant
FTensor::Index< 'i', SPACE_DIM > i
summit Index
OpMixDivTimesScalarImpl(const std::string row_field_name, const std::string col_field_name, ConstantFun alpha_fun, const bool assemble_transpose=false, const bool only_transpose=false)
OpMixDivTimesVecImpl(const std::string row_field_name, const std::string col_field_name, ConstantFun alpha_fun, const bool assemble_transpose, const bool only_transpose=false)
FTensor::Index< 'i', SPACE_DIM > i
summit Index
OpMixDivTimesVecImpl(const std::string row_field_name, const std::string col_field_name, ConstantFun alpha_fun, ScalarFun beta_fun, const bool assemble_transpose, const bool only_transpose=false)
OpMixScalarTimesDivImpl(const std::string row_field_name, const std::string col_field_name, ScalarFun alpha_fun=[](double, double, double) constexpr { return 1;}, const bool assemble_transpose=false, const bool only_transpose=false)
FTensor::Index< 'i', SPACE_DIM > i
summit Index
FTensor::Index< 'i', SPACE_DIM > i
summit Index
OpMixTensorTimesGradImpl(const std::string row_field_name, const std::string col_field_name, ConstantFun alpha_fun, const bool assemble_transpose, const bool only_transpose=false)
OpMixTensorTimesGradImpl(const std::string row_field_name, const std::string col_field_name, ConstantFun alpha_fun, ScalarFun beta_coeff, const bool assemble_transpose, const bool only_transpose=false)
FTensor::Index< 'j', SPACE_DIM > j
summit Index
ConstantFun alphaConstant
FTensor::Index< 'i', SPACE_DIM > i
summit Index
OpMixVectorTimesGradImpl(const std::string row_field_name, const std::string col_field_name, ConstantFun alpha_fun, const bool assemble_transpose=false, const bool only_transpose=false)
FTensor::Index< 'i', SPACE_DIM > i
summit Index
ConstantFun alphaConstant
OpMixVectorTimesGradImpl(const std::string row_field_name, const std::string col_field_name, ConstantFun alpha_fun, const bool assemble_transpose=false, const bool only_transpose=false)