6#ifndef __USER_DATA_OPERATORS_HPP__
7 #define __USER_DATA_OPERATORS_HPP__
22template <
class T,
class A>
28 boost::shared_ptr<ublas::vector<T, A>> data_ptr,
29 const EntityType zero_type = MBVERTEX)
38 boost::shared_ptr<ublas::vector<T, A>> data_ptr,
57 boost::shared_ptr<ublas::vector<T, A>>
dataPtr;
67template <
class T,
class A>
100 if (type ==
zeroType || vec.size() != nb_gauss_pts) {
101 vec.resize(nb_gauss_pts,
false);
114 for (
int i = 0;
i != local_indices.size(); ++
i)
115 if (local_indices[
i] != -1)
123 const size_t nb_base_functions = data.
getN().size2();
126 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
129 for (; bb != nb_dofs; ++bb) {
130 values_at_gauss_pts += field_data * base_function;
135 for (; bb < nb_base_functions; ++bb)
137 ++values_at_gauss_pts;
155template <PetscData::DataContext CTX>
160 const std::string
field_name, boost::shared_ptr<VectorDouble> data_ptr,
161 const EntityType zero_at_type = MBVERTEX)
176 if (type ==
zeroAtType || vec.size() != nb_gauss_pts) {
177 vec.resize(nb_gauss_pts,
false);
182 const size_t nb_dofs = local_indices.size();
187 auto get_array = [&](
const auto ctx,
auto vec) {
193 <<
"In this case field degrees of freedom are read from vector. "
194 "That usually happens when time solver is used, and acces to "
195 "first or second rates is needed. You probably not set ts_u, "
196 "ts_u_t, or ts_u_tt and associated data structure, i.e. "
197 "data_ctx to CTX_SET_X, CTX_SET_X_T, or CTX_SET_X_TT "
202 CHKERR VecGetArrayRead(vec, &array);
206 auto restore_array = [&](
auto vec) {
207 return VecRestoreArrayRead(vec, &array);
222 "That case is not implemented");
225 std::array<double, MAX_DOFS_ON_ENTITY> dot_dofs_vector;
226 for (
int i = 0;
i != local_indices.size(); ++
i)
227 if (local_indices[
i] != -1)
228 dot_dofs_vector[
i] = array[local_indices[
i]];
230 dot_dofs_vector[
i] = 0;
244 "That case is not implemented");
247 const size_t nb_base_functions = data.
getN().size2();
251 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
253 for (; bb != nb_dofs; ++bb) {
254 values_at_gauss_pts += dot_dofs_vector[bb] * base_function;
259 for (; bb < nb_base_functions; ++bb)
261 ++values_at_gauss_pts;
292template <
int Tensor_Dim,
class T,
class L,
class A>
298 boost::shared_ptr<ublas::matrix<T, L, A>> data_ptr,
299 const EntityType zero_type = MBVERTEX)
318 boost::shared_ptr<ublas::matrix<T, L, A>>
dataPtr;
322template <
int Tensor_Dim,
class T,
class L,
class A>
328 "Not implemented for T = %s and dim = %d",
338template <
int Tensor_Dim>
344 boost::shared_ptr<MatrixDouble> data_ptr,
345 const EntityType zero_type = MBVERTEX)
354 boost::shared_ptr<MatrixDouble> data_ptr,
356 const EntityType zero_type = MBVERTEX)
378template <
int Tensor_Dim>
380 Tensor_Dim,
double, ublas::row_major,
385 const size_t nb_gauss_pts = getGaussPts().size2();
386 auto &mat = *dataPtr;
387 if (type == zeroType || mat.size2() != nb_gauss_pts) {
388 mat.resize(Tensor_Dim, nb_gauss_pts,
false);
395 if (dataVec.use_count()) {
396 dotVector.resize(nb_dofs,
false);
398 CHKERR VecGetArrayRead(dataVec, &array);
400 for (
int i = 0;
i != local_indices.size(); ++
i)
401 if (local_indices[
i] != -1)
402 dotVector[
i] = array[local_indices[
i]];
405 CHKERR VecRestoreArrayRead(dataVec, &array);
410 const size_t nb_base_functions = data.
getN().size2();
414 const size_t size = nb_dofs / Tensor_Dim;
415 if (nb_dofs % Tensor_Dim) {
417 "Nb. of DOFs is inconsistent with Tensor_Dim");
419 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
423 if (field_data.l2() != field_data.l2()) {
424 MOFEM_LOG(
"SELF", Sev::error) <<
"field data: " << field_data;
426 "Wrong number in coefficients");
431 for (; bb != size; ++bb) {
435 if (base_function != base_function) {
436 MOFEM_LOG(
"SELF", Sev::error) <<
"base function: " << base_function;
438 "Wrong number number in base functions");
443 values_at_gauss_pts(
I) += field_data(
I) * base_function;
449 for (; bb < nb_base_functions; ++bb)
451 ++values_at_gauss_pts;
455 if (dataVec.use_count()) {
467template <
int Tensor_Dim>
470 Tensor_Dim, double, ublas::row_major, DoubleAllocator> {
473 Tensor_Dim,
double, ublas::row_major,
487template <
int Tensor_Dim, CoordinateTypes COORDINATE_SYSTEM = CARTESIAN>
492 const std::string
field_name, boost::shared_ptr<VectorDouble> data_ptr,
493 const EntityType zero_type = MBVERTEX)
506 if constexpr (COORDINATE_SYSTEM ==
POLAR || COORDINATE_SYSTEM ==
SPHERICAL)
508 "%s coordiante not implemented",
514 vec.resize(nb_gauss_pts,
false);
522 const size_t nb_base_functions = data.
getN().size2();
525 const size_t size = nb_dofs / Tensor_Dim;
527 if (nb_dofs % Tensor_Dim) {
529 "Number of dofs should multiple of dimensions");
534 if constexpr (COORDINATE_SYSTEM ==
CARTESIAN) {
536 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
539 for (; bb != size; ++bb) {
540 values_at_gauss_pts += field_data(
I) * diff_base_function(
I);
542 ++diff_base_function;
546 for (; bb < nb_base_functions; ++bb)
547 ++diff_base_function;
548 ++values_at_gauss_pts;
558 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
561 for (; bb != size; ++bb) {
562 values_at_gauss_pts += field_data(
I) * diff_base_function(
I);
563 values_at_gauss_pts +=
564 base_function * (field_data(0) / t_coords(0));
567 ++diff_base_function;
571 for (; bb < nb_base_functions; ++bb) {
573 ++diff_base_function;
575 ++values_at_gauss_pts;
596template <
int Tensor_Dim, PetscData::DataContext CTX>
601 const std::string
field_name, boost::shared_ptr<MatrixDouble> data_ptr,
602 const EntityType zero_at_type = MBVERTEX,
bool throw_error =
true)
615 const size_t nb_dofs = local_indices.size();
620 mat.resize(Tensor_Dim, nb_gauss_pts,
false);
634 auto get_array = [&](
const auto ctx,
auto vec) {
640 <<
"In this case field degrees of freedom are read from vector. "
641 "That usually happens when time solver is used, and access to "
642 "first or second rates is needed. You probably not set ts_u, "
643 "ts_u_t, or ts_u_tt and associated data structure, i.e. "
644 "data_ctx to CTX_SET_X, CTX_SET_DX, CTX_SET_X_T, or "
645 "CTX_SET_X_TT respectively";
649 CHKERR VecGetArrayRead(vec, &array);
653 auto restore_array = [&](
auto vec) {
654 return VecRestoreArrayRead(vec, &array);
672 "That case is not implemented");
676 for (
int i = 0;
i != local_indices.size(); ++
i)
677 if (local_indices[
i] != -1)
697 "That case is not implemented");
700 const size_t nb_base_functions = data.
getN().size2();
705 const size_t size = nb_dofs / Tensor_Dim;
706 if (nb_dofs % Tensor_Dim) {
709 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
712 for (; bb != size; ++bb) {
713 values_at_gauss_pts(
I) += field_data(
I) * base_function;
719 for (; bb < nb_base_functions; ++bb)
721 ++values_at_gauss_pts;
738template <
int Tensor_Dim>
748template <
int Tensor_Dim>
758template <
int Tensor_Dim>
772template <
int Tensor_Dim0,
int Tensor_Dim1,
class T,
class L,
class A>
778 boost::shared_ptr<ublas::matrix<T, L, A>> data_ptr,
779 const EntityType zero_type = MBVERTEX)
791 boost::shared_ptr<ublas::matrix<T, L, A>>
dataPtr;
795template <
int Tensor_Dim0,
int Tensor_Dim1,
class T,
class L,
class A>
801 "Not implemented for T = %s, dim0 = %d and dim1 = %d",
803 Tensor_Dim0, Tensor_Dim1);
807template <
int Tensor_Dim0,
int Tensor_Dim1>
815 ublas::matrix<double, ublas::row_major, DoubleAllocator>>
817 const EntityType zero_type = MBVERTEX)
828 ublas::matrix<double, ublas::row_major, DoubleAllocator>>
842 boost::shared_ptr<ublas::matrix<double, ublas::row_major, DoubleAllocator>>
849template <
int Tensor_Dim0,
int Tensor_Dim1>
851 Tensor_Dim0, Tensor_Dim1,
double, ublas::row_major,
857 const size_t nb_gauss_pts = data.
getN().size1();
858 if (type == zeroType) {
859 mat.resize(Tensor_Dim0 * Tensor_Dim1, nb_gauss_pts,
false);
865 if (dataVec.use_count()) {
866 dotVector.resize(nb_dofs,
false);
868 CHKERR VecGetArrayRead(dataVec, &array);
870 for (
int i = 0;
i != local_indices.size(); ++
i)
871 if (local_indices[
i] != -1)
872 dotVector[
i] = array[local_indices[
i]];
875 CHKERR VecRestoreArrayRead(dataVec, &array);
880 const size_t nb_base_functions = data.
getN().size2();
882 auto values_at_gauss_pts =
886 const size_t size = nb_dofs / (Tensor_Dim0 * Tensor_Dim1);
887 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
890 for (; bb != size; ++bb) {
891 values_at_gauss_pts(
i,
j) += field_data(
i,
j) * base_function;
895 for (; bb < nb_base_functions; ++bb)
897 ++values_at_gauss_pts;
900 if (dataVec.use_count()) {
912template <
int Tensor_Dim0,
int Tensor_Dim1>
915 Tensor_Dim0, Tensor_Dim1, double, ublas::row_major, DoubleAllocator> {
918 Tensor_Dim0, Tensor_Dim1,
double, ublas::row_major,
927template <
int Tensor_Dim0,
int Tensor_Dim1>
932 boost::shared_ptr<MatrixDouble> data_ptr,
933 const EntityType zero_at_type = MBVERTEX)
948 mat.resize(Tensor_Dim0 * Tensor_Dim1, nb_gauss_pts,
false);
952 const size_t nb_dofs = local_indices.size();
957 for (
size_t i = 0;
i != local_indices.size(); ++
i)
958 if (local_indices[
i] != -1)
964 const size_t nb_base_functions = data.
getN().size2();
967 auto values_at_gauss_pts =
971 const size_t size = nb_dofs / (Tensor_Dim0 * Tensor_Dim1);
972 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
976 for (; bb != size; ++bb) {
977 values_at_gauss_pts(
i,
j) += field_data(
i,
j) * base_function;
981 for (; bb < nb_base_functions; ++bb)
983 ++values_at_gauss_pts;
1002template <
int Tensor_Dim>
1007 const std::string
field_name, boost::shared_ptr<MatrixDouble> data_ptr,
1008 const EntityType zero_type = MBEDGE,
const int zero_side = 0)
1017 const std::string
field_name, boost::shared_ptr<MatrixDouble> data_ptr,
1019 const int zero_side = 0)
1034 mat.resize((Tensor_Dim * (Tensor_Dim + 1)) / 2, nb_gauss_pts,
false);
1043 const double *array;
1046 for (
int i = 0;
i != local_indices.size(); ++
i)
1047 if (local_indices[
i] != -1)
1055 const int nb_base_functions = data.
getN().size2();
1060 const int size = nb_dofs / ((Tensor_Dim * (Tensor_Dim + 1)) / 2);
1061 for (
int gg = 0; gg != nb_gauss_pts; ++gg) {
1064 for (; bb != size; ++bb) {
1065 values_at_gauss_pts(
i,
j) += field_data(
i,
j) * base_function;
1069 for (; bb < nb_base_functions; ++bb)
1071 ++values_at_gauss_pts;
1096template <
int Tensor_Dim>
1101 const std::string
field_name, boost::shared_ptr<MatrixDouble> data_ptr,
1102 const EntityType zero_type = MBEDGE,
const int zero_side = 0)
1115 constexpr auto symm_size = (Tensor_Dim * (Tensor_Dim + 1)) / 2;
1117 mat.resize(symm_size, nb_gauss_pts,
false);
1121 const int nb_dofs = local_indices.size();
1129 <<
"In this case field degrees of freedom are read from vector. "
1130 "That usually happens when time solver is used, and acces to "
1131 "first rates is needed. You probably not set "
1132 "ts_u_t and associated data structure data_ctx to CTX_SET_X_T "
1139 const double *array;
1141 for (
int i = 0;
i != local_indices.size(); ++
i)
1142 if (local_indices[
i] != -1)
1148 const int nb_base_functions = data.
getN().size2();
1154 const int size = nb_dofs / symm_size;
1155 for (
int gg = 0; gg != nb_gauss_pts; ++gg) {
1158 for (; bb != size; ++bb) {
1159 values_at_gauss_pts(
i,
j) += field_data(
i,
j) * base_function;
1163 for (; bb < nb_base_functions; ++bb)
1165 ++values_at_gauss_pts;
1178 static_assert(Dim || !Dim,
"not implemented");
1187 &dotVector[0], &dotVector[1], &dotVector[2], &dotVector[3], &dotVector[4],
1196 &dotVector[0], &dotVector[1], &dotVector[2]);
1210template <
int Tensor_Dim,
class T,
class L,
class A>
1222template <
int Tensor_Dim>
1226 Tensor_Dim, double, ublas::row_major, DoubleAllocator> {
1229 Tensor_Dim,
double, ublas::row_major,
1248template <
int Tensor_Dim>
1250 Tensor_Dim,
double, ublas::row_major,
1255 const size_t nb_gauss_pts = this->getGaussPts().size2();
1256 auto &mat = *this->dataPtr;
1257 if (type == this->zeroType) {
1258 mat.resize(Tensor_Dim, nb_gauss_pts,
false);
1266 const int nb_base_functions = data.
getN().size2();
1271 if (nb_dofs > nb_base_functions)
1273 "Number of base functions inconsistent with number of DOFs "
1275 nb_dofs, nb_base_functions);
1277 if (data.
getDiffN().size2() != nb_base_functions * Tensor_Dim)
1280 "Number of base functions inconsistent with number of derivatives "
1282 data.
getDiffN().size2(), nb_base_functions);
1284 if (data.
getDiffN().size1() != nb_gauss_pts)
1287 "Number of base functions inconsistent with number of integration "
1289 data.
getDiffN().size2(), nb_gauss_pts);
1294 for (
int gg = 0; gg < nb_gauss_pts; ++gg) {
1297 for (; bb != nb_dofs; ++bb) {
1298 gradients_at_gauss_pts(
I) += field_data * diff_base_function(
I);
1300 ++diff_base_function;
1303 for (; bb < nb_base_functions; ++bb)
1304 ++diff_base_function;
1305 ++gradients_at_gauss_pts;
1318template <
int Tensor_Dim>
1321 Tensor_Dim, double, ublas::row_major, DoubleAllocator> {
1323 Tensor_Dim,
double, ublas::row_major,
1331template <
int Tensor_Dim>
1334 Tensor_Dim, double, ublas::row_major, DoubleAllocator> {
1337 Tensor_Dim,
double, ublas::row_major,
1351template <
int Tensor_Dim>
1356 const size_t nb_gauss_pts = this->getGaussPts().size2();
1357 auto &mat = *this->dataPtr;
1358 if (type == this->zeroType) {
1359 mat.resize(Tensor_Dim * Tensor_Dim, nb_gauss_pts,
false);
1367 const int nb_base_functions = data.
getN().size2();
1369 auto &hessian_base = data.
getN(BaseDerivatives::SecondDerivative);
1371 if (hessian_base.size1() != nb_gauss_pts) {
1373 "Wrong number of integration pts (%ld != %d)",
1374 hessian_base.size1(), nb_gauss_pts);
1376 if (hessian_base.size2() != nb_base_functions * Tensor_Dim * Tensor_Dim) {
1378 "Wrong number of base functions (%ld != %d)",
1379 hessian_base.size2() / (Tensor_Dim * Tensor_Dim),
1382 if (hessian_base.size2() < nb_dofs * Tensor_Dim * Tensor_Dim) {
1384 "Wrong number of base functions (%ld < %d)",
1385 hessian_base.size2(), nb_dofs * Tensor_Dim * Tensor_Dim);
1390 &*hessian_base.data().begin());
1392 auto hessian_at_gauss_pts =
1397 for (
int gg = 0; gg != nb_gauss_pts; ++gg) {
1400 for (; bb != nb_dofs; ++bb) {
1401 hessian_at_gauss_pts(
I,
J) +=
1402 field_data * t_diff2_base_function(
I,
J);
1404 ++t_diff2_base_function;
1407 for (; bb < nb_base_functions; ++bb) {
1408 ++t_diff2_base_function;
1411 ++hessian_at_gauss_pts;
1437template <
int Tensor_Dim0,
int Tensor_Dim1,
int S,
class T,
class L,
class A>
1446template <
int Tensor_Dim0,
int Tensor_Dim1,
int S>
1450 Tensor_Dim0, Tensor_Dim1, double, ublas::row_major, DoubleAllocator> {
1453 Tensor_Dim0, Tensor_Dim1,
double, ublas::row_major,
1472template <
int Tensor_Dim0,
int Tensor_Dim1,
int S>
1474 Tensor_Dim0, Tensor_Dim1, S,
double, ublas::row_major,
1480 "Data pointer not allocated");
1482 const size_t nb_gauss_pts = this->getGaussPts().size2();
1483 auto &mat = *this->dataPtr;
1484 if (type == this->zeroType) {
1485 mat.resize(Tensor_Dim0 * Tensor_Dim1, nb_gauss_pts,
false);
1494 if (this->dataVec.use_count()) {
1495 this->dotVector.resize(nb_dofs,
false);
1496 const double *array;
1497 CHKERR VecGetArrayRead(this->dataVec, &array);
1499 for (
int i = 0;
i != local_indices.size(); ++
i)
1500 if (local_indices[
i] != -1)
1501 this->dotVector[
i] = array[local_indices[
i]];
1503 this->dotVector[
i] = 0;
1504 CHKERR VecRestoreArrayRead(this->dataVec, &array);
1508 const int nb_base_functions = data.
getN().size2();
1510 auto gradients_at_gauss_pts =
1514 int size = nb_dofs / Tensor_Dim0;
1515 if (nb_dofs % Tensor_Dim0) {
1517 "Data inconsistency");
1519 for (
int gg = 0; gg < nb_gauss_pts; ++gg) {
1524 if (field_data.l2() != field_data.l2()) {
1525 MOFEM_LOG(
"SELF", Sev::error) <<
"field data " << field_data;
1527 "Wrong number in coefficients");
1532 for (; bb < size; ++bb) {
1535 if (diff_base_function.l2() != diff_base_function.l2()) {
1537 <<
"diff_base_function: " << diff_base_function;
1539 "Wrong number number in base functions");
1544 gradients_at_gauss_pts(
I,
J) += field_data(
I) * diff_base_function(
J);
1546 ++diff_base_function;
1550 for (; bb != nb_base_functions; ++bb)
1551 ++diff_base_function;
1552 ++gradients_at_gauss_pts;
1555 if (this->dataVec.use_count()) {
1572template <
int Tensor_Dim0,
int Tensor_Dim1,
int S = Tensor_Dim0>
1575 double, ublas::row_major,
1579 Tensor_Dim0, Tensor_Dim1, S,
double, ublas::row_major,
1588template <
int Tensor_Dim0,
int Tensor_Dim1>
1593 boost::shared_ptr<MatrixDouble> data_ptr,
1594 const EntityType zero_at_type = MBVERTEX)
1607 const int nb_dofs = local_indices.size();
1608 const int nb_gauss_pts = this->
getGaussPts().size2();
1612 mat.resize(Tensor_Dim0 * Tensor_Dim1, nb_gauss_pts,
false);
1619 const double *array;
1621 for (
int i = 0;
i != local_indices.size(); ++
i)
1622 if (local_indices[
i] != -1)
1628 const int nb_base_functions = data.
getN().size2();
1630 auto gradients_at_gauss_pts =
1634 int size = nb_dofs / Tensor_Dim0;
1635 if (nb_dofs % Tensor_Dim0) {
1639 for (
int gg = 0; gg < nb_gauss_pts; ++gg) {
1642 for (; bb < size; ++bb) {
1643 gradients_at_gauss_pts(
I,
J) += field_data(
I) * diff_base_function(
J);
1645 ++diff_base_function;
1649 for (; bb != nb_base_functions; ++bb)
1650 ++diff_base_function;
1651 ++gradients_at_gauss_pts;
1667template <
int Tensor_Dim0,
int Tensor_Dim1,
class T,
class L,
class A>
1673 const std::string
field_name, boost::shared_ptr<MatrixDouble> data_ptr,
1674 const EntityType zero_type = MBVERTEX)
1680template <
int Tensor_Dim0,
int Tensor_Dim1>
1684 Tensor_Dim0, Tensor_Dim1, double, ublas::row_major, DoubleAllocator> {
1687 const std::string
field_name, boost::shared_ptr<MatrixDouble> data_ptr,
1688 const EntityType zero_type = MBVERTEX)
1710template <
int Tensor_Dim0,
int Tensor_Dim1>
1712 Tensor_Dim0, Tensor_Dim1,
double, ublas::row_major,
1718 "Data pointer not allocated");
1720 const size_t nb_gauss_pts = this->getGaussPts().size2();
1721 constexpr size_t msize = (Tensor_Dim0 * (Tensor_Dim0 + 1)) / 2;
1722 auto &mat = *this->dataPtr;
1723 if (type == this->zeroType) {
1724 mat.resize(msize * Tensor_Dim1, nb_gauss_pts,
false);
1733 const int nb_base_functions = data.
getN().size2();
1735 auto gradients_at_gauss_pts =
1740 int size = nb_dofs / msize;
1741 if (nb_dofs % msize) {
1743 "Data inconsistency");
1745 for (
int gg = 0; gg < nb_gauss_pts; ++gg) {
1748 for (; bb < size; ++bb) {
1749 gradients_at_gauss_pts(
I,
J,
K) +=
1750 field_data(
I,
J) * diff_base_function(
K);
1752 ++diff_base_function;
1756 for (; bb != nb_base_functions; ++bb)
1757 ++diff_base_function;
1758 ++gradients_at_gauss_pts;
1770template <
int Tensor_Dim0,
int Tensor_Dim1>
1773 Tensor_Dim0, Tensor_Dim1, double, ublas::row_major, DoubleAllocator> {
1776 const std::string
field_name, boost::shared_ptr<MatrixDouble> data_ptr,
1777 const EntityType zero_type = MBVERTEX)
1779 Tensor_Dim0, Tensor_Dim1,
double, ublas::row_major,
1783template <
int Tensor_Dim0,
int Tensor_Dim1>
1786 Tensor_Dim0, Tensor_Dim1, double, ublas::row_major, DoubleAllocator> {
1789 Tensor_Dim0, Tensor_Dim1,
double, ublas::row_major,
1808template <
int Tensor_Dim0,
int Tensor_Dim1>
1814 "Data pointer not allocated");
1816 const size_t nb_gauss_pts = this->getGaussPts().size2();
1817 auto &mat = *this->dataPtr;
1818 if (type == this->zeroType) {
1819 mat.resize(Tensor_Dim0 * Tensor_Dim1 * Tensor_Dim1, nb_gauss_pts,
false);
1828 if (this->dataVec.use_count()) {
1829 this->dotVector.resize(nb_dofs,
false);
1830 const double *array;
1831 CHKERR VecGetArrayRead(this->dataVec, &array);
1833 for (
int i = 0;
i != local_indices.size(); ++
i)
1834 if (local_indices[
i] != -1)
1835 this->dotVector[
i] = array[local_indices[
i]];
1837 this->dotVector[
i] = 0;
1838 CHKERR VecRestoreArrayRead(this->dataVec, &array);
1842 const int nb_base_functions = data.
getN().size2();
1844 auto &hessian_base = data.
getN(BaseDerivatives::SecondDerivative);
1846 if (hessian_base.size1() != nb_gauss_pts) {
1848 "Wrong number of integration pts (%d != %d)",
1849 hessian_base.size1(), nb_gauss_pts);
1851 if (hessian_base.size2() !=
1852 nb_base_functions * Tensor_Dim1 * Tensor_Dim1) {
1854 "Wrong number of base functions (%d != %d)",
1855 hessian_base.size2() / (Tensor_Dim1 * Tensor_Dim1),
1858 if (hessian_base.size2() < nb_dofs * Tensor_Dim1 * Tensor_Dim1) {
1860 "Wrong number of base functions (%d < %d)",
1861 hessian_base.size2(), nb_dofs * Tensor_Dim1 * Tensor_Dim1);
1866 &*hessian_base.data().begin());
1868 auto t_hessian_at_gauss_pts =
1875 int size = nb_dofs / Tensor_Dim0;
1877 if (nb_dofs % Tensor_Dim0) {
1879 "Data inconsistency");
1883 for (
int gg = 0; gg < nb_gauss_pts; ++gg) {
1886 for (; bb < size; ++bb) {
1887 t_hessian_at_gauss_pts(
I,
J,
K) +=
1888 field_data(
I) * t_diff2_base_function(
J,
K);
1890 ++t_diff2_base_function;
1894 for (; bb != nb_base_functions; ++bb)
1895 ++t_diff2_base_function;
1896 ++t_hessian_at_gauss_pts;
1899 if (this->dataVec.use_count()) {
1921template <
int DIM_01,
int DIM_23,
int S = 0>
1933 boost::shared_ptr<MatrixDouble> in_mat,
1934 boost::shared_ptr<MatrixDouble> out_mat,
1935 boost::shared_ptr<MatrixDouble> d_mat)
1948 boost::shared_ptr<MatrixDouble> out_mat,
1949 boost::shared_ptr<MatrixDouble> d_mat)
1962 const size_t nb_gauss_pts =
getGaussPts().size2();
1965 outMat->resize((DIM_23 * (DIM_23 + 1)) / 2, nb_gauss_pts,
false);
1967 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
1968 t_out(
i,
j) = t_D(
i,
j,
k,
l) * t_in(
k,
l);
1983 boost::shared_ptr<MatrixDouble>
dMat;
1996 boost::shared_ptr<MatrixDouble> in_mat,
1997 boost::shared_ptr<MatrixDouble> out_mat)
2008 boost::shared_ptr<MatrixDouble> out_mat)
2019 const size_t nb_gauss_pts =
getGaussPts().size2();
2021 outMat->resize((DIM * (DIM + 1)) / 2, nb_gauss_pts,
false);
2023 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
2024 t_out(
i,
j) = (t_in(
i,
j) || t_in(
j,
i)) / 2;
2047 boost::shared_ptr<MatrixDouble> in_mat,
2048 boost::shared_ptr<MatrixDouble> out_mat)
2060 boost::shared_ptr<MatrixDouble> in_mat,
2061 boost::shared_ptr<MatrixDouble> out_mat)
2075 noalias(*
outMat) = (*scalePtr) * (*inMat);
2094template <
int Base_Dim,
int Field_Dim,
class T,
class L,
class A>
2100template <
int Field_Dim>
2106 boost::shared_ptr<MatrixDouble> data_ptr,
2108 const EntityType zero_type = MBEDGE,
2109 const int zero_side = 0)
2112 dataPtr(data_ptr), dataVec(data_vec), zeroType(zero_type),
2113 zeroSide(zero_side) {
2119 boost::shared_ptr<MatrixDouble> data_ptr,
2120 const EntityType zero_type = MBEDGE,
2121 const int zero_side = 0)
2143template <
int Field_Dim>
2145 3, Field_Dim,
double, ublas::row_major,
2149 const size_t nb_integration_points = this->getGaussPts().size2();
2150 if (type == zeroType && side == zeroSide) {
2151 dataPtr->resize(Field_Dim, nb_integration_points,
false);
2158 if (dataVec.use_count()) {
2159 dotVector.resize(nb_dofs,
false);
2160 const double *array;
2161 CHKERR VecGetArrayRead(dataVec, &array);
2163 for (
int i = 0;
i != local_indices.size(); ++
i)
2164 if (local_indices[
i] != -1)
2165 dotVector[
i] = array[local_indices[
i]];
2168 CHKERR VecRestoreArrayRead(dataVec, &array);
2172 const size_t nb_base_functions = data.
getN().size2() / 3;
2176 for (
size_t gg = 0; gg != nb_integration_points; ++gg) {
2179 for (; bb != nb_dofs; ++bb) {
2180 t_data(
i) += t_n_hdiv(
i) * t_dof;
2184 for (; bb != nb_base_functions; ++bb)
2189 if (dataVec.use_count()) {
2199template <
int Base_Dim,
int Field_Dim = Base_Dim>
2202 Base_Dim, Field_Dim, double, ublas::row_major, DoubleAllocator> {
2204 Base_Dim, Field_Dim,
double, ublas::row_major,
2211template <
int Base_Dim,
int Field_Dim = Base_Dim>
2214template <
int Field_Dim>
2219 boost::shared_ptr<MatrixDouble> data_ptr,
2220 const EntityType zero_type = MBEDGE,
2221 const int zero_side = 0)
2224 dataPtr(data_ptr), zeroType(zero_type), zeroSide(zero_side) {
2245template <
int Field_Dim>
2250 const size_t nb_integration_points = this->getGaussPts().size2();
2251 if (type == zeroType && side == zeroSide) {
2252 dataPtr->resize(Field_Dim, nb_integration_points,
false);
2257 const size_t nb_dofs = local_indices.size();
2260 std::array<double, MAX_DOFS_ON_ENTITY> dot_dofs_vector;
2261 const double *array;
2262 CHKERR VecGetArrayRead(getFEMethod()->ts_u_t, &array);
2263 for (
size_t i = 0;
i != nb_dofs; ++
i)
2264 if (local_indices[
i] != -1)
2265 dot_dofs_vector[
i] = array[local_indices[
i]];
2267 dot_dofs_vector[
i] = 0;
2268 CHKERR VecRestoreArrayRead(getFEMethod()->ts_u_t, &array);
2270 const size_t nb_base_functions = data.
getN().size2() / 3;
2274 for (
size_t gg = 0; gg != nb_integration_points; ++gg) {
2276 for (; bb != nb_dofs; ++bb) {
2277 t_data(
i) += t_n_hdiv(
i) * dot_dofs_vector[bb];
2280 for (; bb != nb_base_functions; ++bb)
2296template <
int BASE_DIM,
int SPACE_DIM>
2301 boost::shared_ptr<VectorDouble> data_ptr,
2302 const EntityType zero_type = MBEDGE,
2303 const int zero_side = 0)
2314 const size_t nb_integration_points =
getGaussPts().size2();
2316 dataPtr->resize(nb_integration_points,
false);
2322 const size_t nb_base_functions = data.
getN().size2() /
BASE_DIM;
2327 for (
size_t gg = 0; gg != nb_integration_points; ++gg) {
2330 for (; bb != nb_dofs; ++bb) {
2331 t_data += t_dof * t_n_diff_hdiv(
j,
j);
2335 for (; bb != nb_base_functions; ++bb)
2355template <
int BASE_DIM,
int SPACE_DIM>
2360 boost::shared_ptr<MatrixDouble> data_ptr,
2361 const EntityType zero_type = MBEDGE,
2362 const int zero_side = 0)
2373 const size_t nb_integration_points =
getGaussPts().size2();
2381 const size_t nb_base_functions = data.
getN().size2() /
BASE_DIM;
2386 for (
size_t gg = 0; gg != nb_integration_points; ++gg) {
2389 for (; bb != nb_dofs; ++bb) {
2390 t_data(
i,
j) += t_dof * t_base_diff(
i,
j);
2394 for (; bb != nb_base_functions; ++bb)
2414template <
int BASE_DIM,
int SPACE_DIM>
2419 boost::shared_ptr<MatrixDouble> data_ptr,
2420 const EntityType zero_type = MBEDGE,
2421 const int zero_side = 0)
2432 const size_t nb_integration_points =
getGaussPts().size2();
2442 const int nb_base_functions = data.
getN().size2() /
BASE_DIM;
2445 auto &hessian_base = data.
getN(BaseDerivatives::SecondDerivative);
2446 if (hessian_base.size1() != nb_integration_points) {
2448 "Wrong number of integration pts (%d != %d)",
2449 hessian_base.size1(), nb_integration_points);
2451 if (hessian_base.size2() !=
2454 "Wrong number of base functions (%d != %d)",
2460 "Wrong number of base functions (%d < %d)", hessian_base.size2(),
2472 for (
size_t gg = 0; gg != nb_integration_points; ++gg) {
2475 for (; bb != nb_dofs; ++bb) {
2476 t_data(
i,
j,
k) += t_dof * t_base_diff2(
i,
j,
k);
2481 for (; bb != nb_base_functions; ++bb)
2500template <
int Tensor_Dim1,
int Tensor_Dim2>
2505 boost::shared_ptr<VectorDouble> data_ptr,
2506 const EntityType zero_type = MBEDGE,
2507 const int zero_side = 0)
2518 const size_t nb_integration_points =
getGaussPts().size2();
2520 dataPtr->resize(nb_integration_points,
false);
2525 const int nb_dofs = local_indices.size();
2528 std::array<double, MAX_DOFS_ON_ENTITY> dot_dofs_vector;
2529 const double *array;
2531 for (
size_t i = 0;
i != local_indices.size(); ++
i)
2532 if (local_indices[
i] != -1)
2533 dot_dofs_vector[
i] = array[local_indices[
i]];
2535 dot_dofs_vector[
i] = 0;
2538 const size_t nb_base_functions = data.
getN().size2() / Tensor_Dim1;
2542 for (
size_t gg = 0; gg != nb_integration_points; ++gg) {
2544 for (; bb != nb_dofs; ++bb) {
2546 for (
auto ii = 0; ii != Tensor_Dim2; ++ii)
2547 div += t_n_diff_hdiv(ii, ii);
2548 t_data += dot_dofs_vector[bb] * div;
2551 for (; bb != nb_base_functions; ++bb)
2587 boost::shared_ptr<MatrixDouble> data_ptr,
2588 const EntityType zero_type = MBEDGE,
2589 const int zero_side = 0);
2611 boost::shared_ptr<MatrixDouble> data_ptr,
2612 const EntityType zero_type = MBVERTEX,
2613 const int zero_side = 0);
2636 boost::shared_ptr<MatrixDouble> data_ptr,
2637 const EntityType zero_type = MBVERTEX,
2638 const int zero_side = 0);
2656template <
int Tensor_Dim0,
int Tensor_Dim1>
2661 boost::shared_ptr<MatrixDouble> data_ptr,
2662 boost::shared_ptr<double> scale_ptr,
2664 const EntityType zero_type = MBEDGE,
2665 const int zero_side = 0)
2675 boost::shared_ptr<MatrixDouble> data_ptr,
2676 const EntityType zero_type = MBEDGE,
2677 const int zero_side = 0)
2685 const size_t nb_integration_points =
getGaussPts().size2();
2687 dataPtr->resize(Tensor_Dim0 * Tensor_Dim1, nb_integration_points,
false);
2695 const double *array;
2698 for (
int i = 0;
i != local_indices.size(); ++
i)
2699 if (local_indices[
i] != -1)
2708 const size_t nb_base_functions = data.
getN().size2() / 3;
2713 for (
size_t gg = 0; gg != nb_integration_points; ++gg) {
2716 for (; bb != nb_dofs / Tensor_Dim0; ++bb) {
2717 t_data(
i,
j) += (
scale * t_dof(
i)) * t_n_hvec(
j);
2721 for (; bb < nb_base_functions; ++bb)
2747template <
int Tensor_Dim0,
int Tensor_Dim1>
2754 const std::string
field_name, boost::shared_ptr<MatrixDouble> data_ptr,
2756 boost::shared_ptr<Range> broken_range_ptr =
nullptr,
2757 boost::shared_ptr<double> scale_ptr =
nullptr,
2758 const EntityType zero_type = MBEDGE,
const int zero_side = 0)
2786template <
int Tensor_Dim0,
int Tensor_Dim1>
2791 const size_t nb_integration_points = OP::getGaussPts().size2();
2792 if (type == zeroType) {
2793 dataPtr->resize(Tensor_Dim0 * Tensor_Dim1, nb_integration_points,
false);
2800 if (dataVec.use_count()) {
2801 dotVector.resize(nb_dofs,
false);
2802 const double *array;
2803 CHKERR VecGetArrayRead(dataVec, &array);
2805 for (
int i = 0;
i != local_indices.size(); ++
i)
2806 if (local_indices[
i] != -1)
2807 dotVector[
i] = array[local_indices[
i]];
2810 CHKERR VecRestoreArrayRead(dataVec, &array);
2821 auto get_get_side_face_dofs = [&]() {
2822 auto fe_type = OP::getFEType();
2826 std::vector<int> side_face_dofs;
2827 side_face_dofs.reserve(data.
getIndices().size() / Tensor_Dim0);
2831 auto it = side_dof_map.get<1>().begin();
2832 it != side_dof_map.get<1>().end(); ++it
2835 if ((Tensor_Dim0 * it->dof) >= data.
getIndices().size()) {
2838 if (it->type == brokenType) {
2839 if (brokenRangePtr) {
2840 auto ent = OP::getSideEntity(it->side, brokenType);
2841 if (brokenRangePtr->find(ent) != brokenRangePtr->end()) {
2842 side_face_dofs.push_back(it->dof);
2845 side_face_dofs.push_back(it->dof);
2850 return side_face_dofs;
2853 auto side_face_dofs = get_get_side_face_dofs();
2858 for (
size_t gg = 0; gg != nb_integration_points; ++gg) {
2859 for (
auto b : side_face_dofs) {
2863 t_data(
i,
j) += t_dof(
i) * t_row_base(
j);
2867 *dataPtr *= (scalePtr) ? *scalePtr : 1.0;
2869 if (dataVec.use_count()) {
2883template <
int Tensor_Dim0,
int Tensor_Dim1>
2888 const std::string
field_name, boost::shared_ptr<MatrixDouble> data_ptr,
2889 boost::shared_ptr<double> scale_ptr,
2891 const EntityType zero_type = MBEDGE,
const int zero_side = 0)
2901 boost::shared_ptr<MatrixDouble> data_ptr,
2902 const EntityType zero_type = MBEDGE,
2903 const int zero_side = 0)
2911 const size_t nb_integration_points =
getGaussPts().size2();
2913 dataPtr->resize(Tensor_Dim0 * Tensor_Dim1, nb_integration_points,
false);
2922 const double *array;
2925 for (
int i = 0;
i != local_indices.size(); ++
i)
2926 if (local_indices[
i] != -1)
2935 const size_t nb_base_functions =
2936 data.
getN().size2() / (Tensor_Dim0 * Tensor_Dim1);
2939 auto t_n_hten = data.
getFTensor2N<Tensor_Dim0, Tensor_Dim1>();
2941 for (
size_t gg = 0; gg != nb_integration_points; ++gg) {
2944 for (; bb != nb_dofs; ++bb) {
2945 t_data(
i,
j) += (
scale * t_dof) * t_n_hten(
i,
j);
2949 for (; bb < nb_base_functions; ++bb)
2977template <
int Tensor_Dim0,
int Tensor_Dim1,
2983 boost::shared_ptr<MatrixDouble> data_ptr,
2985 const EntityType zero_type = MBEDGE,
2986 const int zero_side = 0)
2996 boost::shared_ptr<MatrixDouble> data_ptr,
2997 const EntityType zero_type = MBEDGE,
2998 const int zero_side = 0)
3005 const size_t nb_integration_points =
getGaussPts().size2();
3007 dataPtr->resize(Tensor_Dim0, nb_integration_points,
false);
3015 const double *array;
3018 for (
int i = 0;
i != local_indices.size(); ++
i)
3019 if (local_indices[
i] != -1)
3027 const size_t nb_base_functions = data.
getN().size2() / 3;
3034 for (
size_t gg = 0; gg != nb_integration_points; ++gg) {
3037 for (; bb != nb_dofs / Tensor_Dim0; ++bb) {
3038 double div = t_n_diff_hvec(
j,
j);
3039 t_data(
i) += t_dof(
i) * div;
3041 t_data(
i) += t_base(0) * (t_dof(
i) / t_coords(0));
3047 for (; bb < nb_base_functions; ++bb) {
3080template <
int Tensor_Dim0,
int Tensor_Dim1,
3088 const std::string
field_name, boost::shared_ptr<MatrixDouble> data_ptr,
3090 boost::shared_ptr<Range> broken_range_ptr =
nullptr,
3091 boost::shared_ptr<double> scale_ptr =
nullptr,
3092 const EntityType zero_type = MBEDGE)
3104 const size_t nb_integration_points =
getGaussPts().size2();
3105 if (type ==
zeroType && side == 0) {
3106 dataPtr->resize(Tensor_Dim0, nb_integration_points,
false);
3115 const double *array;
3118 for (
int i = 0;
i != local_indices.size(); ++
i)
3119 if (local_indices[
i] != -1)
3134 auto get_get_side_face_dofs = [&]() {
3139 std::vector<int> side_face_dofs;
3140 side_face_dofs.reserve(data.
getIndices().size() / Tensor_Dim0);
3144 auto it = side_dof_map.get<1>().begin();
3145 it != side_dof_map.get<1>().end(); ++it
3148 if ((Tensor_Dim0 * it->dof) >= data.
getIndices().size()) {
3155 side_face_dofs.push_back(it->dof);
3158 side_face_dofs.push_back(it->dof);
3163 return side_face_dofs;
3166 auto side_face_dofs = get_get_side_face_dofs();
3172 for (
size_t gg = 0; gg != nb_integration_points; ++gg) {
3173 for (
auto b : side_face_dofs) {
3178 double div = t_diff_base(
j,
j);
3179 t_data(
i) += t_dof(
i) * div;
3181 t_data(
i) += t_base(0) * (t_dof(
i) / t_coords(0));
3212template <
int Tensor_Dim,
typename OpBase>
3216 boost::shared_ptr<MatrixDouble> data_ptr,
3217 boost::shared_ptr<double> scale_ptr,
3218 const EntityType zero_type = MBEDGE,
3219 const int zero_side = 0)
3227 boost::shared_ptr<MatrixDouble> data_ptr,
3228 const EntityType zero_type = MBEDGE,
3229 const int zero_side = 0)
3236 const size_t nb_integration_points = OpBase::getGaussPts().size2();
3237 if (type ==
zeroType && side == 0) {
3238 dataPtr->resize(Tensor_Dim, nb_integration_points,
false);
3244 auto t_normal = OpBase::getFTensor1NormalsAtGaussPts();
3245 const size_t nb_base_functions = data.
getN().size2() / 3;
3248 for (
size_t gg = 0; gg != nb_integration_points; ++gg) {
3250 t_normalized_normal(
j) = t_normal(
j);
3254 for (; bb != nb_dofs / Tensor_Dim; ++bb) {
3256 (scale_val * t_dof(
i)) * (t_base(
j) * t_normalized_normal(
j));
3260 for (; bb < nb_base_functions; ++bb) {
3303 boost::shared_ptr<MatrixDouble> inv_jac_ptr);
3317 template <
int D1,
int D2,
int J1,
int J2>
3321 static_assert(D2 == J2,
"Dimension of jacobian and dimension of <out> "
3322 "directive does not match");
3324 size_t nb_functions = diff_n.size2() / D1;
3326 size_t nb_gauss_pts = diff_n.size1();
3331 "Wrong number of Gauss Pts");
3332 if (diff_n.size2() % D1)
3334 "Number of directives of base functions and D1 dimension does "
3338 diffNinvJac.resize(diff_n.size1(), D2 * nb_functions,
false);
3345 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg, ++t_inv_jac) {
3346 for (
size_t dd = 0; dd != nb_functions; ++dd) {
3347 t_diff_n(
i) = t_inv_jac(
k,
i) * t_diff_n_ref(
k);
3392template <
int DERIVARIVE = 1>
3399template <
int DERIVARIVE = 1>
3406template <
int DERIVARIVE = 1>
3410 boost::shared_ptr<MatrixDouble> inv_jac_ptr)
3414template <
int DERIVARIVE = 1>
3418 boost::shared_ptr<MatrixDouble> inv_jac_ptr)
3436 invJacPtr(inv_jac_ptr) {}
3490 boost::shared_ptr<MatrixDouble> inv_jac_ptr);
3527 boost::shared_ptr<MatrixDouble> jac_ptr)
3544 2>::OpSetContravariantPiolaTransformOnFace2DImpl;
3727 boost::shared_ptr<VectorDouble> det_ptr,
3728 boost::shared_ptr<MatrixDouble> out_ptr)
3748 "Pointer for inPtr matrix not allocated");
3751 "Pointer for detPtr matrix not allocated");
3753 const auto nb_rows = inPtr->size1();
3754 const auto nb_integration_pts = inPtr->size2();
3758 detPtr->resize(nb_integration_pts,
false);
3761 for (
size_t gg = 0; gg != nb_integration_pts; ++gg) {
3770 outPtr->resize(nb_rows, nb_integration_pts,
false);
3774 for (
size_t gg = 0; gg != nb_integration_pts; ++gg) {
3797 boost::shared_ptr<VectorDouble> out_ptr)
3818 "Pointer for inPtr matrix not allocated");
3820 const auto nb_integration_pts = inPtr->size2();
3823 outPtr->resize(nb_integration_pts,
false);
3827 for (
size_t gg = 0; gg != nb_integration_pts; ++gg) {
3850 boost::shared_ptr<VectorDouble> out_ptr)
3871 "Pointer for inPtr matrix not allocated");
3873 const auto nb_integration_pts = inPtr->size2();
3876 outPtr->resize(nb_integration_pts,
false);
3880 for (
size_t gg = 0; gg != nb_integration_pts; ++gg) {
ForcesAndSourcesCore::UserDataOperator UserDataOperator
Tensor1< T, Tensor_Dim > normalize()
static const char *const CoordinateTypesNames[]
Coordinate system names.
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
FieldSpace
approximation spaces
@ L2
field with C-1 continuity
@ HCURL
field with continuous tangents
#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 CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
#define MOFEM_LOG(channel, severity)
Log.
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'J', DIM1 > J
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
VecAllocator< double > DoubleAllocator
implementation of Data Operators for Forces and Sources
static FTensor::Tensor3< FTensor::PackPtr< T *, 1 >, Tensor_Dim0, Tensor_Dim1, Tensor_Dim2 > getFTensor3FromMat(ublas::matrix< T, L, A > &data)
Get tensor rank 3 (non symmetries) form data matrix.
static FTensor::Dg< FTensor::PackPtr< T *, 1 >, Tensor_Dim01, Tensor_Dim2 > getFTensor3DgFromMat(ublas::matrix< T, L, A > &data)
Get symmetric tensor rank 3 on the first two indices from form data matrix.
auto getFTensor2FromPtr(double *ptr)
Make Tensor2 from pointer.
static FTensor::Ddg< FTensor::PackPtr< T *, 1 >, Tensor_Dim01, Tensor_Dim23 > getFTensor4DdgFromMat(ublas::matrix< T, L, A > &data)
Get symmetric tensor rank 4 on first two and last indices from form data matrix.
FTensor::Tensor1< FTensor::PackPtr< T *, S >, Tensor_Dim > getFTensor1FromMat(ublas::matrix< T, L, A > &data)
Get tensor rank 1 (vector) form data matrix.
static MoFEMErrorCode invertTensor(FTensor::Tensor2< T1, DIM, DIM > &t, T2 &det, FTensor::Tensor2< T3, DIM, DIM > &inv_t)
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.
FTensor::Tensor1< FTensor::PackPtr< double *, S >, DIM > getFTensor1FromPtr(double *ptr)
Make Tensor1 from pointer.
static auto determinantTensor(FTensor::Tensor2< T, DIM, DIM > &t)
Calculate the determinant of a tensor of rank DIM.
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
OpCalculateScalarFieldValuesFromPetscVecImpl< PetscData::CTX_SET_X_T > OpCalculateScalarFieldValuesDot
auto getFTensor1FromArray(VectorDouble &data)
Get FTensor1 from array.
constexpr IntegrationType I
constexpr auto field_name
multi_index_container< DofsSideMapData, indexed_by< ordered_non_unique< tag< TypeSide_mi_tag >, composite_key< DofsSideMapData, member< DofsSideMapData, EntityType, &DofsSideMapData::type >, member< DofsSideMapData, int, &DofsSideMapData::side > > >, ordered_unique< tag< EntDofIdx_mi_tag >, member< DofsSideMapData, int, &DofsSideMapData::dof > > > > DofsSideMap
Map entity stype and side to element/entity dof index.
std::array< bool, MBMAXTYPE > doEntities
If true operator is executed for entity.
default operator for EDGE element
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::Tensor2< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 >, Tensor_Dim0, Tensor_Dim1 > getFTensor2N(FieldApproximationBase base)
Get base functions for Hdiv/Hcurl spaces.
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.
MatrixDouble & getDiffN(const FieldApproximationBase base)
get derivatives of base functions
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....
const VectorDouble & getFieldData() const
get dofs values
FTensor::Tensor3< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 *Tensor_Dim2 >, Tensor_Dim0, Tensor_Dim1, Tensor_Dim2 > getFTensor3Diff2N(FieldApproximationBase base)
Get second derivatives of base functions for Hvec space.
const VectorInt & getLocalIndices() const
get local indices of dofs on entity
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1FieldData()
Return FTensor of rank 1, i.e. vector from field data coefficients.
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0FieldData()
Resturn scalar files as a FTensor of rank 0.
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1N(FieldApproximationBase base)
Get base functions for Hdiv/Hcurl spaces.
FTensor::Tensor2_symmetric< FTensor::PackPtr< double *,(Tensor_Dim *(Tensor_Dim+1))/2 >, Tensor_Dim > getFTensor2SymmetricFieldData()
Return symmetric FTensor rank 2, i.e. matrix from field data coefficients.
FTensor::Tensor2< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 >, Tensor_Dim0, Tensor_Dim1 > getFTensor2FieldData()
Return FTensor rank 2, i.e. matrix from field data coefficients.
const VectorInt & getIndices() const
Get global indices of dofs on entity.
default operator for TRI element
default operator for Flat Prism element
default operator for Flat Prism element
FlatPrism finite element.
EntityType getFEType() const
Get dimension of finite element.
auto getFTensor1CoordsAtGaussPts()
Get coordinates at integration points assuming linear geometry.
EntityHandle getSideEntity(const int side_number, const EntityType type)
Get the side entity.
@ OPROW
operator doWork function is executed on FE rows
@ OPSPACE
operator do Work is execute on space data
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
MatrixDouble & getGaussPts()
matrix of integration (Gauss) points for Volume Element
structure to get information form mofem into EntitiesFieldData
friend class UserDataOperator
Calculate divergence of tonsorial field using vectorial base.
boost::shared_ptr< double > scalePtr
boost::shared_ptr< Range > brokenRangePtr
SmartPetscObj< Vec > dataVec
OpCalculateBrokenHVecTensorDivergence(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, SmartPetscObj< Vec > data_vec, EntityType broken_type, boost::shared_ptr< Range > broken_range_ptr=nullptr, boost::shared_ptr< double > scale_ptr=nullptr, const EntityType zero_type=MBEDGE)
boost::shared_ptr< MatrixDouble > dataPtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
const EntityHandle zeroType
Get tensor field for H-div approximation.
boost::shared_ptr< MatrixDouble > dataPtr
boost::shared_ptr< Range > brokenRangePtr
boost::shared_ptr< double > scalePtr
OpCalculateBrokenHVecTensorField(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, SmartPetscObj< Vec > data_vec, EntityType broken_type, boost::shared_ptr< Range > broken_range_ptr=nullptr, boost::shared_ptr< double > scale_ptr=nullptr, const EntityType zero_type=MBEDGE, const int zero_side=0)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Calculate values of vector field at integration points.
const EntityHandle zeroType
SmartPetscObj< Vec > dataVec
Calculate field values (template specialization) for tensor field rank 1, i.e. vector field.
boost::shared_ptr< VectorDouble > dataPtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpCalculateDivergenceVectorFieldValues(const std::string field_name, boost::shared_ptr< VectorDouble > data_ptr, const EntityType zero_type=MBVERTEX)
const EntityHandle zeroType
Calculate tenor field using tensor base, i.e. Hdiv/Hcurl.
VectorDouble dotVector
Keeps temporary values of time derivatives.
OpCalculateHTensorTensorField(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, boost::shared_ptr< double > scale_ptr, SmartPetscObj< Vec > data_vec=SmartPetscObj< Vec >(), const EntityType zero_type=MBEDGE, const int zero_side=0)
SmartPetscObj< Vec > dataVec
boost::shared_ptr< MatrixDouble > dataPtr
const EntityHandle zeroType
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpCalculateHTensorTensorField(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, const EntityType zero_type=MBEDGE, const int zero_side=0)
boost::shared_ptr< double > scalePtr
Calculate divergence of tonsorial field using vectorial base.
boost::shared_ptr< MatrixDouble > dataPtr
SmartPetscObj< Vec > dataVec
VectorDouble dotVector
Keeps temporary values of time derivatives.
const EntityHandle zeroType
OpCalculateHVecTensorDivergence(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, SmartPetscObj< Vec > data_vec, const EntityType zero_type=MBEDGE, const int zero_side=0)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpCalculateHVecTensorDivergence(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, const EntityType zero_type=MBEDGE, const int zero_side=0)
Calculate tenor field using vectorial base, i.e. Hdiv/Hcurl.
OpCalculateHVecTensorField(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, boost::shared_ptr< double > scale_ptr, SmartPetscObj< Vec > data_vec=SmartPetscObj< Vec >(), const EntityType zero_type=MBEDGE, const int zero_side=0)
const EntityHandle zeroType
boost::shared_ptr< MatrixDouble > dataPtr
OpCalculateHVecTensorField(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, const EntityType zero_type=MBEDGE, const int zero_side=0)
SmartPetscObj< Vec > dataVec
VectorDouble dotVector
Keeps temporary values of time derivatives.
boost::shared_ptr< double > scalePtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Calculate trace of vector (Hdiv/Hcurl) space.
const EntityHandle zeroType
boost::shared_ptr< MatrixDouble > dataPtr
boost::shared_ptr< double > scalePtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
OpCalculateHVecTensorTrace(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, const EntityType zero_type=MBEDGE, const int zero_side=0)
FTensor::Index< 'j', Tensor_Dim > j
FTensor::Index< 'i', Tensor_Dim > i
OpCalculateHVecTensorTrace(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, boost::shared_ptr< double > scale_ptr, const EntityType zero_type=MBEDGE, const int zero_side=0)
const EntityHandle zeroType
boost::shared_ptr< MatrixDouble > dataPtr
OpCalculateHVecVectorFieldDot(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, const EntityType zero_type=MBEDGE, const int zero_side=0)
Get vector field for H-div approximation.
OpCalculateHVecVectorField_General(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, SmartPetscObj< Vec > data_vec, const EntityType zero_type=MBEDGE, const int zero_side=0)
boost::shared_ptr< MatrixDouble > dataPtr
OpCalculateHVecVectorField_General(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, const EntityType zero_type=MBEDGE, const int zero_side=0)
SmartPetscObj< Vec > dataVec
const EntityHandle zeroType
Get vector field for H-div approximation.
Get vector field for H-div approximation.
Calculate gradient of vector field.
boost::shared_ptr< MatrixDouble > dataPtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
const EntityHandle zeroType
OpCalculateHVecVectorGradient(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, const EntityType zero_type=MBEDGE, const int zero_side=0)
Calculate gradient of vector field.
boost::shared_ptr< MatrixDouble > dataPtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
const EntityHandle zeroType
OpCalculateHVecVectorHessian(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, const EntityType zero_type=MBEDGE, const int zero_side=0)
boost::shared_ptr< MatrixDouble > dataPtr
const EntityHandle zeroType
const EntityHandle zeroType
boost::shared_ptr< MatrixDouble > dataPtr
boost::shared_ptr< MatrixDouble > dataPtr
const EntityHandle zeroType
Calculate curl of vector field.
Calculate divergence of vector field dot.
boost::shared_ptr< VectorDouble > dataPtr
OpCalculateHdivVectorDivergenceDot(const std::string field_name, boost::shared_ptr< VectorDouble > data_ptr, const EntityType zero_type=MBEDGE, const int zero_side=0)
const EntityHandle zeroType
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Calculate divergence of vector field.
boost::shared_ptr< VectorDouble > dataPtr
OpCalculateHdivVectorDivergence(const std::string field_name, boost::shared_ptr< VectorDouble > data_ptr, const EntityType zero_type=MBEDGE, const int zero_side=0)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
const EntityHandle zeroType
Calculate inverse of jacobian for face element.
OpCalculateInvJacForFatPrism(MatrixDouble &inv_jac)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpCalculateInvJacForFatPrism(boost::shared_ptr< MatrixDouble > inv_jac_ptr)
const boost::shared_ptr< MatrixDouble > invJacPtr
Calculate inverse of jacobian for face element.
OpCalculateInvJacForFlatPrism(MatrixDouble &inv_jac_f3)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
calculate gradient values of scalar field at integration points
Evaluate field gradient values for scalar field, i.e. gradient is tensor rank 1 (vector)
Get field gradients at integration pts for scalar field rank 0, i.e. vector field.
Evaluate field gradient values for scalar field, i.e. gradient is tensor rank 1 (vector),...
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
calculate gradient values of scalar field at integration points
Get rate of scalar field at integration points.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
boost::shared_ptr< VectorDouble > dataPtr
OpCalculateScalarFieldValuesFromPetscVecImpl(const std::string field_name, boost::shared_ptr< VectorDouble > data_ptr, const EntityType zero_at_type=MBVERTEX)
const EntityHandle zeroAtType
Scalar field values at integration points.
OpCalculateScalarFieldValues_General(const std::string field_name, boost::shared_ptr< ublas::vector< T, A > > data_ptr, SmartPetscObj< Vec > data_vec, const EntityType zero_type=MBVERTEX)
boost::shared_ptr< ublas::vector< T, A > > dataPtr
SmartPetscObj< Vec > dataVec
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
calculate values of scalar field at integration points
const EntityHandle zeroType
OpCalculateScalarFieldValues_General(const std::string field_name, boost::shared_ptr< ublas::vector< T, A > > data_ptr, const EntityType zero_type=MBVERTEX)
Get value at integration points for scalar field.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
calculate values of scalar field at integration points
Get time direvarive values at integration pts for tensor field rank 2, i.e. matrix field.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpCalculateTensor2FieldValuesDot(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, const EntityType zero_at_type=MBVERTEX)
EntityType zeroAtType
Zero values at Gauss point at this type.
VectorDouble dotVector
Keeps temporary values of time derivatives.
boost::shared_ptr< MatrixDouble > dataPtr
Data computed into this matrix.
SmartPetscObj< Vec > dataVec
OpCalculateTensor2FieldValues_General(const std::string field_name, boost::shared_ptr< ublas::matrix< double, ublas::row_major, DoubleAllocator > > data_ptr, SmartPetscObj< Vec > data_vec, const EntityType zero_type=MBVERTEX)
OpCalculateTensor2FieldValues_General(const std::string field_name, boost::shared_ptr< ublas::matrix< double, ublas::row_major, DoubleAllocator > > data_ptr, const EntityType zero_type=MBVERTEX)
boost::shared_ptr< ublas::matrix< double, ublas::row_major, DoubleAllocator > > dataPtr
const EntityHandle zeroType
Calculate field values for tenor field rank 2.
const EntityHandle zeroType
OpCalculateTensor2FieldValues_General(const std::string field_name, boost::shared_ptr< ublas::matrix< T, L, A > > data_ptr, const EntityType zero_type=MBVERTEX)
boost::shared_ptr< ublas::matrix< T, L, A > > dataPtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Get values at integration pts for tensor field rank 2, i.e. matrix field.
OpCalculateTensor2SymmetricFieldGradient_General(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, const EntityType zero_type=MBVERTEX)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
calculate values of vector field at integration points
Evaluate field gradient values for symmetric 2nd order tensor field, i.e. gradient is tensor rank 3.
OpCalculateTensor2SymmetricFieldGradient_General(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, const EntityType zero_type=MBVERTEX)
Get field gradients at integration pts for symmetric tensorial field rank 2.
OpCalculateTensor2SymmetricFieldGradient(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, const EntityType zero_type=MBVERTEX)
Calculate symmetric tensor field rates ant integratio pts.
const EntityHandle zeroType
boost::shared_ptr< MatrixDouble > dataPtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpCalculateTensor2SymmetricFieldValuesDot(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, const EntityType zero_type=MBEDGE, const int zero_side=0)
Calculate symmetric tensor field values at integration pts.
OpCalculateTensor2SymmetricFieldValues(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, SmartPetscObj< Vec > data_vec, const EntityType zero_type=MBEDGE, const int zero_side=0)
OpCalculateTensor2SymmetricFieldValues(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, const EntityType zero_type=MBEDGE, const int zero_side=0)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
SmartPetscObj< Vec > dataVec
boost::shared_ptr< MatrixDouble > dataPtr
const EntityHandle zeroType
Calculates the trace of an input matrix.
boost::shared_ptr< MatrixDouble > inPtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
FTensor::Index< 'i', DIM > i
boost::shared_ptr< VectorDouble > outPtr
OpCalculateTraceFromMat(boost::shared_ptr< MatrixDouble > in_ptr, boost::shared_ptr< VectorDouble > out_ptr)
Calculates the trace of an input matrix.
FTensor::Index< 'i', DIM > i
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
boost::shared_ptr< VectorDouble > outPtr
boost::shared_ptr< MatrixDouble > inPtr
OpCalculateTraceFromSymmMat(boost::shared_ptr< MatrixDouble > in_ptr, boost::shared_ptr< VectorDouble > out_ptr)
Get field gradients time derivative at integration pts for scalar field rank 0, i....
boost::shared_ptr< MatrixDouble > dataPtr
Data computed into this matrix.
OpCalculateVectorFieldGradientDot(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, const EntityType zero_at_type=MBVERTEX)
EntityType zeroAtType
Zero values at Gauss point at this type.
VectorDouble dotVector
Keeps temporary values of time derivatives.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
calculate values of vector field at integration points
Evaluate field gradient values for vector field, i.e. gradient is tensor rank 2.
Get field gradients at integration pts for scalar field rank 0, i.e. vector field.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
calculate values of vector field at integration points
Approximate field values for given petsc vector.
boost::shared_ptr< MatrixDouble > dataPtr
OpCalculateVectorFieldValuesFromPetscVecImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, const EntityType zero_at_type=MBVERTEX, bool throw_error=true)
const EntityHandle zeroAtType
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
boost::shared_ptr< MatrixDouble > dataPtr
OpCalculateVectorFieldValues_General(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, const EntityType zero_type=MBVERTEX)
OpCalculateVectorFieldValues_General(const std::string field_name, boost::shared_ptr< MatrixDouble > data_ptr, SmartPetscObj< Vec > data_vec, const EntityType zero_type=MBVERTEX)
SmartPetscObj< Vec > dataVec
const EntityHandle zeroType
Calculate field values for tenor field rank 1, i.e. vector field.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
calculate values of vector field at integration points
OpCalculateVectorFieldValues_General(const std::string field_name, boost::shared_ptr< ublas::matrix< T, L, A > > data_ptr, const EntityType zero_type=MBVERTEX)
const EntityHandle zeroType
boost::shared_ptr< ublas::matrix< T, L, A > > dataPtr
Get values at integration pts for tensor field rank 1, i.e. vector field.
boost::shared_ptr< VectorDouble > detPtr
OpInvertMatrix(boost::shared_ptr< MatrixDouble > in_ptr, boost::shared_ptr< VectorDouble > det_ptr, boost::shared_ptr< MatrixDouble > out_ptr)
boost::shared_ptr< MatrixDouble > outPtr
boost::shared_ptr< MatrixDouble > inPtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Make Hdiv space from Hcurl space in 2d.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Operator for fat prism element updating integration weights in the volume.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpMultiplyDeterminantOfJacobianAndWeightsForFatPrisms()
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
boost::shared_ptr< MatrixDouble > outMat
DEPRECATED OpScaleMatrix(const std::string field_name, const double scale, boost::shared_ptr< MatrixDouble > in_mat, boost::shared_ptr< MatrixDouble > out_mat)
boost::shared_ptr< double > scalePtr
OpScaleMatrix(boost::shared_ptr< double > scale_ptr, boost::shared_ptr< MatrixDouble > in_mat, boost::shared_ptr< MatrixDouble > out_mat)
boost::shared_ptr< MatrixDouble > inMat
OpSetInvJacH1ForFaceEmbeddedIn3DSpace(boost::shared_ptr< MatrixDouble > inv_jac_ptr)
OpSetInvJacH1ForFace(boost::shared_ptr< MatrixDouble > inv_jac_ptr)
Transform local reference derivatives of shape functions to global derivatives.
const boost::shared_ptr< MatrixDouble > invJacPtr
OpSetInvJacH1ForFatPrism(MatrixDouble &inv_jac)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpSetInvJacH1ForFatPrism(boost::shared_ptr< MatrixDouble > inv_jac_ptr)
Transform local reference derivatives of shape functions to global derivatives.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpSetInvJacH1ForFlatPrism(MatrixDouble &inv_jac_f3)
boost::shared_ptr< MatrixDouble > invJacPtr
MatrixDouble diffHcurlInvJac
OpSetInvJacHcurlFaceImpl(boost::shared_ptr< MatrixDouble > inv_jac_ptr)
Transform local reference derivatives of shape function to global derivatives for face.
OpSetInvJacL2ForFaceEmbeddedIn3DSpace(boost::shared_ptr< MatrixDouble > inv_jac_ptr)
OpSetInvJacL2ForFace(boost::shared_ptr< MatrixDouble > inv_jac_ptr)
Transform local reference derivatives of shape functions to global derivatives.
MoFEMErrorCode applyTransform(MatrixDouble &diff_n)
Apply transformation to the input matrix.
OpSetInvJacToScalarBasesBasic(FieldSpace space, boost::shared_ptr< MatrixDouble > inv_jac_ptr)
boost::shared_ptr< MatrixDouble > invJacPtr
FTensor::Index< 'i', DIM > i
boost::shared_ptr< MatrixDouble > outMat
OpSymmetrizeTensor(boost::shared_ptr< MatrixDouble > in_mat, boost::shared_ptr< MatrixDouble > out_mat)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
FTensor::Index< 'j', DIM > j
DEPRECATED OpSymmetrizeTensor(const std::string field_name, boost::shared_ptr< MatrixDouble > in_mat, boost::shared_ptr< MatrixDouble > out_mat)
boost::shared_ptr< MatrixDouble > inMat
FTensor::Index< 'i', DIM_01 > i
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
FTensor::Index< 'k', DIM_23 > k
boost::shared_ptr< MatrixDouble > dMat
OpTensorTimesSymmetricTensor(boost::shared_ptr< MatrixDouble > in_mat, boost::shared_ptr< MatrixDouble > out_mat, boost::shared_ptr< MatrixDouble > d_mat)
FTensor::Index< 'j', DIM_01 > j
boost::shared_ptr< MatrixDouble > inMat
FTensor::Index< 'l', DIM_23 > l
DEPRECATED OpTensorTimesSymmetricTensor(const std::string field_name, boost::shared_ptr< MatrixDouble > in_mat, boost::shared_ptr< MatrixDouble > out_mat, boost::shared_ptr< MatrixDouble > d_mat)
boost::shared_ptr< MatrixDouble > outMat
static constexpr Switches CtxSetX
static constexpr Switches CtxSetX_TT
std::bitset< 8 > Switches
static constexpr Switches CtxSetX_T
intrusive_ptr for managing petsc objects