G+Smo  25.01.0
Geometry + Simulation Modules
 
Loading...
Searching...
No Matches
gsMaterialMatrixIntegrate.h
Go to the documentation of this file.
1
16#pragma once
17
18#include <gsCore/gsFunction.h>
20
21namespace gismo
22{
23
24template <class T, enum MaterialOutput out> class gsMaterialMatrixIntegrateSingle;
25
26template <class T, enum MaterialOutput out>
27class gsMaterialMatrixIntegrate : public gsFunction<T>
28{
29public:
31 gsMaterialMatrixIntegrate( const gsMaterialMatrixContainer<T> & materialMatrices,
32 const gsFunctionSet<T> * deformed);
33
35 gsMaterialMatrixIntegrate( gsMaterialMatrixBase<T> * materialMatrix,
36 const gsFunctionSet<T> * deformed);
37
39 gsMaterialMatrixIntegrate( typename gsMaterialMatrixBase<T>::uPtr & materialMatrix,
40 const gsFunctionSet<T> * deformed);
41
43 gsMaterialMatrixIntegrate( const gsMaterialMatrixContainer<T> & materialMatrices,
44 const gsFunctionSet<T> * undeformed,
45 const gsFunctionSet<T> * deformed);
46
48 gsMaterialMatrixIntegrate( gsMaterialMatrixBase<T> * materialMatrix,
49 const gsFunctionSet<T> * undeformed,
50 const gsFunctionSet<T> * deformed);
51
53 gsMaterialMatrixIntegrate( typename gsMaterialMatrixBase<T>::uPtr & materialMatrix,
54 const gsFunctionSet<T> * undeformed,
55 const gsFunctionSet<T> * deformed);
56
58 ~gsMaterialMatrixIntegrate();
59
61 short_t domainDim() const {return 2;}
62
70 short_t targetDim() const { return this->piece(0).targetDim(); }
71
73 const gsFunction<T> & piece(const index_t p) const
74 {
75 return *m_pieces[p];
76 }
77
79 void eval_into(const gsMatrix<T>& /*u*/, gsMatrix<T>& /*result*/) const
81
82protected:
83 void _makePieces(const gsFunctionSet<T> * deformed);
84
85 void _makePieces(const gsFunctionSet<T> * undeformed, const gsFunctionSet<T> * deformed);
86
87protected:
88 gsMaterialMatrixContainer<T> m_materialMatrices;
89 const gsFunctionSet<T> * m_deformed;
90 mutable std::vector<gsMaterialMatrixIntegrateSingle<T,out> *> m_pieces;
91};
92
101template <class T, enum MaterialOutput out>
103{
104public:
105
108 gsMaterialMatrixBase<T> * materialMatrix,
109 const gsFunctionSet<T> * deformed);
110
113 gsMaterialMatrixBase<T> * materialMatrix,
114 const gsFunctionSet<T> * undeformed,
115 const gsFunctionSet<T> * deformed);
116
118
120 short_t domainDim() const {return 2;}
121
129 short_t targetDim() const { return targetDim_impl<out>(); }
130
131
132private:
134 template<enum MaterialOutput _out>
135 typename std::enable_if<_out==MaterialOutput::Density , short_t>::type targetDim_impl() const { return 1; };
136
138 template<enum MaterialOutput _out>
139 typename std::enable_if<_out==MaterialOutput::VectorN ||
140 _out==MaterialOutput::CauchyVectorN ||
141 _out==MaterialOutput::VectorM ||
142 _out==MaterialOutput::CauchyVectorM , short_t>::type targetDim_impl() const { return 3; };
143
145 template<enum MaterialOutput _out>
146 typename std::enable_if<_out==MaterialOutput::MatrixA ||
147 _out==MaterialOutput::MatrixB ||
148 _out==MaterialOutput::MatrixC ||
149 _out==MaterialOutput::MatrixD , short_t>::type targetDim_impl() const { return 9; };
150
152 template<enum MaterialOutput _out>
153 typename std::enable_if<_out==MaterialOutput::PStressN ||
154 _out==MaterialOutput::PStressM ||
155 _out==MaterialOutput::PStrainN ||
156 _out==MaterialOutput::PStrainM , short_t>::type targetDim_impl() const { return 2; };
157
159 template<enum MaterialOutput _out>
160 typename std::enable_if<_out==MaterialOutput::Stretch , short_t>::type targetDim_impl() const { return 3; };
161
163 template<enum MaterialOutput _out>
164 typename std::enable_if<_out==MaterialOutput::StretchDir, short_t>::type targetDim_impl() const { return 9; };
165
166protected:
168 void setPatch(index_t p) {m_pIndex = p; }
169
170public:
172 void eval_into(const gsMatrix<T>& u, gsMatrix<T>& result) const;
173
174private:
176 template<enum MaterialOutput _out>
177 typename std::enable_if<_out==MaterialOutput::Density , void>::type eval_into_impl(const gsMatrix<T>& u, gsMatrix<T>& result) const;
178
180 template<enum MaterialOutput _out>
181 typename std::enable_if<_out==MaterialOutput::VectorN ||
182 _out==MaterialOutput::CauchyVectorN , void>::type eval_into_impl(const gsMatrix<T>& u, gsMatrix<T>& result) const;
183
185 template<enum MaterialOutput _out>
186 typename std::enable_if<_out==MaterialOutput::VectorM ||
187 _out==MaterialOutput::CauchyVectorM , void>::type eval_into_impl(const gsMatrix<T>& u, gsMatrix<T>& result) const;
188
190 template<enum MaterialOutput _out>
191 typename std::enable_if<_out==MaterialOutput::MatrixA ||
192 _out==MaterialOutput::MatrixB ||
193 _out==MaterialOutput::MatrixC ||
194 _out==MaterialOutput::MatrixD , void>::type eval_into_impl(const gsMatrix<T>& u, gsMatrix<T>& result) const;
195
197 template<enum MaterialOutput _out>
198 typename std::enable_if<_out==MaterialOutput::PStressN ||
199 _out==MaterialOutput::PStressM , void>::type eval_into_impl(const gsMatrix<T>& u, gsMatrix<T>& result) const;
200
202 template<enum MaterialOutput _out>
203 typename std::enable_if<_out==MaterialOutput::PStrainN ||
204 _out==MaterialOutput::PStrainM , void>::type eval_into_impl(const gsMatrix<T>& u, gsMatrix<T>& result) const;
205
207 template<enum MaterialOutput _out>
208 typename std::enable_if<_out==MaterialOutput::Stretch , void>::type eval_into_impl(const gsMatrix<T>& u, gsMatrix<T>& result) const;
209
211 template<enum MaterialOutput _out>
212 typename std::enable_if<_out==MaterialOutput::StretchDir, void>::type eval_into_impl(const gsMatrix<T>& u, gsMatrix<T>& result) const;
213
214protected:
215
231 T getMoment() const { return getMoment_impl<out>(); }
232
233private:
235 template<enum MaterialOutput _out>
236 typename std::enable_if<_out==MaterialOutput::VectorN ||
237 _out==MaterialOutput::PStressN||
238 _out==MaterialOutput::CauchyVectorN, T>::type getMoment_impl() const { return 0; };
239
241 template<enum MaterialOutput _out>
242 typename std::enable_if<_out==MaterialOutput::VectorM ||
243 _out==MaterialOutput::PStressM||
244 _out==MaterialOutput::CauchyVectorM, T>::type getMoment_impl() const { return 1; };
245
247 template<enum MaterialOutput _out>
248 typename std::enable_if<_out==MaterialOutput::MatrixA, T>::type getMoment_impl() const { return 0; };
249
251 template<enum MaterialOutput _out>
252 typename std::enable_if<_out==MaterialOutput::MatrixB, T>::type getMoment_impl() const { return 1; };
253
255 template<enum MaterialOutput _out>
256 typename std::enable_if<_out==MaterialOutput::MatrixC, T>::type getMoment_impl() const { return 1; }; // must be 1
257
259 template<enum MaterialOutput _out>
260 typename std::enable_if<_out==MaterialOutput::MatrixD, T>::type getMoment_impl() const { return 2; }; // must be 2
261
263 template<enum MaterialOutput _out>
264 typename std::enable_if<!(_out==MaterialOutput::VectorN || _out==MaterialOutput::VectorM ||
265 _out==MaterialOutput::CauchyVectorN || _out==MaterialOutput::CauchyVectorM ||
266 _out==MaterialOutput::MatrixA || _out==MaterialOutput::MatrixB ||
267 _out==MaterialOutput::MatrixC || _out==MaterialOutput::MatrixD ||
268 _out==MaterialOutput::PStressN || _out==MaterialOutput::PStressM )
270
271protected:
272
280 void integrateZ_into(const gsMatrix<T> & u, const index_t moment, gsMatrix<T> & result) const;
281
305 void multiplyLinZ_into(const gsMatrix<T> & u, const index_t moment, gsMatrix<T> & result) const;
306
321 void multiplyZ_into(const gsMatrix<T> & u, index_t moment, gsMatrix<T> & result) const;
322
323protected:
324
333 gsMatrix<T> _eval3D(const gsMatrix<T>& u, const gsMatrix<T>& Z) const;
334
335
345 gsMatrix<T> _eval (const gsMatrix<T>& u) const;
346
347private:
349 template<enum MaterialOutput _out>
350 typename std::enable_if<_out==MaterialOutput::VectorN ||
351 _out==MaterialOutput::VectorM , gsMatrix<T>>::type eval3D_impl(const gsMatrix<T>& u, const gsMatrix<T>& Z) const;
352
354 template<enum MaterialOutput _out>
355 typename std::enable_if<_out==MaterialOutput::CauchyVectorN ||
356 _out==MaterialOutput::CauchyVectorM , gsMatrix<T>>::type eval3D_impl(const gsMatrix<T>& u, const gsMatrix<T>& Z) const;
357
358
360 template<enum MaterialOutput _out>
361 typename std::enable_if<_out==MaterialOutput::MatrixA ||
362 _out==MaterialOutput::MatrixB ||
363 _out==MaterialOutput::MatrixC ||
364 _out==MaterialOutput::MatrixD , gsMatrix<T>>::type eval3D_impl(const gsMatrix<T>& u, const gsMatrix<T>& Z) const;
365
367 template<enum MaterialOutput _out>
368 typename std::enable_if<_out==MaterialOutput::PStress ||
369 _out==MaterialOutput::PStressN ||
370 _out==MaterialOutput::PStressM , gsMatrix<T>>::type eval3D_impl(const gsMatrix<T>& u, const gsMatrix<T>& Z) const;
371
373 template<enum MaterialOutput _out>
374 typename std::enable_if<!(_out==MaterialOutput::VectorN || _out==MaterialOutput::VectorM ||
375 _out==MaterialOutput::CauchyVectorN || _out==MaterialOutput::CauchyVectorM ||
376 _out==MaterialOutput::MatrixA || _out==MaterialOutput::MatrixB ||
377 _out==MaterialOutput::MatrixC || _out==MaterialOutput::MatrixD ||
378 _out==MaterialOutput::PStressN || _out==MaterialOutput::PStressM )
379 , gsMatrix<T>>::type eval3D_impl(const gsMatrix<T>& /*u*/, const gsMatrix<T>& /*Z*/) const
381
382protected:
383 index_t m_pIndex;
384 gsMaterialMatrixBase<T> * m_materialMat;
385 gsMatrix<T> m_z;
386};
387
388} // namespace
389
390
391#ifndef GISMO_BUILD_LIB
392#include GISMO_HPP_HEADER(gsMaterialMatrixIntegrate.hpp)
393#endif
394
Interface for the set of functions defined on a domain (the total number of functions in the set equa...
Definition gsFunctionSet.h:219
A function from a n-dimensional domain to an m-dimensional image.
Definition gsFunction.h:60
This class defines the base class for material matrices.
Definition gsMaterialMatrixBase.h:33
memory::unique_ptr< gsMaterialMatrixBase > uPtr
Unique pointer for gsGeometry.
Definition gsMaterialMatrixBase.h:44
This class serves as the integrator of material matrices, based on gsMaterialMatrixBase.
Definition gsMaterialMatrixIntegrate.h:103
std::enable_if< _out==MaterialOutput::MatrixB, T >::type getMoment_impl() const
Implementation of getMoment for MaterialOutput::MatrixB; the moment is 1.
Definition gsMaterialMatrixIntegrate.h:252
std::enable_if<!(_out==MaterialOutput::VectorN||_out==MaterialOutput::VectorM||_out==MaterialOutput::CauchyVectorN||_out==MaterialOutput::CauchyVectorM||_out==MaterialOutput::MatrixA||_out==MaterialOutput::MatrixB||_out==MaterialOutput::MatrixC||_out==MaterialOutput::MatrixD||_out==MaterialOutput::PStressN||_out==MaterialOutput::PStressM), index_t >::type getMoment_impl() const
Implementation of getMoment for MaterialOutput other than VectorN, VectorM, MatrixA,...
Definition gsMaterialMatrixIntegrate.h:269
std::enable_if< _out==MaterialOutput::CauchyVectorN||_out==MaterialOutput::CauchyVectorM, gsMatrix< T > >::type eval3D_impl(const gsMatrix< T > &u, const gsMatrix< T > &Z) const
Specialisation of eval3D for vectors.
std::enable_if< _out==MaterialOutput::VectorN||_out==MaterialOutput::VectorM, gsMatrix< T > >::type eval3D_impl(const gsMatrix< T > &u, const gsMatrix< T > &Z) const
Specialisation of eval3D for vectors.
Definition gsMaterialMatrixIntegrate.hpp:410
std::enable_if< _out==MaterialOutput::StretchDir, short_t >::type targetDim_impl() const
Implementation of targetDim for principal stress directions.
Definition gsMaterialMatrixIntegrate.h:164
std::enable_if< _out==MaterialOutput::Density, void >::type eval_into_impl(const gsMatrix< T > &u, gsMatrix< T > &result) const
Specialisation of eval_into for densities.
Definition gsMaterialMatrixIntegrate.hpp:170
std::enable_if< _out==MaterialOutput::VectorN||_out==MaterialOutput::CauchyVectorN||_out==MaterialOutput::VectorM||_out==MaterialOutput::CauchyVectorM, short_t >::type targetDim_impl() const
Implementation of targetDim for stress tensors.
Definition gsMaterialMatrixIntegrate.h:142
std::enable_if<!(_out==MaterialOutput::VectorN||_out==MaterialOutput::VectorM||_out==MaterialOutput::CauchyVectorN||_out==MaterialOutput::CauchyVectorM||_out==MaterialOutput::MatrixA||_out==MaterialOutput::MatrixB||_out==MaterialOutput::MatrixC||_out==MaterialOutput::MatrixD||_out==MaterialOutput::PStressN||_out==MaterialOutput::PStressM), gsMatrix< T > >::type eval3D_impl(const gsMatrix< T > &, const gsMatrix< T > &) const
Specialisation of eval3D for other types.
Definition gsMaterialMatrixIntegrate.h:379
std::enable_if< _out==MaterialOutput::Stretch, short_t >::type targetDim_impl() const
Implementation of targetDim for principal stretch fields.
Definition gsMaterialMatrixIntegrate.h:160
std::enable_if< _out==MaterialOutput::MatrixA, T >::type getMoment_impl() const
Implementation of getMoment for MaterialOutput::MatrixA; the moment is 0.
Definition gsMaterialMatrixIntegrate.h:248
std::enable_if< _out==MaterialOutput::MatrixD, T >::type getMoment_impl() const
Implementation of getMoment for MaterialOutput::MatrixD; the moment is 2.
Definition gsMaterialMatrixIntegrate.h:260
std::enable_if< _out==MaterialOutput::MatrixA||_out==MaterialOutput::MatrixB||_out==MaterialOutput::MatrixC||_out==MaterialOutput::MatrixD, short_t >::type targetDim_impl() const
Implementation of targetDim for material tensors.
Definition gsMaterialMatrixIntegrate.h:149
void multiplyZ_into(const gsMatrix< T > &u, index_t moment, gsMatrix< T > &result) const
Multiplies the evaluation at z=0 with 2.0/(moment+1) * Thalf^(moment + 1)
Definition gsMaterialMatrixIntegrate.hpp:367
std::enable_if< _out==MaterialOutput::MatrixA||_out==MaterialOutput::MatrixB||_out==MaterialOutput::MatrixC||_out==MaterialOutput::MatrixD, void >::type eval_into_impl(const gsMatrix< T > &u, gsMatrix< T > &result) const
Specialisation of eval_into for the moments of the material matrices.
std::enable_if< _out==MaterialOutput::VectorN||_out==MaterialOutput::PStressN||_out==MaterialOutput::CauchyVectorN, T >::type getMoment_impl() const
Implementation of getMoment for MaterialOutput::VectorN; the moment is 0.
Definition gsMaterialMatrixIntegrate.h:238
void setPatch(index_t p)
Sets the patch index.
Definition gsMaterialMatrixIntegrate.h:168
std::enable_if< _out==MaterialOutput::Stretch, void >::type eval_into_impl(const gsMatrix< T > &u, gsMatrix< T > &result) const
Specialisation of eval_into for the stretches.
std::enable_if< _out==MaterialOutput::PStressN||_out==MaterialOutput::PStressM, void >::type eval_into_impl(const gsMatrix< T > &u, gsMatrix< T > &result) const
Specialisation of eval_into for the membrane and flexural principle stresses.
short_t targetDim() const
Target dimension.
Definition gsMaterialMatrixIntegrate.h:129
std::enable_if< _out==MaterialOutput::VectorM||_out==MaterialOutput::CauchyVectorM, void >::type eval_into_impl(const gsMatrix< T > &u, gsMatrix< T > &result) const
Specialisation of eval_into for the flexural stress tensor M.
void eval_into(const gsMatrix< T > &u, gsMatrix< T > &result) const
Implementation of eval_into, see gsFunction.
Definition gsMaterialMatrixIntegrate.hpp:160
std::enable_if< _out==MaterialOutput::MatrixC, T >::type getMoment_impl() const
Implementation of getMoment for MaterialOutput::MatrixC; the moment is 1.
Definition gsMaterialMatrixIntegrate.h:256
std::enable_if< _out==MaterialOutput::PStressN||_out==MaterialOutput::PStressM||_out==MaterialOutput::PStrainN||_out==MaterialOutput::PStrainM, short_t >::type targetDim_impl() const
Implementation of targetDim for principal stress fields.
Definition gsMaterialMatrixIntegrate.h:156
std::enable_if< _out==MaterialOutput::VectorM||_out==MaterialOutput::PStressM||_out==MaterialOutput::CauchyVectorM, T >::type getMoment_impl() const
Implementation of getMoment for MaterialOutput::VectorM; the moment is 1.
Definition gsMaterialMatrixIntegrate.h:244
std::enable_if< _out==MaterialOutput::VectorN||_out==MaterialOutput::CauchyVectorN, void >::type eval_into_impl(const gsMatrix< T > &u, gsMatrix< T > &result) const
Specialisation of eval_into for the membrane stress tensor N.
short_t domainDim() const
Domain dimension, always 2 for shells.
Definition gsMaterialMatrixIntegrate.h:120
void integrateZ_into(const gsMatrix< T > &u, const index_t moment, gsMatrix< T > &result) const
Integrates through-thickness using Gauss integration.
Definition gsMaterialMatrixIntegrate.hpp:280
void multiplyLinZ_into(const gsMatrix< T > &u, const index_t moment, gsMatrix< T > &result) const
Uses the top and bottom parts of the thickness to compute the integral exactly.
Definition gsMaterialMatrixIntegrate.hpp:333
std::enable_if< _out==MaterialOutput::Density, short_t >::type targetDim_impl() const
Implementation of targetDim for densities.
Definition gsMaterialMatrixIntegrate.h:135
gsMatrix< T > _eval3D(const gsMatrix< T > &u, const gsMatrix< T > &Z) const
Integrateuates the base class in 3D.
Definition gsMaterialMatrixIntegrate.hpp:401
std::enable_if< _out==MaterialOutput::MatrixA||_out==MaterialOutput::MatrixB||_out==MaterialOutput::MatrixC||_out==MaterialOutput::MatrixD, gsMatrix< T > >::type eval3D_impl(const gsMatrix< T > &u, const gsMatrix< T > &Z) const
Specialisation of eval3D for matrix.
std::enable_if< _out==MaterialOutput::PStrainN||_out==MaterialOutput::PStrainM, void >::type eval_into_impl(const gsMatrix< T > &u, gsMatrix< T > &result) const
Specialisation of eval_into for the membrane and flexural principle stresses.
T getMoment() const
Gets the moment based on the output type.
Definition gsMaterialMatrixIntegrate.h:231
gsMatrix< T > _eval(const gsMatrix< T > &u) const
Integrateuates the base class in 3D, but on Z=0.
Definition gsMaterialMatrixIntegrate.hpp:393
A matrix with arbitrary coefficient type and fixed or dynamic size.
Definition gsMatrix.h:41
#define short_t
Definition gsConfig.h:35
#define index_t
Definition gsConfig.h:32
#define GISMO_NO_IMPLEMENTATION
Definition gsDebug.h:129
Provides declaration of Function abstract interface.
Provides a container for material matrices for thin shells.
The G+Smo namespace, containing all definitions for the library.