8#ifndef __HENCKY_OPS_HPP__
9#define __HENCKY_OPS_HPP__
13constexpr double eps = std::numeric_limits<float>::epsilon();
15auto f = [](
double v) {
return 0.5 * std::log(
v); };
16auto d_f = [](
double v) {
return 0.5 /
v; };
17auto dd_f = [](
double v) {
return -0.5 / (
v *
v); };
20 static inline auto check(
const double &
a,
const double &b) {
28inline auto is_eq(
const double &
a,
const double &b) {
33 std::array<double, DIM> tmp;
34 std::copy(ptr, &ptr[DIM], tmp.begin());
35 std::sort(tmp.begin(), tmp.end());
36 isEq::absMax = std::max(std::abs(tmp[0]), std::abs(tmp[DIM - 1]));
37 return std::distance(tmp.begin(), std::unique(tmp.begin(), tmp.end(),
is_eq));
45 std::max(std::max(std::abs(eig(0)), std::abs(eig(1))), std::abs(eig(2)));
47 int i = 0,
j = 1,
k = 2;
49 if (
is_eq(eig(0), eig(1))) {
53 }
else if (
is_eq(eig(0), eig(2))) {
57 }
else if (
is_eq(eig(1), eig(2))) {
64 eigen_vec(
i, 0), eigen_vec(
i, 1), eigen_vec(
i, 2),
66 eigen_vec(
j, 0), eigen_vec(
j, 1), eigen_vec(
j, 2),
68 eigen_vec(
k, 0), eigen_vec(
k, 1), eigen_vec(
k, 2)};
76 eigen_vec(
i,
j) = eigen_vec_c(
i,
j);
80struct CommonData :
public boost::enable_shared_from_this<CommonData> {
95 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
100 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
105 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &
matLogC);
109 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &
matTangent);
113template <
int DIM, IntegrationType I,
typename DomainEleOp>
116template <
int DIM, IntegrationType I,
typename DomainEleOp>
119template <
int DIM, IntegrationType I,
typename DomainEleOp>
122template <
int DIM, IntegrationType I,
typename DomainEleOp,
int S>
125template <
int DIM, IntegrationType I,
typename DomainEleOp,
int S>
128template <
int DIM, IntegrationType I,
typename DomainEleOp,
int S>
131template <
int DIM, IntegrationType I,
typename DomainEleOp,
int S>
134template <
int DIM, IntegrationType I,
typename DomainEleOp,
int S>
137template <
int DIM, IntegrationType I,
typename DomainEleOp,
int S>
140template <
int DIM, IntegrationType I,
typename DomainEleOp,
int S>
143template <
int DIM,
typename DomainEleOp>
147 boost::shared_ptr<CommonData> common_data)
149 commonDataPtr(common_data) {
150 std::fill(&DomainEleOp::doEntities[MBEDGE],
151 &DomainEleOp::doEntities[MBMAXTYPE],
false);
163 const size_t nb_gauss_pts = DomainEleOp::getGaussPts().size2();
165 auto t_grad = getFTensor2FromMat<DIM, DIM>(*(commonDataPtr->matGradPtr));
167 commonDataPtr->matEigVal.resize(DIM, nb_gauss_pts,
false);
168 commonDataPtr->matEigVec.resize(DIM * DIM, nb_gauss_pts,
false);
169 auto t_eig_val = getFTensor1FromMat<DIM>(commonDataPtr->matEigVal);
170 auto t_eig_vec = getFTensor2FromMat<DIM, DIM>(commonDataPtr->matEigVec);
172 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
182 for (
int ii = 0; ii != DIM; ii++)
183 for (
int jj = 0; jj != DIM; jj++)
184 eigen_vec(ii, jj) = C(ii, jj);
186 CHKERR computeEigenValuesSymmetric(eigen_vec, eig);
187 for (
auto ii = 0; ii != DIM; ++ii)
188 eig(ii) = std::max(std::numeric_limits<double>::epsilon(), eig(ii));
192 if constexpr (DIM == 3) {
198 t_eig_val(
i) = eig(
i);
199 t_eig_vec(
i,
j) = eigen_vec(
i,
j);
213template <
int DIM,
typename DomainEleOp>
217 boost::shared_ptr<CommonData> common_data)
219 commonDataPtr(common_data) {
220 std::fill(&DomainEleOp::doEntities[MBEDGE],
221 &DomainEleOp::doEntities[MBMAXTYPE],
false);
231 const size_t nb_gauss_pts = DomainEleOp::getGaussPts().size2();
232 constexpr auto size_symm = (DIM * (DIM + 1)) / 2;
233 commonDataPtr->matLogC.resize(
size_symm, nb_gauss_pts,
false);
235 auto t_eig_val = getFTensor1FromMat<DIM>(commonDataPtr->matEigVal);
236 auto t_eig_vec = getFTensor2FromMat<DIM, DIM>(commonDataPtr->matEigVec);
238 auto t_logC = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matLogC);
240 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
254template <
int DIM,
typename DomainEleOp>
258 boost::shared_ptr<CommonData> common_data)
260 commonDataPtr(common_data) {
261 std::fill(&DomainEleOp::doEntities[MBEDGE],
262 &DomainEleOp::doEntities[MBMAXTYPE],
false);
273 const size_t nb_gauss_pts = DomainEleOp::getGaussPts().size2();
274 constexpr auto size_symm = (DIM * (DIM + 1)) / 2;
276 auto t_logC_dC = getFTensor4DdgFromMat<DIM, DIM>(commonDataPtr->matLogCdC);
277 auto t_eig_val = getFTensor1FromMat<DIM>(commonDataPtr->matEigVal);
278 auto t_eig_vec = getFTensor2FromMat<DIM, DIM>(commonDataPtr->matEigVec);
280 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
283 t_logC_dC(
i,
j,
k,
l) =
285 nb_uniq)(
i,
j,
k,
l);
299template <
int DIM,
typename DomainEleOp,
int S>
304 boost::shared_ptr<CommonData> common_data)
306 commonDataPtr(common_data) {
307 std::fill(&DomainEleOp::doEntities[MBEDGE],
308 &DomainEleOp::doEntities[MBMAXTYPE],
false);
320 const size_t nb_gauss_pts = DomainEleOp::getGaussPts().size2();
321 auto t_D = getFTensor4DdgFromMat<DIM, DIM, S>(*commonDataPtr->matDPtr);
322 auto t_logC = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matLogC);
323 constexpr auto size_symm = (DIM * (DIM + 1)) / 2;
324 commonDataPtr->matHenckyStress.resize(
size_symm, nb_gauss_pts,
false);
325 auto t_T = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matHenckyStress);
327 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
328 t_T(
i,
j) = t_D(
i,
j,
k,
l) * t_logC(
k,
l);
341template <
int DIM,
typename DomainEleOp,
int S>
346 const std::string
field_name, boost::shared_ptr<VectorDouble> temperature,
347 boost::shared_ptr<CommonData> common_data,
348 boost::shared_ptr<VectorDouble> coeff_expansion_ptr,
349 boost::shared_ptr<double> ref_temp_ptr)
351 commonDataPtr(common_data), coeffExpansionPtr(coeff_expansion_ptr),
352 refTempPtr(ref_temp_ptr) {
353 std::fill(&DomainEleOp::doEntities[MBEDGE],
354 &DomainEleOp::doEntities[MBMAXTYPE],
false);
368 const size_t nb_gauss_pts = DomainEleOp::getGaussPts().size2();
369 auto t_D = getFTensor4DdgFromMat<DIM, DIM, S>(*commonDataPtr->matDPtr);
370 auto t_logC = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matLogC);
371 auto t_logC_dC = getFTensor4DdgFromMat<DIM, DIM>(commonDataPtr->matLogCdC);
372 constexpr auto size_symm = (DIM * (DIM + 1)) / 2;
373 commonDataPtr->matHenckyStress.resize(
size_symm, nb_gauss_pts,
false);
374 commonDataPtr->matFirstPiolaStress.resize(DIM * DIM, nb_gauss_pts,
false);
375 commonDataPtr->matSecondPiolaStress.resize(
size_symm, nb_gauss_pts,
false);
376 auto t_T = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matHenckyStress);
377 auto t_P = getFTensor2FromMat<DIM, DIM>(commonDataPtr->matFirstPiolaStress);
379 getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matSecondPiolaStress);
380 auto t_grad = getFTensor2FromMat<DIM, DIM>(*(commonDataPtr->matGradPtr));
381 auto t_temp = getFTensor0FromVec(*tempPtr);
384 t_coeff_exp(
i,
j) = 0;
386 t_coeff_exp(d, d) = (*coeffExpansionPtr)[d];
389 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
390#ifdef HENCKY_SMALL_STRAIN
391 t_P(
i,
j) = t_D(
i,
j,
k,
l) *
392 (t_grad(
k,
l) - t_coeff_exp(
k,
l) * (t_temp - (*refTempPtr)));
394 t_T(
i,
j) = t_D(
i,
j,
k,
l) *
395 (t_logC(
k,
l) - t_coeff_exp(
k,
l) * (t_temp - (*refTempPtr)));
398 t_S(
k,
l) = t_T(
i,
j) * t_logC_dC(
i,
j,
k,
l);
399 t_P(
i,
l) = t_F(
i,
k) * t_S(
k,
l);
421template <
int DIM,
typename DomainEleOp,
int S>
426 boost::shared_ptr<CommonData> common_data,
427 boost::shared_ptr<MatrixDouble> mat_D_ptr,
428 const double scale = 1)
430 scaleStress(
scale), matDPtr(mat_D_ptr) {
431 std::fill(&DomainEleOp::doEntities[MBEDGE],
432 &DomainEleOp::doEntities[MBMAXTYPE],
false);
434 matLogCPlastic = commonDataPtr->matLogCPlastic;
446 const size_t nb_gauss_pts = DomainEleOp::getGaussPts().size2();
447 auto t_D = getFTensor4DdgFromMat<DIM, DIM, S>(*matDPtr);
448 auto t_logC = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matLogC);
449 auto t_logCPlastic = getFTensor2SymmetricFromMat<DIM>(*matLogCPlastic);
450 constexpr auto size_symm = (DIM * (DIM + 1)) / 2;
451 commonDataPtr->matHenckyStress.resize(
size_symm, nb_gauss_pts,
false);
452 auto t_T = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matHenckyStress);
454 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
455 t_T(
i,
j) = t_D(
i,
j,
k,
l) * (t_logC(
k,
l) - t_logCPlastic(
k,
l));
456 t_T(
i,
j) /= scaleStress;
473template <
int DIM,
typename DomainEleOp,
int S>
478 boost::shared_ptr<CommonData> common_data)
480 commonDataPtr(common_data) {
481 std::fill(&DomainEleOp::doEntities[MBEDGE],
482 &DomainEleOp::doEntities[MBMAXTYPE],
false);
496 const size_t nb_gauss_pts = DomainEleOp::getGaussPts().size2();
497#ifdef HENCKY_SMALL_STRAIN
498 auto t_D = getFTensor4DdgFromMat<DIM, DIM, S>(*commonDataPtr->matDPtr);
500 auto t_logC = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matLogC);
501 auto t_logC_dC = getFTensor4DdgFromMat<DIM, DIM>(commonDataPtr->matLogCdC);
502 constexpr auto size_symm = (DIM * (DIM + 1)) / 2;
503 commonDataPtr->matFirstPiolaStress.resize(DIM * DIM, nb_gauss_pts,
false);
504 commonDataPtr->matSecondPiolaStress.resize(
size_symm, nb_gauss_pts,
false);
505 auto t_P = getFTensor2FromMat<DIM, DIM>(commonDataPtr->matFirstPiolaStress);
506 auto t_T = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matHenckyStress);
508 getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matSecondPiolaStress);
509 auto t_grad = getFTensor2FromMat<DIM, DIM>(*(commonDataPtr->matGradPtr));
511 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
513#ifdef HENCKY_SMALL_STRAIN
514 t_P(
i,
j) = t_D(
i,
j,
k,
l) * t_grad(
k,
l);
518 t_S(
k,
l) = t_T(
i,
j) * t_logC_dC(
i,
j,
k,
l);
519 t_P(
i,
l) = t_F(
i,
k) * t_S(
k,
l);
528#ifdef HENCKY_SMALL_STRAIN
540template <
int DIM,
typename DomainEleOp,
int S>
543 boost::shared_ptr<CommonData> common_data,
544 boost::shared_ptr<MatrixDouble> mat_D_ptr =
nullptr)
546 commonDataPtr(common_data) {
547 std::fill(&DomainEleOp::doEntities[MBEDGE],
548 &DomainEleOp::doEntities[MBMAXTYPE],
false);
552 matDPtr = commonDataPtr->matDPtr;
569 const size_t nb_gauss_pts = DomainEleOp::getGaussPts().size2();
570 commonDataPtr->matTangent.resize(DIM * DIM * DIM * DIM, nb_gauss_pts);
572 getFTensor4FromMat<DIM, DIM, DIM, DIM, 1>(commonDataPtr->matTangent);
574 auto t_D = getFTensor4DdgFromMat<DIM, DIM, S>(*matDPtr);
575 auto t_eig_val = getFTensor1FromMat<DIM>(commonDataPtr->matEigVal);
576 auto t_eig_vec = getFTensor2FromMat<DIM, DIM>(commonDataPtr->matEigVec);
577 auto t_T = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matHenckyStress);
579 getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matSecondPiolaStress);
580 auto t_grad = getFTensor2FromMat<DIM, DIM>(*(commonDataPtr->matGradPtr));
581 auto t_logC_dC = getFTensor4DdgFromMat<DIM, DIM>(commonDataPtr->matLogCdC);
583 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
585#ifdef HENCKY_SMALL_STRAIN
586 dP_dF(
i,
j,
k,
l) = t_D(
i,
j,
k,
l);
601 P_D_P_plus_TL(
i,
j,
k,
l) =
603 (t_logC_dC(
i,
j, o, p) * t_D(o, p,
m,
n)) * t_logC_dC(
m,
n,
k,
l);
604 P_D_P_plus_TL(
i,
j,
k,
l) *= 0.5;
607 t_F(
i,
k) * (P_D_P_plus_TL(
k,
j, o, p) * dC_dF(o, p,
m,
n));
630template <
int DIM,
typename AssemblyDomainEleOp,
int S>
634 const std::string row_field_name,
const std::string col_field_name,
635 boost::shared_ptr<CommonData> elastic_common_data_ptr,
636 boost::shared_ptr<VectorDouble> coeff_expansion_ptr);
638 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
639 EntitiesFieldData::EntData &col_data);
646template <
int DIM,
typename AssemblyDomainEleOp,
int S>
649 const std::string row_field_name,
const std::string col_field_name,
650 boost::shared_ptr<CommonData> elastic_common_data_ptr,
651 boost::shared_ptr<VectorDouble> coeff_expansion_ptr)
654 elasticCommonDataPtr(elastic_common_data_ptr),
655 coeffExpansionPtr(coeff_expansion_ptr) {
659template <
int DIM,
typename AssemblyDomainEleOp,
int S>
662 iNtegrate(EntitiesFieldData::EntData &row_data,
663 EntitiesFieldData::EntData &col_data) {
666 auto &locMat = AssemblyDomainEleOp::locMat;
668 const auto nb_integration_pts = row_data.getN().size1();
669 const auto nb_row_base_functions = row_data.getN().size2();
670 auto t_w = this->getFTensor0IntegrationWeight();
673 auto t_row_diff_base = row_data.getFTensor1DiffN<DIM>();
674 auto t_D = getFTensor4DdgFromMat<DIM, DIM, 0>(*elasticCommonDataPtr->matDPtr);
676 getFTensor2FromMat<DIM, DIM>(*(elasticCommonDataPtr->matGradPtr));
678 getFTensor4DdgFromMat<DIM, DIM>(elasticCommonDataPtr->matLogCdC);
691 t_coeff_exp(
i,
j) = 0;
693 t_coeff_exp(d, d) = (*coeffExpansionPtr)[d];
696 t_eigen_strain(
i,
j) = (t_D(
i,
j,
k,
l) * t_coeff_exp(
k,
l));
698 for (
auto gg = 0; gg != nb_integration_pts; ++gg) {
702 double alpha = this->getMeasure() * t_w;
704 for (; rr != AssemblyDomainEleOp::nbRows / DIM; ++rr) {
705 auto t_mat = getFTensor1FromMat<DIM, 1>(locMat, rr * DIM);
706 auto t_col_base = col_data.getFTensor0N(gg, 0);
707 for (
auto cc = 0; cc != AssemblyDomainEleOp::nbCols; cc++) {
708#ifdef HENCKY_SMALL_STRAIN
710 (t_row_diff_base(
j) * t_eigen_strain(
i,
j)) * (t_col_base * alpha);
712 t_mat(
i) -= (t_row_diff_base(
j) *
713 (t_F(
i, o) * ((t_D(
m,
n,
k,
l) * t_coeff_exp(
k,
l)) *
714 t_logC_dC(
m,
n, o,
j)))) *
715 (t_col_base * alpha);
724 for (; rr != nb_row_base_functions; ++rr)
737 template <
int DIM, IntegrationType I>
740 template <
int DIM, IntegrationType I>
743 template <
int DIM, IntegrationType I>
746 template <
int DIM, IntegrationType I,
int S>
750 template <
int DIM, IntegrationType I,
int S>
754 template <
int DIM, IntegrationType I,
int S>
758 template <
int DIM, IntegrationType I,
int S>
762 template <
int DIM, IntegrationType I,
int S>
765 template <
int DIM, IntegrationType I,
typename AssemblyDomainEleOp,
int S>
773 boost::ptr_deque<ForcesAndSourcesCore::UserDataOperator> &pip,
774 std::string block_name,
775 boost::shared_ptr<MatrixDouble> mat_D_Ptr, Sev sev,
779 PetscBool plane_strain_flag = PETSC_FALSE;
780 CHKERR PetscOptionsGetBool(PETSC_NULLPTR,
"",
"-plane_strain",
781 &plane_strain_flag, PETSC_NULLPTR);
786 std::vector<const CubitMeshSets *> meshset_vec_ptr,
787 double scale, PetscBool plane_strain_flag)
791 planeStrainFlag(plane_strain_flag) {
793 "Can not get data from block");
796 MoFEMErrorCode doWork(
int side, EntityType type,
797 EntitiesFieldData::EntData &data) {
800 for (
auto &b : blockData) {
802 if (b.blockEnts.find(getFEEntityHandle()) != b.blockEnts.end()) {
803 CHKERR getMatDPtr(matDPtr, b.bulkModulusK * scaleYoungModulus,
804 b.shearModulusG * scaleYoungModulus,
810 CHKERR getMatDPtr(matDPtr, bulkModulusKDefault * scaleYoungModulus,
811 shearModulusGDefault * scaleYoungModulus,
817 boost::shared_ptr<MatrixDouble> matDPtr;
818 const double scaleYoungModulus;
819 const PetscBool planeStrainFlag;
823 double shearModulusG;
827 double bulkModulusKDefault;
828 double shearModulusGDefault;
829 std::vector<BlockData> blockData;
833 std::vector<const CubitMeshSets *> meshset_vec_ptr,
837 for (
auto m : meshset_vec_ptr) {
839 std::vector<double> block_data;
840 CHKERR m->getAttributes(block_data);
841 if (block_data.size() != 2) {
843 "Expected that block has two attribute");
845 auto get_block_ents = [&]() {
848 m_field.
get_moab().get_entities_by_handle(
m->meshset, ents,
true);
867 MoFEMErrorCode getMatDPtr(boost::shared_ptr<MatrixDouble> mat_D_ptr,
872 auto set_material_stiffness = [&]() {
882 auto t_D = getFTensor4DdgFromMat<DIM, DIM, 0>(*mat_D_ptr);
889 constexpr auto size_symm = (DIM * (DIM + 1)) / 2;
891 set_material_stiffness();
899 PetscOptionsBegin(PETSC_COMM_WORLD,
"",
"",
"none");
900 CHKERR PetscOptionsScalar(
"-young_modulus",
"Young modulus",
"",
E, &
E,
902 CHKERR PetscOptionsScalar(
"-poisson_ratio",
"poisson ratio",
"", nu, &nu,
908 pip.push_back(
new OpMatBlocks(
912 m_field.
getInterface<MeshsetsManager>()->getCubitMeshsetPtr(std::regex(
914 (boost::format(
"%s(.*)") % block_name).str()
917 scale, plane_strain_flag
924template <
int DIM, IntegrationType I,
typename DomainEleOp>
927 boost::ptr_deque<ForcesAndSourcesCore::UserDataOperator> &pip,
928 std::string
field_name, std::string block_name, Sev sev,
double scale = 1) {
930 auto common_ptr = boost::make_shared<HenckyOps::CommonData>();
931 common_ptr->matDPtr = boost::make_shared<MatrixDouble>();
932 common_ptr->matGradPtr = boost::make_shared<MatrixDouble>();
935 common_ptr->matDPtr, sev,
scale),
940 pip.push_back(
new OpCalculateVectorFieldGradient<DIM, DIM>(
942 pip.push_back(
new typename H::template OpCalculateEigenVals<DIM, I>(
945 new typename H::template OpCalculateLogC<DIM, I>(
field_name, common_ptr));
946 pip.push_back(
new typename H::template OpCalculateLogC_dC<DIM, I>(
949 pip.push_back(
new typename H::template OpCalculateHenckyStress<DIM, I, 0>(
951 pip.push_back(
new typename H::template OpCalculatePiolaStress<DIM, I, 0>(
957template <
int DIM, AssemblyType A, IntegrationType I,
typename DomainEleOp>
960 boost::ptr_deque<ForcesAndSourcesCore::UserDataOperator> &pip,
961 std::string
field_name, boost::shared_ptr<HenckyOps::CommonData> common_ptr,
965 using B =
typename FormsIntegrators<DomainEleOp>::template Assembly<
975template <
int DIM, AssemblyType A, IntegrationType I,
typename DomainEleOp>
978 boost::ptr_deque<ForcesAndSourcesCore::UserDataOperator> &pip,
979 std::string
field_name, std::string block_name, Sev sev,
double scale = 1) {
990template <
int DIM, AssemblyType A, IntegrationType I,
typename DomainEleOp>
993 boost::ptr_deque<ForcesAndSourcesCore::UserDataOperator> &pip,
994 std::string
field_name, boost::shared_ptr<HenckyOps::CommonData> common_ptr,
998 using B =
typename FormsIntegrators<DomainEleOp>::template Assembly<
1000 using OpKPiola =
typename B::template OpGradTensorGrad<1, DIM, DIM, 1>;
1004 pip.push_back(
new typename H::template OpHenckyTangent<DIM, I, 0>(
1012template <
int DIM, AssemblyType A, IntegrationType I,
typename DomainEleOp>
1015 boost::ptr_deque<ForcesAndSourcesCore::UserDataOperator> &pip,
1016 std::string
field_name, std::string block_name, Sev sev,
double scale = 1) {
#define MOFEM_TAG_AND_LOG(channel, severity, tag)
Tag and log in channel.
DomainEle::UserDataOperator DomainEleOp
Finire element operator type.
Kronecker Delta class symmetric.
#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()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_DATA_INCONSISTENCY
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
FTensor::Index< 'i', SPACE_DIM > i
const double v
phase velocity of light in medium (cm/ns)
const double n
refractive index of diffusive medium
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
auto getMat(A &&t_val, B &&t_vec, Fun< double > f)
Get the Mat object.
auto getDiffMat(A &&t_val, B &&t_vec, Fun< double > f, Fun< double > d_f, const int nb)
Get the Diff Mat object.
auto getDiffDiffMat(A &&t_val, B &&t_vec, Fun< double > f, Fun< double > d_f, Fun< double > dd_f, C &&t_S, const int nb)
Get the Diff Diff Mat object.
auto get_uniq_nb(double *ptr)
auto sort_eigen_vals(FTensor::Tensor1< double, DIM > &eig, FTensor::Tensor2< double, DIM, DIM > &eigen_vec)
MoFEMErrorCode opFactoryDomainLhs(MoFEM::Interface &m_field, boost::ptr_deque< ForcesAndSourcesCore::UserDataOperator > &pip, std::string field_name, boost::shared_ptr< HenckyOps::CommonData > common_ptr, Sev sev)
MoFEMErrorCode opFactoryDomainRhs(MoFEM::Interface &m_field, boost::ptr_deque< ForcesAndSourcesCore::UserDataOperator > &pip, std::string field_name, boost::shared_ptr< HenckyOps::CommonData > common_ptr, Sev sev)
MoFEMErrorCode addMatBlockOps(MoFEM::Interface &m_field, boost::ptr_deque< ForcesAndSourcesCore::UserDataOperator > &pip, std::string block_name, boost::shared_ptr< MatrixDouble > mat_D_Ptr, Sev sev, double scale=1)
auto commonDataFactory(MoFEM::Interface &m_field, boost::ptr_deque< ForcesAndSourcesCore::UserDataOperator > &pip, std::string field_name, std::string block_name, Sev sev, double scale=1)
auto is_eq(const double &a, const double &b)
FormsIntegrators< DomainEleOp >::Assembly< A >::LinearForm< I >::OpGradTimesTensor< 1, FIELD_DIM, SPACE_DIM > OpGradTimesTensor
double young_modulus
Young modulus.
double poisson_ratio
Poisson ratio.
constexpr auto field_name
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpGradTensorGrad< 1, SPACE_DIM, SPACE_DIM, 1 > OpKPiola
[Only used for dynamics]
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpGradTimesTensor< 1, SPACE_DIM, SPACE_DIM > OpInternalForcePiola
FTensor::Index< 'm', 3 > m
auto getMatHenckyStress()
auto getMatFirstPiolaStress()
boost::shared_ptr< MatrixDouble > matLogCPlastic
MatrixDouble matFirstPiolaStress
boost::shared_ptr< MatrixDouble > matDPtr
MatrixDouble matSecondPiolaStress
MatrixDouble matHenckyStress
boost::shared_ptr< MatrixDouble > matGradPtr
OpCalculateEigenValsImpl(const std::string field_name, boost::shared_ptr< CommonData > common_data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
boost::shared_ptr< CommonData > commonDataPtr
boost::shared_ptr< MatrixDouble > matDPtr
boost::shared_ptr< MatrixDouble > matLogCPlastic
boost::shared_ptr< CommonData > commonDataPtr
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
OpCalculateHenckyPlasticStressImpl(const std::string field_name, boost::shared_ptr< CommonData > common_data, boost::shared_ptr< MatrixDouble > mat_D_ptr, const double scale=1)
boost::shared_ptr< CommonData > commonDataPtr
OpCalculateHenckyStressImpl(const std::string field_name, boost::shared_ptr< CommonData > common_data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
boost::shared_ptr< VectorDouble > coeffExpansionPtr
OpCalculateHenckyThermalStressImpl(const std::string field_name, boost::shared_ptr< VectorDouble > temperature, boost::shared_ptr< CommonData > common_data, boost::shared_ptr< VectorDouble > coeff_expansion_ptr, boost::shared_ptr< double > ref_temp_ptr)
boost::shared_ptr< double > refTempPtr
boost::shared_ptr< CommonData > commonDataPtr
boost::shared_ptr< VectorDouble > tempPtr
boost::shared_ptr< CommonData > elasticCommonDataPtr
boost::shared_ptr< VectorDouble > coeffExpansionPtr
boost::shared_ptr< CommonData > commonDataPtr
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
OpCalculateLogCImpl(const std::string field_name, boost::shared_ptr< CommonData > common_data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
OpCalculateLogC_dCImpl(const std::string field_name, boost::shared_ptr< CommonData > common_data)
boost::shared_ptr< CommonData > commonDataPtr
boost::shared_ptr< CommonData > commonDataPtr
OpCalculatePiolaStressImpl(const std::string field_name, boost::shared_ptr< CommonData > common_data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
OpHenckyTangentImpl(const std::string field_name, boost::shared_ptr< CommonData > common_data, boost::shared_ptr< MatrixDouble > mat_D_ptr=nullptr)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
boost::shared_ptr< CommonData > commonDataPtr
boost::shared_ptr< MatrixDouble > matDPtr
static auto check(const double &a, const double &b)
virtual moab::Interface & get_moab()=0
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.