Detailed Description

This module is about solving elasticity problems on thin shells using the Kirchhoff-Love shell formulations. References for this implementation include the original work by Kiendl et. al. (2009) [I], the PhD thesis of Josef Kiendl (2011) [II], the PhD thesis of Anmol Goyal (2015) [III], the MSc thesis of Hugo Verhelst (2019) [IV]. Anmol contributed to the first Kirchhoff-Love implementation of G+Smo and Hugo to the current one.

The variational formulation for the Kirchhoff-Love shell is:

\begin{eqnarray*} \int_\Omega \mathbf{f}\cdot\mathbf{v}\ + p\mathbf{\hat{n}}:\text{d}\Omega - \int_\Omega \mathbf{n}(\mathbf{u}):\mathbf{\varepsilon}^\prime(\mathbf{u},\mathbf{v}) + \mathbf{m}(\mathbf{u}):\mathbf{\kappa}^\prime(\mathbf{u},\mathbf{v}) \text{d}\Omega = 0 \end{eqnarray*}

The left-hand side of this formulation in fact is the residual of the problem, being the balance between external and internal forces. In this equation, \(\Omega\) is the domain, \( \mathbf{f} \) is a vector with a distributed load acting on the mid-plane of the shell, \( \mathbf{v} \) is the three-dimensional test function, \( p\) is a follower pressure acting on the shell, \( \mathbf{\hat{n}}\) is the shell normal, \( n \) is the normal-force tensor, \( \mathbf{\varepsilon}^\prime \) is the variation of the membrane strain tensor, \( \mathbf{m} \) is the bending moment tensor, \( \mathbf{\kappa}^\prime \) is the variation of the bending strain tensor and \(\mathbf{u} \) is the displacement of the mid-plane of the shell; thus the solution to the problem.

In order to solve the variational formulation for \( \mathbf{u} \), the second variation is derived, such that a lineararized form can be found or such that the Jacobian and for Newton iterations can be found. The second variation is:

\begin{eqnarray*} j(\mathbf{u},\mathbf{v},\mathbf{w}) = \int_\Omega \mathbf{n}(\mathbf{u}):\mathbf{\varepsilon}^{\prime\prime}(\mathbf{u},\mathbf{v},\mathbf{v}) + \mathbf{n}^\prime(\mathbf{u},\mathbf{v}):\mathbf{\varepsilon}^{\prime}(\mathbf{u},\mathbf{v}) + \mathbf{m}(\mathbf{u}):\mathbf{\kappa}^{\prime\prime}(\mathbf{u},\mathbf{v},\mathbf{w}) + \mathbf{m}^\prime(\mathbf{u},\mathbf{v}):\mathbf{\kappa}^\prime(\mathbf{u},\mathbf{v}) \text{d}\Omega - \int_\Omega p \mathbf{v}\cdot\mathbf{\hat{n}}^\prime\text{d}\Omega \end{eqnarray*}

In this equation, the primed ( \(\cdot^\prime\)) expressions again denote the first variation and the double-primed ( \(\cdot^{\prime\prime}\)) expressions denote second variations.

From the above equations, a system of equations can be assembled, which will later be done by the gsExprAssembler . In case of a lineararized system, only the single-primed expressions are non-zero, together with the term containing the external load.

The expressions for the strains, normal foce and bending moment are as follows (see sec. 3.2 and 3.3 of [III]):

