Nonlinear static analysis of an axisymmetric shell storage container in spherical polar coordinates with constraint volume

Nonlinear static analysis of an axisymmetric shell storage container in spherical polar coordinates with constraint volume

Engineering Structures 68 (2014) 111–120 Contents lists available at ScienceDirect Engineering Structures journal homepage: www.elsevier.com/locate/...

2MB Sizes 2 Downloads 98 Views

Engineering Structures 68 (2014) 111–120

Contents lists available at ScienceDirect

Engineering Structures journal homepage: www.elsevier.com/locate/engstruct

Nonlinear static analysis of an axisymmetric shell storage container in spherical polar coordinates with constraint volume Weeraphan Jiammeepreecha a, Somchai Chucheepsakul a,⇑, Tseng Huang b a b

Department of Civil Engineering, Faculty of Engineering, King Mongkut’s University of Technology Thonburi, Bangkok 10140, Thailand Department of Civil Engineering, University of Texas at Arlington, Arlington, TX 76019, USA

a r t i c l e

i n f o

Article history: Received 17 January 2013 Revised 21 February 2014 Accepted 23 February 2014 Available online 27 March 2014 Keywords: Axisymmetric membrane shell Incompressible fluid Half drop shell Hydrostatic pressure Spherical polar coordinate

a b s t r a c t An axisymmetric membrane shell fully filled with an incompressible fluid is investigated and it is modeled as a half drop shell storage container under hydrostatic pressure subjected to the volume constraint conditions of the shell and contained fluid. The shell geometry is simulated using one-dimensional beam elements described in spherical polar coordinates. Energy functional of the shell expressed in the appropriate forms is derived from the variational principle in terms of displacements and the obtained nonlinear equation can be solved by an iterative procedure. Numerical results of the shell displacements with various water depths, thicknesses, and internal pressures are demonstrated. Ó 2014 Elsevier Ltd. All rights reserved.

1. Introduction Axisymmetric shells as structural elements are widely used in many engineering fields, such as structural, mechanical, aeronautical and offshore engineering [1–8]. Examples of shell structures in structural and mechanical engineering are concrete arch domes, and liquid tanks, including the development of pressure vessels and underwater storage containers. Aircrafts, ship hulls, and submarines are examples of the usage of shells in aeronautical and offshore engineering. In the field of biomechanics, axisymmetric shells are commonly used in the cornea and crystalline lens of a human eye model [9–14]. These structures are very effective in resisting external loading [15]. The geometrically nonlinear analysis of thin elastic shells of revolution subjected to arbitrary load has been presented by Ball [16,17]. The analysis was based on Sanders nonlinear thin shell theory and solved by finite difference formulation. Based on this method, the nonlinear terms were treated in the form of pseudoloads. The advantage of this method is the considerable reduction in execution time. However, many researchers have presented the use of finite element formulations for analysis of axisymmetric shells [18–24]. Delpak and Peshkam [20] developed the variational method in order to study the geometrically nonlinear behavior of a ⇑ Corresponding author. Tel.: +66 2 470 9146; fax: +66 2 427 9063. E-mail address: [email protected] (S. Chucheepsakul). http://dx.doi.org/10.1016/j.engstruct.2014.02.014 0141-0296/Ó 2014 Elsevier Ltd. All rights reserved.

rotational shell. The formulation is based on the second variation of the total potential energy equation and implemented using finite element method, as in the research of Peshkam and Delpak [21]. New finite formulations for analysis of shells of revolution have been presented by Delpak [18], Teng and Rotter [19], and Hong and Teng [22]. In addition, Polat and Calayir [23] proposed formulations for investigating shells of revolution based on the total Lagrangian approach, and where the material behavior was assumed to be linearly elastic. The numerical solutions were obtained by using Newmark integration technique coupled with Newton–Raphson iteration procedure. Wu [24] presented a vector form intrinsic finite element (VFIFE) for the dynamic nonlinear analysis of shell structures. This method was based on the theory of vector form analysis. The numerical results are accurate and efficient for solving the shell problems includes geometrical and material nonlinearity. Moreover, Sekhon and Bhatia [25] presented a method for generating stiffness coefficients and fixed edge forces for a spherical shell element. The formulation was based on an approximate analytical solution of the differential equations. The numerical results are accurate and computing time is reduced considerably because of the use of fewer elements. Lang et al. [26] proposed the nonlinear static analysis of shells of revolution using a ring element. The displacement field in the circumferential direction of the ring element is defined by Fourier series. For analyzing a shell is fully filled with an incompressible fluid, Sharma et al. [27] presented the effect of internal fluid height level on

112

W. Jiammeepreecha et al. / Engineering Structures 68 (2014) 111–120

