v0.14.0
Loading...
Searching...
No Matches
Public Member Functions | Protected Attributes | List of all members
EshelbianPlasticity::OpTractionBc Struct Reference

#include <users_modules/eshelbian_plasticty/src/EshelbianPlasticity.hpp>

Inheritance diagram for EshelbianPlasticity::OpTractionBc:
[legend]
Collaboration diagram for EshelbianPlasticity::OpTractionBc:
[legend]

Public Member Functions

 OpTractionBc (std::string row_field, FeTractionBc &bc_fe)
 
MoFEMErrorCode doWork (int side, EntityType type, EntData &data)
 

Protected Attributes

FeTractionBcbcFe
 
MatrixDouble matPP
 
MatrixDouble vecPv
 

Detailed Description

Definition at line 582 of file EshelbianPlasticity.hpp.

Constructor & Destructor Documentation

◆ OpTractionBc()

EshelbianPlasticity::OpTractionBc::OpTractionBc ( std::string  row_field,
FeTractionBc bc_fe 
)
inline

Definition at line 583 of file EshelbianPlasticity.hpp.

584 : FaceUserDataOperator(row_field, FaceUserDataOperator::OPROW),
585 bcFe(bc_fe) {}
FaceElementForcesAndSourcesCore::UserDataOperator FaceUserDataOperator

Member Function Documentation

◆ doWork()

MoFEMErrorCode EshelbianPlasticity::OpTractionBc::doWork ( int  side,
EntityType  type,
EntData data 
)
Examples
EshelbianOperators.cpp.

Definition at line 944 of file EshelbianOperators.cpp.