\begin{align*} \varepsilon_{\alpha\beta} &= \frac{\partial\mathbf{c}}{\partial\theta_\alpha}\cdot\frac{\partial\mathbf{c}}{\partial\theta_\beta} - \frac{\partial\mathbf{C}}{\partial\theta_\alpha}\cdot\frac{\partial\mathbf{C}}{\partial\theta_\beta}\\ \varepsilon^\prime_{\alpha\beta} &= \frac{1}{2}\frac{\partial\mathbf{v}}{\partial\theta_\alpha}\cdot\frac{\partial\mathbf{c}}{\partial\theta_\beta} + \frac{1}{2}\frac{\partial\mathbf{c}}{\partial\theta_\alpha}\cdot\frac{\partial\mathbf{v}}{\partial\theta_\beta} \\ \kappa_{\alpha\beta} &= \frac{\partial^2\mathbf{C}}{\partial\theta_\alpha\partial\theta_\beta}\cdot\mathbf{\hat{N}} - \frac{\partial^2\mathbf{c}}{\partial\theta_\alpha\partial\theta_\beta}\cdot\mathbf{\hat{n}}\\ \kappa^\prime_{\alpha\beta} &= \frac{\partial^2\mathbf{v}}{\partial\theta_\alpha\partial\theta_\beta}\cdot\mathbf{\hat{n}} + \frac{\partial^2\mathbf{c}}{\partial\theta_\alpha\partial\theta_\beta}\cdot\mathbf{\hat{n}}^\prime(\mathbf{u})\\ \mathbf{\hat{n}}^\prime(\mathbf{v}) &= \mathbf{m}_{\mathbf{v}} - (\mathbf{\hat{n}}\cdot \mathbf{m}_{\mathbf{v}})\mathbf{\hat{n}} \\ \mathbf{m}_{\mathbf{v}} &= \frac{1}{\vert \frac{\partial\mathbf{c}}{\partial\theta_1}\times \frac{\partial\mathbf{c}}{\partial\theta_2}\vert } \left( \frac{\partial\mathbf{v}}{\partial\theta_1}\times \frac{\partial\mathbf{c}}{\partial\theta_2} + \frac{\partial\mathbf{c}}{\partial\theta_1} \times \frac{\partial\mathbf{v}}{\partial\theta_2} \right)\\ \varepsilon^{\prime\prime}_{\alpha\beta} &= \frac{\partial\mathbf{v}}{\partial\theta_\alpha}\cdot\frac{\partial\mathbf{w}}{\partial\theta_\beta}\\ \kappa^{\prime\prime}_{\alpha\beta} &= \frac{\partial^2\mathbf{v}}{\partial\theta_\alpha\partial\theta_\beta}\cdot\mathbf{\hat{n}}^\prime(\mathbf{w}) + \frac{\partial^2\mathbf{w}}{\partial\theta_\alpha\partial\theta_\beta}\cdot\mathbf{\hat{n}}^\prime(\mathbf{v}) + \frac{\partial^2\mathbf{c}}{\partial\theta_\alpha\partial\theta_\beta}\cdot\mathbf{\hat{n}}^{\prime\prime}(\mathbf{v},\mathbf{w})\\ \mathbf{\hat{n}}^{\prime\prime}(\mathbf{v},\mathbf{w}) &= \mathbf{m}_{\mathbf{v}}^\prime - ( \mathbf{m}_{\mathbf{v}}\cdot\mathbf{\hat{n}}^\prime(\mathbf{w}) + \mathbf{\hat{n}}\cdot\mathbf{m}_{\mathbf{v}}^\prime(\mathbf{w}) )\mathbf{\hat{n}} - ( \mathbf{\hat{n}}\cdot \mathbf{m}_{\mathbf{v}} )\mathbf{\hat{n}}^\prime(\mathbf{w})\\ \mathbf{m}_{\mathbf{v}}^{\prime}(\mathbf{w}) &= \mathbf{m}_{\mathbf{vw}}-(\mathbf{n}\cdot\mathbf{m}_{\mathbf{w}})\mathbf{m}_{\mathbf{w}}\\ \mathbf{m}_{\mathbf{vw}} &= \frac{1}{\vert \frac{\partial\mathbf{c}}{\partial\theta_1}\times \frac{\partial\mathbf{c}}{\partial\theta_2}\vert } \left( \frac{\partial\mathbf{v}}{\partial\theta_1}\times \frac{\partial\mathbf{w}}{\partial\theta_2} + \frac{\partial\mathbf{w}}{\partial\theta_1} \times \frac{\partial\mathbf{v}}{\partial\theta_2} \right)\\ n_{\alpha\beta} &= t \mathcal{C}^{\alpha\beta\gamma\delta} \varepsilon_{\alpha\beta}\\ n_{\alpha\beta}^\prime &= \mathcal{C}^{\alpha\beta\gamma\delta} \varepsilon_{\alpha\beta}^{\prime}\\ m_{\alpha\beta} &= t \mathcal{C}^{\alpha\beta\gamma\delta} \kappa{\alpha\beta}\\ m_{\alpha\beta}^\prime &= \mathcal{C}^{\alpha\beta\gamma\delta} \kappa_{\alpha\beta}^{\prime} \end{align*}