the natural frequency of a vertical clamped-free cylindrical shell. It was found that the natural frequency of the cylindrical shell decreased when increasing the internal fluid height level. In the field of offshore work, the applications of an axisymmetric shell for investigating the behavior of an underwater spherical shell has received much attention from many researchers [5,6,28–31]. Yasuzawa [28] proposed the static and dynamic responses of an underwater half drop shell by the theory of thin shells of revolution and finite element method. It was found that the displacement and membrane stress distribution are uniform along the meridian except at the bottom part. Furthermore, the optimal dome shape of submerged spherical domes has been presented by Vo et al. [29] and Wang et al. [30]. They proposed the membrane analysis and minimum weight design of submerged spherical domes. The numerical solutions were obtained using shooting-optimization technique. It was found that the variation of the shell thickness of spherical domes can be accurately defined by the first-nine terms in the power series for practical applications. Recently, Jiammeepreecha et al. [31] adopted the finite element technique used by Goan [11] in finding the static equilibrium configurations of a deep water half drop shell with constraint volume. In order to solve this problem and obtain the correct solution, the meridian line should be divided into many finite elements in two sub-regions. At the junction of two adjacent sub-regions, the function values, displacements and slopes are made all continuous. The problem can be solved by using the Lagrange multiplier technique; in the process four Lagrange multipliers are used. To alleviate this difficulty, this paper presents a new technique of discretization by using spherical polar coordinates. The advantage of this discretization technique is that the meridian curve is not divided into many regions. The meridian curve is divided into number of sub-regions with equal arc length. Thus, the Lagrange multiplier technique is no longer used in this analysis. The purpose of this paper is to present the nonlinear static analysis of an axisymmetric membrane shell storage container subjected to hydrostatic pressure and constraint volume of the shell. Since the geometry of the shell is always axisymmetric, any meridional curve may be considered as the generating curve. Therefore, the shell is simulated using one-dimensional beam elements and is described in spherical polar coordinates. Third-order shape functions are used in the finite element formulation. In this study, large displacements and rotations are considered. Thus, the strain energy density is expressed as a quadratic function of Lagrangian strains, and the material behavior is assumed to be linearly elastic. By using the strain–displacement relations, the strain energy density is given in terms of displacements and their derivatives. This problem is formulated in a variational form by using shell theory [15], and written in the appropriate forms [32] which are introduced in order to reduce the computation time. The principle of virtual work [33] and finite element method are used to solve the problem, and then the nonlinear equilibrium equations are derived. These equations are solved by an iterative process. The final deformed configuration of the shell is to be determined.

Fig. 1. Three states of the shell.

are assumed to be zero. Finally, this shell is subjected to several external loadings such as linearly hydrostatic pressure, uniform force or imposed displacement. The final deformed configuration is to be determined. This state is referred to as the deformed state or equilibrium state 2 (ES2). 2.1. Shell geometry at the reference state (ES1) A portion of the shell at reference state (ES1) is shown in Fig. 2. ^ be the unit Let (X, Y, Z) be the rectangular coordinates and ð^i; ^j; kÞ vectors along the coordinate axes. A surface may be defined by parametric parameters (h, /); that is X = X(h, /), Y = Y(h, /), and Z = Z(h, /) where (h, /) are the two surface parameters defining the position of a point on the meridian and longitude, respectively. These two surface parameters specify the orthogonal curvilinear

2. Analytical model Consider the three states of a shell, shown in Fig. 1. At an undeformed state, the empty shell without any strains is designated as the initial unstrained state (IUS). When the initially unstrained axisymmetric shell is fully filled with an incompressible fluid and the internal pressure is assumed to be constant, this state is called the reference state or equilibrium state 1 (ES1), in which the geometric configuration is known. The initial strains and displaced configuration at the reference state are small and can be determined by traditional shell analysis [34]. It is noted that the initial unstrained state and the reference state are the same when the initial strains

Fig. 2. Shell reference surface.

113

W. Jiammeepreecha et al. / Engineering Structures 68 (2014) 111–120

coordinate lines on the surface. Also, X = X(h, /), Y = Y(h, /), and Z = Z(h, /) are single valued, continuous, and differentiable functions. Furthermore, there shall be a one-to-one correspondence between pairs of two surface parameter (h, /) values and points on the surface P. Due to symmetry, the shell surface may be generated by rotating a plane curve about the axis of symmetry. For the case of a shell of revolution, any meridional curve may be considered as the generating curve on the r–Z plane, as shown in Fig. 2. Therefore, only the change shape of the generating curve need be considered. Let r be the position vector of any point P. Therefore, the position vector r can be defined by using parallel circle radius r as

^ r ¼ r cos /^i þ r sin /^j þ Z k

ð1Þ

where r = r(h) and Z = Z(h). The total differential line element of r is given by dr ¼ r h dh þ r / d/, where subscripts (h, /) denote partial derivatives along the shell coordinates. The first fundamental form of the reference surface (S) is defined by

dr  dr ¼ E dh2 þ 2F dh d/ þ G d/2

ð2Þ

where E, F, and G are the metric tensor components of the reference surface (S), and are given by

2.2. Displacements and deformed surface When the shell is deformed, its reference surface (S) transfers to a new surface (S). The position vector R of a point on the surface at the deformed state (ES2) is established. It is expressed in terms of displacements measured from the configuration at reference state (ES1) as follows:

r r  ¼ r þ phffiffiffi u þ p/ffiffiffiffi v þ n ^w R ¼ r þ q E G

ð9Þ

where u, v, and w are the displacement components along meridian, longitudinal, and normal directions, respectively. Since an axial pffiffiffiffi shell problem is considered herein, the term ðr / = GÞv is equal to q ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffi pffiffiffiffi zero. Let A ¼ E ¼ r 2h þ Z 2h and B ¼ G ¼ r. Then Rh and R/ can be written as

  e  r h  e ^ Rh ¼ A þ uh  w þ u þ wh n A A A   r/ Bh g R/ ¼ B þ u  w B A B

ð10aÞ

ð10bÞ

The metric tensor components of the deformed surface (S) can be written as

ð3aÞ

 2 e 2  e u þ wh E ¼ Rh  Rh ¼ A þ uh  w þ A A

ð11aÞ

F ¼ r h  r/ ¼ 0

ð3bÞ

F  ¼ Rh  R/ ¼ 0

ð11bÞ

G ¼ r/  r / ¼ r 2

ð3cÞ

G ¼ R/  R/ ¼

E ¼ rh  r h ¼ r2h þ Z 2h

The unit vector normal to the shell surface at point P can be determined by

^ r  r / rZ h cos /^i  rZ h sin /^j þ rr h k ^¼ h n ¼   D jr h  r / j qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi EG  F 2 ¼ r r2h þ Z 2h

2.3. Strain–displacement relations

