1 /** @file gsElasticityFunctions.h
2 
3     @brief Provides useful classes derived from gsFunction which can be used
4     for visualization or coupling.
5 
6     This file is part of the G+Smo library.
7 
8     This Source Code Form is subject to the terms of the Mozilla Public
9     License, v. 2.0. If a copy of the MPL was not distributed with this
10     file, You can obtain one at http://mozilla.org/MPL/2.0/.
11 
12     Author(s):
13         A.Shamanskiy (2016 - ...., TU Kaiserslautern)
14 */
15 
16 #pragma once
17 
18 #include <gsCore/gsMultiPatch.h>
19 #include <gsElasticity/gsBaseUtils.h>
20 #include <gsIO/gsOptionList.h>
21 
22 namespace gismo
23 {
24 
25 /** @brief Compute Cauchy stresses for a previously computed/defined displacement field.
26  *         Can be pushed into gsPiecewiseFunction to construct gsField for visualization in Paraview.
27 */
28 template <class T>
29 class gsCauchyStressFunction : public gsFunction<T>
30 {
31 public:
32 
33     gsCauchyStressFunction(index_t patch, stress_components::components comp,
34                            const gsOptionList & options,
35                            const gsMultiPatch<T> * geometry,
36                            const gsMultiPatch<T> * displacement,
37                            const gsMultiPatch<T> * pressure = nullptr,
38                            const gsMultiPatch<T> * velocity = nullptr)
m_geometry(geometry)39         : m_geometry(geometry),
40           m_displacement(displacement),
41           m_pressure(pressure),
42           m_velocity(velocity),
43           m_patch(patch),
44           m_dim(m_geometry->patch(m_patch).parDim()),
45           m_options(options),
46           m_type(comp)
47     {}
48 
49 
50 
domainDim()51     virtual short_t domainDim() const
52     {
53         return m_geometry->patch(m_patch).parDim();
54     }
55 
targetDim()56     virtual short_t targetDim() const
57     {
58         switch(m_type)
59         {
60         case stress_components::von_mises: return 1;
61         case stress_components::all_2D_vector: return 3;
62         case stress_components::all_2D_matrix: return 2;
63         case stress_components::normal_3D_vector: return 3;
64         case stress_components::shear_3D_vector: return 3;
65         case stress_components::all_3D_matrix: return 3;
66         default: return 0;
67         };
68     }
69 
70     /** @brief Each column of the input matrix (u) corresponds to one evaluation point.
71      *         Columns of the output matrix (result) correspond to a set of stress components for vector cases,
72      *         or consist of one number in case of stress_type::von_mises,
73      *         or form square matrices concatinated in the col-direction.
74      */
eval_into(const gsMatrix<T> & u,gsMatrix<T> & result)75     virtual void eval_into(const gsMatrix<T> & u, gsMatrix<T> & result) const
76     {
77         switch (material_law::law(m_options.getInt("MaterialLaw")))
78         {
79         case material_law::hooke : linearElastic(u,result); return;
80         case material_law::saint_venant_kirchhoff : nonLinearElastic(u,result); return;
81         case material_law::neo_hooke_ln : nonLinearElastic(u,result); return;
82         case material_law::neo_hooke_quad : nonLinearElastic(u,result); return;
83         case material_law::mixed_hooke : mixedLinearElastic(u,result); return;
84         case material_law::mixed_neo_hooke_ln : mixedNonLinearElastic(u,result); return;
85         case material_law::muscle : mixedNonLinearElastic(u,result); return;
86         default: return;
87         }
88     }
89 
90 protected:
91 
92     /// size of the output matrix according to the m_type
outputCols(index_t inputCols)93     index_t outputCols(index_t inputCols) const
94     {
95         switch (m_type)
96         {
97         case stress_components::von_mises: return inputCols;
98         case stress_components::all_2D_vector: return inputCols;
99         case stress_components::all_2D_matrix: return 2*inputCols;
100         case stress_components::normal_3D_vector: return inputCols;
101         case stress_components::shear_3D_vector: return inputCols;
102         case stress_components::all_3D_matrix: return 3*inputCols;
103         default: return 0;
104         }
105     }
106     /// save components of the stress tensor to the output matrix according to the m_type
107     void saveStress(const gsMatrix<T> & S, gsMatrix<T> & result, index_t q) const;
108 
109     /// computation routines for different material laws
110     void linearElastic(const gsMatrix<T> & u, gsMatrix<T> & result) const;
111     void nonLinearElastic(const gsMatrix<T> & u, gsMatrix<T> & result) const;
112     void mixedLinearElastic(const gsMatrix<T> & u, gsMatrix<T> & result) const;
113     void mixedNonLinearElastic(const gsMatrix<T> & u, gsMatrix<T> & result) const;
114 
115 protected:
116     const gsMultiPatch<T> * m_geometry;
117     const gsMultiPatch<T> * m_displacement;
118     const gsMultiPatch<T> * m_pressure;
119     const gsMultiPatch<T> * m_velocity;
120     index_t m_patch;
121     short_t m_dim;
122     const gsOptionList & m_options;
123     stress_components::components m_type;
124 
125 }; // class definition ends
126 
127 
128 
129 /** @brief Compute jacobian determinant of the geometry mapping.
130  *         Can be pushed into gsPiecewiseFunction to construct gsField for visualization in Paraview.
131 */
132 template <class T>
133 class gsDetFunction : public gsFunction<T>
134 {
135 public:
136 
gsDetFunction(const gsMultiPatch<T> & geo,index_t patch)137     gsDetFunction(const gsMultiPatch<T> & geo, index_t patch)
138         : m_geo(geo),
139           m_patch(patch)
140     {}
141 
domainDim()142     virtual short_t domainDim() const { return m_geo.domainDim(); }
143 
targetDim()144     virtual short_t targetDim() const { return 1; }
145 
146     /** @brief Each column of the input matrix (u) corresponds to one evaluation point.
147      *         Each column of the output matrix is the jacobian determinant of the mapping at this point.
148      */
149     virtual void eval_into(const gsMatrix<T> & u, gsMatrix<T> & result) const;
150 
151 protected:
152 
153     gsMultiPatch<T> const & m_geo;
154     index_t m_patch;
155 
156 }; // class definition ends
157 
158 /** @brief Loading function to transfer fluid action to the solid.
159  * Used in Fluid-Structure Interaction simulation.
160  * Different parametrizations can be used for the geometry+ALE and velocity+pressure
161 */
162 template <class T>
163 class gsFsiLoad : public gsFunction<T>
164 {
165 public:
166 
gsFsiLoad(const gsMultiPatch<T> & geoRef,const gsMultiPatch<T> & ALEdisplacement,index_t patchGeo,boxSide sideGeo,const gsMultiPatch<T> & velocity,const gsMultiPatch<T> & pressure,index_t patchVelPres,T viscosity,T density)167     gsFsiLoad(const gsMultiPatch<T> & geoRef, const gsMultiPatch<T> & ALEdisplacement,
168               index_t patchGeo, boxSide sideGeo,
169               const gsMultiPatch<T> & velocity, const gsMultiPatch<T> & pressure,
170               index_t patchVelPres, T viscosity, T density)
171         : m_geo(geoRef),
172           m_ale(ALEdisplacement),
173           m_patchGeo(patchGeo),
174           m_sideGeo(sideGeo),
175           m_vel(velocity),
176           m_pres(pressure),
177           m_patchVP(patchVelPres),
178           m_viscosity(viscosity),
179           m_density(density)
180     {}
181 
domainDim()182     virtual short_t domainDim() const { return m_geo.domainDim(); }
183 
targetDim()184     virtual short_t targetDim() const { return m_geo.domainDim(); }
185 
186     /** @brief Each column of the input matrix (u) corresponds to one evaluation point.
187      *         Each column of the output matrix is the jacobian determinant of the mapping at this point.
188      */
189     virtual void eval_into(const gsMatrix<T> & u, gsMatrix<T> & result) const;
190 
191 protected:
192 
193     gsMultiPatch<T> const & m_geo;
194     gsMultiPatch<T> const & m_ale;
195     index_t m_patchGeo;
196     boxSide m_sideGeo;
197     gsMultiPatch<T> const & m_vel;
198     gsMultiPatch<T> const & m_pres;
199     index_t m_patchVP;
200     T m_viscosity;
201     T m_density;
202 
203 }; // class definition ends
204 
205 } // namespace ends
206 
207 
208 #ifndef GISMO_BUILD_LIB
209 #include GISMO_HPP_HEADER(gsElasticityFunctions.hpp)
210 #endif
211