v0.14.0
Loading...
Searching...
No Matches
Public Member Functions | Public Attributes | List of all members
SmallStrainTranverslyIsotropic< TYPE > Struct Template Reference

Hook equation. More...

#include <users_modules/basic_finite_elements/nonlinear_elastic_materials/src/SmallTransverselyIsotropic.hpp>

Inheritance diagram for SmallStrainTranverslyIsotropic< TYPE >:
[legend]
Collaboration diagram for SmallStrainTranverslyIsotropic< TYPE >:
[legend]

Public Member Functions

 SmallStrainTranverslyIsotropic ()
 
MoFEMErrorCode calculateStrain ()
 
MoFEMErrorCode calculateLocalStiffnesMatrix ()
 
MoFEMErrorCode calculateAxisAngleRotationalMatrix ()
 Function to Calculate the Rotation Matrix at a given axis and angle of rotation. More...
 
MoFEMErrorCode stressTransformation ()
 Function to Calculate Stress Transformation Matrix This function computes the stress transformation Matrix at a give axis and angle of rotation
One can also output the axis/angle rotational Matrix. More...
 
MoFEMErrorCode strainTransformation ()
 Function to Calculate Strain Transformation Matrix
This function computes the strain transformation Matrix at a give axis and angle of rotation
One can also output the axis/angle rotational Matrix. More...
 
MoFEMErrorCode calculateGlobalStiffnesMatrix ()
 
virtual MoFEMErrorCode calculateAngles ()
 