ð5Þ

eL ¼

 2

ð6Þ

where e, f, and g are the curvature tensor components of the reference surface (S), and given by

r Z  r hh Z h ^ ¼ hqhhffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi e ¼ r hh  n r 2h þ Z 2h

ð7aÞ

1 ðds Þ  ðds0 Þ 2 2 ðds0 Þ

^¼0 f ¼ r h/  n

ð7bÞ

rZ h ^ ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi g ¼ r //  n r 2h þ Z 2h

ð7cÞ

Accordingly, the curvature j of a normal section of the surface can be determined by

e dh2 þ 2f dh d/ þ g d/2 E dh2 þ 2F dh d/ þ G d/2

2

eL ¼ e0

ð12Þ

ds

2 ds0

2

þe

ds

ð13Þ

2

ds0

in which

e0 ¼

 2

2

1 ðdsÞ  ðds0 Þ ; 2 2 ðdsÞ



1 ðds Þ  ðdsÞ 2 2 ðdsÞ

2

ð14a-bÞ

It is apparent that e0 is the initial Eulerian strains component at the reference state (ES1) and e is the added strains component associated with the surface deformation from the reference state (ES1) to the deformed state (ES2). For the case of a symmetrical shell, the shearing strains c0h/ = 0. Hence, the initial Eulerian strains e0 can be expressed in terms of the metric tensor components as follows:

ð8Þ

In the case of the axisymmetric shell, the lines of principal curvature coincide with the coordinate lines. That means F = f = 0. Therefore, the principal curvatures can be expressed as j1 = e/E and j2 = g/G.

2

Since all the displacements are measured from the reference state (ES1), it is necessary to express the strain component in terms of arc length ds at the reference state (ES1) and separate it into two parts, as follows:

2



ð11cÞ

Consider an infinitesimal line element of length ds0 in the initial unstrained state (IUS) and ds in the deformed state (ES2); the total Lagrangian strains can be expressed as follows:

^¼n ^ h dh þ n ^ / d/, the second fundamental form of the Since dn reference surface (S) is determined by

^ ¼ e dh2 þ 2f dh d/ þ gd/2 dr  dn

2 Bh g u w B A



ð4Þ

in which

D ¼ jr h  r / j ¼



e0h ¼

  1 E0 ; 1 2 E

e0/ ¼

  1 G0 1 2 G

ð15a-bÞ

where E0, F0, and G0 are the metric tensor components at the initial unstrained state (IUS). It should be noted that F0 is zero due to the symmetrical shell. Thus, Eq. (5) becomes

114

W. Jiammeepreecha et al. / Engineering Structures 68 (2014) 111–120

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffi D0 ¼ E0 G0 ¼ D ð1  2e0h Þð1  2e0/ Þ

ð16Þ

Similarly, the added strains e are related to the metric tensor components by

eh ¼

  1 E 1 ; 2 E

e/ ¼

  1 G 1 2 G

The first term in the integrand of Eq. (26) is irrelevant, so it is dropped hereafter. Substituting Eq. (21) into Eq. (26), the strain energy can be expressed in the appropriate forms [32]; their first and second derivatives with respect to gi can be obtained by only changing coefficients without re-calculating the matrices.

ð17a-bÞ U¼

Z

h2



c0k g k þ

h1

Finally, the total Lagrangian strains can be expressed in a matrix form

in which

feL g ¼ ½Tðfe0 g þ fegÞ

c0k ¼ C ij ei0 Ljk ;

ð18Þ

in which

" ½T ¼

#

1 12e0h

0

0

1 12e0/

ð19Þ

where [T] is the diagonal material-element matrix. It is formed in order to transform the strain components from the reference state (ES1) to the initial unstrained state (IUS). By substituting Eq. (11) into Eq. (17), the added strains can be expressed in terms of displacements as follows:

eh

 2  2 1 e 1 1 e 1 e 1 uh  2 w þ ¼ uh  2 w þ u þ wh 2 A 2 A 2 A A A A

ð20aÞ

e/

 2 Bh g 1 Bh g ¼ u 2wþ u 2w 2 AB AB B B

ð20bÞ

N i

ei ¼ e þ e ¼ Lik

Lik g k

1 þ Hikl g k g l 2

Hikl

  1 n2kn ¼ C ij Hikl Hjmn þ Hink Hjlm g l g m 2

ð28eÞ

According to Eq. (28), the matrices c1, k, n1, and n2 are symmetric and the variation of strain energy dU can be obtained as follows:

dU ¼



h1

DV ¼

  1 m E0 1  m2 m 1

DV ¼



h2

h1

1 ðfe0 gT þ fegT Þ½Cðfe0 g þ fegÞdh 2

ð22Þ

ð24Þ

ð25Þ

Since [T] and [C0 ] are symmetric metric, [C] is also symmetric and Cij = Cji. Using the index notations, the strain energy of the shell can be written as



Z

h2

h1

  1 1 C ij ei0 ej0 þ C ij ei0 ej þ C ij ei ej dh 2 2

ð29Þ

1 3

Z

h2

h1

Z 2p

Rh  R/  R  r h  r /  r d/ dh

ð30Þ

0

1 3

Z

h2

h1

Z 2p

ðv 1 þ v 2 þ v 3 Þd/ dh

ð31Þ

0

    Be Ag Be ^ Þ  2 ðr  rh Þ u þ  ðr  n ^ Þ  ðr  n ^ Þ þ AB w Bh ðr  n B A A   B ^ ÞÞuh þ  ðr  r h Þ wh ð32aÞ þ ðBðr  n A

    Bh e Be 2 Bh e eg ^Þ þ 2 ðr  r h Þ þ Bh uw u þ  2 ðr  n  3 ðr  r h Þ  A A A A B     Bh Bh  ^   ðr  nÞ uuh þ  2 ðr  r h Þ  B uwh þ A A    g  eg Ag Be ^Þ  ^ Þ þ B wuh w2 þ  ðr  n þ ðr  n  AB B A B g  þ ðr  r h Þ wwh ð32bÞ AB      eg  g  Bh e Bh 2 u wh þ ¼  2 u3 þ u2 w þ  uwwh AB B A A      g  eg  Bh Bh e uuh w þ  uh w2 þ  2 uw2 þ þ ð32cÞ w3 B AB A A

v2 ¼

ð23Þ

in which

½C ¼ 2p½TT ½C 0 ½TtD0

    1 1 dg k c0k þ c1kn þ kkn þ n1kn þ n2kn g n dh 2 3

in which

v1 ¼

where [C0 ] is the shell material property matrix, t is the shell thickness, E0 is the Young’s modulus, and m is the Poisson’s ratio. Substitution of Eq. (18) into Eq. (22) yields

Z

h2

Substituting Eqs. (1), (9) and (10) into Eq. (30) yields

in which

½C 0  ¼

Z

The volume change of the shell from reference state (ES1) to deformed state (ES2) can be expressed in terms of displacements as follows:

The shell is assumed to be a linearly elastic material of constant thickness. Then the strain energy of the shell can be expressed as follows:

Z 2p 1 L T 0 L fe g ½C fe gtD0 d/ dh 2 0

ð28a-cÞ ð28dÞ

ð21Þ

2.4. Strain energy of the shell

h2

kkn ¼ C ij Lik Ljn

2.5. Volume change of the shell

where and are column and symmetric matrices, respectively. These matrices depend on the reference surface (S) characteristics and are identified by Eq. (20). The strains e0, eLi , and eNi are constant, linear, and nonlinear in terms of displacements, respectively.

Z

c1kn ¼ C ij ei0 Hjkn ;

ð27Þ

  n1kn ¼ C ij Lik Hjmn þ Lim Hjnk þ Lin Hjkm g m

h1

Let fggT ¼ bu w uh wh c. Then the added strains can be separated into two parts, a linear and nonlinear, and written in the following index form: L i

   1 1 1 1 1 2 ckn þ kkn þ n1kn þ nkn g k g n dh 2 2 6 12

v3

where v1, v2, and v3 contain the linear, quadratic, and cubic terms of the displacements in the volume change DV, respectively. Using the 1 , a 2 , and a 3 ; the unit index notations, designate Rh , R/ , and R by a ^ can be denoted by ^i1 , ^i2 , and ^i3 , respectively. vectors r h =A, r / =B, and n 1 , a 2 , and a 3 can be written as Then the vectors a

  i ¼ aij^ij ¼ aij þ bijk g k ^ij a

ð33Þ

ð26Þ Consequently, using the permutation symbols (eijk) as follows:

115

W. Jiammeepreecha et al. / Engineering Structures 68 (2014) 111–120

Rh  R/  R ¼ eijk a1i a2j a3k

ð34Þ

Therefore, the volume change can be expressed in the appropriate forms [35] as follows:

Z

DV ¼



h2

v ck g k þ

h1

   1 K 1 v kn þ v Nkn g k g n dh 2 6

ð35Þ

in which

v ck ¼ v

K kn

 2p  1 2 3 eijl ai aj blk þ a1i b2jk a3l þ b1ik a2j a3l 3

ð36aÞ

 2p  1 2 3 ¼ eijl ai bjk bln þ b1ik a2j b3ln þ b1ik b2jn a3l þ a1i b2jn b3lk þ b1in a2j b3lk þ b1in b2jk a3l 3 ð36bÞ

2p  1 2 3 eijl bik bjm bln þ b1im b2jn b3lk þ b1in b2jk b3lm þ b1in b2jm b3lk 3 

v Nkn ¼

þb1ik b2jn b3lm þ b1im b2jk b3ln g m

Z



h2

dg k



1 2

v ck þ v Kkn þ v Nkn

h1

  g n dh

ð36cÞ

Assume that the internal pressure at reference state (ES1) is linearly proportional to the volumetric strain by the relation

ð38Þ

~ is the bulk modulus of the fluid, V is the unstrained fluid where k 0W volume, and DV0W is the volume change of fluid from the unstrained fluid state to the reference state (ES1). Then the strain energy in the enclosed fluid is determined by the relation

1 2

C ¼ k~

  DV 0W þ DV 2 V 0W V 0W

ð39Þ

where DV is the volume change of fluid and shell which are unchangeable. Associated with the variation of volume change, the variation of strain energy due to internal fluid dC can be written as

  ~ DV 0W þ DV dðDVÞ dC ¼ k V 0W

ð40Þ

ð41Þ

Physically, k, which represents the change of pressure from the reference state (ES1) to the deformed state (ES2) is defined as

~ DV k ¼ k V 0W

ð42Þ

According to Eqs. (38) and (42), the constraint equation is considered by the relation

DV þ

V k¼0 ~ k  p0

ð44Þ

where qw is the density of the sea water, g is the specific gravity, and Zw is the vertical distance from the sea water level. The virtual work done by linearly hydrostatic pressure dX can be expressed as follows:

Z

dX ¼

h2

h1

Z 2p 0

pw fdwgD d/ dh ¼ 2p

Z

h2

pw fdwgD dh

ð45Þ

h1

Based on the principle of virtual work, the equilibrium equation of the shell can be obtained by setting the total virtual work of the shell to be zero, as follows:

dp ¼ dU þ dC þ dX ¼ 0

ð46Þ