944 {
946
949
950 auto t_normal = getFTensor1Normal();
951 const double nrm2 = sqrt(t_normal(i) * t_normal(i));
952 FTensor::Tensor1<double, 3> t_unit_normal;
953 t_unit_normal(i) = t_normal(i) / nrm2;
954 int nb_dofs = data.getFieldData().size();
955 int nb_integration_pts = data.getN().size1();
956 int nb_base_functions = data.getN().size2() / 3;
957 double ts_t = getFEMethod()->ts_t;
958
959 auto integrate_matrix = [&]() {
961
962 auto t_w = getFTensor0IntegrationWeight();
963 matPP.resize(nb_dofs / 3, nb_dofs / 3, false);
964 matPP.clear();
965
966 auto t_row_base_fun = data.getFTensor1N<3>();
967 for (int gg = 0; gg != nb_integration_pts; ++gg) {
968
969 int rr = 0;
970 for (; rr != nb_dofs / 3; ++rr) {
971 const double a = t_row_base_fun(i) * t_unit_normal(i);
972 auto t_col_base_fun = data.getFTensor1N<3>(gg, 0);
973 for (int cc = 0; cc != nb_dofs / 3; ++cc) {
974 const double b = t_col_base_fun(i) * t_unit_normal(i);
975 matPP(rr, cc) += t_w * a * b;
976 ++t_col_base_fun;
977 }
978 ++t_row_base_fun;
979 }
980
981 for (; rr != nb_base_functions; ++rr)
982 ++t_row_base_fun;
983
984 ++t_w;
985 }
986
988 };
989
990 auto integrate_rhs = [&](auto &bc) {
992
993 auto t_w = getFTensor0IntegrationWeight();
994 vecPv.resize(3, nb_dofs / 3, false);
995 vecPv.clear();
996
997 auto t_row_base_fun = data.getFTensor1N<3>();
998 double ts_t = getFEMethod()->ts_t;
999
1000 for (int gg = 0; gg != nb_integration_pts; ++gg) {
1001 int rr = 0;
1002 for (; rr != nb_dofs / 3; ++rr) {
1003 const double t = ts_t * t_w * t_row_base_fun(i) * t_unit_normal(i);
1004 for (int dd = 0; dd != 3; ++dd)
1005 if (bc.flags[dd])
1006 vecPv(dd, rr) += t * bc.vals[dd];
1007 ++t_row_base_fun;
1008 }
1009 for (; rr != nb_base_functions; ++rr)
1010 ++t_row_base_fun;
1011 ++t_w;
1012 }
1014 };
1015
1016 auto integrate_rhs_cook = [&](auto &bc) {
1018
1019 vecPv.resize(3, nb_dofs / 3, false);
1020 vecPv.clear();
1021
1022 auto t_w = getFTensor0IntegrationWeight();
1023 auto t_row_base_fun = data.getFTensor1N<3>();
1024 auto t_coords = getFTensor1CoordsAtGaussPts();
1025
1026 for (int gg = 0; gg != nb_integration_pts; ++gg) {
1027
1028 auto calc_tau = [](double y) {
1029 y -= 44;
1030 y /= (60 - 44);
1031 return -y * (y - 1) / 0.25;
1032 };
1033
1034 const double tau = calc_tau(t_coords(1));
1035
1036 int rr = 0;
1037 for (; rr != nb_dofs / 3; ++rr) {
1038 const double t = ts_t * t_w * t_row_base_fun(i) * t_unit_normal(i);
1039 for (int dd = 0; dd != 3; ++dd)
1040 if (bc.flags[dd])
1041 vecPv(dd, rr) += t * tau * bc.vals[dd];
1042 ++t_row_base_fun;
1043 }
1044
1045 for (; rr != nb_base_functions; ++rr)
1046 ++t_row_base_fun;
1047 ++t_w;
1048 ++t_coords;
1049 }
1051 };
1052
1053 // get entity of face
1054 EntityHandle fe_ent = getFEEntityHandle();
1055 for (auto &bc : *(bcFe.bcData)) {
1056 if (bc.faces.find(fe_ent) != bc.faces.end()) {
1057
1058 int nb_dofs = data.getFieldData().size();
1059 if (nb_dofs) {
1060
1061 CHKERR integrate_matrix();
1062 if (std::regex_match(bc.blockName, std::regex(".*COOK.*")))
1063 CHKERR integrate_rhs_cook(bc);
1064 else
1065 CHKERR integrate_rhs(bc);
1066
1067 const auto info =
1068 lapack_dposv('L', nb_dofs / 3, 3, &*matPP.data().begin(),
1069 nb_dofs / 3, &*vecPv.data().begin(), nb_dofs / 3);
1070 if (info > 0)
1071 SETERRQ1(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
1072 "The leading minor of order %i is "
1073 "not positive; definite;\nthe "
1074 "solution could not be computed",
1075 info);
1076
1077 for (int dd = 0; dd != 3; ++dd)
1078 if (bc.flags[dd])
1079 for (int rr = 0; rr != nb_dofs / 3; ++rr)
1080 data.getFieldDofs()[3 * rr + dd]->getFieldData() = vecPv(dd, rr);
1081 }
1082 }
1083 }
1084
1086}
constexpr double a
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ MOFEM_OPERATION_UNSUCCESSFUL
Definition: definitions.h:34
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
FTensor::Index< 'i', SPACE_DIM > i
static __CLPK_integer lapack_dposv(char uplo, __CLPK_integer n, __CLPK_integer nrhs, __CLPK_doublereal *a, __CLPK_integer lda, __CLPK_doublereal *b, __CLPK_integer ldb)
Definition: lapack_wrap.h:211
FTensor::Index< 'j', 3 > j
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
constexpr double t
plate stiffness
Definition: plate.cpp:59
boost::shared_ptr< TractionBcVec > bcData
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::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1N(FieldApproximationBase base)
Get base functions for Hdiv/Hcurl spaces.
const VectorDofs & getFieldDofs() const
get dofs data stature FEDofEntity

Member Data Documentation

◆ bcFe

FeTractionBc& EshelbianPlasticity::OpTractionBc::bcFe
protected
Examples
EshelbianOperators.cpp.

Definition at line 589 of file EshelbianPlasticity.hpp.

◆ matPP

MatrixDouble EshelbianPlasticity::OpTractionBc::matPP
protected
Examples
EshelbianOperators.cpp.

Definition at line 590 of file EshelbianPlasticity.hpp.

◆ vecPv

MatrixDouble EshelbianPlasticity::OpTractionBc::vecPv
protected
Examples
EshelbianOperators.cpp.

Definition at line 591 of file EshelbianPlasticity.hpp.


The documentation for this struct was generated from the following files: