A new technique for large deflection analysis of non-prismatic cantilever beams

MECHANICS RESEARCH COMMUNICATIONS

Mechanics Research Communications 32 (2005) 692–703 www.elsevier.com/locate/mechrescom

A new technique for large deflection analysis of non-prismatic cantilever beams Mohammad Dado a, Samir Al-Sadder a

b,*

Department of Mechanical Engineering, Faculty of Engineering, University of Jordan, Amman 11943, Jordan b Department of Civil Engineering, Faculty of Engineering, Hashemite University, Zarqa 13115, Jordan Available online 10 February 2005

Abstract This paper studies the very large deflection behavior of prismatic and non-prismatic cantilever beams subjected to various types of loadings. The formulation is based on representing the angle of rotation of the beam by a polynomial on the position variable along the deflected beam axis. The coefficients of the polynomial are obtained by minimizing the integral of the residual error of the governing differential equation and by applying the beamÕs boundary conditions. Several numerical examples are presented covering prismatic and non-prismatic cantilever beams subjected to uniform, non-uniform distributed loads and tip concentrated loadings in vertical and horizontal directions. The loads considered in this study are restricted to the non-follower type loads. Cases with different loadings and geometries are compared with MSC/NASTRAN computer package. However, for some very large deflection case, the MSC/NASTRAN failed to predict the deflected shape due to divergence problems.
2005 Elsevier Ltd. All rights reserved. Keywords: MSC/NASTRAN; Non-prismatic cantilever; Slender cantilever beams; Very large deflection

1. Introduction Closed form solutions of large deflection problems for cantilever beams with general loading conditions using elliptic integrals are not available. If a non-prismatic beam is considered, the complexity of the problem becomes much greater. In this situation a numerical scheme is the only approach available. Available numerical approaches such as finite elements, finite differences and numerical integration are cumbersome due to cost, difficulties in setting up the proper model and the divergence of solution. This paper presents a general numerical scheme that yields accurate results for prismatic and non-prismatic cantilever beams *

Corresponding author. Tel.: +962 795048993. E-mail address: [email protected] (S. Al-Sadder).

0093-6413/$ - see front matter
2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.mechrescom.2005.01.004

M. Dado, S. Al-Sadder / Mechanics Research Communications 32 (2005) 692–703

693

subjected to various types of loadings. The present scheme is based on the integrated least square error of the non-linear governing differential equation in which the angle of rotation is presented by a polynomial. Most of the research that deals with large deflection problems has used four different approaches. The first approach is based on elliptic integral formulation (e.g. Barten, 1945; Bisshop and Drucker, 1945; Timoshenko and Gere, 1961; Lau, 1974; Mattiasson, 1981; Chucheepsakul et al., 1994; Bona and Zelenika, 1997; Wang et al., 1997; Chucheepsakul et al., 1999; Coffin and Bloom, 1999; Ohtsuki and Ellyin, 2001). This approach is tedious and only suitable for simple loading cases. For example, it cannot solve any non-prismatic or prismatic beam with simple uniformly distributed load in the vertical or horizontal directions. The second approach uses numerical integration with iterative shooting techniques (e.g. Freeman, 1946; Conway, 1947; Holden, 1972; Wang and Watson, 1980; Wang, 1981; Watson and Wang, 1981; Wang and Watson, 1982; Watson and Wang, 1983; Mau, 1990; Wang and Kitipornchai, 1992; Lee and Oh, 2000; Lee, 2001; Magnusson et al., 2001). However, it is suitable only for beams subjected to loads producing moderate deflections and it fails in cases incorporating very large deflection. The third approach utilizes incremental finite element method in connection with Newton–Rhapson iteration techniques for solving elastic problems (e.g. Schmidt, 1977; Golley, 1984; Kooi, 1985; Golley, 1997). This method requires the use of expensive commercial packages and the generation of a very fine mesh requiring huge computational time. In addition, the method may experience divergence problems in very large deflection cases and it requires special types of numerical techniques such as arc-length method. Moreover, an experience user of this program is needed to setup the model and the solution method in the proper form. In the fourth approach, the incremental finite differences method in connection with Newton–Rhapson iteration techniques is used (e.g. Kooi and Kuipers, 1984; Saje and Srpcic, 1985; Srpcic and Saje, 1986). This method requires a very large number of nodes for accurate results and it is prone to divergence in very large deflection cases. This paper studies the very large deflection behavior of prismatic and non-prismatic cantilever beams subjected to various types of loadings. The formulation is based on representing the angle of rotation by a polynomial and substituting it into the governing differential equation and boundary conditions. The governing differential equation is set in the residual form and the integrated residual error over the deflected beam axis is minimized. Representative numerical examples are studied covering prismatic and non-prismatic cantilever beams subjected to uniform, non-uniform distributed loads and tip concentrated loadings in vertical and horizontal directions. Some of the results were compared with MSC/NASTRAN for Windows 95 (1995) computer package. However, for very large deflection cases the comparisons cannot be made due to the divergence of NASTRAN commercial package.

2. Formulation of basic differential equations Consider an inextensible non-prismatic slender cantilever beam of length L and a variable flexural stiffness EI(s) subjected to variable distributed loads in global x, y-axes as shown in Fig. 1. Where s is the curved coordinate along the deflected axis of the beam. qx(s) and qy(s) are the horizontal and vertical distributed loads and Fx, Fy and Me are concentrated tip loads and moment. By definition, these loads are of non-follower type. Consider the right segment of the beam as a free body diagram as shown in Fig. 2 where the length of the segment is equal to (L
s). Horizontal and vertical equilibriums of this segment gives: Z s¼L H ðsÞ ¼ qx ðsÞ ds þ F x ð1Þ s¼s

694

M. Dado, S. Al-Sadder / Mechanics Research Communications 32 (2005) 692–703

y

Fy Me

Fx qy (s)

L

qx (s) s x Fig. 1. Non-prismatic cantilever beam subjected to general loading conditions.

Fy

y

Me Fx qy (s) qx (s) M (s) H (s) V (s)

s

M (s)

V (s)

H (s) x

Fig. 2. Free body diagram of a beam segment.

V ðsÞ ¼

Z

s¼L

qy ðsÞ ds þ F y

ð2Þ

s¼s

Consider an infinitesimal element of length ds of the beam as shown in Fig. 3. The internal forces that act on this element are the horizontal force H, the vertical force V and the bending moment M. x and y are the coordinates of a point on this beam. Angle h represents the angle of rotation of the beam with respect to the positive x-axis as shown in Fig. 3. Moment equilibrium equation of the infinitesimal element yields: dM ds þ V dx
H dy ¼ 0 ds

ð3Þ

dM dy dx ¼H
V ds ds ds

ð4Þ

or

But from geometrical relationship of the infinitesimal element (Fig. 3): dy ¼ sin h; ds

dx ¼ cos h ds

ð5Þ

M. Dado, S. Al-Sadder / Mechanics Research Communications 32 (2005) 692–703 y

V+

dV ds ds M+

ds

qy(s)

695

dM ds ds dH H+ ds ds dy

M θ

qx(s)

H V

dx

x

Fig. 3. Internal and external forces acting on an infinitesimal element (ds) in global coordinate.

Then dM ¼ H sin h
V cos h ds Substituting Eqs. (1) and (2) into Eq. (6), results in Z s¼L  Z s¼L  dM ¼
qy ðsÞ ds þ F y cos h þ qx ðsÞ ds þ F x sin h ds s¼s s¼s

ð6Þ

ð7Þ

From Euler–Bernoulli law, which states that local bending moment is proportional to local curvature; MðsÞ ¼ EIðsÞ

dh ds

ð8Þ

Differentiating Eq. (8) once with respect to s, leads to dM d2 h dIðsÞ dh ¼ EIðsÞ 2 þ E ds ds ds ds Substituting Eq. (9) into Eq. (7), results in Z s¼L  Z s¼L  d2 h dIðsÞ dh þ EIðsÞ 2 þ E qy ðsÞ ds þ F y cos h
qx ðsÞ ds þ F x sin h ¼ 0 ds ds ds s¼s s¼s

ð9Þ

ð10Þ

Eq. (10) is the differential equation governing the large deflection behavior of non-prismatic cantilever beams subjected to various types of loadings.

3. Method of solution The numerical solution scheme for Eq. (10) is based on approximating the angle of rotation (h) along the deflected beam axis by a polynomial on the non-dimensional parameter s
s ¼ ð11Þ L The polynomial approximation is defined as hð
sÞ ¼