Here, \( t\) is the thickness of the shell and \( \mathcal{C}^{\alpha\beta\gamma\delta}\) is the material tensor (see [III] eq. 3.15).

The Kirchhoff-Love shell is implemented in the gsThinShellAssembler class and the constitutive relations are included in separate classes based on gsMaterialMatrixBase.

The gsThinShellAssembler is a class that implements the kirchhoff-Love_example.cpp and adds advanced functionalities. That is, the class uses the gsExprAssembler to assemble matrices and vectors for shell modelling. Key features of the class are:

Geometric nonlinearities, loading nonlinearities and material nonlinearities (via the gsMaterialMatrixBase),
Distributed loads, point loads and follower loads (nonlinear)
Multipatch modelling by accepting the gsMappedBasis and gsMappedSpline

The expressions for the assembly of the system of equations are partially coming from gsExpressions (e.g. the jac_expr, the grad_expr) and shell-specific expressions are implemented in the gsThinShellUtils.h. Furthermore, gsThinShellFunctions.h contains a helper-class for plotting the stress field within the shell.

The material relations are handled in separate classes, derived from gsMaterialMatrixBase. The material matrix generally behaves as a gsFunctionSet with implemented evaluation functions. The quantities that can be computed using material matrices are the thickness, the density, the (membrane/bending) stress, the force/moment through-thickness and principal stretches or stresses. Helper classes are employed to evaluate the stress or material tensor at a point or to find the integral over the thickness. The constitutive relations currently implemented are:

Linear materials (Saint-Venant Kirchhoff) in gsMaterialMatrixLinear
Linear laminates in gsMaterialMatrixComposite
Nonlinear materials models (Neo-Hookean, Mooney-Rivlin and Ogden) via gsMaterialMatrix

To evaluate (the integral of) any material quantity, the following classes can be used:

gsMaterialMatrixEval evaluates any quantity given a through-thickness coordinate,
gsMaterialMatrixIntegral evaluates the through-thickness integral of quantity.

Helper functions for defining material matrices are provided in getMaterialMatrix.h, which constructs a pointer to a material matrix based on a gsOptionsList

Typically, the aforementioned classes are used by first defining a gsMaterialMatrix and passing this to a gsThinShellAssembler.

// Define the material matrix options
options.addInt("Material","Material model: (0): SvK | (1): NH | (2): NH_ext | (3): MR | (4): Ogden",0);
options.addInt("Implementation","Implementation: (0): Composites | (1): Analytical | (2): Generalized | (3): Spectral",1);
// Construct a material matrix pointer, for a matrix with geometric dimension DIM and T as a double type
gsMaterialMatrixBase<T> * materialMatrix = getMaterialMatrix<DIM,T>(mp,t,parameters,rho,options);
// Construct a gsThinShellAssembler class with geometric dimension DIM, T as double type and the flag BENDING to specify whether bending stiffness should be taken into account. The
gsThinShellAssemblerBase<real_t>* assembler;
assembler = new gsThinShellAssembler<DIM,T, BENDING>(
        geom,
        basis,
        boundary_conditions,
        force,
        materialMatrix);
// Call assembly routines:
assembler->assemble();
    // Stiffness matrix
gsSparseMatrix<T> K = assember->matrix();
    // External force vector
gsVector<T>       vector = assember->rhs();
    // Mass matrix
assembler->assembleMass();
gsSparseMatrix<T> M = assember->matrix();
    // Mass matrix lumped into a vector
assembler->assembleMass(true);
gsSparseMatrix<T> M = assember->rhs();
    // Jacobian
assembler->assembleMatrix(deformed_geometry);
gsSparseMatrix<T> K_NL = assember->matrix();
    // Residual
assembler->assembleVector(deformed_geometry);
gsVector<T>       R = assember->rhs();

In the folder gsKLShell/tutorials, some tutorials are provided, explaining the use of the gsKLShell module. The tutorials can be compiled using the target gsKLShell-tutorials. The available tutorials are:

Few examples are provided in this module

