v0.15.0
Loading...
Searching...
No Matches
SurfacePressureComplexForLazy.hpp
Go to the documentation of this file.
1/* \file SurfacePressureComplexForLazy.hpp
2 */
3
4
5
6#ifndef __COMPLEX_FOR_LAZY_NEUMANN_FORCES_HPP__
7#define __COMPLEX_FOR_LAZY_NEUMANN_FORCES_HPP__
8
9/** \brief NonLinear surface pressure element (obsolete implementation)
10 * \ingroup nonlinear_elastic_elem
11
12 \todo This is old implementation, need to be reimplemented, using
13 auto-differentiation. It is well tested and works. well.
14
15 */
17
21
23 AuxMethodSpatial(const string &field_name, MyTriangleSpatialFE *my_ptr,
24 const char type);
25 MoFEMErrorCode doWork(int side, EntityType type,
26 EntitiesFieldData::EntData &data);
27 };
28
31
34 const char type);
35 MoFEMErrorCode doWork(int side, EntityType type,
36 EntitiesFieldData::EntData &data);
37 };
38
40
41 double *sCaleLhs;
42 double *sCaleRhs;
44
46 const double eps;
47 bool uSeF;
49
50 Mat Aij;
51 Vec F;
52
53 MyTriangleSpatialFE(MoFEM::Interface &_mField, Mat _Aij, Vec &_F,
54 double *scale_lhs, double *scale_rhs,
55 std::string spatial_field_name = "SPATIAL_POSITION",
56 std::string mat_field_name = "MESH_NODE_POSITIONS");
57
58 int getRule(int order) { return max(1, order); };
59
60 double *N;
61 double *N_face;
62 double *N_edge[3];
63 double *diffN;
64 double *diffN_face;
65 double *diffN_edge[3];
66
68 int order_edge[3];
69 double *dofs_x;
70 double *dofs_x_edge[3];
71 double *dofs_x_face;
72 double *idofs_x;
73 double *idofs_x_edge[3];
74 double *idofs_x_face;
78
81 double *dofs_X;
82 double *dofs_X_edge[3];
83 double *dofs_X_face;
84 double *idofs_X;
85 double *idofs_X_edge[3];
86 double *idofs_X_face;
87
89
90 VectorDouble tLoc, tGlob;
91 MatrixDouble tLocNodal, tGlobNodal;
92 double *t_loc;
93
95 ublas::vector<ublas::vector<int>> dOfs_x_edge_indices;
97 ublas::vector<ublas::vector<int>> dOfs_X_edge_indices;
98
99 VectorDouble dOfs_x, dOfs_x_face;
100 ublas::vector<VectorDouble> dOfs_x_edge;
101 VectorDouble dOfs_X, dOfs_X_face;
102 ublas::vector<VectorDouble> dOfs_X_edge;
103
104 VectorDouble fExtNode, fExtFace;
105 ublas::vector<VectorDouble> fExtEdge;
106 double *Fext_edge[3];
107
109 ublas::vector<MatrixDouble> kExtEdgeNode;
110 double *Kext_edge_node[3];
111
113 ublas::vector<MatrixDouble> kExtEdgeFace;
114 double *Kext_edge_face[3];
115
116 ublas::vector<MatrixDouble> kExtFaceEdge, kExtNodeEdge;
117 ublas::matrix<MatrixDouble> kExtEdgeEdge;
118 double *Kext_node_edge[3];
119 double *Kext_face_edge[3];
120 double *Kext_edge_edge[3][3];
121
122 virtual MoFEMErrorCode calcTraction();
123 virtual MoFEMErrorCode rHs();
124 virtual MoFEMErrorCode lHs();
125
126 MoFEMErrorCode preProcess();
127 MoFEMErrorCode operator()();
128
129 MoFEMErrorCode addForce(int ms_id);
130 MoFEMErrorCode addPressure(int ms_id);
131
132 DEPRECATED MoFEMErrorCode addPreassure(int ms_id) {
133 return addPressure(ms_id);
134 }
135
136 struct bCForce {
137 ForceCubitBcData data;
139 };
140 map<int, bCForce> mapForce;
141 struct bCPressure {
142 PressureCubitBcData data;
144 };
145 map<int, bCPressure> mapPressure;
146 MoFEMErrorCode reBaseToFaceLoocalCoordSystem(MatrixDouble &t_glob_nodal);
147
148 boost::ptr_vector<MethodForForceScaling> methodsOp;
149 };
150
153
155
156 double *sCale;
157 MoFEMErrorCode setForceScale(double scale) {
159 *sCale = scale;
161 }
162
164
166 MoFEM::Interface &m_field, Mat _Aij, Vec _F, double *scale_lhs,
167 double *scale_rhs, std::string spatial_field_name = "SPATIAL_POSITION",
168 std::string material_field_name = "MESH_NODE_POSITIONS");
170 MoFEM::Interface &m_field, Mat _Aij, Vec _F,
171 std::string spatial_field_name = "SPATIAL_POSITION",
172 std::string material_field_name = "MESH_NODE_POSITIONS");
173
174private:
175 const std::string spatialField;
176 const std::string materialField;
177};
178
179/**
180 * \depracted do not use name with spelling mistake.
181 */
184
185#endif //__COMPLEX_FOR_LAZY_NEUMANN_FORCES_HPP__
DEPRECATED typedef NeumannForcesSurfaceComplexForLazy NeummanForcesSurfaceComplexForLazy
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define DEPRECATED
Definition definitions.h:17
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
constexpr int order
double scale
Definition plastic.cpp:123
constexpr auto field_name
#define _F(n)
Definition quad.c:25
Deprecated interface functions.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
AuxMethodMaterial(const string &field_name, MyTriangleSpatialFE *my_ptr, const char type)
AuxMethodSpatial(const string &field_name, MyTriangleSpatialFE *my_ptr, const char type)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
MoFEMErrorCode reBaseToFaceLoocalCoordSystem(MatrixDouble &t_glob_nodal)
MoFEMErrorCode operator()()
function is run for every finite element
MoFEMErrorCode preProcess()
function is run at the beginning of loop
MyTriangleSpatialFE(MoFEM::Interface &_mField, Mat _Aij, Vec &_F, double *scale_lhs, double *scale_rhs, std::string spatial_field_name="SPATIAL_POSITION", std::string mat_field_name="MESH_NODE_POSITIONS")
NonLinear surface pressure element (obsolete implementation)
NeumannForcesSurfaceComplexForLazy(MoFEM::Interface &m_field, Mat _Aij, Vec _F, double *scale_lhs, double *scale_rhs, std::string spatial_field_name="SPATIAL_POSITION", std::string material_field_name="MESH_NODE_POSITIONS")