i¼N X i¼0

ai
si

ð12Þ

696

M. Dado, S. Al-Sadder / Mechanics Research Communications 32 (2005) 692–703

where ai(i = 0, 1, 2, . . . N) are unknown coefficients to be determined by the solution procedure and N is the order of the polynomial approximating the angle of rotation. Substituting the approximation of h(s) and its derivates into Eq. (10) leads to a residual error e(s) defined as ! Z 1  N N N X X dIð
sÞ X i
2 i
1 i eð
sÞ ¼ EIð
sÞ iði
1Þai
s þ E iai
s þ qy ð
sÞ d
s þ F y cos ai
s d
s i¼0
s i¼0 i¼0 ! Z 1  N X i
qx ð
sÞ d
s þ F x sin ai
s ð13Þ
s

i¼0

where eð
sÞ is the error at position
s along the deflected axis of the beam. The aim of the procedure is to find the optimum values of ai that minimize the integrated error over the deflected beam axis. The integrated or total residual error is defined as TRE ¼

Z

1

ðeð
sÞÞ2 d~ s

ð14Þ

0

¼ 0, The polynomial coefficients are obtained by minimizing the total residual error TRE by setting oTRE oai (i = 0, 1, 2, . . . , N) which in turn leads to the following set of simultaneous equations on ai: Z

1

eð
sÞ 0

oeð
sÞ d
s ¼ 0 oai

i ¼ 0; 1; 2; . . . ; N

ð15Þ

The boundary conditions for the cantilever beam are as follows: (1) Zero slope at the fixed end; hð0Þ ¼ 0

ð16Þ

(2) Zero deflection at the fixed end; yð0Þ ¼ 0

ð17Þ

(3) End moment at the free moment is equal to the applied external moment; Mð1Þ ¼

 EIð1Þ dh ¼ Me L d
s 
s¼1

ð18Þ

Since the internal shear equation is used in expressing the governing differential equation, the shear boundary conditions should not be included here. Writing these conditions in terms of the approximated angle of rotation, leads to the following relations; ao ¼ 0 N X i¼0

iai
si
1 ¼

ð19Þ M eL EIð1Þ

ð20Þ

where Eqs. (19) and (20) correspond to Eqs. (16) and (18), respectively. Eq. (17) is clearly satisfied by observing Eq. (5).

M. Dado, S. Al-Sadder / Mechanics Research Communications 32 (2005) 692–703

697

The different load components are define in a non-dimensional form as,
qx ð
sÞ ¼

qx ð
sÞL3 EIð0Þ


qy ð
sÞ ¼

qy ð
sÞL3 EIð0Þ

Fx ¼

F x L2 EIð0Þ

ð21Þ

F y L2 EIð0Þ M eL Me ¼ EIð0Þ Fy ¼

Eqs. (15), (19) and (20) are solved simultaneously using Newton–Raphson iterative method for the unknown parameters ai (i = 0, 1, 2, . . . N). The form of these equation allows the exact formulation for the residual vector and the Jacobian matrix so that numerical differentiation is not used. However, a numerical integration scheme using 16 points Gaussian Quadrature is utilized in performing the required integrals.

4. Large displacement finite element analysis using MSC/NASTRAN As a check for the results of the present study, a large displacement finite element analysis using the multi-purpose computer program MSC/NASTRAN is performed by dividing the cantilever beam into sufficient number of two-dimensional beam elements depending on beam geometry and loading conditions. The following fundamental parameters are introduced to the program: (1) (2) (3) (4) (5) (6)

Number of loading increments = 50, Maximum number of iterations during each increment = 25, Number of iterations before updating tangent stiffness matrix = 15, Method of updating tangent stiffness matrix = iterative, Tolerance for displacement convergence criterion = 1*10
4, Solution strategy = conventional Newton–Raphson iteration method.

5. Numerical examples and discussions Several representative examples are presented in this section. These examples demonstrate the robustness of the present technique in dealing with very large deflection problems. As stated earlier, these examples were also solved using MSC/NASTRAN for comparison. Unfortunately, MSC/NASTRAN had failed in predicting cases with very large deflections. The running time for MSC/NASTRAN examples was more than one hour on the average while the present technique required only few seconds using the same computing facility, namely, P4-2.4 GHz-PC. Example 1. In this example, a prismatic slender cantilever beam is subjected to uniformly distributed load in the vertical direction as shown in Fig. 4. Several intensity values of the load are considered. Fig. 5 shows the deflection shapes for the loads
qy ð
sÞ ¼
4,
10,
20,
40,
100. This example was solved using MSC/