example_shell2D.cpp : simple 2D elasticity including nonlinear materials
example_shell3D.cpp : 3D shell examples including nonlinear materials
gsMaterialMatrix_test.cpp : simple file to demonstrate evaluation of material matrices

[I] J. Kiendl, K.-U. Bletzinger, J. Linhard, and R. Wüchner, "Isogeometric shell analysis with Kirchhoff–Love elements," Comput. Methods Appl. Mech. Eng., vol. 198, no. 49–52, pp. 3902–3914, Nov. 2009.

[II] J. Kiendl, "Isogeometric analysis and shape optimal design of shell structures," Technische Universität München, 2011.

[III] A. Goyal, "Isogeometric Shell Discretizations for Flexible Multibody Dynamics," Technische Universität Kaiserslautern, 2015.

[IV] H. M. Verhelst, "Modelling Wrinkling Behaviour of Large Floating Thin Offshore Structures: An application of Isogeometric Structural Analysis for Post-Buckling Analyses," Delft University of Technology, 2019.

[8] H. M. Verhelst, "Isogeometric Analysis of Wrinkling.", PhD Thesis, TU Delft, 2024

[7] H.M. Verhelst, A. Mantzaflaris, M. Möller, J.H. Den Besten, "Goal-Adaptive Meshing of Isogeometric Kirchhoff-Love Shells", Engineering with Computers, 2024

[6] H.M. Verhelst, J.H. Den Besten, M. Möller, "An adaptive parallel arc-length method", Computers & Structures, 2024

[5] H.M. Verhelst, P. Weinmüller, A. Mantzaflaris, T. Takacs, D. Toshniwal, "A comparison of smooth basis constructions for isogeometric analysis", Computer Methods in Applied Mechanics and Engineering, 2024

[4] A. Farahat, H.M. Verhelst, J. Kiendl, M. Kapl, "Isogeometric analysis for multi-patch structured Kirchhoff–Love shells", Computer Methods in Applied Mechanics and Engineering, 2023

[3] H.M. Verhelst, M. Möller, J.H. Den Besten, A. Mantzaflaris, M.L. Kaminski, "Stretch-based hyperelastic material formulations for isogeometric Kirchhoff–Love shells with application to wrinkling", Computer-Aided Design. 2021 Oct 1;139:103075.

[2] H.M. Verhelst, M. Möller, J.H. Den Besten, F.J. Vermolen, M.L. Kaminski, "Equilibrium Path Analysis Including Bifurcations with an Arc-Length Method Avoiding A Priori Perturbations", Numerical Mathematics and Advanced Applications ENUMATH 2019: European Conference, Egmond aan Zee, The Netherlands, September 30-October 4. Cham: Springer International Publishing, 2020.

[1] H.M. Verhelst, "Modelling Wrinkling Behaviour of Large Floating Thin Offshore Structures: An application of Isogeometric Structural Analysis for Post-Buckling Analyses.", MSc. Thesis, TU Delft, 2019.

Author: Hugo Verhelst – h.m.v.nosp@m.erhe.nosp@m.lst@t.nosp@m.udel.nosp@m.ft.nl

Classes
struct	decodeMat_id< id >
	Decodes the material model and implementation. More...

struct	encodeMat_id< material, implementation >
	Encodes the material model and implementation. More...

class	gsMaterialMatrixBase< T >
	This class defines the base class for material matrices. More...

class	gsMaterialMatrixBaseDim< dim, T >
	This class defines the base class for material matrices. More...

class	gsMaterialMatrixComposite< dim, T >
	This class defines a linear material laminate. More...

class	gsMaterialMatrixContainer< T >
	This class serves as the evaluator of material matrices, based on gsMaterialMatrixBase. More...

class	gsMaterialMatrixEvalSingle< T, out >
	This class serves as the evaluator of material matrices, based on gsMaterialMatrixBase. More...

class	gsMaterialMatrixIntegrateSingle< T, out >
	This class serves as the integrator of material matrices, based on gsMaterialMatrixBase. More...

class	gsMaterialMatrixNonlinear< dim, T, matId, comp, mat, imp >
	This class defines hyperelastic material matrices. More...