Substituting Eqs. (29), (37), (41), and (45) into Eq. (46) gives

   1 1 c0k þ c1kn þ kkn þ n1kn þ n2kn g n 2 3 h1     Z h2 1 pw fdwgDdh ¼ 0 ð47Þ ðp0 þ kÞ v ck þ v Kkn þ v Nkn g n dh þ 2p 2 h1

Z



h2

dg k

Two highly nonlinear differential equations in terms of u(h) and w(h) are embedded in the above Euler’s equation. Before the hydrostatic pressure is applied, the shell is in equilibrium at the reference state (ES1). Thus, setting pw, k, and gn to be zero, Eq. (47) is also valid. This requires

Z

h2

dg k ðc0k  p0 v ck Þdh ¼ 0

ð48Þ

h1

This equation should be satisfied everywhere, and it can be used to predict the value of initial strains e0. 3.1. Constraint equation ~ approaches infinity and the Since the fluid is incompressible, k last term in Eq. (43) becomes zero. Thus, k in Eq. (43) may be interpreted as a Lagrange multiplier associated with the constraint volume (DV = 0). Finally, the constraint equation can be written as

Z

h2

h1



v ck g k þ

   1 K 1 v kn þ v Nkn g k g n dh ¼ 0 2 6

ð49Þ

4. Finite element method

Substituting Eq. (38) into Eq. (40) yields

dC ¼ ðp0 þ kÞdðDVÞ

pw ¼ qw gZ w

ð37Þ

2.6. Strain energy due to internal fluid

~ DV 0W p0 ¼ k V 0W

The linearly hydrostatic pressure acting on the normal surface of the shell is given by

3. Equilibrium equation requirement

where v ck , v Kkn , and v Nkn are linear, quadratic, and cubic in terms of displacement gradients g, respectively. However, the values of ai and bij can be derived by Eq. (32). According to Eqs. (36b), (36c), the matrices v Kkn and v Nkn are symmetric and directly obtained from the variation of volume change d(DV) as follows:

dðDVÞ ¼

2.7. Virtual work done by linearly hydrostatic pressure

ð43Þ

where V is the shell volume at reference state (ES1). It is noted that the numerical value of k is an unknown, which will be obtained by solving the entire problem as will be presented.

To solve the problem by using the finite element method, the shell is divided along the h coordinate into many finite ring elements. Consider a general single element with the local coordinate u, the shell global coordinate h and the angle a = h2  h1, as shown in Fig. 3. The local coordinate u is related to the global coordinate h by u = h  h1, and the derivatives of any quantity with respect to u and h are equal. Therefore, using the C1 continuity in finite element method [36], the displacements u(u) and w(u) within each element are approximated by a third-order polynomial of the local coordinate u

uðuÞ ¼ b1 þ b2 u þ b3 u2 þ b4 u3

ð50aÞ

wðuÞ ¼ b5 þ b6 u þ b7 u2 þ b8 u3

ð50bÞ

116

W. Jiammeepreecha et al. / Engineering Structures 68 (2014) 111–120

T

fddg

Z

h2

½wT ðfc0 g  ðp0 þ kÞfv c gÞdh

h1

 1 1 ½wT ½c1  þ ½k þ ½n1  þ ½n2  2 3 h1    1 ðpo þ kÞ ½v K  þ ½v N  ½wdhfdg þ ff g ¼ 0 2

þ

Z

h2

ð55Þ

in which

ff g ¼ 2pfdwgT

Z

h2

 pw fwgD dh

ð56Þ

h1

Since the global degree of freedom {Q} is the same as the local degree of freedom {d}, the global equilibrium equation can be obtained by assembly process using Eq. (55). The results are

fC 0 g  ðp0 þ kÞfVCg    1 1 1 þ ½C 1  þ ½K þ ½N1  þ ½N2   ðp0 þ kÞ ½VK þ ½VN fQ g 2 3 2

Fig. 3. Deep water axisymmetric half drop shell.

where bi ði ¼ 1; 2; . . . ; 8Þ are unknown coefficients. Their first derivatives with respect to u or h are as follows:

uh ðuÞ ¼ b2 þ 2b3 u þ 3b4 u

2

ð51aÞ

wh ðuÞ ¼ b6 þ 2b7 u þ 3b8 u2

ð51bÞ

Consider the eight-unknown coefficient (bi) in Eq. (50). In finite element formulation the displacements u and w are expressed in terms of element nodal degrees of freedom {d} via the cubic polynomial shape function. Therefore, the displacement gradient vector {g} can be expressed as

fgg ¼ ½wfdg

þ fFg ¼ f0g

ð57Þ

Similarly, the constraint equation, Eq. (49), becomes

    1 1 fQ gT fVCg þ ½VK þ ½VN fQg ¼ 0 2 6

ð58Þ

Finally, the equilibrium equation, Eq. (57), and the constraint equation, Eq. (58), are combined into a symmetrical matrix form, as follows:

ð52Þ

in which

fggT ¼ b uðuÞ wðuÞ uh ðuÞ wh ðuÞ c

ð53aÞ

T

fdg ¼ b uð0Þ wð0Þ uh ð0Þ wh ð0Þ uðaÞ wðaÞ uh ðaÞ wh ðaÞ c ð53bÞ 2

N1

6 6 0 ½w ¼ 6 6N 4 1;u 0

0

N2

0

N3

0

N4

N1

0

N2

0

N3

0

0

N2;u

0

N 3;u

0

N4;u

N1;u

0

N2;u

0

N3;u

0

0

3

7 N4 7 7 0 7 5

ð59Þ ð53cÞ

N4;u

u ¼ 0;

where [w] is the cubic polynomial shape function. This function and its derivatives can be expressed as follows:



u2 u3 6 u u2 N1 ¼ 1  3 2 þ 2 3 ; N 1;u ¼  þ 2 a a a a a N2 ¼ u  2