698

M. Dado, S. Al-Sadder / Mechanics Research Communications 32 (2005) 692–703

y

qy(s) x

constant EI

s

L

Fig. 4. Prismatic cantilever beam subjected to uniformly distributed load in y-direction.

0.1 0 Prsent Study MSC/NASTRAN

-0.1 -0.2

qy = −4

y/L

-0.3 -0.4 -0.5 -0.6 -0.7

-10

-0.8 -20

-0.9

-40 -100

-1 -0.1

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

x/L

Fig. 5. Deflection shapes for Example 1.

NASTRAN with 100 two-dimensional beam elements. Both results up to q
y ð
sÞ ¼
20 are in excellent agreement. However, for
qy ð
sÞ ¼
40 and
100, the MSC/NASTRAN failed in predicting results due to divergence problem. The results shown in Fig. 5 emphasize the stability of the present technique especially for the case of
qy ð
sÞ ¼
100 where the beam is almost aligned with the load for a large portion of the beam length. The order of the polynomial for the case of
qy ð
sÞ ¼
100 is (N = 15) and the total residual error TRE as defined in Eq. (14) is equal to (2 · 10
4). Example 2. This example is similar to Example (1) except that the beam is non-prismatic but tapered as shown in Fig. 6. Fig. 7 shows the deflection shapes for the loads
qy ð
sÞ ¼
4;
6;
8;
10;
12;
14;
16;
20;
30;
40;
100. This example was solved using MSC/NASTRAN with 200 two-dimensional beam elements. Both results up to
qy ð
sÞ ¼
12 are in excellent agreement. However, for
qy ð
sÞ ¼
14 to
100, the MSC/NASTRAN failed in predicting results due to divergence problem. Even with the presence of tapered cross-section the results as shown in Fig. 7 are stable and reach the physical limit of the deflection shape. The order of the polynomial for the case of
qy ð
sÞ ¼
100 is (N = 15) and the total residual error is (1.5 · 10
4).

M. Dado, S. Al-Sadder / Mechanics Research Communications 32 (2005) 692–703

699

y qy(s) A

x Variable EI

A

s

L

h

h1

h2 b Section A-A

Fig. 6. Tapered cantilever beam subjected to uniformly distributed load in the y-direction.

0 Present Study MSC/NASTRAN

-0.1 -0.2 -0.3

y/L

-0.4 qy = −4

-0.5 -0.6 -6

-0.7

-8

-0.8 -20

-0.9

-40

-14 -16

-12

-10

-30

-100

-1 -0.1

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

x/L

Fig. 7. Deflection shapes for Example 2.

Example 3. The loading case for this example is a linearly distributed load q
y ð
sÞ in the vertical direction and a uniformly distributed load
qx ð
sÞ in the horizontal direction and the cantilever beam is a prismatic as shown in Fig. 8. The maximum intensity of the linearly distributed vertical load
qy ð
sÞ was held at a constant value of
5, while the uniformly distributed horizontal load
qx ð
sÞ was given the values
1,
5,
10,
15,
20.
30,
50,
80,
100,
120. Fig. 9 shows the deflection shapes for the different loading combinations. This example was solved using MSC/NASTRAN with 100 two-dimensional beam elements. Both results up to
qx ð
sÞ ¼
20 are in excellent agreement. However, for
qx ð
sÞ ¼
30 to
120, the MSC/NASTRAN failed in predicting results due to divergence problem. For the case of
qx ð
sÞ ¼
120, the order of the polynomial is (N = 15) and the total residual error equals to (1.4 · 10
3).

700

M. Dado, S. Al-Sadder / Mechanics Research Communications 32 (2005) 692–703

qy(s)

y

x s Constant EI

qx(s)

L

Fig. 8. Prismatic cantilever beam subjected to linearly distributed load in the y-direction and uniformly distributed load in the x-direction.

0

qy = −5

-0.4

y/L

qx = −1

-120 -100 -80

-0.2

-50

-0.6

-5

-30 -20

-0.8

Present Study MSC/NASTRAN

