An enhanced strain 3D element for large

Computational Mechanics 34 (2004) 38–52 Ó Springer-Verlag 2004
DOI 10.1007/s00466-004-0551-7
An enhanced strain 3D element for large deformation elastoplastic
thin-shell applications
R. A. Fontes Valente, R. J. Alves de Sousa, R. M. Natal Jorge
38
Abstract In this work a previously proposed solid-shell
finite element, entirely based on the Enhanced Assumed
Strain (EAS) formulation, is extended in order to account
for large deformation elastoplastic thin-shell problems. An
optimal number of 12 enhanced (internal) variables is
employed, leading to a computationally efficient performance when compared to other 3D or solid-shell enhanced
elements. This low number of enhanced variables is sufficient to (directly) eliminate either volumetric and
transverse shear lockings, the first one arising, for instance, in the fully plastic range, whilst the last appears for
small thickness’ values. The enhanced formulation comprises an additive split of the Green-Lagrange material
strain tensor, turning the inclusion of nonlinear kinematics a straightforward task. Finally, some shell-type
numerical benchmarks are carried out with the present
formulation, and good results are obtained, compared to
well-established formulations in the literature.
Keywords Solid-shell elements, Enhanced strains,
Volumetric and transverse shear lockings, Geometric and
material nonlinearities, Thin shells
1
Introduction
Finite element analysis of shell structures goes back in
time until the onset of the so-called degenerated approach
in works of Ahmad et al. [1] and Zienkiewicz et al. [101],
as well as in early papers of Ramm [78], and afterwards
with Hughes and Liu [55] and Hughes and Carnoy [57],
among others. Soon it was verified that brick elements
were prone to the appearance of volumetric and transverse
Received: 25 August 2003 / Accepted: 13 January 2004
Published online: 27 February 2004
R. A. F. Valente (&), R. J. A. de Sousa
Department of Mechanical Engineering,
University of Aveiro Campus de Santiago,
3810-193 Aveiro, Portugal
E-mail: [email protected]
R. A. F. Valente, R. M. N. Jorge
IDMEC, Faculty of Engineering,
University of Porto, Porto, Portugal
Funding by Ministe´rio da Cieˆncia e do Ensino Superior (FCT
and FSE) (Portugal) under grant PRAXIS XXI/ BD/21662/99; as
well as the funding by FEDER, under grant POCTI/EME/47289/
2002, are gratefully acknowledged.
shear locking effects. The first one is characteristic of
common metal plasticity models, where plastic deformation is taken to be isochoric or, in other words, incompressible [19]. The second one comes from the analysis of
thin shells, were the limit between ‘‘thick’’ and ‘‘thin’’
geometries is somewhat difficult to establish, with the
occurrence of locking not only strictly relying on thickness/length ratios, as demonstrated by Chapelle, Bathe and
co-workers [11, 12, 29, 30, 31].
In order to circumvent these parasitic phenomena,
selective reduced integration (or, equivalently, u/p formulation, mean-dilatation technique and B-bar methods)
– for volumetric locking – and the ‘‘mixed interpolation of
tensorial components’’/assumed strain method – for
transverse shear locking - had arisen as possible and
successful techniques. In the literature see, for instance,
references [4, 45, 53, 54, 61, 68, 69, 74, 85, 86, 96] for the
grounds of computational treatment of incompressibility,
in elastic and elastoplastic finite element cases, and also [9,
38, 56, 67] for earlier works dealing with transverse shear
locking.
In the specific case of shell elements, original planestress assumptions were enough to avoid or postpone
incompressibility issues in the nonlinear material range
[5, 47, 55, 78, 88], although at the expense of a rotation
tensor inclusion. As more generality was needed, higher
order theories including thickness change via extensible
director fields and/or ‘‘layerwise’’ approaches were
developed, including (or not) rotational variables, as in
[6, 7, 14–18, 22–25, 39, 40, 41, 52, 84, 89, 93], to name but
a few.
Despite the good results obtained by these formulations
in thick and thin shell problems, interest in trilinear
brick-type elements, resting just on translation-type degrees-of-freedom, has been increasing over the last decade.
A relative advantage gained with this kind of formulation
would then be the avoidance of a specific treatment for
rotation variables. On the other side, for this kind of elements, locking pathologies must be appropriately treated
while keeping its scope of application independent of
thickness values. Such an hexahedral solid element should
also naturally incorporate kinematical formulations typical
of shell approaches with, at the same time, the automatic
account for thickness variations.
According to Wriggers et al. [99], reliable threedimensional elements for shell-type applications with
finite strains can be obtained using the Enhanced Assumed
Strain (EAS) method of Simo and co-workers [87, 90, 92].
Representative lines of research in this field are, for
example, the intensive work of Schweizerhof et al. [37]
Freischlager and Schweizerhof [44], Harnau and
Schweizerhof [48], Hauptmann and Schweizerhof [49],
Hauptmann et al. [50], Klinkel and Wagner [62], Klinkel
et al. [63], Wagner et al. [98], Miehe [73] and recently VuQuoc and Tan [97] and Legay and Combescure [65]. All
these works have the common feature that enhanced
assumed strain, assumed strain method and/or selective
integration procedures have been combined in order to
obtain a wide class of solid-shell elements with good
performances. For typical shell problems, solid-shell elements can then represent an alternative with, as stated
before, a simpler formulation when compared to shell
elements, although more advantages can be specified. In
metal forming simulations involving two-sided contact
along the thickness direction (presence of blank-holder)
and in composites delamination problems (with a more
accurate evaluation of interlaminar shear and normal
stresses), numerical simulations can be effectively carried
out with this class of finite elements.
The grounds of the present work rely on the recent
paper of Alves de Sousa et al. [2], where a new class of
three-dimensional EAS elements for incompressible cases
was introduced. Starting with a sound analysis of the
deformation subspace granting the incompressibility
condition ðdiv u ¼ 0Þ, an enhanced strain field was
developed and introduced into the functional of the classical displacement-based solid element. It was then shown
that the inclusion of 6 enhanced variables, acting on the
volumetric components of the strain field, was sufficient to
avoid the volumetric locking phenomenon. A first proposal for a new 3D element, characterized by a total of 18
internal variables has proved to be effective in solving
general three-dimensional problems (HCiS18 solid element). The adopted EAS approach avoids the direct use of
classical selective reduced integration, which is consistent
only for material models with decoupled isochoric and
volumetric behavior. Another important feature was that
the element has proved to be reliable in thin shell problems. However, for the specific case of shell structures, the
authors have verified, also in the last reference, that the use
of only 12 enhanced parameters (leading to more computational efficiency) was enough for the obtention of
sound results. In this case, 6 enhanced variables are
responsible for the elimination of transverse shear locking
effects, without resorting to assumed strain methods, and
following previous works of the authors [28, 43]. This last
element, then coined HCiS12 solid-shell element, is now
extended and applied in large deformation elastoplastic
shell problems.
The distinguishing characteristic of the present
formulation can be summarized by the fact that only the
enhanced assumed strain method is used to simultaneously treat volumetric and transverse shear locking in
classical thin-shell problems. This point contrast with
the generalized use of the assumed natural strain
approach (for the transverse shear locking) and/or the
selective reduced integration technique (for nearincompressibility constraints) in well-established solidshell formulations in the literature. As a first step in the
formulation, linear benchmarks were provided in
reference [2], with the extension of the methodology to
account for nonlinear geometric as well as elastoplastic
problems being carried out in the present work. Besides
leading to an unified and ‘‘neat’’ formulation for the
solid-shell element as a whole, the present proposal
(coming from a subspace analysis detailed in [28] and
[2]) relies upon an enhanced strain field based on the
derivatives of a three-dimensional ‘‘bubble-function’’.
This specific choice of functions is grounded on improved results obtained for distorted meshes in incompressibility conditions, and also in bending-dominated
situations for two-dimensional problems, as can be
inferred from previous works of the authors (references
[26, 27], respectively).
This paper is organized as follows. In Sect. 2, kinematic
aspects of the HCiS12 solid-shell element are revisited,
along with expressions for stress and strain tensors in the
convective frame. Section 3 deals with the enhanced strain
variational formulation, establishing the grounds for the
finite element implementation described in Sect. 4. In the
latter, an overview of the algorithmic aspects related to
nonlinear material and geometric behaviors is performed,
focusing on the advantages of both the adoption of a
convective frame for tensorial representation and the
additive strain enhancement employed. Finally, results are
given in Sect. 5 for a class of benchmarks in shell and
solid-shell elements evaluation, and comparisons are
carried out with distinct formulations well-established in
the literature.
2
Kinematics of the solid-shell element
During deformation, it is theoretically useful to establish
some configurations related to whom each particle in the
analyzed body can be referred to. In this sense, it is possible to consider a reference (or material) and current (or
spatial) configurations M R3 and S R3 , respectively.
A third configuration employed is the parametric
configuration P R3 , defined by a set of curvilinear
(convective) coordinates
n ¼ ðn1 ; n2 ; n3 Þ 2 ( ½1; 1 ½1; 1 ½1; 1
ð1Þ
Without loss of generality, and within the incrementaliterative context, the reference configuration can be related
to a converged state ðnÞ (last increment) whereas the
current configuration points to the unknown configuration
ðn þ 1Þ (corresponding to the next load increment). For
the solid-shell topology treated in this work, any point in
the reference configuration can be defined by a position
vector as
n
xðnÞ ¼
n 1
n 1
1 þ n3 xu n1 ; n2 þ 1 n3 xl n1 ; n2
2
2
ð2Þ
The corresponding position after deformation (current
configuration) can be defined by an analogous expression,
now referred to state ðn þ 1Þ. In Eq. (2) it is worth noting
the use of position vectors xu ð; Þ and xl ð; Þ, of auxiliary
points corresponding to upper ðn3 ¼ þ1Þ and lower
surfaces ðn3 ¼ 1Þ, respectively, and referred to a fixed
orthonormal global frame ðe1 ; e2 ; e3 Þ. From expression (2)
39
it is straightforward to arrive at different formulations
commonly used in the literature ([73], [48], [97], among
others)
n
40
1 n 1 2 n 1 2 xu n ; n þ xl n ; n
2
1 þ n3 n xu n1 ; n2 n xl n1 ; n2
ð3aÞ
2
1 n ð3bÞ
¼ n xm n1 ; n2 þ n3 n a n1 ; n2 v n1 ; n2
2
xðnÞ ¼
The material Green-Lagrange strain tensor will be additively enhanced with incompatible (element-wise defined)
strain terms and, in conjunction with the stress tensor S,
will form the variational structure of the proposed element, as shown in the next section. In a departure from the
majority of works in the literature, only the enhanced
assumed strain (EAS) approach will be employed, in an
unified way, to solve either the transverse shear and volumetric lockings, with a minimum number of variables
and following previous works of the authors ([2, 26–28,
43).
The right-hand side of Eq. (3b) is characteristic of shell
formulations, with the introduction of the mid-surface
position vector xm, the shell
thickness a ¼ a n1 ; n2 and
the unit director v n1 ; n2 . In all the equations before, it is
implicit the definition of a preferred through-thickness
orientation, common in solid-shell approaches (see, for
instance, references [37, 48, 49, 50, 62, 63, 73, 97, 98]). In
this context, ðn1 Þ and ðn2 Þ represent inplane curvilinear
coordinate axes while ðn3 Þ denotes the through-thickness
orientation.
The displacement field for any point from a converged
state until the current position can be given by the relation
3
EAS variational formulation
The starting point of the present formulation is the
Enhanced Assumed Strain method, in its linear version as
originally presented by Simo and Rifai [87]. Even dealing
with nonlinearities, the proposed approach keeps the
original frame of additive enhancement over the displacement-based Green-Lagrange strain tensor, in a way
successfully advocated at first by Andelfinger and Ramm
[3] (linear cases), Bischoff and Ramm [18] (nonlinear
cases), after that by Klinkel and Wagner [62] and Klinkel
ðiÞ
nþ1
¼ nþ1 xðnÞðiÞ n xðnÞ
ð4Þ et al. [63] and, more recently, by Fontes Valente et al. [43]
n uðnÞ
for a given iteration ðiÞ. Also from the position vectors it is and Vu-Quoc and Tan [97]. As showed in these references,
this approach is indeed computationally simpler (and
possible to define the covariant base vectors from the
partial derivatives (dropping the iteration index for brevity leading to virtually the same results) than the one
advocated by Simo and Armero [90] and subsequently
reasons)
successfully used, for example, by Miehe [73].
onþ1 xðnÞ
In general terms, the starting point is the Hu-Washizunþ1
gk ðnÞ ¼
ð5aÞ
de Veubeke 3-field functional for static cases [18]
onk
Z
HWV
on xðnÞ
n
P
ðu;
E;
SÞ
¼
Ws ðEÞdV
gk ðnÞ ¼
ð5bÞ
onk
V
Z
1 T
each one directly related to its contravariant counterparts
F F I2 E dV Pext
ð9aÞ
þ S:
nþ1 k
g and n gk , respectively, with ðk ¼ 1; 2; 3Þ [8].
2
V
The relative deformation gradient between configuraZ
Z
tions ðnÞ and ðn þ 1Þ in the above defined convective basis
ext
u bqdV
þ
u tdS
ð9bÞ
P ¼
then takes the usual form
nþ1
nF
V
onþ1 x
¼ n ¼ nþ1 gk n gk
o x
ð6Þ
It is now possible to define the (displacement-based)
Green-Lagrange strain tensor Eu (as well as its
components), between states ðnÞ and ðn þ 1Þ as
1 nþ1 T nþ1
u n k n l
nF
n F I2 ¼ Ekl g g
2
1 n
onþ1n u onþ1n u n
nþ1 u
E
¼
g
þ
gl
k
n kl
2
onl
onk
onþ1n u onþ1n u
þ
onk
onl
nþ1 u
nE
¼
ð7aÞ
ð7bÞ
and including the second-order identity tensor I2 . The
strain tensor in Eqs. (7) is a material entity work-conjugated to the 2nd Piola-Kirchhoff stress tensor S
S ¼ Skl n gk n gl
ð8Þ
Sr
where the displacement ðuÞ, the Green-Lagrange strain ðEÞ
and the 2nd Piola-Kirchhoff stress ðSÞ are independent
variables. In the previous expressions, it is worth noting
the displacement-driven strain energy ðWs Þ, the traction
and volume force vectors ðtÞ and ðqbÞ, respectively, altogether with the corresponding prescribed fields ðtÞ and
over control area Sr and volume V. It is worth noting
ðbÞ,
that all variables are referred to the reference configuration, while the boundary conditions for the displacement
field were omitted in Eq. (9).
The total strain field is then assumed to be composed of
a compatible (displacement-based) (7) and incompatible
(element-wise) parts, in the form [87]
E ¼ Eu þ Ea
ð10Þ
where the left indexes relating to the configuration were
omitted for the sake of simplicity and the right indexes
reports to the driven field of strain (displacement – u or
enhanced variables – a). The substitution of Eq. (10) into
(9a), together with (9b), and the orthogonality condition
[87, 90, 92]
Once the major application field of the present formulation relies on shell structures, in each of the
Z
ð2 2 2Þ Gauss points within an element, local frames
S : Ea d V ¼ 0
ð11Þ are constructed, based on the inherent thickness direction and following the guidelines often taken for
V
degenerated shell elements [10]. This local frame is upreduces the number of independent variables in the oridated, from the last converged configuration, using the
ginal functional to just two. The weak form of this modiincremental rotation tensor ðnþ1n RÞ, extracted via a polar
fied functional is obtained with the Gateaux or directional
decomposition from its respective incremental deformaderivative, leading to the total variation [97]
tion gradient (6). The strain-displacement operator (15)
dPðu; Ea Þ ¼ dPint dPext
ð12aÞ is defined in this orthonormal frame and, as a conseZ
quence, the strain tensor assumes automatically a corooWs ðEu þ Ea Þ
int
tational character, as well as its work-conjugated stress
dP ¼ ðdEu þ dEa Þ :
dV
ð12bÞ
oðEu þ Ea Þ
tensor. This approach is responsible to provide a
V
straightforward constitutive update in the nonlinear
Z
Z
ext
du b qdV þ
du t dS
ð12cÞ range (for further details, see references [71, 72] and,
dP ¼
particularly, [43]).
V
Sr
The weak form can be expanded via a truncated Taylor
series
about
the solution (fixed point) at the kth state
ujk ; Ea jk [18]
4.2
Interpolation of the enhanced strain tensor
Along with the strain tensor ðEu Þ, it is necessary to
dP ujkþ1 ; Ea jkþ1 dP ujk ; Ea jk
define its enhanced counterpart ðEa Þ, responsible for the
a
a
attenuation
of volumetric and transverse shear locking
þ D½dP ujk ; E jk ðDu; DE Þ
inherent to a conventional displacement formulation. In
ð13Þ formal terms, the enhanced strain field is analogous to
where, in the present context, the ðDÞ operator relates to a the displacement-based one, being possible defined the
finite variation between ðkÞ and ðk þ 1Þ states. The finite relations
element interpolation for the displacement-based and
Ea ¼ Ma ðnÞa; dEa ¼ Ma ðnÞda; DEa ¼ Ma ðnÞDa
enhanced strain fields is described next, along with the
ð17Þ
explicit expression for the D½dP operator and the main
a
advantages of including the additive approach as in (10). where the enhanced operator ðM Þ involves a set of
internal variables ðaÞ, discontinuous between elements and
4
eliminated by a static condensation procedure at each
Finite element implementation including nonlinearities
element level. The enhanced strain field, in components
form, is equivalent to its displacement-based counterpart,
4.1
as presented in Eq. (16).
General aspects
The grounds of the present formulation points to the
In each element’s domain, the displacement field (with the recent work of Alves de Sousa et al. [2]. In this work, the
corresponding variation and increment) is interpolated in original (displacement-based) trilinear 3D element was
the form
analyzed in detail, with the description of its subspace
basis components able to enforce the incompressibility
u uh ¼ NðnÞd; du duh ¼ NðnÞdd;
condition ðdivðuÞ ¼ 0Þ. These linearly independent comð14Þ
Du Duh ¼ NðnÞDd
ponents are responsible for the range of possible deforwhere matrix N encopreses the usual isoparametric com- mation modes within an element while, on the other hand,
patible shape functions for a 8-node 3D element, relating the absence of specific modes leads invariably to the onset
the continuum displacement field and the corresponding of volumetric locking.
In order to suppress this parasitic phenomenon, the
vector of elemental degrees-of-freedom ðdÞ (in the present
missing
deformation modes are directly introduced into
case with a total of 24 translational components) [10].
the
formulation.
The enhanced assumed strain approach
The displacement-based Green-Lagrange strain tensor
was
then
chosen
to achieve this goal, and a first mixeddefined in Eqs. (7) can then be related (as well as its
formulated
element
with 18 internal variables was introvariations) to vector ðdÞ at the element level, i.e.,
duced in reference [2]. The HCiS18 element represented
Eu ¼ Mu ðnÞd; dEu ¼ Mu ðnÞdd; DEu ¼ Mu ðnÞDd ð15Þ then a general approach, with reliable results either in
where Mu is the conventional strain-displacement matrix, three-dimensional and thin-shell cases.
However, for shell applications, an approach involving
for the present case with 6 24 components, and the
lesser
variables then the latter was also introduced by
specific arrangement of the convective strain tensor
Alves
de
Sousa et al. [2] for linear applications. From the
components is according to
use
of
the
enhanced assumed strain method to treat
n
oT
u
u
u
u
u
u
u
transverse
shear locking in thin shells, adopted initially for
E ¼ En1 n1 En2 n2 En3 n3 En1 n2 En1 n3 En2 n3
ð16Þ
linear cases by Ce´sar de Sa´ et al. [28] and further
41
42
developed for nonlinear cases by Fontes Valente et al. [43],
it was possible to devise a solid-shell element with a
relative low computational cost (due to the use of 12
enhanced variables), with no rotations and good performances irrespective of thickness/length ratios. It is worth
emphasize that the referred 12 internal variables are
responsible for the treatment of volumetric and transverse
shear locking altogether, with no further enhancement or
assumed strain method being necessary. This HCiS12
solid-shell element can be characterized by an interpolation matrix as follows
2 oN
a
1
6 on
6
6 0
6
6
6 0
a M( HCiS12 ¼ 6
6
6 0
6
6 0
6
4
0
4.3
Algorithmic aspect of the linearized discrete weak
form and the constitutive update procedure
After the description of the interpolation functions and
variables for both the displacement-based and enhanced
strain fields, the second member of the right-hand side of
the linearized weak form (13) can be stated (dropping the
kth iteration indices) as [18, 62, 63, 97]
0
0
o2 Na
on1 on2
o2 Na
on1 on3
o2 Na
on2 on3
0
0
0
0
0
oNa
on2
0
o2 Na
on1 on2
o2 Na
on1 on3
o2 Na
on2 on3
0
0
0
0
0
0
oNa
on3
o2 Na
on1 on2
o2 Na
on1 on3
o2 Na
on2 on3
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
oNa
on1
oNa
on2
o2 Na
on1 on2
0
0
0
0
0
0
0
0
0
0
oNa
on1
oNa
on2
and resorting to the three-dimensional bubble function
1
Na ðnÞ ¼ ð1 n1 n1 Þð1 n2 n2 Þð1 n3 n3 Þ
ð19Þ
2
About the specific form of the enhanced operator (18), it is
noticeable the distinctions between the present formulation and classical works in the enhanced strain field for
three-dimensional elements [3, 24], namely in the last 9
columns of ðMa( jHCiS12 Þ matrix, and specifically in the low
number of internal variables involved in the present case.
Another point of interest in the present formulation is that
an enhancement on the in-plane displacement-based
strain terms is not strictly necessary. In fact, no substantial
gain were obtained if – instead of using the present
approach – a full 3D enhanced element (including 18
enhanced variables to improve each strain component
coming from the displacement formulation) is adopted for
membrane-dominated linear shell problems (as showed by
the authors in reference [2]).
After the adoption of the enhanced interpolation matrix
(18), defined in the convective frame ðn1 ; n2 ; n3 Þ, its
transformation onto the local corotational frame is carried
out with the mapping
Ma jHCiS12 ¼
instance, in the works of Andelfinger and Ramm [3] and
Klinkel et al. [63].
detJð0Þ
T0 Ma( jHCiS12
detJðnÞ
ð20Þ
3
0
7
7
0 7
7
7
0 7
7
7
0 7
7
0 7
7
5
ð18Þ
o2 Na
on1 on2
D½dPðd; aÞ ðDd; DaÞ ¼
oðdPint dPext Þ
ðDd; DaÞ
oðd; aÞ
ð21Þ
After including the interpolation functions (14), (15) and
(17), the variations in (21) take the matrix form
int
dP ðd; aÞ ¼ dd
ðTÞ
Z
V
þ da
ðTÞ
dPext ðdÞ ¼ ddðTÞ
ðMu ÞðTÞ S dV
Z
Z
V
V
Z
ðMa ÞðTÞ S dV
ð22aÞ
qdV
NðTÞ b
þ ddðTÞ
NðTÞt dS
ð22bÞ
Sr
a
where the ðM Þ matrix refers to the enhanced strain
operator in the local frame and presented in (20).
Focusing on the variation of the internal part (22a) of
the whole potential, it is possible to state that
oðdPint Þ
oðdPint Þ
Dd þ
Da
od
oa
nlg
¼ ddT ðKlg
uu þ Kuu ÞDd þ Kua Da
D½dPint ðDd; DaÞ ¼
involving the Jacobian operators directly relating the
convective and local frames, and evaluated at the center
of the element and in each Gauss point (Jð0Þ and JðnÞ
þ daT ½Kau Dd þ Kaa Da
respectively). Furthermore, transformation matrix ðT0 Þ
involves the components of the Jacobian inverse, being
ð8dd; 8daÞ
ð23Þ
also evaluated at each element’s center. This scaling
procedure is common in literature dealing with the
The linear and nonlinear geometric (initial stress) stiffness
lg
nlg
enhanced assumed strain approach, and the formal
matrices (Kuu and Kuu , respectively) are defined as in a
structure for the operator ðT0 Þ, can be found, for
fully displacement-based formulation [10]. The main
result of the inclusion of the enhanced parameters into the
variational formulation is the appearance of the coupling
stiffness matrices Kau and Kua , as well as the introduction
of the fully-enhanced stiffness operator Kaa , all of them
possessing the same structure as in the linear formulation
of Simo and Rifai [87]. In fact, the adopted additive approach (10) leads to a straightforward algorithmic extension from the linear case, with no inclusion of nonlinear
geometric stiffness matrices associated with the enhanced
variables, as in [90] and [73], and generating a final system
of equations, on matrix form, with the structure [18, 62,
63, 97]
"
lg
nlg
Kuu þ Kuu
Kua
#
Dd
Da
Kau
Kaa
8R T
9
R T
R
u T
>
< N bqdV þ N tdS ðM Þ S dV>
=
V
V
Sr
¼
R
>
>
ðMa ÞT S dV
:
;
V
ð24Þ
and elastoplastic counterpart. The implementation steps
followed, to some extent, those detailed for the specific
case of shell elements in the classical work of Brank et al.
[21].
5
Numerical examples
In the following, numerical benchmark problems are
considered in order to evaluate the performance of the
HCiS12 solid-shell element, as well as the validity of the
presented kinematical and constitutive formulation. Relative convergence tolerances for forces and displacements’
norms in Newton-Raphson procedure were set to
1:0 105 . Both convergence indicators were treated
simultaneously in each analysis. Problems including nonsmooth (unstable) load-deflection paths were solved
resorting to the ‘‘cylindrical’’ standard arc-length method
of Crisfield [33, 34] with the refinements introduced by de
Souza Neto and Feng [35], and following the algorithmic
guidelines detailed by the authors in [43].
The internal forceRvectors related toRdisplacement and
5.1
enhanced fields, ð ðMu ÞT S dV and ðMa ÞT S dV, respec- Thick-wall sphere problem with geometric nonlinearity
V
V
Enhanced strain methods are known to provide nontively) come from the discrete form of Eq. (12b).
Besides these nonlinear geometric aspects, the constit- stable response in the nonlinear range for large homoutive update of the stress tensor needs a specific treatment geneous compressive strain states. This pathology has
been identified by a number of authors, with a sound
due to nonlinearities coming from the account of
analysis being carried out initially by Wriggers and Reese
plasticity.
[100] and, subsequently, in references [64, 79, 80–82], to
The spatial incremental Cauchy stress tensor (in the
name but a few. The common point in all these works is
local frame) can be directly related to the (corotational)
2nd Piola-Kirchhoff stress tensor [13]. This approach leads the focusing on plane-strain and full three-dimensional
problems.
to the adoption of a hypoelastic constitutive model repAs the formulation for the present case, on the other
resentative of a Green-Naghdi objective stress rate, which
will derive in an additive constitutive update of the stress side, is devoted to simulation of typical shell problems, the
tensor also in the nonlinear material range. Doing so, and analysis of a (free) thick-wall sphere subjected to an
with the last converged increment representing the refer- internal pressure field was chosen in order to attest the
level of occurrence of numerical instabilities in compresence configuration (updated Lagrangian approach) it is
sion loading cases.
possible to define an evolution equation for the spatial
The problem is accounted for following the guidelines
stress field in terms of the material (mixed) strain field, in
of Kasper and Taylor [60], who considered a geometrically
the form
linear, nearly incompressible, problem. The authors have
nþ1
r ¼ n r þ nþ1n r
also previously analysed this fully-linear case in reference
[2], with encouraging results. In the present work, mate¼ n r þ nþ1n S
rial, boundary conditions and geometric data are kept the
ð25Þ same (elastic modulus E ¼ 250, inner radius R ¼ 7:5,
¼ n r þ C nþ1n E
i
external
radius
R
¼
10:0
and
a
varying
Poisson’s
ratio m).
e
In order to obtain a second order accurate procedure, the
However,
a
nonlinear
geometric
behavior
is
now
also
mid-point configuration between states ðnÞ and ðn þ 1Þ is
taken into account, altogether with a higher pressure load
used for the evaluation of the displacement-based and
(internal pressure p ¼ 2:0). Mesh topology also follows the
enhanced strain components
one presented in reference [60], being reproduced in
nþ12 u nþ1
nþ12 a nþ1
nþ1
ð26Þ Fig. 1, where a symmetric (undeformed) one-eight of the
nM
nM
nE ¼
nd þ
na
total volume is depicted.
The evolution of the radial displacements with the
Doing so, the continuous update of the local frame with
Poisson’s ratio, for node points located at both the internal
resort to the relative rotation tensor is theoretically
and external radius, is shown in Table 1. Apart from the
equivalent to successive push-forward and pull-back
variation of the results as incompressibility condition is
operations over the spatial Cauchy stress tensor [36].
In algorithmic terms, a conventional backward-Euler progressively achieved, it is worth reporting that, for valprocedure is employed in the update of the stress tensor. ues of m 0:49999999, convergence is completely lost and
hourglass patterns appears (even for the relatively low
In this sense, the incremental Green-Lagrange strain
tensor is consider to be additively composed of an elastic displacement values involved).
43
44
Fig. 1. Nonlinear geometric thick-wall sphere – Finite elements’
mesh adopted, with a total of 2100 solid-shell elements
Fig. 2. Membrane (in-plane) bending benchmark – Initial configurations for 3 different meshes
Table 1. Nonlinear geometric thick-wall sphere – Evolution of
radial displacements for an internal pressure level p ¼ 2:0
Poisson’s
ratio (m)
radial displacement ð102 Þ
Ri ¼ 7:5
0.3
0.49
0.499
0.4999
0.49999
0.499999
0.4999999
8.750
8.049
8.014
8.000
8.000
8.000
7.989
Re ¼ 10:0
6.305
4.589
4.509
4.510
4.540
4.583
4.598
5.2
Elastic large deflections (membrane) bending problem
This example relates to the analysis of a clamped beam
which is loaded by a transverse force F ¼ 1000 on its free
edge (see, for instance, references [88, 89]). The resultant
in-plane bending deformation is reproduced using a mesh
of 10 solid-shell elements, in both regular and distorted
patterns [15, 73]. The elastic properties refer to a bulk
modulus of j ¼ 83:33 105 and a Lame’s parameter of
l ¼ 38:46 105 , whilst the geometry is characterized by
the relationship height/width/length of h=w=l ¼
0:1=0:1=1:0 consistent units.
In order to infer the effect of mesh distortion in the
nonlinear geometric range, three meshes are considered as
represented in Fig. 2. The first two meshes are defined
following the previous references, the skewed pattern of
the second mesh obtained, in the present work, with the
translation of nodes (0.05 unit) along the beam axis. A
third mesh is taken into account (coming from the 90
rotation of the second one) in order to evaluate a real
tridimensional mesh distortion level. The load-deflection
curves for the displacement of nodes A and B is presented
in Fig. 3 and compared to the solution presented by Simo,
Fox and Rifai [88, 89]. For the shown meshes the results
are almost the same, being in good agreement with those
presented in the last references. It is worth noting that
results with HCiS12 solid-shell element also correspond to
Fig. 3. Membrane (in-plane) bending benchmark – Evolution of
displacements (in the load direction) with load-level for points A
and B
those obtained by Miehe [73] (with a solid-shell enhancedþassumed strain formulation) and by Betsch et al.
[15] (using a bilinear shell formulation incorporating
extensibility of the director field and also enhanced
strains).
5.3
Elastic large deflections (out-of-plane) bending problem
Consider now an elastic cantilever beam with length
L ¼ 10 and rectangular section with constant width w ¼ 1,
clamped in one end and subjected to (out-of-plane) point
loads on the opposite (free) end. This example has been
treated by a number of authors, either with extensibledirector shell or solid-shell formulations, such as in [24],
[41], [73], [76] and [89], among others.
Following the last two references, the thickness of the
beam is taken as a ¼ 0:1, and the material properties are
defined as E ¼ 107 and m ¼ 0:3. The external load considered has a constant value of F ¼ 40 k, where for the
present case k is a geometrical factor, function of the
thickness, and defined as ðk ¼ 103 a3 Þ. In the present
work three mapped meshes were considered, with 10, 16
and 20 HCiS12 elements along the length direction.
For the small deformation theory, a solution of 16:0
consistent unities for the tip displacement can be advanced
[41], according to linear beam theory. In this case, HCiS12
element gives results of 15:820 (10 elements’ mesh), 15:880
(16 elements’ mesh) and 15:890 (20 elements’ mesh).
In case of large rotations and displacements, the same
load level as before is now applied in 10 equal steps, as
proposed by el-Abbasi and Meguid [41]. Following this
reference, solutions are compared to a theoretical one
relying on inextensional elastica, and coming from the
work of Frisch-Fay [46]. The analysis of the present results
for the three meshes and the theoretical one is presented in
Fig. 4. For the three meshes, it is noticeable the performance of the proposed element, with solutions in good
agreement with the reference one. Still referring to the
work of el-Abbasi and Meguid [41], no Poisson’s locking
appears with HCiS12 element in this test case.
Focusing on the 16 1 1 mesh, and based on the
proposal of Hauptmann et al. [50], a set of numerical
analysis with different Poisson’s ratio is carried out. The
load level is the same as before, being likewise applied in
10 equal steps. The results presented in Fig. 5 clearly show
a virtually insensitivity of the load-deflection curve for the
various Poisson’s coefficients presented.
Elastic constitutive parameters are Young modulus
E ¼ 2:0685 107 and Poisson coefficient m ¼ 0:3. The
length of the cylinder is L ¼ 3:048, with mean radius of
R ¼ 1:016 and thickness a ¼ 0:03. Nominal load in Fig. 6
is Ftot ¼ 1600 k, where the load factor employed is
supposed to vary between k ¼ 0:0 to k ¼ 1:0. Due to
symmetry, only 1=4 of the structure needs to be meshed.
The results obtained with HCiS12 element for the
deflection of the point under the concentrated load is
represented in Fig. 7, and compared to those of Brank et al.
[20]. The value of the cylinder radius is highlighted in the
picture, establishing the physical limit of the deformation.
From the last picture it is clear that the refined mesh of
20 20 1 is necessary in order the results can be in good
agreement with those from the reference. However, the
5.4
Nonlinear geometric pinching of a clamped cylinder
In this test problem an elastic cylindrical shell, fully
clamped at one end, is subjected to a pair of concentrated
loads at its free end. Following references dealing with
shell elements [20, 43, 58, 94] a regular mesh of 16 16
Fig. 5. Out-of-plane bending benchmark – Influence of Poisson’s
elements is employed (Fig. 6). An additional mesh topol- coefficient on the deflection of 16 1 1 mesh
ogy with 20 20 solid-shell elements over the midsurface
of the cylinder is also considered. In both cases, only 1
element along the radial direction is employed.
Fig. 4. Out-of-plane bending benchmark – Evolution of displacements with load level for different meshes and m ¼ 0:3
Fig. 6. Clamped cylinder problem – Mesh, loading and boundary
conditions
45
46
Fig. 7. Clamped cylinder problem – Deflection curve for loaded
points
Fig. 9. Shallow roof problem – Geometric model with load and
boundary conditions
Fig. 10. Shallow roof problem – Results for points A and B
Fig. 8. Clamped cylinder problem – Sequence of deformed configurations for displacements of: a 0.26 R, b 0.58 R, c 0.72 R, d
0.97 R
overall predictive capability of the proposed formulation
in an example traditionally analyzed with shell elements is
worth noting. In Fig. 8 successive deformed equilibrium
configurations at different load stages are shown for the
coarser mesh over one half of the cylinder, until a value of
displacement near its radius.
5.5
Unstable behavior of a shallow roof structure
In this classical shell example, the snap-through and snapback load-displacement path of a cylindrical shell is analyzed (see references [32, 33, 42, 43, 51, 70, 75, 83, 95], to
name but a few) . The structure, schematically represented
in Fig. 9, is mapped with 5 5 HCiS12 elements over onequarter of its surface, along with 2 elements in the thickness direction. The imposition of these two elements is
related to the proper reproduction of the hinged support
over the straight edges. About the input data for this
problem, linear dimensions are L1 ¼ 508:0 and
L2 ¼ 507:15, with a nominal radius R ¼ 2540:0 and
thickness value of a ¼ 6:35. Material parameters are
E ¼ 3102:75 and m ¼ 0:3. The maximum load level attained
is equal to Ftot ¼ 1000.
The displacement along OZ direction for points A and
B is reproduced in Fig. 10, plotted against the reference
load level and compared to the solution advanced by
Horrigmoe and Bergan [51], coming from a shell formulation. It is noticeable the agreement between both
solutions, with the proposed approach spanning the
whole nonlinear range in a total of 39 arc-length controlled steps.
5.6
Geometric- and material-nonlinear analysis
of a pinched hemispherical shell
The well-known nonlinear geometric hemispherical shell
test case, introduced by Simo et al. [88], is now considered
with the inclusion of elastoplastic effects. This combined
nonlinear behavior has been previously investigated by
Masud and Tham [72], based on a class of reduced
integrated solid elements [66, 71].
According to reference [72], geometric and elastic
parameters, as well as restraint conditions, are kept the
same as in the work of Simo et al. [88], while a new set of
plastic properties in coherent unities (initial yield stress
r0 ¼ 6:825 105 ; isotropic linear hardening factor
Hiso ¼ 6:825 106 ) are now introduced.
In [72], a maximum load level of F ¼ 400:0 is proposed
(in opposition to the value of 200:0 in [88]), along with
mapped meshes of 16 16 2 and 18 18 2 elements.
For the sake of comparison, the results using the last
topology are included in this work.
In the present simulation, a coarser mesh of 16 16 1
HCiS12 solid-shell elements was adopted. The results
obtained for the displacement along the OX and OY directions (traction and compression external loads, respectively) are represented in Fig. 11. It is noticeable the good
correspondence between the present and reference results,
even with the lower number of elements in the earlier case. It
is also worth noting that the complete deformation path is
obtained in 20 steps, 5 times less than the number of
increments adopted by [72]. The deformed configuration
for the maximum load level is shown in Fig. 12.
the hole shell is showed. Previous publications dealing
with this example include, among others, the works of
Peng and Crisfield [77], Sansour and Bufler [83], Jiang and
Chernuka [59], Brank et al. [20], Masud et al. [71], Ibrahimbegovic´ et al. [58] and Fontes Valente et al. [43]. In all
these cases, only geometric nonlinearities were considered,
whereas the nonlinear material analysis have previously
been considered by Masud and Tham [72].
The starting point relies on an initially cylindrical
geometry characterized by a length of L ¼ 10:35, radius of
R ¼ 4:953 and a constant thickness value a ¼ 0:094.
Material properties are: Young modulus E ¼ 10:5 106 ,
and Poisson’s ratio m ¼ 0:3125. No boundary conditions
are applied to the free ends of the shell, being the load pair
responsible for the equilibrium of the cylinder. Plastic
parameters are the initial yield stress r0 ¼ 1:05 105 , and
a linear isotropic hardening coefficient of
Hiso ¼ 10:5 105 , as adopted in reference [72].
Mesh topologies are analogous to those employed in
[43], that is: 12 8 1 HCiS12 elements, regular and
distorted, and also 16 8 1 elements (regular only). The
numbers relate to elements along the periphery, the semilength and radial directions, respectively. In Fig. 13, the
distorted pattern with 12 8 1 elements is schematically
5.7
Elastic and elastoplastic stretching
of a short cylinder with free ends
A cylindrical shell is submitted to a pair of concentrated
forces, inducing large displacements and rotations to the
elements. A schematic representation is presented in
Fig. 13 where, due to symmetry reasons, only one octant of
Fig. 12. Elastoplastic pinched hemispherical shell – Orthogonal
views of deformed geometry
Fig. 11. Elastoplastic pinched hemispherical shell – Loaddeflection diagram
Fig. 13. Stretching of a cylinder – Schematic view (one octant)
with 12 8 1 elements in a distorted mesh
47
48
represented. The main goal is to evaluate the effects of
mesh distortions in the final solution of the proposed
solid-shell element. For this distorted mesh, the element
most far away from the load point is 10 times larger than
the smaller one.
The maximum load level, both for the elastic and elastoplastic cases, is Ftot ¼ 40000 k, where k ranges from
0:0 until 1:0. Results for the deflection of points A and B,
as a function of the k parameter, are given respectively in
Figs. 14 and 15 (for the purely nonlinear geometric case)
and in Figs. 16 and 17 (for the elastoplastic nonlinear
geometric case). For the present HCiS12 element, all curves
were obtained with the arc-length procedure, as stated
before, in a total of 26 automatic load steps.
In both the elastic and elastoplastic cases, the present
solution tends to follow the one obtained in [71] and [72],
but with less elements. For the coarser meshes, in the
elastic case, there is a deviation of results when compared
to the earlier reference. It is also interesting to note, for
this example, the increase in displacements when distorted
elements are employed, in opposition to what happened
with the class of fully enhanced shell elements, recently
presented by the authors [43]. For the refined mesh
employed, the results are acceptable either in the elastic
and elastoplastic regimens, still keeping a lower number of
increments to achieve the full deformation path.
5.8
Elastoplastic analysis of a simply supported plate
with pressure loads
In this test case the inflation of a square plate is analyzed.
This example has been treated before in a number of
references, including Miehe [73], Betsch and Stein [17],
Eberlein and Wriggers [40], Doll et al. [37], Hauptmann
et al. [50] and Harnau and Schweizerhof [48]. Along these
references, a variety of mesh topologies is employed, while
Fig. 16. Stretching of a cylinder – Results for point A, elastoFig. 14. Stretching of a cylinder – Results for point A, elastic case plastic case
Fig. 15. Stretching of a cylinder – Results for point B, elastic
cases
Fig. 17. Stretching of a cylinder – Results for point B, elastoplastic case
strong shape modifications. Earlier works dealing with this
problem includes the contributions of Simo and Kennedy
[91], Wriggers et al. [99], Hauptmann and Schweizerhof
[49], Miehe [73], Eberlein and Wriggers [40] and Wagner
et al. [98], among others.
For the present case, comparisons will be carried out
with the results presented in references [73] and [99].
Following these guidelines, the cylinder geometry is characterized by a relation radius/thickness/length of
300=3=300 consistent unities, respectively. The boundary
conditions are such that the circular shape of the cylinder’s
end is preserved although free deformation in longitudinal
direction is allowed. Once each element has upper and
lower nodes, a ‘‘hard support’’, in the sense of Wagner et al.
[98], is considered. A von Mises yield criterion is assumed,
with yield stress r0 ¼ 24:3 and a linear isotropic hardening
parameter Hiso ¼ 300. Elastic parameters are j ¼ 2500 and
l ¼ 1154 [73]. Two discretizations were accounted for: a
coarse mesh of 16 16 1 elements [98, 99] and a finer
one with 32 32 1 elements [40, 73, 99], applied over
one eighth of the cylinder, benefiting from symmetric
geometry. The simulation is performed in a load-controlled
way and the whole path is covered, with the arc-length
method, in 65 steps for both meshes. The obtained results
for the displacement value in the external load direction
(against this same load) are presented in Fig. 20. A good
agreement between the present and reference results is
obtained, with the coarser mesh leading to a load-deflection path similar to the obtained with the refined mesh. For
this last case, the HCiS12 element agrees quite well with the
results presented in the work of Miehe [73]. It is also worth
noting that the apparent ‘‘snap-through’’ behavior noticeable in the 16 16 1 mesh almost disappear with the
adoption of the 32 32 1 mesh, reproducing a mesh5.9
dependent behavior already pointed out by Hauptmann
Elastoplastic nonlinear geometric response
and Schweizerhof [49]. In Fig. 21 a sequence of deformed
of a pinched cylinder
This example deals with the elastoplastic deformation of a configurations for the adopted meshes is presented, startthin-walled cylinder, submitted to a pair of concentrated ing from the undeformed point and ranging until the
forces. It is a classical test to analyze the behavior of finite physically acceptable displacement value of 300 consistent
unities (equal to the cylinder radius).
element in problems involving localized plasticity and
in references [37, 48] a higher number of Gauss points (6)
is employed for the integration in thickness direction.
The geometric properties for this test are defined by the
relations length/width/thickness of 508=508=2:54 consistent unities [17]. The plate is submitted to an uniformly
distributed load of p ¼ 60 p0 , where p0 ¼ 102 . Material
properties are given as E ¼ 6:9 104 , with a perfectly
plastic law characterized by an initial yield stress of
r0 ¼ 248 (Hiso ¼ 0). Boundary conditions only restrain the
displacements in the direction normal to the plate, being
applied to just the lower nodes of the mesh (defined over
one quarter of the plate). Do to this fact, it is valid the
occurrence of a sort of ‘‘edge-rotations’’ and, as the pressure
value increases, the plate assumes a ‘‘pillow-type’’ deformation mode, changing from a bending dominated deformation (in the beginning) to a membrane dominated one.
From the meshes available for comparison, a
15 15 1 topology, refined in the corners (inspired by
Eberlein and Wriggers [40]) and a second one with
24 24 1 mapped elements (following Betsch and Stein
[17]) are adopted. The so-called ‘‘6-parameter’’ formulation on reference [40] is the one chosen for comparison
purposes. For both meshes, the maximum load level is
attained in 45 load steps, and the resulting out-of-plane
displacement curves are shown in Fig. 18, where the central node of the plate was monitored. It can be seen that
the present results are in agreement with both reference
solutions, although the predictive capability of few HCiS12
elements approaches the values obtained with the refined
mesh, even with a two-point integration rule across
thickness direction. The deformed configurations for the
full load and both meshes are represented in Fig. 19.
Fig. 18. Simply supported plate – Results for the out of plane
deflection of the center node
Fig. 19. Simply supported plate – Deformed configurations for
the analyzed meshes at full load
49
problems typically solved by shell elements. The adopted
additive enhanced strain field, involving a minimum
number of internal variables, is responsible for the treatment of transverse shear locking (for low limiting thickness values) and volumetric locking (arising on the full
plastic range) in an unified way, without resort to other
techniques, such as the reduced integration approach or
the assumed natural strain method. The presented
numerical simulations show high predictability characteristics when compared to other well-established (shell
and solid-shell) formulations in the literature.
50
References
Fig. 20. Pinched cylinder problem – Load deflection curve for the
center point in the middle of the cylinder
Fig. 21. Pinched cylinder problem – Sequence of deformed
configurations for both meshes. a Initial configurations;
Deformed meshes at b w 100, c w 200, d w 300
6
Conclusions
In this work a recently proposed formulation for a new
solid-shell element, entirely based on the Enhanced
Assumed Strain approach, was extended and tested in
large deformation elastoplastic and nonlinear geometric
1. Ahmad S, Irons BM, Zienkiewicz OC (1970) Analysis of
thick and thin shell structures by curved finite elements. Int
J Numer. Meth. Eng. 2: 419–451
2. Alves de Sousa RJ, Natal Jorge RM, Fontes Valente RA,
Ce´sar de Sa´ JMA (2003) A new volumetric and shear locking-free 3D enhanced strain element. Eng. Comput. 20: 896–
925 DOI 10.1108/02644400310502036
3. Andelfinger U, Ramm E (1993) EAS-Elements for twodimensional, three-dimensional, plate and shells and their
equivalence to HR-elements. Int. J. Numer. Meth. Eng. 36:
1413–1449
4. Argyris J, Dunne P, Angelopoulus T, Bichat B (1974) Large
natural strains and some special difficulties due to nonlinearity and incompressibility in finite elements. Comput.
Meth. Appl Mech Eng 4: 219–278
5. Basar Y, Ding Y, Kra¨tzig WB (1992) Finite-rotation shell
elements via mixed formulation. Comput. Mech. 10: 289–306
6. Basar Y, Ding Y (1997) Shear deformation models for largestrain shell analysis. Int J Solids. Struct. 34: 1687–1708
7. Basar Y, Itskov M (1999) Constitutive model and finite
element formulation for large strain elasto-plastic analysis of
shells. Comput. Mech. 23: 466–481
8. Basar Y, Weicher D (2000) Nonlinear Continuum Mechanics
of Solids. Springer-Verlag, Berlin Heidelberg New York
9. Bathe KJ, Dvorkin E (1986) A formulation of general shell
elements - the use of mixed interpolation of tensorial
components. Int. J. Numer. Meth. Eng. 22: 697–722
10. Bathe KJ (1996) Finite Element Procedures. 2nd edn.
Prentice-Hall, New Jersey
11. Bathe KJ, Iosilevich A, Chapelle D (2000) An evaluation of
the MITC shell elements. Comput. Struct. 75: 1–30
12. Bathe KJ, Chapelle D, Lee PS (2003) A shell problem ‘‘highly
sensitive’’ to thickness changes. Int. J. Numer. Meth. Eng. 57:
1039–1052
13. Belytschko T, Liu WK, Moran B (2000) Nonlinear Finite
Elements for Continua and Structures. John Wiley & Sons,
West Sussex, England
14. Betsch P, Stein E (1995) An assumed strain approach
avoiding artificial thickness straining for a nonlinear 4-node
shell element. Comm. Numer. Meth. Eng. 11: 899–909
15. Betsch P, Gruttmann F, Stein E (1996) A 4-node finite shell
element for the implementation of general hyperelastic 3Delasticity at finite strains. Comput. Meth. Appl. Mech. Eng.
130: 57–79
16. Betsch P, Stein E (1996) A nonlinear extensible 4-node shell
element based on continuum theory and assumed strain
interpolations. J. Nonlinear. Sci. 6: 169–199
17. Betsch P, Stein E (1999) Numerical implementation of
multiplicative elasto-plasticity into assumed strain elements
with application to shells at large strains. Comput. Meth.
Appl. Mech. Eng. 179: 215–245
18. Bischoff M, Ramm E (1997) Shear deformable shell elements
for large strains and rotations. Int. J. Numer. Meth. Eng. 40:
4427–4449
19. Boe¨r CR, Rebelo N, Rydstad H, Schro¨der G (1986) Process
Modelling of Metal Forming and Thermomechanical Treatment. Springer-Verlag, Berlin Heidelberg New York
20. Brank B, Peric D, Damjanic FB (1995) On implementation
of a nonlinear four node shell finite element for thin multilayered elastic shells. Comput. Mech. 16: 341–359
21. Brank B, Peric D, Damjanic FB (1997) On large deformations of thin elasto-plastic shells: implementation of a finite
rotation model for a quadrilateral shell element. Int.
J. Numer. Meth. Eng 40: 689–726
22. Brank B, Korelc J, Ibrahimbegovic´ A (2002) Nonlinear shell
problem formulation accounting for through-the-thickness
stretching and its finite element interpolation. Comput.
Struct. 80: 699–717
23. Braun M, Bischoff M, Ramm E (1994) Nonlinear shell formulations for complete three-dimensional constitutive laws
including composites and laminates. Comput. Mech.
15: 1–18
24. Bu¨chter N, Ramm E, Roehl D (1994) Three-dimensional
extension of nonlinear shell formulation based on the enhanced assumed strain concept. Int. J. Numer. Meth. Eng.
37: 2551–2568
25. Cardoso RPR, Yoon JW, Gra´cio JJA, Barlat F, Ce´sar de Sa´
JMA (2002) Development of a one point quadrature shell
element for nonlinear applications with contact and
anisotropy. Comput. Meth. Appl. Mech. Eng. 191: 5177–5206
26. Ce´sar de Sa´ JMA, Natal Jorge RM (1999) New enhanced
strain elements for incompressible problems. Int. J. Numer.
Meth. Eng. 44: 229–248
27. Ce´sar de Sa´ JMA, Areias PMA, Natal Jorge RM (2001)
Quadrilateral elements for the solution of elasto-plastic finite
strain problems. Int. J. Numer. Meth. Eng. 51: 883–917
28. Ce´sar de Sa´ JMA, Natal Jorge RM, Fontes Valente RA,
Areias PMA (2002) Development of shear locking-free shell
elements using an enhanced assumed strain formulation. Int
J Numer Meth Eng 53: 1721–1750
29. Chapelle D, Bathe KJ (1998) Fundamental considerations for
the finite element analysis of shell structures. Comput.
Struct. 66: 19–36
30. Chapelle D, Bathe KJ (2000) The mathematical shell model
underlying general shell elements. Int. J. Numer. Meth. Eng.
48: 289–313
31. Chapelle D, Bathe KJ (2003) The Finite Element Analysis of
Shells - Fundamentals. Springer-Verlag Berlin Heidelberg
New York
32. Cho C, Park HC, Lee SW (1998) Stability analysis using a
geometrically nonlinear assumed strain solid-shell element
model. Finite. El. Anal. Design. 29: 121–135
33. Crisfield MA (1981) A fast incremental/iterative solution
procedure that handles ‘‘snap-through’’. Comput. Struct .13:
55–62
34. Crisfield MA (1983) An arc-length method including line
searches and accelerations. Int. J. Numer. Meth. Eng. 19:
1269–1289
35. de Souza Neto EA, Feng YT (1999) On the determination of
the path direction for arc-length methods in the presence of
bifurcations and ‘‘snap-backs’’. Comput. Meth. Appl. Mech.
Eng. 179: 81–89
36. Doghri I (2000) Mechanics of Deformable Solids - Linear,
Nonlinear, Analytical and Computational Aspects. SpringerVerlag, Berlin Heidelberg New York
37. Doll S, Schweizerhof K, Hauptmann R, Freischla¨ger C
(2000) On volumetric locking of low-order solid and solidshell elements for finite elastoviscoplastic deformations and
selective reduced integration. Eng. Comput. 17: 874–902
38. Dvorkin E, Bathe KJ (1984) A continuum mechanics based
four-node shell element for general nonlinear analysis. Eng.
Comput. 1: 77–88
39. Dvorkin EN, Pantuso D, Repetto EA (1995) A formulation of
the MITC4 shell element for finite strain elastoplastic analysis. Comput. Meth. Appl. Mech. Eng. 125: 17–40
40. Eberlein R, Wriggers P (1999) Finite element concepts for
finite elastoplastic strains and isotropic stress response in
shells: theoretical and computational analysis. Comput.
Meth. Appl. Mech. Eng. 171: 243–279
41. el-Abbasi N, Meguid SA (2000) A new shell element
accounting for through-thickness deformation. Comput.
Meth. Appl. Mech. Eng. 189: 841–862
42. Eriksson A, Pacoste C (2002) Element formulation and
numerical techniques for stability problems in shells.
Comput. Meth. Appl. Mech. Eng. 191: 3775–3810
43. Fontes Valente RA, Natal Jorge RM, Cardoso RPR, Ce´sar de
Sa´ JMA, Gra´cio JJ (2003) On the use of an enhanced transverse shear strain shell element for problems involving large
rotations. Comput. Mech. 30: 286–296 DOI 10.1007/s00466002-0388-x
44. Freischlager C, Schweizerhof K (1996) On a systematic
development of trilinear three-dimensional solid elements
based on Simo’s enhanced strain formulation. Int. J. Solids.
Struct. 33: 2993–3017
45. Fried I (1974) Finite element analysis of incompressible
material by residual energy balancing. Int. J. Solids. Struct.
10: 993–1002
46. Frisch-Fay R (1962) Flexible Bars, Butterworths, London
47. Gruttmann F, Wagner W, Wriggers P (1992) A nonlinear
quadrilateral shell element with drilling degrees of freedom.
Arch. Appl. Mech. 62: 474–486
48. Harnau M, Schweizerhof K (2002) About linear and quadratic ‘‘solid-shell’’ elements at large deformations. Comput.
Struct. 80: 805–817
49. Hauptmann R, Schweizerhof K (1998) A systematic development of ‘‘solid-shell’’ element formulations for linear and
non-linear analyses employing only displacement degrees of
freedom. Int. J. Numer. Meth. Eng. 42: 49–69
50. Hauptmann R, Schweizerhof K, Doll S (2000) Extension of
the solid-shell concept for application to large elastic and
large elastoplastic deformations. Int. J. Numer. Meth. Eng.
49: 1121–1141
51. Horrigmoe G, Bergan PG (1978) Nonlinear analysis of freeform shells by flat finite elements. Comput. Meth. Appl.
Mech. Eng. 16: 11–35
52. Huettel C, Matzenmiller A (1999) Consistent discretization
of thickness strains in thin shells including 3D-material
models. Comm. Numer. Meth. Eng. 15: 283–293
53. Hughes TJR (1977) Equivalence of finite elements for nearly
incompressible elasticity. J. Appl. Mechanics 44: 181–183
54. Hughes TJR (1980) Generalisation of selective integration
procedures to anisotropic and nonlinear media. Int.
J. Numer. Meth. Eng. 15: 1413–1418
55. Hughes TJR, Liu WK (1981) Nonlinear finite element analysis of shells: Part I. Three-dimensional shells. Comput.
Meth. Appl. Mech. Eng. 26: 331–362
56. Hughes TJR, Tezduyar TE (1981) Finite elements based
upon Mindlin plate theory with particular reference to the
four-node bilinear isoparametric element. J. Appl. Mechanics 48: 587–596
57. Hughes TJR, Carnoy E (1983) Nonlinear finite shell element
formulation accounting for large membrane strains.
Comput. Meth. Appl. Mech. Eng. 39: 69–82
58. Ibrahimbegovic´ A, Brank B, Courtois P (2001) Stress
resultant geometrically exact form of classical shell model
and vector-like parameterization of constrained finite rotations. Int. J. Numer. Meth. Eng. 52: 1235–1252
59. Jiang L, Chernuka MW (1994) A simple four-noded corotational shell element for arbitrarily large rotations. Comput.
Struct. 53: 1123–1132
51
52
82. Reese S (2003) On a consistent hourglass stabilization
60. Kasper EP, Taylor RL (2000) A mixed-enhanced strain
technique to treat large inelastic deformations and
method. Part I: Geometrically linear problems. Comput.
thermo-mechanical coupling in plane strain problems. Int.
Struct .75: 237–250
J. Numer. Meth. Eng. 57: 1095–1127
61. Key SW (1969) A variational principle for incompressible
83. Sansour C, Bufler H (1992) An exact finite rotation shell
and nearly incompressible anisotropic elasticity. Int. J. Soltheory, its mixed variational formulation and its finite eleids. Struct. 5: 951–964
ment implementation. Int. J. Numer. Meth. Eng. 34: 73–115
62. Klinkel S, Wagner W (1997) A geometrical nonlinear brick
element based on the eas-method. Int. J. Numer. Meth. Eng. 84. Sansour C (1995) A theory and finite element formulation of
shells at finite deformations including thickness change:
40: 4529–4545
Circunventing the use of a rotation tensor. Arch. Appl.
63. Klinkel S, Gruttman F, Wagner W (1999) A continuum
Mech. 65: 194–216
based 3D-shell element for laminated structures. Comput.
85. Simo JC, Taylor RL (1982) Penalty-function formulations for
Struct. 71: 43–62
incompressible non-linear elastostatics. Comput. Meth.
64. Korelc J, Wriggers P (1996) Consistent gradient formulation
Appl. Mech. Eng. 35: 107–118
for a stable enhanced strain method for large deformations.
86. Simo JC, Taylor RL, Pister KS (1985) Variational and proEng. Comput. 13: 103–123
jection methods for the volume constraint in finite defor65. Legay A, Combescure A (2003) Elastoplastic stability analmation elasto-plasticity. Comput. Meth. Appl. Mech. Eng.
ysis of shells using the physically stabilized finite element
51: 177–208
SHB8PS. Int. J. Numer. Meth. Eng. 57: 1299–1322
87. Simo JC, Rifai MS (1990) A class of mixed assumed strain
66. Liu WK, Guo Y, Tang S, Belytschko T (1998) A multiplemethods and the method of incompatible modes. Int.
quadrature eight-node hexahedral finite element for large
J. Numer. Meth. Eng. 29: 1595–1638
deformation elasto-plastic analysis. Comput. Meth. Appl.
88. Simo JC, Fox DD, Rifai MS (1990) On a stress resultant
Mech. Eng. 154: 69–132
geometrically exact shell model. Part III: Computational
67. MacNeal RH (1982) Derivation of element stiffness matrices
aspects of the nonlinear theory. Comput. Meth. Appl. Mech.
by assumed strain distribution. Nucl. Eng. Des. 70: 3–12
Eng. 79: 21–70
68. Malkus DS (1976) A finite element displacement model valid
for any value of the compressibility. Int. J. Solids. Struct. 11: 89. Simo JC, Fox DD, Rifai MS (1990) On a stress resultant
geometrically exact shell model. Part IV: Variable thickness
731–738
shells with through-the-thickness stretching. Comput. Meth.
69. Malkus DS, Hughes TJR (1978) Mixed finite-element methAppl. Mech. Eng. 81: 91–126
ods - reduced and selective integration techniques - a unification of concepts. Comput. Meth. Appl. Mech. Eng. 15: 63– 90. Simo JC, Armero F (1992) Geometrically non-linear enhanced strain mixed methods and the method of incom81
patible modes. Int. J. Numer. Meth. Eng. 33: 1413–1449
70. Massin P, al Mikdad M (2002) Nine-node and seven-node
thick shell elements with large displacements and rotations. 91. Simo JC, Kennedy JG (1992) On a stress resultant geometrically exact shell model. Part V: Nonlinear plasticity - forComput. Struct. 80: 835–847
mulation and integration algorithms. Comput. Meth. Appl.
71. Masud A, Tham CL, Liu WK (2000) A stabilized 3D coMech. Eng. 96: 133–171
rotational formulation for geometrically nonlinear analysis
of multi-layered composite shells. Comput. Mech. 26(1): 1–12 92. Simo JC, Armero F, Taylor RL (1993) Improved versions of
assumed enhanced strain tri-linear elements for 3D finite
72. Masud A, Tham CL (2000) Three-dimensional corotational
deformation problems. Comput. Meth. Appl. Mech. Eng.
framework for finite deformation elasto-plastic analysis of
110: 359–386
multilayered composite shells. AIAA Journal 38(9): 1–8
93. Soric´ J, Montag U, Kra¨tzig WB (1997) An efficient formu73. Miehe C (1998) A theoretical and computational model for
lation of integration algorithms for elastoplastic shell analisotropic elastoplastic stress analysis in shells at large
ysis based on layered finite element approach. Comput.
strains. Comput. Meth. Appl. Mech. Eng. 155: 193–233
Meth. Appl. Mech. Eng. 148: 315–328
74. Nagtegaal JC, Parks DM, Rice JR (1974) On numerically
94. Stander N, Matzenmiller A, Ramm E (1989) An assessment
accurate finite element solutions in the fully plastic range.
of assumed strain methods in finite rotation shell analysis.
Comput. Meth. Appl. Mech. Eng. 4: 153–177
Eng. Comput. 6: 58–66
75. Oliver J, On˜ate E (1984) A total lagrangian formulation for
95. Surana KS (1983) Geometrically non-linear formulation for
geometrically nonlinear analysis of structures using finite
the curved shell elements. Int. J. Numer. Meth. Eng. 19: 581–
elements - Part I: Two-dimensional problems: shell and plate
615
structures. Int. J. Numer. Meth. Eng. 20: 2253–2281
76. Parisch H (1995) A continuum based shell element for non- 96. Taylor RL, Pister KS, Herrmann LR (1968) On a variational
theorem for incompressible and nearly incompressible orlinear applications. Int. J. Numer. Meth. Eng. 38: 1855–1883
thotropic elasticity. Int. J. Solids. Struct. 4: 875–883
77. Peng X, Crisfield MA (1992) A consistent corotational for97. Vu-Quoc L, Tan XG (2003) Optimal solid shells for nonmulation for shells using the constant stress/constant molinear analysis of multilayer composites. Part I: Statics.
ment triangle. Int. J. Numer. Meth. Eng. 35: 1829–1847
Comput. Meth. Appl. Mech. Eng. 192: 975–1016
78. Ramm E (1977) A plate/shell element for large deflections
and rotations. In: K.J. Bathe, J.T. Oden, W. Wunderlich (eds) 98. Wagner W, Klinkel S, Gruttmann F (2002) Elastic and
plastic analysis of thin-walled structures using improved
Formulations and Computational Algorithms in Finite
hexahedral elements. Comput. Struct. 80: 857–869
Element Analysis, M.I.T. Press
99. Wriggers P, Eberlein R, Reese S (1996) A comparison of
79. Reese S, Ku¨ssner M, Reddy BD (1999) A new stabilization
three-dimensional continuum and shell elements for finite
technique for finite elements in non-linear elasticity. Int.
plasticity. Int. J. Solids. Struct. 33: 3309–3326
J. Numer. Meth. Eng. 44: 1617–1652
100. Wriggers P, Reese S (1996) A note on enhanced strain
80. Reese S, Wriggers P (2000) A stabilization technique to
methods for large deformations. Comput. Meth. Appl. Mech.
avoid hourglassing in finite elasticity. Int. J. Numer. Meth.
Eng. 135: 201–209
Eng. 48: 79–109
101. Zienkiewicz OC, Taylor RL, Too JM (1979) Reduced inte81. Reese S (2002) On the equivalence of mixed element forgration techniques in general analysis of plates and shells.
mulations and the concept of reduced integration in large
Int. J. Numer. Meth. Eng. 3: 275–290
deformation problems. Int. J. Nonlinear Sci. and Num.
Simul. 3: 1–33