u2 u3 6 u u2 N3 ¼ 3 2  2 3 ; N3;u ¼  a a a a a2 N4 ¼ 

ð54a-bÞ

wh ¼ 0

ð60Þ

The supported condition is considered to be fully fixed at the sea bed. Therefore

u ¼ 0;



u2 u3 u u2 þ 2 ; N2;u ¼ 1  4 þ 3 2 a a a a

Since an axially symmetrical shell is considered, the boundary conditions at the top are

w ¼ 0;

uh ¼ 0;

wh ¼ 0

ð61Þ

The system of nonlinear equations in Eq. (59), which is constrained by both boundary conditions Eqs. (60) and (61), can be solved numerically by an iterative procedure.

ð54c-dÞ 5. Numerical example and results



u2 u3 u u2 þ 2 ; N4;u ¼ 2 þ 3 2 a a a a

ð54e-fÞ

ð54g-hÞ

Substituting Eq. (52) into the matrices c0k , c1kn , kkn , n1kn , n2kn , v ck , v Kkn , and v Nkn in Eq. (47) yields

In order to present the finite element formulation of the membrane shell theory, one has to study the behaviors of the axisymmetric half drop shell storage container installed in deep water, as shown in Fig. 3. A computer program developed by Goan [11] is modified to solve the problem, and the independent variable to h coordinate is used. This independent variable is generally for a spherical shell having a constant Gaussian curvature. In the case

117

W. Jiammeepreecha et al. / Engineering Structures 68 (2014) 111–120

of discretization by using h coordinate, the meridian curve is not divided into many regions. Therefore, the size of the global matrix is reduced in parts of the four Lagrange multipliers when compared with the previous work of Jiammeepreecha et al. [31]. To validate the accuracy of the present solutions, consider a half drop shell, as shown in Fig. 3, submerged at a water depth of H = 1745 mm. The shell geometry and material are: a = 220 mm, t = 2.5 mm, E0 = 757 kgf/mm2, and m = 0.36, and the specific weight of external fluid is cw = 1.0  106 kgf/mm3. The present results show the tangential and normal displacements for a linearly distributed hydrostatic pressure along the sea depth and a constant hydrostatic pressure at the sea bed. As shown in Fig. 4, it can be seen that the results of constant hydrostatic pressure are in close agreement with Yasuzawa’s results [28], except the normal displacement near the support. In this study, the hydrostatic pressure is varied along the sea depth, while there is no information on the hydrostatic pressure in Yasuzawa’s [28] work. However, the numerical results from this study were verified with Roark’s formula for a spherical subjected to uniform external pressure [37], and were found to be conformable. The input parameters employed in this analysis are tabulated in Table 1.

Table 2 Convergence of deflection of half drop shell at the apex. Number of elements

wapex (103 m)

8 12 16 20 24 28

0.019425 0.020044 0.020356 0.020544 0.020669 0.020759

5.1. Half drop shell behavior subjected to hydrostatic pressure Table 2 shows the convergence of the apex displacement for half drop shell with constraint volume. It can be seen that the higher mesh models gives the more accurate result. However, the difference of apex movement is less than 0.50% between the model with 24 and 28 elements. In this paper the model with 24 elements is assumed to be sufficient for accurate results. Based on the results of the first part of the study, the shell at the deformed state (ES2) subjected to hydrostatic pressure is shown in Fig. 5. The present results show very good agreement with the previous work of Jiammeepreecha et al. [31].

Fig. 5. Configuration of the half drop shell at deformed state (ES2).

5.2. Effects of hydrostatic pressure on half drop shell Linearly varying hydrostatic pressure has no effect on the displacement response of the half drop shell, as shown in Figs. 6 and 7. Fig. 8 describes the values of k versus the sea water level. It can be seen that the change of pressure from the reference state (ES1) to the deformed state (ES2) is linearly proportional to the hydrostatic pressure; that is, the value of k increases under large hydrostatic pressure and decreases when the hydrostatic pressure becomes small. 5.3. Effects of radius-to-thickness ratio on half drop shell Using the main data in Table 1 and varying the radius-to-thickness ratio (a/t ratio) from 25 to 200, the tangential and normal displacements of the half drop shell are shown in Figs. 9 and 10,

Fig. 4. Comparison of displacement responses with Yasuzawa’s results [28].

Table 1 Input parameter data. Parameter

Value

Young’s modulus, E0 (N/m2) Poisson’s ratio, m Sea water level, H (m) Radius of shell, a (m) Thickness of shell, t (m) Initial internal pressure, p0 (N/m2) Density of sea water, qw (kg/m3)

2.04  1011 0.30 40 5 0.20 50  103 1025

Fig. 6. Effects of linearly varying hydrostatic pressure on tangential displacement of half drop shell.

118

W. Jiammeepreecha et al. / Engineering Structures 68 (2014) 111–120

Fig. 7. Effects of linearly varying hydrostatic pressure on normal displacement of half drop shell.

Fig. 10. Effects of thickness variation on normal displacement of half drop shell.

Fig. 11. Effects of thickness variation on the change of pressure k.

Fig. 8. Effects of linearly varying hydrostatic pressure on the change of pressure k.

respectively. It can be seen that the radius-to-thickness ratio has a significant effect on the displacements on a deep water half drop shell. On the contrary, changing of the radius-to-thickness ratios may have little effect on the values of k, as shown in Fig. 11. Furthermore, the results show that the point of intersection is on the same location, as shown in Fig. 10. This intersection point location is independent of the radius-to-thickness ratio. 5.4. Effects of initial internal pressure on half drop shell Using the main data in Table 1 and varying the initial internal pressure (p0) from 50  103 to 50  106 N/m2, the tangential and normal displacements of the half drop shell are shown in Figs. 12

Fig. 12. Effects of initial internal pressure on tangential displacement of half drop shell.

