v0.16.3
Loading...
Searching...
No Matches
Public Member Functions | Public Attributes | List of all members
ContactOps::OpConstrainBoundaryLhs_dUImpl< DIM, GAUSS, AssemblyBoundaryEleOp > Struct Template Reference

#include "tutorials/adv-1_contact/src/ContactOps.hpp"

Inheritance diagram for ContactOps::OpConstrainBoundaryLhs_dUImpl< DIM, GAUSS, AssemblyBoundaryEleOp >:
[legend]
Collaboration diagram for ContactOps::OpConstrainBoundaryLhs_dUImpl< DIM, GAUSS, AssemblyBoundaryEleOp >:
[legend]

Public Member Functions

 OpConstrainBoundaryLhs_dUImpl (const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< CommonData > common_data_ptr, bool is_axisymmetric=false)
 
MoFEMErrorCode iNtegrate (EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
 

Public Attributes

SurfaceDistanceFunction surfaceDistanceFunction = surface_distance_function
 
GradSurfaceDistanceFunction gradSurfaceDistanceFunction
 
HessSurfaceDistanceFunction hessSurfaceDistanceFunction
 
boost::shared_ptr< CommonData > commonDataPtr
 
bool isAxisymmetric
 

Detailed Description

template<int DIM, typename AssemblyBoundaryEleOp>
struct ContactOps::OpConstrainBoundaryLhs_dUImpl< DIM, GAUSS, AssemblyBoundaryEleOp >

Definition at line 556 of file ContactOps.hpp.

Constructor & Destructor Documentation

◆ OpConstrainBoundaryLhs_dUImpl()

template<int DIM, typename AssemblyBoundaryEleOp >
ContactOps::OpConstrainBoundaryLhs_dUImpl< DIM, GAUSS, AssemblyBoundaryEleOp >::OpConstrainBoundaryLhs_dUImpl ( const std::string  row_field_name,
const std::string  col_field_name,
boost::shared_ptr< CommonData >  common_data_ptr,
bool  is_axisymmetric = false 
)

Definition at line 988 of file ContactOps.hpp.

993 : AssemblyBoundaryEleOp(row_field_name, col_field_name,
994 AssemblyBoundaryEleOp::OPROWCOL),
995 commonDataPtr(common_data_ptr), isAxisymmetric(is_axisymmetric) {
996 AssemblyBoundaryEleOp::sYmm = false;
997}
PetscBool is_axisymmetric
Definition contact.cpp:91
FormsIntegrators< BoundaryEleOp >::Assembly< A >::OpBase AssemblyBoundaryEleOp

Member Function Documentation

◆ iNtegrate()

template<int DIM, typename AssemblyBoundaryEleOp >
MoFEMErrorCode ContactOps::OpConstrainBoundaryLhs_dUImpl< DIM, GAUSS, AssemblyBoundaryEleOp >::iNtegrate ( EntitiesFieldData::EntData &  row_data,
EntitiesFieldData::EntData &  col_data 
)

Definition at line 1001 of file ContactOps.hpp.

1003 {
1005
1006 FTensor::Index<'i', DIM> i;
1007 FTensor::Index<'j', DIM> j;
1008 FTensor::Index<'k', DIM> k;
1009
1010 const size_t nb_gauss_pts = AssemblyBoundaryEleOp::getGaussPts().size2();
1011 auto &locMat = AssemblyBoundaryEleOp::locMat;
1012
1013 auto t_traction = getFTensor1FromMat<DIM>(commonDataPtr->contactTraction);
1014 auto t_coords = AssemblyBoundaryEleOp::getFTensor1CoordsAtGaussPts();
1015
1016 auto t_w = AssemblyBoundaryEleOp::getFTensor0IntegrationWeight();
1017 auto t_row_base = row_data.getFTensor1N<3>();
1018 size_t nb_face_functions = row_data.getN().size2() / 3;
1019
1020 auto m_spatial_coords = get_spatial_coords(
1021 BoundaryEleOp::getFTensor1CoordsAtGaussPts(),
1022 getFTensor1FromMat<DIM>(commonDataPtr->contactDisp), nb_gauss_pts);
1023 auto m_normals_at_pts = get_normalize_normals(
1024 BoundaryEleOp::getFTensor1NormalsAtGaussPts(), nb_gauss_pts);
1025
1026 auto t_normal = getFTensor1FromMat<3>(m_normals_at_pts);
1027
1028 auto ts_time = AssemblyBoundaryEleOp::getTStime();
1029 auto ts_time_step = AssemblyBoundaryEleOp::getTStimeStep();
1030
1031 // placeholder to pass boundary block id to python
1032 int block_id = 0;
1033
1034 auto v_sdf =
1035 surfaceDistanceFunction(ts_time_step, ts_time, nb_gauss_pts,
1036 m_spatial_coords, m_normals_at_pts, block_id);
1037
1038 auto m_grad_sdf =
1039 gradSurfaceDistanceFunction(ts_time_step, ts_time, nb_gauss_pts,
1040 m_spatial_coords, m_normals_at_pts, block_id);
1041
1042 auto m_hess_sdf =
1043 hessSurfaceDistanceFunction(ts_time_step, ts_time, nb_gauss_pts,
1044 m_spatial_coords, m_normals_at_pts, block_id);
1045
1046 auto t_sdf = getFTensor0FromVec(v_sdf);
1047 auto t_grad_sdf = getFTensor1FromMat<3>(m_grad_sdf);
1048 auto t_hess_sdf = getFTensor2SymmetricFromMat<3>(m_hess_sdf);
1049
1050 for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
1051
1052 double jacobian = 1.;
1053 if (isAxisymmetric) {
1054 jacobian = 2. * M_PI * t_coords(0);
1055 }
1056 const double alpha = t_w * jacobian * AssemblyBoundaryEleOp::getMeasure();
1057
1058 auto tn = -t_traction(i) * t_grad_sdf(i);
1059 auto c = constrain(t_sdf, tn);
1060
1062 t_cP(i, j) = (c * t_grad_sdf(i)) * t_grad_sdf(j);
1064 t_cQ(i, j) = kronecker_delta(i, j) - t_cP(i, j);
1065
1067 t_res_dU(i, j) = kronecker_delta(i, j) + t_cP(i, j);
1068
1069 if (c > 0) {
1070 t_res_dU(i, j) +=
1071 (c * cn_contact) *
1072 (t_hess_sdf(i, j) * (t_grad_sdf(k) * t_traction(k)) +
1073 t_grad_sdf(i) * t_hess_sdf(k, j) * t_traction(k)) +
1074 c * t_sdf * t_hess_sdf(i, j);
1075 }
1076
1077 size_t rr = 0;
1078 for (; rr != AssemblyBoundaryEleOp::nbRows / DIM; ++rr) {
1079
1080 auto t_mat = getFTensor2FromArray<DIM, DIM, DIM>(locMat, DIM * rr);
1081
1082 const double row_base = t_row_base(i) * t_normal(i);
1083
1084 auto t_col_base = col_data.getFTensor0N(gg, 0);
1085 for (size_t cc = 0; cc != AssemblyBoundaryEleOp::nbCols / DIM; ++cc) {
1086 const double beta = alpha * row_base * t_col_base;
1087
1088 t_mat(i, j) -= beta * t_res_dU(i, j);
1089
1090 ++t_col_base;
1091 ++t_mat;
1092 }
1093
1094 ++t_row_base;
1095 }
1096 for (; rr < nb_face_functions; ++rr)
1097 ++t_row_base;
1098
1099 ++t_traction;
1100 ++t_coords;
1101 ++t_w;
1102 ++t_normal;
1103 ++t_sdf;
1104 ++t_grad_sdf;
1105 ++t_hess_sdf;
1106 }
1107
1109}
SmartPetscObj< Mat > jacobian
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
FTensor::Index< 'i', SPACE_DIM > i
const double c
speed of light (cm/ns)
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
auto get_normalize_normals(FTensor::Tensor1< T1, DIM1 > &&t_normal_at_pts, size_t nb_gauss_pts)
double cn_contact
Definition contact.cpp:95
auto get_spatial_coords(FTensor::Tensor1< T1, DIM1 > &&t_coords, FTensor::Tensor1< T2, DIM2 > &&t_disp, size_t nb_gauss_pts)
double constrain(double sdf, double tn)
constrain function
Tensor2_Expr< Kronecker_Delta< T >, T, Dim0, Dim1, i, j > kronecker_delta(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
Rank 2.
int nbRows
number of dofs on rows
MatrixDouble locMat
local entity block matrix
int nbCols
number if dof on column

Member Data Documentation

◆ commonDataPtr

template<int DIM, typename AssemblyBoundaryEleOp >
boost::shared_ptr<CommonData> ContactOps::OpConstrainBoundaryLhs_dUImpl< DIM, GAUSS, AssemblyBoundaryEleOp >::commonDataPtr

Definition at line 571 of file ContactOps.hpp.

◆ gradSurfaceDistanceFunction

template<int DIM, typename AssemblyBoundaryEleOp >
GradSurfaceDistanceFunction ContactOps::OpConstrainBoundaryLhs_dUImpl< DIM, GAUSS, AssemblyBoundaryEleOp >::gradSurfaceDistanceFunction
Initial value:
=
MatrixDouble grad_surface_distance_function(double delta_t, double t, int nb_gauss_pts, MatrixDouble &m_spatial_coords, MatrixDouble &m_normals_at_pts, int block_id)

Definition at line 566 of file ContactOps.hpp.

◆ hessSurfaceDistanceFunction

template<int DIM, typename AssemblyBoundaryEleOp >
HessSurfaceDistanceFunction ContactOps::OpConstrainBoundaryLhs_dUImpl< DIM, GAUSS, AssemblyBoundaryEleOp >::hessSurfaceDistanceFunction
Initial value:
=
MatrixDouble hess_surface_distance_function(double delta_t, double t, int nb_gauss_pts, MatrixDouble &m_spatial_coords, MatrixDouble &m_normals_at_pts, int block_id)

Definition at line 568 of file ContactOps.hpp.

◆ isAxisymmetric

template<int DIM, typename AssemblyBoundaryEleOp >
bool ContactOps::OpConstrainBoundaryLhs_dUImpl< DIM, GAUSS, AssemblyBoundaryEleOp >::isAxisymmetric

Definition at line 572 of file ContactOps.hpp.

◆ surfaceDistanceFunction

Definition at line 565 of file ContactOps.hpp.


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