class	gsShellStressFunction< T >
	Compute Cauchy stresses for a previously computed/defined displacement field. Can be pushed into gsPiecewiseFunction to construct gsField for visualization in Paraview. More...

class	gsThinShellAssembler< d, T, bending >
	Assembles the system matrix and vectors for 2D and 3D shell problems, including geometric nonlinearities and loading nonlinearities. The material nonlinearities are handled by the gsMaterialMatrixIntegrate class. More...

class	gsThinShellAssemblerBase< T >
	Base class for the gsThinShellAssembler. More...

class	gsThinShellAssemblerDWRBase< T >
	Base class for the gsThinShellAssembler. More...

struct	stress_type
	Specifies the type of stresses to compute. More...

Enumerations
enum	Implementation : short_t
	This class describes the way material models are implemented. More...

enum	Material : short_t
	This class describes a material model.

enum	MaterialOutput : short_t
	This class describes the output type. More...

enum	MatIntegration : short_t
	This class describes if an object is integrated through-thickness or not. More...

Functions
template<short_t d, class T >
gsMaterialMatrixBase< T > *	getMaterialMatrix (const gsMultiPatch< T > &mp, const gsFunctionSet< T > &thickness, const std::vector< gsFunctionSet< T > * > &parameters, const gsFunctionSet< T > &rho, const gsOptionList &options)
	Gets a material matrix based on options. More...

template<short_t d, class T >
gsMaterialMatrixBase< T > *	getMaterialMatrix (const gsMultiPatch< T > &mp, const gsFunctionSet< T > &thickness, const std::vector< gsFunctionSet< T > * > &parameters, const gsOptionList &options)
	Gets a material matrix based on options. More...

Enumeration Type Documentation

enum Implementation : short_t

strong

This class describes the way material models are implemented.

Composite: laminate material model Analytical: The expressions for Cijkl and Sij have to be provided in closed form Generalized: Uses a generalized way, where only the derivatives of psi have to be implemented Spectral: Implementation based on derivatives of psi w.r.t. principal stretches

enum MaterialOutput : short_t

strong

This class describes the output type.

Generic: unspecified output type Density VectorN: Membrane forces VectorM: Bending moments MatrixA: 0th order moment of the differential of the material matrix w.r.t. membrane strains MatrixB: 1st order moment of the differential of the material matrix w.r.t. membrane strains MatrixC: 0th order moment of the differential of the material matrix w.r.t. bending strains MatrixD: 1st order moment of the differential of the material matrix w.r.t. bending strains PStressN: membrane principal stress PStressN: bending principal stress Stretch: Principal stretch StretchDir: Principal stretch directions

enum MatIntegration : short_t

strong

This class describes if an object is integrated through-thickness or not.

NotIntegrated: The object has to be integrated Integrated: The object is integrated Constant: The object is constant through thickness, but is not integrated Linear: The object is linear through thickness, but is not integrated

Function Documentation

gsMaterialMatrixBase<T>* gismo::getMaterialMatrix	(	const gsMultiPatch< T > &	mp,
		const gsFunctionSet< T > &	thickness,
		const std::vector< gsFunctionSet< T > * > &	parameters,
		const gsFunctionSet< T > &	rho,
		const gsOptionList &	options
	)

Gets a material matrix based on options.

Parameters

[in]	mp	The undeformed geometry
[in]	thickness	The thickness
[in]	parameters	The parameters
[in]	rho	The density
[in]	options	The option list

Template Parameters

d	The dimension of the problem (2 = planar, 3 = surface)
T	Real type

Returns: The material matrix.

gsMaterialMatrixBase<T>* gismo::getMaterialMatrix	(	const gsMultiPatch< T > &	mp,
		const gsFunctionSet< T > &	thickness,
		const std::vector< gsFunctionSet< T > * > &	parameters,
		const gsOptionList &	options
	)

Gets a material matrix based on options.

Parameters

[in]	mp	The undeformed geometry
[in]	thickness	The thickness
[in]	parameters	The parameters
[in]	options	The option list

Template Parameters

d	The dimension of the problem (2 = planar, 3 = surface)
T	Real type

Returns: The material matrix.

Detailed Description

Classes

Enumerations

Functions

Enumeration Type Documentation

Function Documentation