Fig. 9. Effects of thickness variation on tangential displacement of half drop shell.

Fig. 13. Effects of initial internal pressure on normal displacement of half drop shell.

119

W. Jiammeepreecha et al. / Engineering Structures 68 (2014) 111–120

ee ¼

ds  ds0 ds0

ðA:2Þ

In the present study, the initial engineering strains ee are assumed to be small, and the small displacement theory is used to calculate the initial strains. Furthermore, by neglecting the quadratic term in Eq. (A.1), eL0  ee . The initial engineering strains ee can be computed by using the membrane theory equilibrium equation of the shell; that is

Nh N/ þ ¼ p0 r1 r2

Fig. 14. Effects of initial internal pressure on the change of pressure k.

and 13, respectively. When a low value of initial internal pressure is applied to the shell, the shell has a large effect versus the high initial internal pressure. Fig. 14 shows the effects of the increase of initial internal pressure on the values of k. It can be seen that the values of k decrease when the initial internal pressure becomes large. However, the values of k are unchangeable under high initial internal pressure.

ðA:3Þ

where p0 is the internal pressure. For the case of a reference surface (S) of a spherical shell having a constant Gaussian curvature, the total tension forces are Nh = N/ = N and the principal curvatures are 1/ r1 = 1/r2 = 1/a. Therefore, the tension force N is given by



1 p a 2 0

ðA:4Þ

The initial engineering strains ee can be determined by

r N p a ee ¼ 0 ð1  mÞ ¼ 0 ð1  mÞ ¼ 0 0 ð1  mÞ Et

E

2E t

ðA:5Þ

Let a0 be the radius of a spherical shell at the initial unstrained state (IUS). The initial Lagrangian strain eL0 becomes

6. Conclusions The nonlinear static responses of a deep-water axisymmetric half drop shell storage container with constraint volume condition by using membrane theory are presented in this paper. The problem is formulated by using the variational principle and the finite element method in terms of displacements, which are expressed in the appropriate forms. The change of pressure from the reference state to the deformed state may be explained as a Lagrange multiplier. In the present study, small displacement theory and traditional shell analysis are used to calculate the initial strains and displaced configuration of the half drop shell, respectively. The numerical results indicate that the Lagrange multiplier represents the parameter for adjusting the internal pressure in order to sustain the shell in equilibrium position under the constraint volume condition. The numerical results show that the effect of radius-tothickness ratio has a major impact on the displacements on the shell, whereas changing the hydrostatic pressure has no effect. However, by varying linearly hydrostatic pressure, the change of pressure from the reference state to the deformed state is linearly proportional to the hydrostatic pressure. For a large value of initial internal pressure, the change of pressure is unchangeable.

eL0 ¼ ee ¼

1 a2  a20 2 a20

Finally, the initial Eulerian strain e0 component at the reference state (ES1) can be determined by 2

e0 ¼

2

1 ðdsÞ  ðds0 Þ 1 a2  a20 ¼ 2 2 2 a2 ðdsÞ

The first and second authors gratefully acknowledge financial support by the Thailand Research Fund (TRF) and King Mongkut’s University of Technology Thonburi (KMUTT) through the Royal Golden Jubilee Ph.D. program (Grant No. PHD/0134/2552).

Referring to the position vector in Eq. (1), the reference surface (S) of a spherical shell having a radius a can be defined by

^ r ¼ a sin h cos /^i þ a sin h sin /^j þ a cos hk

ðB:1Þ

in which

r ¼ a sin h;

r h ¼ a cos h;

Z ¼ a cos h;

Z h ¼ a sin h;

r hh ¼ a sin h Z hh ¼ a cos h

ðB:2a-cÞ ðB:3a-cÞ

The metric tensor components of the reference surface (S) are

F ¼ 0;

G ¼ r2 ¼ a2 sin 2 h

ðB:4a-cÞ

Also,

D¼r

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi r 2h þ Z 2h ¼ a2 sin h

ðB:5Þ

pffiffiffiffi pffiffiffi Let A ¼ E and B ¼ G, then

A ¼ a;

Appendix A. Derivation of initial Eulerian strain

ðA:7Þ

Appendix B. Characteristic quantities of the reference surface

E ¼ r 2h þ Z 2h ¼ a2 ;

Acknowledgements

ðA:6Þ

B ¼ a sin h;

Bh ¼ a cos h

ðB:6a-cÞ

Therefore, the unit vector normal is given by

eL0

The initial engineering strains ee and initial Lagrangian strains are related by [11] 2

eL0 ¼

2

1 ðdsÞ  ðds0 Þ 1 ¼ ee þ e2e 2 2 2 ðds0 Þ

in which

ðA:1Þ

^ rZ h cos /^i  rZ h sin /^j þ rr h k D ^ ¼ sin h cos /^i þ sin h sin /^j þ cos hk

^¼ n

ðB:7Þ

The curvature tensor components of the reference surface (S) are

120

W. Jiammeepreecha et al. / Engineering Structures 68 (2014) 111–120

r h Z hh  r hh Z h e ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ a; r2h þ Z 2h

f ¼ 0;

rZ h 2 g ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ a sin h 2 2 rh þ Zh ðB:8a-cÞ

Since F = f = 0, the coordinate lines (h, /) are also lines of principal curvature. Accordingly, the principal curvatures can be determined by

e E

1 a

g G

j1 ¼ ¼  ; j2 ¼ ¼ 

1 a

ðB:9a-bÞ

^ are given From Eqs. (B.1) and (B.7), the quantities r  r h and r  n by

r  r h ¼ 0;

r  n ^¼a

ðB:10a-bÞ