-15

-10

-1 -1

-0.8

-0.6

-0.4

-0.2

0

0.2

0.4

0.6

0.8

1

x/L

Fig. 9. Deflection shapes for Example 3.

Example 4. A tip concentrated loads in the vertical and the horizontal directions: F y , F x along with a tip concentrated moment M e were applied on a prismatic cantilever beam as shown in Fig. 10. Several different combination of these loads were examined. These loads were selected in such away to produce deflection shape with an inflection point. Fig. 11 shows these load combinations and the resulting deflection shapes. As in the previous examples the MSC/NASTRAN was limited in its deflection shapes prediction. For the case of F x ¼
7, F y ¼ 14 and M e ¼
7, the order of the polynomial is (N = 12) and the total residual error equals to (6.6 · 10
5). Example 5. A tapered cantilever beam with several taper ratios is considered. The load consist of a uniformly distributed vertical load
qy ð
sÞ and a tip concentrated moment M e as shown in Fig. 12. The loads were held constant while the taper ratio was given the values 0.2, 0.3, 0.4,. 0.6, 0.8 and 1. Fig. 13 illustrates the large deflection beam behavior for these different taper ratios. This example was solved using MSC/ NASTRAN with 200 two-dimensional beam elements. Both results were in excellent agreement for taper ratios from 0.4 to 1. However, for taper ratios 0.2 and 0.3, the MSC/NASTRAN failed in predicting results

y Fy Fx

s

Constant EI

L

x

Me

Fig. 10. Prismatic cantilever beam subjected to tip concentrated forces and tip moment.

M. Dado, S. Al-Sadder / Mechanics Research Communications 32 (2005) 692–703

701

1 0.9

Present Study MSC/NASTRAN

0.8 0.7 0.6 0.5

y/L

-7 14 -7

0.4

-4 8 -4

-6 12 -6

-2 4 -2

0.3 0.2

Fx = – 1

0.1

Fy = 2 0

Me = – 1

-0.1 -0.1

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

x/L

Fig. 11. Deflection shapes for Example 4.

y qy(s) A

s

Me

x

A Variable EI

L

h

h1

h2 b Section A-A

Fig. 12. Tapered cantilever beam subjected to uniformly distributed load in the y-direction and to tip concentrated moment.

due to divergence problem. It is observed that for a taper ratio of 0.3 the beam made half a loop in the tip region and two loops for the case of taper ratio 0.2. This is due to low flexural stiffness at the tip. For the case of taper ratio = 0.2, the order of the polynomial is (N = 12) and the total residual error equals to (8 · 10
6). The accurate prediction of this complicated deflection behavior demonstrates the strength of the present method.

702

M. Dado, S. Al-Sadder / Mechanics Research Communications 32 (2005) 692–703 0.1 0 Present Study MSC/NASTRAN

-0.1 -0.2 -0.3

y/L

-0.4 -0.5

h2 /h1 =0.2 -0.6

0.3 q = –50

-0.7

0.4

y

1.0 0.6

Me = 1

-0.8

0.8

-0.9 -1 -0.1

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

x/L

Fig. 13. Deflection shapes for Example 5.

6. Conclusion A robust and stable numerical scheme for the solution of very large deflection problem for prismatic and non-prismatic slender cantilever beams is presented in this paper. The method is based on approximating the angle of rotation of the deflected beam axis by an Nth order polynomial. The governing differential equation along with the beams boundary conditions are used to evaluate the polynomial coefficients. Several examples ranging from moderately large deflection to very large deflection were considered. It was shown that the presented scheme handles extreme cases with high accuracy and solution stability while the finite element method as depicted by MSC/NASTRAN failed in predicting these extreme cases. A recommended extension for this work is the application of this scheme on beams with different boundary conditions such as: pinned-pinned, fixed-pinned, and fixed-fixed. In addition, continuous beams may also be considered using the same scheme.