virtual MoFEMErrorCode calculateP_PiolaKirchhoffI (const NonlinearElasticElement::BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
 Calculate global stress. More...
 
virtual MoFEMErrorCode calculateElasticEnergy (const NonlinearElasticElement::BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
 calculate density of strain energy More...
 
MoFEMErrorCode calculateFibreAngles ()
 
- Public Member Functions inherited from NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >
 FunctionsToCalculatePiolaKirchhoffI ()
 
virtual ~FunctionsToCalculatePiolaKirchhoffI ()=default
 
MoFEMErrorCode calculateC_CauchyDeformationTensor ()
 
MoFEMErrorCode calculateE_GreenStrain ()
 
MoFEMErrorCode calculateS_PiolaKirchhoffII ()
 
virtual MoFEMErrorCode calculateP_PiolaKirchhoffI (const BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
 Function overload to implement user material. More...
 
virtual MoFEMErrorCode calculateCauchyStress (const BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
 Function overload to implement user material. More...
 
virtual MoFEMErrorCode setUserActiveVariables (int &nb_active_variables)
 add additional active variables More...
 
virtual MoFEMErrorCode setUserActiveVariables (VectorDouble &activeVariables)
 Add additional independent variables More complex physical models depend on gradient of defamation and some additional variables. For example can depend on temperature. This function adds additional independent variables to the model. More...
 
virtual MoFEMErrorCode calculateElasticEnergy (const BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
 Calculate elastic energy density. More...
 
virtual MoFEMErrorCode calculatesIGma_EshelbyStress (const BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
 Calculate Eshelby stress. More...
 
virtual MoFEMErrorCode getDataOnPostProcessor (std::map< std::string, std::vector< VectorDouble > > &field_map, std::map< std::string, std::vector< MatrixDouble > > &grad_map)
 Do operations when pre-process. More...
 

Public Attributes

ublas::matrix< TYPE > sTrain
 
ublas::vector< TYPE > voightStrain
 
TYPE tR
 
double nu_p
 
double nu_pz
 
double E_p
 
double E_z
 
double G_zp
 
ublas::symmetric_matrix< TYPE, ublas::upper > localStiffnessMatrix
 
ublas::matrix< TYPE > aARotMat
 
ublas::vector< TYPE > axVector
 
TYPE axAngle
 
ublas::matrix< TYPE > stressRotMat
 
ublas::matrix< TYPE > strainRotMat
 
ublas::matrix< TYPE > dR
 
ublas::matrix< TYPE > globalStiffnessMatrix
 
ublas::vector< TYPE > voigtStress
 
VectorDouble normalizedPhi
 
VectorDouble axVectorDouble
 
double axAngleDouble
 
- Public Attributes inherited from NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >
FTensor::Index< 'i', 3 > i
 
FTensor::Index< 'j', 3 > j
 
FTensor::Index< 'k', 3 > k
 
double lambda
 
double mu
 
MatrixBoundedArray< TYPE, 9 > F
 
MatrixBoundedArray< TYPE, 9 > C
 
MatrixBoundedArray< TYPE, 9 > E
 
MatrixBoundedArray< TYPE, 9 > S
 
MatrixBoundedArray< TYPE, 9 > invF
 
MatrixBoundedArray< TYPE, 9 > P
 
MatrixBoundedArray< TYPE, 9 > sIGma
 
MatrixBoundedArray< TYPE, 9 > h
 
MatrixBoundedArray< TYPE, 9 > H
 
MatrixBoundedArray< TYPE, 9 > invH
 
MatrixBoundedArray< TYPE, 9 > sigmaCauchy
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_F
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_C
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_E
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_S
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_invF
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_P
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_sIGma
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_h
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_H
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_invH
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_sigmaCauchy
 
TYPE J
 
TYPE eNergy
 
TYPE detH
 
TYPE detF
 
int gG
 Gauss point number. More...
 
CommonDatacommonDataPtr
 
MoFEM::VolumeElementForcesAndSourcesCore::UserDataOperatoropPtr
 pointer to finite element tetrahedral operator More...
 

Additional Inherited Members

- Static Protected Member Functions inherited from NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >
static auto resizeAndSet (MatrixBoundedArray< TYPE, 9 > &m)
 

Detailed Description

template<typename TYPE>
struct SmallStrainTranverslyIsotropic< TYPE >

Hook equation.

Definition at line 19 of file SmallTransverselyIsotropic.hpp.

Constructor & Destructor Documentation

◆ SmallStrainTranverslyIsotropic()

template<typename TYPE >
SmallStrainTranverslyIsotropic< TYPE >::SmallStrainTranverslyIsotropic ( )
inline

Member Function Documentation

◆ calculateAngles()

template<typename TYPE >
virtual MoFEMErrorCode SmallStrainTranverslyIsotropic< TYPE >::calculateAngles ( )
inlinevirtual

Reimplemented in SmallStrainTranverslyIsotropicDouble.

Definition at line 236 of file SmallTransverselyIsotropic.hpp.

236 {
239 }
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440

◆ calculateAxisAngleRotationalMatrix()

template<typename TYPE >
MoFEMErrorCode SmallStrainTranverslyIsotropic< TYPE >::calculateAxisAngleRotationalMatrix ( )
inline

Function to Calculate the Rotation Matrix at a given axis and angle of rotation.

This function computes the rotational matrix for a given axis of rotation and angle of rotation about that angle

Parameters
axVectorA vector representing the axis of rotation
axAngleAngle of rotation along the axis (in radians)

Definition at line 82 of file SmallTransverselyIsotropic.hpp.

82 {
84
85 aARotMat.resize(3,3,false);
86 aARotMat.clear();
87
88 TYPE norm = sqrt(pow(axVector[0],2) + pow(axVector[1],2) + pow(axVector[2],2));
89
90 aARotMat(0,0) = 1-((1-cos(axAngle))*(pow(axVector[1],2)+pow(axVector[2],2))/pow(norm,2));
91 aARotMat(1,1) = 1-((1-cos(axAngle))*(pow(axVector[0],2)+pow(axVector[2],2))/pow(norm,2));
92 aARotMat(2,2) = 1-((1-cos(axAngle))*(pow(axVector[0],2)+pow(axVector[1],2))/pow(norm,2));
93
94 aARotMat(0,1) = ((1-cos(axAngle))*axVector[0]*axVector[1]-norm*axVector[2]*sin(axAngle))/pow(norm,2);
95 aARotMat(1,0) = ((1-cos(axAngle))*axVector[0]*axVector[1]+norm*axVector[2]*sin(axAngle))/pow(norm,2);
96
97 aARotMat(0,2) = ((1-cos(axAngle))*axVector[0]*axVector[2]+norm*axVector[1]*sin(axAngle))/pow(norm,2);
98 aARotMat(2,0) = ((1-cos(axAngle))*axVector[0]*axVector[2]-norm*axVector[1]*sin(axAngle))/pow(norm,2);
99
100 aARotMat(1,2) = ((1-cos(axAngle))*axVector[1]*axVector[2]-norm*axVector[0]*sin(axAngle))/pow(norm,2);
101 aARotMat(2,1) = ((1-cos(axAngle))*axVector[1]*axVector[2]+norm*axVector[0]*sin(axAngle))/pow(norm,2);
102
104 }

◆ calculateElasticEnergy()

template<typename TYPE >
virtual MoFEMErrorCode SmallStrainTranverslyIsotropic< TYPE >::calculateElasticEnergy ( const NonlinearElasticElement::BlockData  block_data,
boost::shared_ptr< const NumeredEntFiniteElement >  fe_ptr 
)
inlinevirtual

calculate density of strain energy

Reimplemented from NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >.

Definition at line 284 of file SmallTransverselyIsotropic.hpp.

287 {
289
295 axAngle = -axAngle;
299
300 voigtStress.resize(6,false);
302 this->eNergy = 0.5*inner_prod(voigtStress,voightStrain);
304 }
static PetscErrorCode ierr
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:483
MoFEMErrorCode stressTransformation()
Function to Calculate Stress Transformation Matrix This function computes the stress transformation M...
MoFEMErrorCode strainTransformation()
Function to Calculate Strain Transformation Matrix This function computes the strain transformation ...
MoFEMErrorCode calculateAxisAngleRotationalMatrix()
Function to Calculate the Rotation Matrix at a given axis and angle of rotation.

◆ calculateFibreAngles()

template<typename TYPE >
MoFEMErrorCode SmallStrainTranverslyIsotropic< TYPE >::calculateFibreAngles ( )
inline

Definition at line 310 of file SmallTransverselyIsotropic.hpp.

310 {
312
313 try {
314
315 int gg = this->gG; // number of integration point
316 MatrixDouble &phi = (this->commonDataPtr->gradAtGaussPts["POTENTIAL_FIELD"][gg]);
317 normalizedPhi.resize(3,false);
318 double nrm2_phi = sqrt(pow(phi(0,0),2)+pow(phi(0,1),2)+pow(phi(0,2),2));
319 for(int ii = 0;ii<3;ii++) {
320 normalizedPhi[ii] = -phi(0,ii)/nrm2_phi;
321 }
322
323 axVectorDouble.resize(3,false);
324 const double zVec[3]={ 0.0,0.0,1.0 };
325 axVectorDouble[0] = normalizedPhi[1]*zVec[2]-normalizedPhi[2]*zVec[1];
326 axVectorDouble[1] = normalizedPhi[2]*zVec[0]-normalizedPhi[0]*zVec[2];
327 axVectorDouble[2] = normalizedPhi[0]*zVec[1]-normalizedPhi[1]*zVec[0];
328 double nrm2_ax_vector = norm_2(axVectorDouble);
329 const double eps = 1e-12;
330 if(nrm2_ax_vector<eps) {
331 axVectorDouble[0] = 0;
332 axVectorDouble[1] = 0;
333 axVectorDouble[2] = 1;
334 nrm2_ax_vector = 1;
335 }
336 axAngleDouble = asin(nrm2_ax_vector);
337
338 } catch (const std::exception& ex) {
339 std::ostringstream ss;
340 ss << "throw in method: " << ex.what() << std::endl;
341 SETERRQ(PETSC_COMM_SELF,1,ss.str().c_str());
342 }
343
345 }
static const double eps
static double phi
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
std::map< std::string, std::vector< MatrixDouble > > gradAtGaussPts

◆ calculateGlobalStiffnesMatrix()

template<typename TYPE >
MoFEMErrorCode SmallStrainTranverslyIsotropic< TYPE >::calculateGlobalStiffnesMatrix ( )
inline

Definition at line 225 of file SmallTransverselyIsotropic.hpp.

225 {
227
228 dR.resize(6,6,false);
229 noalias(dR) = prod(localStiffnessMatrix,strainRotMat);
230 globalStiffnessMatrix.resize(6,6,false);
231 noalias(globalStiffnessMatrix) = prod(stressRotMat,dR);
232
234 }
ublas::symmetric_matrix< TYPE, ublas::upper > localStiffnessMatrix

◆ calculateLocalStiffnesMatrix()

template<typename TYPE >
MoFEMErrorCode SmallStrainTranverslyIsotropic< TYPE >::calculateLocalStiffnesMatrix ( )
inline

◆ calculateP_PiolaKirchhoffI()

template<typename TYPE >
virtual MoFEMErrorCode SmallStrainTranverslyIsotropic< TYPE >::calculateP_PiolaKirchhoffI ( const NonlinearElasticElement::BlockData  block_data,
boost::shared_ptr< const NumeredEntFiniteElement >  fe_ptr 
)
inlinevirtual

Calculate global stress.

This is small strain approach, i.e. Piola stress is like a Cauchy stress, since configurations are notation distinguished.

Reimplemented from NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >.

Definition at line 250 of file SmallTransverselyIsotropic.hpp.

253 {
260 axAngle = -axAngle;
264
265 voigtStress.resize(6,false);
267 this->P.resize(3,3,false);
268 this->P(0,0) = voigtStress[0];
269 this->P(1,1) = voigtStress[1];
270 this->P(2,2) = voigtStress[2];
271 this->P(0,1) = voigtStress[3];
272 this->P(1,2) = voigtStress[4];
273 this->P(0,2) = voigtStress[5];
274 this->P(1,0) = this->P(0,1);
275 this->P(2,1) = this->P(1,2);
276 this->P(2,0) = this->P(0,2);
277
279 }

◆ calculateStrain()

template<typename TYPE >
MoFEMErrorCode SmallStrainTranverslyIsotropic< TYPE >::calculateStrain ( )
inline

Definition at line 29 of file SmallTransverselyIsotropic.hpp.

29 {
31 sTrain.resize(3,3,false);
32 noalias(sTrain) = this->F;
33 for(int dd = 0;dd<3;dd++) {
34 sTrain(dd,dd) -= 1;
35 }
36 sTrain += trans(sTrain);
37 sTrain *= 0.5;
38 voightStrain.resize(6,false);
39 voightStrain[0] = sTrain(0,0);
40 voightStrain[1] = sTrain(1,1);
41 voightStrain[2] = sTrain(2,2);
42 voightStrain[3] = 2*sTrain(0,1);
43 voightStrain[4] = 2*sTrain(1,2);
44 voightStrain[5] = 2*sTrain(2,0);
46 }
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

◆ strainTransformation()

template<typename TYPE >
MoFEMErrorCode SmallStrainTranverslyIsotropic< TYPE >::strainTransformation ( )
inline

Function to Calculate Strain Transformation Matrix
This function computes the strain transformation Matrix at a give axis and angle of rotation
One can also output the axis/angle rotational Matrix.

Definition at line 171 of file SmallTransverselyIsotropic.hpp.

171 {
173
174 strainRotMat.resize(6,6,false);
175 strainRotMat.clear();
176
177 strainRotMat(0, 0) = aARotMat(0,0) * aARotMat(0,0);
178 strainRotMat(0, 1) = aARotMat(1,0) * aARotMat(1,0);
179 strainRotMat(0, 2) = aARotMat(2,0) * aARotMat(2,0);
180 strainRotMat(0, 3) = aARotMat(1,0) * aARotMat(0,0);
181 strainRotMat(0, 4) = aARotMat(2,0) * aARotMat(1,0);
182 strainRotMat(0, 5) = aARotMat(0,0) * aARotMat(2,0);
183
184 strainRotMat(1, 0) = aARotMat(0,1) * aARotMat(0,1);
185 strainRotMat(1, 1) = aARotMat(1,1) * aARotMat(1,1);
186 strainRotMat(1, 2) = aARotMat(2,1) * aARotMat(2,1);
187 strainRotMat(1, 3) = aARotMat(1,1) * aARotMat(0,1);
188 strainRotMat(1, 4) = aARotMat(2,1) * aARotMat(1,1);
189 strainRotMat(1, 5) = aARotMat(0,1) * aARotMat(2,1);
190
191 strainRotMat(2, 0) = aARotMat(0,2) * aARotMat(0,2);
192 strainRotMat(2, 1) = aARotMat(1,2) * aARotMat(1,2);
193 strainRotMat(2, 2) = aARotMat(2,2) * aARotMat(2,2);
194 strainRotMat(2, 3) = aARotMat(1,2) * aARotMat(0,2);
195 strainRotMat(2, 4) = aARotMat(2,2) * aARotMat(1,2);
196 strainRotMat(2, 5) = aARotMat(0,2) * aARotMat(2,2);
197
198 strainRotMat(3, 0) = 2.0 * aARotMat(0,1) * aARotMat(0,0);
199 strainRotMat(3, 1) = 2.0 * aARotMat(1,1) * aARotMat(1,0);
200 strainRotMat(3, 2) = 2.0 * aARotMat(2,1) * aARotMat(2,0);
201 strainRotMat(3, 3) = ( aARotMat(1,1) * aARotMat(0,0) + aARotMat(0,1) * aARotMat(1,0) );
202 strainRotMat(3, 4) = ( aARotMat(2,1) * aARotMat(1,0) + aARotMat(1,1) * aARotMat(2,0) );
203 strainRotMat(3, 5) = ( aARotMat(0,1) * aARotMat(2,0) + aARotMat(2,1) * aARotMat(0,0) );
204
205 strainRotMat(4, 0) = 2.0 * aARotMat(0,2) * aARotMat(0,1);
206 strainRotMat(4, 1) = 2.0 * aARotMat(1,2) * aARotMat(1,1);
207 strainRotMat(4, 2) = 2.0 * aARotMat(2,2) * aARotMat(2,1);
208 strainRotMat(4, 3) = ( aARotMat(1,2) * aARotMat(0,1) + aARotMat(0,2) * aARotMat(1,1) );
209 strainRotMat(4, 4) = ( aARotMat(2,2) * aARotMat(1,1) + aARotMat(1,2) * aARotMat(2,1) );
210 strainRotMat(4, 5) = ( aARotMat(0,2) * aARotMat(2,1) + aARotMat(2,2) * aARotMat(0,1) );
211
212 strainRotMat(5, 0) = 2.0 * aARotMat(0,0) * aARotMat(0,2);
213 strainRotMat(5, 1) = 2.0 * aARotMat(1,0) * aARotMat(1,2);
214 strainRotMat(5, 2) = 2.0 * aARotMat(2,0) * aARotMat(2,2);
215 strainRotMat(5, 3) = ( aARotMat(1,0) * aARotMat(0,2) + aARotMat(0,0) * aARotMat(1,2) );
216 strainRotMat(5, 4) = ( aARotMat(2,0) * aARotMat(1,2) + aARotMat(1,0) * aARotMat(2,2) );
217 strainRotMat(5, 5) = ( aARotMat(0,0) * aARotMat(2,2) + aARotMat(2,0) * aARotMat(0,2) );
218
220 }

◆ stressTransformation()

template<typename TYPE >
MoFEMErrorCode SmallStrainTranverslyIsotropic< TYPE >::stressTransformation ( )
inline

Function to Calculate Stress Transformation Matrix This function computes the stress transformation Matrix at a give axis and angle of rotation
One can also output the axis/angle rotational Matrix.

Definition at line 113 of file SmallTransverselyIsotropic.hpp.

113 {
115
116 stressRotMat.resize(6,6,false);
117 stressRotMat.clear();
118
119 stressRotMat(0, 0) = aARotMat(0,0) * aARotMat(0,0);
120 stressRotMat(0, 1) = aARotMat(1,0) * aARotMat(1,0);
121 stressRotMat(0, 2) = aARotMat(2,0) * aARotMat(2,0);
122 stressRotMat(0, 3) = 2.0 * aARotMat(1,0) * aARotMat(0,0);
123 stressRotMat(0, 4) = 2.0 * aARotMat(2,0) * aARotMat(1,0);
124 stressRotMat(0, 5) = 2.0 * aARotMat(0,0) * aARotMat(2,0);
125
126 stressRotMat(1, 0) = aARotMat(0,1) * aARotMat(0,1);
127 stressRotMat(1, 1) = aARotMat(1,1) * aARotMat(1,1);
128 stressRotMat(1, 2) = aARotMat(2,1) * aARotMat(2,1);
129 stressRotMat(1, 3) = 2.0 * aARotMat(1,1) * aARotMat(0,1);
130 stressRotMat(1, 4) = 2.0 * aARotMat(2,1) * aARotMat(1,1);
131 stressRotMat(1, 5) = 2.0 * aARotMat(0,1) * aARotMat(2,1);
132
133 stressRotMat(2, 0) = aARotMat(0,2) * aARotMat(0,2);
134 stressRotMat(2, 1) = aARotMat(1,2) * aARotMat(1,2);
135 stressRotMat(2, 2) = aARotMat(2,2) * aARotMat(2,2);
136 stressRotMat(2, 3) = 2.0 * aARotMat(1,2) * aARotMat(0,2);
137 stressRotMat(2, 4) = 2.0 * aARotMat(2,2) * aARotMat(1,2);
138 stressRotMat(2, 5) = 2.0 * aARotMat(0,2) * aARotMat(2,2);
139
140 stressRotMat(3, 0) = aARotMat(0,1) * aARotMat(0,0);
141 stressRotMat(3, 1) = aARotMat(1,1) * aARotMat(1,0);
142 stressRotMat(3, 2) = aARotMat(2,1) * aARotMat(2,0);
143 stressRotMat(3, 3) = ( aARotMat(1,1) * aARotMat(0,0) + aARotMat(0,1) * aARotMat(1,0) );
144 stressRotMat(3, 4) = ( aARotMat(2,1) * aARotMat(1,0) + aARotMat(1,1) * aARotMat(2,0) );
145 stressRotMat(3, 5) = ( aARotMat(0,1) * aARotMat(2,0) + aARotMat(2,1) * aARotMat(0,0) );
146
147 stressRotMat(4, 0) = aARotMat(0,2) * aARotMat(0,1);
148 stressRotMat(4, 1) = aARotMat(1,2) * aARotMat(1,1);
149 stressRotMat(4, 2) = aARotMat(2,2) * aARotMat(2,1);
150 stressRotMat(4, 3) = ( aARotMat(1,2) * aARotMat(0,1) + aARotMat(0,2) * aARotMat(1,1) );
151 stressRotMat(4, 4) = ( aARotMat(2,2) * aARotMat(1,1) + aARotMat(1,2) * aARotMat(2,1) );
152 stressRotMat(4, 5) = ( aARotMat(0,2) * aARotMat(2,1) + aARotMat(2,2) * aARotMat(0,1) );
153
154 stressRotMat(5, 0) = aARotMat(0,0) * aARotMat(0,2);
155 stressRotMat(5, 1) = aARotMat(1,0) * aARotMat(1,2);
156 stressRotMat(5, 2) = aARotMat(2,0) * aARotMat(2,2);
157 stressRotMat(5, 3) = ( aARotMat(1,0) * aARotMat(0,2) + aARotMat(0,0) * aARotMat(1,2) );
158 stressRotMat(5, 4) = ( aARotMat(2,0) * aARotMat(1,2) + aARotMat(1,0) * aARotMat(2,2) );
159 stressRotMat(5, 5) = ( aARotMat(0,0) * aARotMat(2,2) + aARotMat(2,0) * aARotMat(0,2) );
160
162 }

Member Data Documentation

◆ aARotMat

template<typename TYPE >
ublas::matrix<TYPE> SmallStrainTranverslyIsotropic< TYPE >::aARotMat

Definition at line 69 of file SmallTransverselyIsotropic.hpp.

◆ axAngle

template<typename TYPE >
TYPE SmallStrainTranverslyIsotropic< TYPE >::axAngle

Definition at line 71 of file SmallTransverselyIsotropic.hpp.

◆ axAngleDouble

template<typename TYPE >
double SmallStrainTranverslyIsotropic< TYPE >::axAngleDouble

Definition at line 308 of file SmallTransverselyIsotropic.hpp.

◆ axVector

template<typename TYPE >
ublas::vector<TYPE> SmallStrainTranverslyIsotropic< TYPE >::axVector

Definition at line 70 of file SmallTransverselyIsotropic.hpp.

◆ axVectorDouble

template<typename TYPE >
VectorDouble SmallStrainTranverslyIsotropic< TYPE >::axVectorDouble

Definition at line 307 of file SmallTransverselyIsotropic.hpp.

◆ dR

template<typename TYPE >
ublas::matrix<TYPE> SmallStrainTranverslyIsotropic< TYPE >::dR

Definition at line 222 of file SmallTransverselyIsotropic.hpp.

◆ E_p

template<typename TYPE >
double SmallStrainTranverslyIsotropic< TYPE >::E_p

Definition at line 48 of file SmallTransverselyIsotropic.hpp.

◆ E_z

template<typename TYPE >
double SmallStrainTranverslyIsotropic< TYPE >::E_z

Definition at line 48 of file SmallTransverselyIsotropic.hpp.

◆ G_zp

template<typename TYPE >
double SmallStrainTranverslyIsotropic< TYPE >::G_zp

Definition at line 48 of file SmallTransverselyIsotropic.hpp.

◆ globalStiffnessMatrix

template<typename TYPE >
ublas::matrix<TYPE> SmallStrainTranverslyIsotropic< TYPE >::globalStiffnessMatrix

Definition at line 223 of file SmallTransverselyIsotropic.hpp.

◆ localStiffnessMatrix

template<typename TYPE >
ublas::symmetric_matrix<TYPE,ublas::upper> SmallStrainTranverslyIsotropic< TYPE >::localStiffnessMatrix

Definition at line 49 of file SmallTransverselyIsotropic.hpp.

◆ normalizedPhi

template<typename TYPE >
VectorDouble SmallStrainTranverslyIsotropic< TYPE >::normalizedPhi

Definition at line 306 of file SmallTransverselyIsotropic.hpp.

◆ nu_p

template<typename TYPE >
double SmallStrainTranverslyIsotropic< TYPE >::nu_p

Definition at line 48 of file SmallTransverselyIsotropic.hpp.

◆ nu_pz

template<typename TYPE >
double SmallStrainTranverslyIsotropic< TYPE >::nu_pz

Definition at line 48 of file SmallTransverselyIsotropic.hpp.

◆ sTrain

template<typename TYPE >
ublas::matrix<TYPE> SmallStrainTranverslyIsotropic< TYPE >::sTrain

Definition at line 25 of file SmallTransverselyIsotropic.hpp.

◆ strainRotMat

template<typename TYPE >
ublas::matrix<TYPE> SmallStrainTranverslyIsotropic< TYPE >::strainRotMat

Definition at line 164 of file SmallTransverselyIsotropic.hpp.

◆ stressRotMat

template<typename TYPE >
ublas::matrix<TYPE> SmallStrainTranverslyIsotropic< TYPE >::stressRotMat

Definition at line 106 of file SmallTransverselyIsotropic.hpp.

◆ tR

template<typename TYPE >
TYPE SmallStrainTranverslyIsotropic< TYPE >::tR

Definition at line 27 of file SmallTransverselyIsotropic.hpp.

◆ voightStrain

template<typename TYPE >
ublas::vector<TYPE> SmallStrainTranverslyIsotropic< TYPE >::voightStrain

Definition at line 26 of file SmallTransverselyIsotropic.hpp.

◆ voigtStress

template<typename TYPE >
ublas::vector<TYPE> SmallStrainTranverslyIsotropic< TYPE >::voigtStress

Definition at line 241 of file SmallTransverselyIsotropic.hpp.


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