References [1] Yang T, Kapania R. Shell elements for cooling tower analysis. J Eng Mech ASCE 1983;109:1270–89. [2] Jianping P, Harik IE. Axisymmetric general shells and jointed shells of revolution. J Struct Eng ASCE 1992;118:3186–202. [3] Lay KS. Seismic coupled modeling of axisymmetric tanks containing liquid. J Eng Mech ASCE 1993;119:1747–61. [4] Bucalem ML, Bathe KJ. Finite element analysis of shell structures. Arch Comput Meth Eng 1997;4:3–61. [5] Toyota K, Yasuzawa Y, Kagawa K. Hydroelastic response analysis of a large underwater shell of revolution. In: Proceedings of the 12th international offshore and polar engineering conference, Kitakyushu; 2002. p. 456–63. [6] Huang T. A concept of deep water axisymmetric shell storage container equatorially anchored. In: Proceedings of the 12th international offshore and polar engineering conference, Kitakyushu; 2002. [7] Grigolyuk EI, Lopanitsyn YA. The axisymmetric postbuckling behaviour of shallow spherical domes. J Appl Math Mech 2002;66:605–16. [8] Hou C, Yin Y, Wang C. Axisymmetric nonlinear stability of a shallow conical shell with a spherical cap of arbitrary variable shell thickness. J Eng Mech ASCE 2006;132:1146–9. [9] Huang T, Bisarnsin T, Schachar RA, Black TD. Corneal curvature change due to structural alternation by radial keratotomy. J Biomech Eng ASME 1988;110:249–53. [10] Yeh HL, Huang T, Schachar RA. A closed shell structured eyeball model with application to radial keratotomy. J Biomech Eng ASME 2000;122:504–10. [11] Goan LA. An analysis of an axisymmetrical closed shell subjected to equatorial pull with application to accommodation of the crystalline lens. PhD thesis. University of Texas at Arlington; 2000. [12] Shung WV. An analysis of a crystalline lens subjected to equatorial periodic pulls. PhD thesis. University of Texas at Arlington; 2002. [13] Shung WV, Goan LA, Huang T. Analysis of an axisymmetric shell equatorially anchored with constrained volume. In: Proceedings of 16th ASCE engineering mechanics conference, Seattle; 2003.

[14] Chien CM, Huang T, Schachar RA. Analysis of human crystalline lens accommodation. J Biomech 2006;39:672–80. [15] Langhaar HL. Foundations of practical shell analysis. Urbana-Champaign, IL: University of Illinois; 1964. [16] Ball RE. A geometrically non-linear analysis of arbitrarily loaded shells of revolution. NASA contract report, CR-909. National Aeronautics and Space Administration; 1968. [17] Ball RE. A program for the nonlinear static and dynamic analysis of arbitrarily loaded shells of revolution. Comput Struct 1972;2:141–62. [18] Delpak R. Static analysis of thin rotational shells. Comput Struct 1980;11:305–25. [19] Teng JG, Rotter JM. Elastic–plastic large deflection analysis of axisymmetric shells. Comput Struct 1989;31:211–33. [20] Delpak R, Peshkam V. A variational approach to geometrically non-linear analysis of asymmetrically loaded rotational shells – I. Theory and formulation. Comput Struct 1991;39:317–26. [21] Peshkam V, Delpak R. A variational approach to geometrically non-linear analysis of asymmetrically loaded rotational shells – II. Finite element application. Comput Struct 1993;46:1–11. [22] Hong T, Teng JG. Non-linear analysis of shells of revolution under arbitrary loads. Comput Struct 2002;80:1547–68. [23] Polat C, Calayir Y. Nonlinear static and dynamic analysis of shells of revolution. Mech Res Commun 2010;37:205–9. [24] Wu TY. Dynamic nonlinear analysis of shell structures using a vector form intrinsic finite element. Eng Struct 2013;56:2028–40. [25] Sekhon GS, Bhatia RS. Generation of exact stiffness matrix for a spherical shell element. Comput Struct 2000;74:335–49. [26] Lang C, Meiswinkel R, Filippou FC. Non-linear analysis of shells of revolution with ring elements. Eng Struct 2002;24:163–77. [27] Sharma CB, Darvizeh M, Darvizeh A. Natural frequency response of vertical cantilever composite shells containing fluid. Eng Struct 1998;20:732–7. [28] Yasuzawa Y. Structural response of underwater half drop shaped shell. In: Proceedings of the 3rd international offshore and polar engineering conference, Singapore; 1993. p. 475–81. [29] Vo KK, Wang CM, Chai YH. Membrane analysis and optimization of submerged domes with allowance for selfweight and skin cover load. Arch Appl Mech 2006;75:235–47. [30] Wang CM, Vo KK, Chai YH. Membrane analysis and minimum weight design of submerged spherical domes. J Struct Eng ASCE 2006;132:253–9. [31] Jiammeepreecha W, Chucheepsakul S, Huang T. Nonlinear static analysis of deep water axisymmetric half drop shell storage container with constrained volume. In: Proceedings of the 22nd international offshore and polar engineering conference, Rhodes; 2012. p. 863–71. [32] Rajasekaran S, Murray DW. Incremental finite element matrices. J Struct Div ASCE 1973;99:2423–38. [33] Langhaar HL. Energy methods in applied mechanics. New York: John Wiley & Sons; 1962. [34] Flügge W. Stresses in shells. 2nd ed. Berlin: Springer-Verlag; 1973. [35] Chen JS, Huang T. Appropriate forms in nonlinear analysis. J Eng Mech ASCE 1985;111:1215–26. [36] Cook RD, Malkus DS, Plesha ME, Witt RJ. Concepts and applications of finite element analysis. 4th ed. New York: John Wiley & Sons; 2002. [37] Young WC, Budynas RG. Roark’s formulas for stress and strain. 7th ed. New York: McGraw-Hill; 2002.