References Barten, H.J., 1945. On the deflection of a cantilever beam. Quart. J. Appl. Math. 3, 275. Bisshop, K.E., Drucker, D.C., 1945. Large deflections of cantilever beams. Quart. J. Appl. Math. 3, 272–275. Bona, F., Zelenika, S., 1997. A generalized elastica-type approach to the analysis of large displacements of spring-strips. Proc. Instn. Mech. Engrs. Part C 21, 509–517. Chucheepsakul, S., Buncharoen, S., Wang, C.M., 1994. Large deflection of beams under moment gradient. ASCE J. Eng. Mech. 120, 1848. Chucheepsakul, S., Wang, C.M., He, X.Q., Monprapussorn, T., 1999. Double curvature bending of variable-arc-length elastica. J. Appl. Mech. 66, 87–94. Coffin, D.W., Bloom, F., 1999. Elastica solution for the hygrothermal buckling of a beam. Int. J. Non-lin. Mech. 34, 935. Conway, H.D., 1947. Large deflection of simply supported beams. Philos. Mag. 38, 905. Freeman, J.G., 1946. Mathematical theory of deflection of beam. Philos. Mag. 37, 551. Golley, B.W., 1984. The finite element solution of a class of elastica problems. Comput. Meth. Appl. Mech. Eng. 46, 159–168. Golley, B.W., 1997. The solution of open and closed elasticas using intrinsic coordinate finite elements. Comput. Meth. Appl. Mech. Eng. 146, 127–134.

M. Dado, S. Al-Sadder / Mechanics Research Communications 32 (2005) 692–703

703

Holden, J.T., 1972. On the finite deflections of thin beams. Int. J. Solids Struct. 8, 1051–1055. Kooi, B.W., Kuipers, M., 1984. A unilateral contact problem with the heavy elastica. Int. J. Non-lin. Mech. 19, 309–321. Kooi, B.W., 1985. A unilateral contact problem with the heavy elastica solve by use of finite element. Int. J. Non-lin. Mech. 21, 95–103. Lau, J.H., 1974. Large deflections of beams with combined loads. ASCE J. Eng. Mech. Div. 12, 140. Lee, B.K., Oh, S.J., 2000. Elastica and buckling load of simple tapered columns with constant volume. Int. J. Solid Struct. 37, 2507– 2518. Lee, K., 2001. Post-buckling of uniform cantilever column under a combined load. Int. J. Mech. Non-lin. Mech. 36, 813–816. Magnusson, A., Ristinmaa, M., Ljung, C., 2001. Behavior of the extensible elastica column. Int. J. Solid Struct. 38, 8441–8457. Mattiasson, K., 1981. Numerical results from large deflection beam and frame problems analysis by means of elliptic integrals. Int. J. Num. Meth. Eng. 16, 145. Mau, S.T., 1990. Elastica solution of braced struts. J. Eng. Mech. 116, 688. MSC/NASTRAN for Windows 95, 1995. McNeal-Schwendler Corporation, CA, USA. Ohtsuki, A., Ellyin, F., 2001. Analytical approach to large deformation problems of frame structures (in case of a square frame with rigid joints). JSME Int. J. 44, 89–93. Saje, M., Srpcic, S., 1985. Large deformations of in-plane beams. Int. J. Solids Struct. 21, 1181. Schmidt, W.F., 1977. Nonlinear bending of beams using the finite element method. Int. J. Comput. Struct. 8, 153. Srpcic, S., Saje, M., 1986. Large deformations of thin curved plane beam of constant initial curvature. Int. J. Mech. Sci. 28, 275. Timoshenko, S.P., Gere, J.M., 1961. Theory of Elastic Stability. McGraw-Hill, NY. Wang, C.M., Kitipornchai, S., 1992. Shooting-optimization technique for large deflection analysis of structural members. Eng. Struct. 14, 231–240. Wang, C.Y., Watson, L.T., 1980. On large deformations of C-shaped springs. Int. J. Mech. Sci. 22, 395–400. Wang, C.Y., 1981. Large deformations of a heavy cantilever. Quart. J. Appl. Math. 39, 261–273. Wang, C.Y., Watson, L.T., 1982. The elastic catenary. Int. J. Mech. Sci. 24, 349–357. Wang, C.M., Lam, K.Y., He, X.Q., Chucheepsakul, S., 1997. Large deformation of an end supported beam subjected to a point load. Int. J. Mech. Non-lin. Mech. 32, 63–72. Watson, L.T., Wang, C.Y., 1981. Hanging an elastic ring. Int. J. Mech. Sci. 23, 161–167. Watson, L.T., Wang, C.Y., 1983. Periodically supported heavy elastica sheet. ASCE J. Eng. Mech. Div. 109, 811.