Large deflection analysis of curved beam problem with varying curvature and moving boundaries

Engineering Science and Technology, an International Journal 21 (2018) 408–420

Contents lists available at ScienceDirect

Engineering Science and Technology, an International Journal journal homepage: www.elsevier.com/locate/jestch

Full Length Article

Large deflection analysis of curved beam problem with varying curvature and moving boundaries Sushanta Ghuku, Kashi Nath Saha ⇑ Department of Mechanical Engineering, Jadavpur University, Kolkata 700032, India

a r t i c l e

i n f o

Article history: Received 4 December 2017 Revised 21 March 2018 Accepted 11 April 2018 Available online 26 April 2018 Keywords: Initially curved beam Large deflection Geometric nonlinearity Combined bending-stretching-shear Moving boundary Updated Lagrangian approach

a b s t r a c t The paper presents experimental and theoretical large deflection analysis of non-uniformly curved beam with moving boundaries under static loading within elastic domain. A master leaf spring is considered as physical model of the curved beam problem and its load–deflection behaviour is studied experimentally in a specially designed testing rig. Beside direct deflection measurement at some discrete points within the specimen domain, image processing technique is also used to obtain complete deflection profiles under loaded conditions. The indirect deflection measurement through post processing photographs of loaded master leaf is implemented manually in AutoCADÒ. Deflection behaviour of the physical system involves strong geometric nonlinearity coming from non-uniform initial curvature, moving boundaries, nonlinear kinematics due to coupling between bending, stretching, shear deformation and large deflection, asymmetry in beam geometry and eccentricity in load application point with respect to geometric centre. All of these complicating effects are considered in the mathematical model of the physical system. As large deflection involves a large rigid body motion and the induced strain coming from deformation displacement is rather small, present analysis is carried out within elastic regime where material constitutive relation remains linear. Hence system governing equation is derived within the framework of geometric nonlinearity and small strain assumption, using energy principle based variational method. The nonlinear governing equation, in association with complicated moving boundary conditions, is solved iteratively through incremental loading using an updated Lagrangian approach. After each incremental load step, kinetic relation is also satisfied through shear force balance. Numerical results are generated for the same loading conditions of the experimental work and comparisons between theoretical and experimental results are quite good. However, the comparison study leads to identification of several geometric parameters of the physical system, incorporation of which may provide more realistic simulation. Ó 2018 Karabuk University. Publishing services by Elsevier B.V. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).

1. Introduction Various structures and machine elements in civil, aerospace and mechanical engineering disciplines are generally modeled as beam. Most of such practical members are initially curved and show nonlinearity in their deformation behaviour. Hence precise design of such structural elements calls for nonlinear analysis. Nonlinearities in beam bending problem are generally manifested through nonlinear kinematic and material constitutive models and they are known as geometric and material nonlinearity. Large deformation of flexible members induces a large rigid body motion and small strain. Hence linear material modeling is generally used for large ⇑ Corresponding author. E-mail address: [email protected] (K.N. Saha).

deflection analysis of such slender structures within elastic limit. The present large deflection problem only focuses on nonlinearities associated with nonlinear kinematic and kinetic relations, and material nonlinearity is out of scope of the present paper. Thus relevant research papers regarding geometric nonlinear analysis of beam and equivalent structures like leaf spring, arch, compliant mechanisms, etc., are critically reviewed and several observations are presented in the following paragraphs. Trivial solution of beam bending problem generally linearizes Euler-Bernoulli moment-curvature relation assuming small deflection and hence cannot be used for beam undergoing large deflection. This produces the basic difference between small and large deflection analyses. One of the classical approaches, widely used to large deflection beam bending problem, is elliptic integral approach which provides solutions for nonlinear displacement field in terms of elliptic integrals [1–3]. This analytical approach

Peer review under responsibility of Karabuk University. https://doi.org/10.1016/j.jestch.2018.04.007 2215-0986/Ó 2018 Karabuk University. Publishing services by Elsevier B.V. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).

S. Ghuku, K.N. Saha / Engineering Science and Technology, an International Journal 21 (2018) 408–420

becomes inconsistent when generalized loading condition is considered, even for an initially straight uniform beam [4,5]. In addition, spatial variation of beam stiffness along length enforces difficulty. In such cases, classical mechanics based system governing equation is solved through series approximation for unknown field [6–9] or using some iterative shooting process [3,10–12], depending on the strength of nonlinearity present in governing equation. For inextensible beam, iterative process is generally converged with the constant beam length. Constancy in beam length remains effective for deformation under follower loading [13,14]. However, conservative loading results some centre line stretching and in addition, presence of multiple solutions makes iterative approaches unsuitable to undertake more complicated problems [15]. Considerable stretching effect is generally encountered for beams with high slenderness ratio, whereas, shear deformation becomes significant for stub beams. Shear deformable beam is generally analyzed in the frame work of Timoshenko beam model and higher order shear deformation theories incorporating warping of beam cross-section [16,17]. Generally, combined effect of stretching and shear deformation is not considered in large deflection analysis. However, actual deformation process of a flexible structural member is vastly complicated as it includes combined effects of bending, stretching, shear deformation, torsion, warping of beam cross-section, and many such other complicated phenomena. In addition, most of the natural and manmade structures are arbitrarily curved in space which increases complexity of a problem [18]. Rigorous analysis of such three dimensional deformation characteristics of generalized naturally curved beam like structures tends towards geometrically exact beam theory, popularly known as Simo-Reissner beam theory [19]. An extensive review on development of geometric nonlinear beam theory is well documented in a recently published review paper [20]. Complicating effects, generally encountered in large deflection beam problem, arise from several problem parameters like loading condition, boundary condition, initial geometry, etc. Complications associated with loading pattern generally include change in intensity of distributed load in the course of large deformation and presence of singularity point within problem domain [4,5,11,21]. Nonlinear effects, predominantly coming from boundary condition, include development of friction force at simple supports [2] and elastic restraints of boundaries [3,22]. The friction force generated at supports and in addition, initial curvature of beam geometry and coupled transverse-in-plane displacement field add more complexity in nonlinear system response [6,11,23,24]. When such complicating effects are considered, even to some extent, classical mechanics based approach becomes inappropriate and analysis is carried out in the framework of variational mechanics [5,21,25–30]. In this approach, governing equation is derived through minimization of error arising from force [5,25,26,30] or energy [21,27–30] balance. Solution of variational equation, considering whole problem domain leads to semi analytical method [5,21], whereas solution through domain decomposition leads to a pure numerical method, e.g., finite element method [23,27–29]. Due to strong dependency of system response on deformation state, change in geometry greatly influences deformation characteristics of beam with complicated geometry, undergoing very large rotation and translation. Geometry updation is sometimes implemented considering initial configuration as reference, which leads to total Lagrangian approach [31–33]. Whereas, solution with last calculated geometry during incremental loading as reference provides more realistic prediction of system behaviour and known as updated Lagrangian approach [33,34]. In spite of numerous theoretical research works on large deflection of beam, experimental works are rarely reported. Most of the reported experimental works are performed with several beam [31,35,36] and equivalent structures like leaf spring [11,37,38],

409

arch [39], etc., to validate theoretical models. In most of the experimental work, deflection is only measured at some specified points within physical domain, generally at load application point, using several precise instruments [31,35] and displacement sensors [9,21,37,40]. Whereas measurement of complete deflection profile is also reported in a research paper [11], where deflection measurement involves image processing technique. Large deflection analysis of slender beam having initial straight profile, with classical boundary condition has been reported in a large number. Analysis of slender beams with initial uniform curvature is also reported in a moderate number, by modeling most of them as circular segment. Whereas geometric nonlinear analysis of extensible and shear deformable beam undergoing large deflection with generalized non-uniform initial curvature, asymmetric geometry, eccentric loading and moving boundaries is rare. Hence large deflection behaviour of a master leaf spring is simulated experimentally and theoretically as it involves all the mentioned complicating effects. Experiment is performed in a specially designed three point bending set-up. Several physical parameters of the system are identified and incorporated in the mathematical model and analysis is carried out through energy principle based geometry updation technique. 2. Experiment As mentioned earlier, the present work simulates loaddeflection behaviour of asymmetric non-uniformly curved beam with moving boundaries under eccentric static loading. A master leaf spring is considered as specimen of the present experimental work whose geometry under no-load condition is shown in Fig. 1 along with some major dimensions. Complete profile of the specimen is obtained in a Cartesian frame (X 0e ; Y 0e ) by measuring coordinates of twenty-one points, equally spaced along X 0e axis (marked as 2 to 22 in Fig. 1). The two end points of measurement domain are marked as points 1 and 23. These end points of curved beam domain (A1 ; B1 ) are defined as points of tangential intersections of master leaf centre line with pitch circles of the eye ends. Midpoint of the straight line joining the eye centres (span) is considered as origin of Cartesian coordinate system (X 0e ; Y 0e ) in the present geometry measurement and complete profile of the master leaf is presented numerically in Table 1. It is obvious from Fig. 1 and the tabulated values in Table 1 that no-load profile of the master leaf is asymmetric and have non-uniform curvature throughout the domain. In addition, centre line of hole present in the specimen (refer Supplementary Fig. S1 for enlarged view of the hole) has an eccentricity of 1.6 mm with respect to the origin of (X 0e ; Y 0e ), which will produce eccentric loading as described later on in the following paragraphs. To obtain deflection characteristics of the above specified master leaf spring under static load, experiment is carried out on a specially designed experimental set-up. Photograph and schematic diagram of the set-up are shown in Fig. 2(a) and (b) respectively. Main components of the set-up fall into two categories and they are support structures (items 3, 5, 6, 10, 11) and load imparting components. Load imparting components mainly consist of load connector (item 7), vertical guide rod (item 1) and bush (item 2). Detailed drawings of the components are not presented separately to maintain compactness of the paper. Eye ends of the master leaf spring are assembled with roller support sub-assemblies (item 9), and this roller ended master leaf spring is placed on the C.I. bed (item 10) of the test rig underneath the load connector. However before placing the master leaf in the test rig, twenty-one equally spaced points are marked on the master leaf spring along its centre line, using prick punch (refer Supplementary Fig. S1 for clear understanding of the markings). These marked points on the

410

S. Ghuku, K.N. Saha / Engineering Science and Technology, an International Journal 21 (2018) 408–420

Fig. 1. Profile of the master leaf spring in its free state (all dimensions in mm).

Table 1 Numerical representation of the master leaf profile in its free state: (X 0e ; Y 0e ) coordinates of the 23 marked points in mm, with respect to origin at mid-span. Loc.

1

2

3

4

5

6

7

8

9

10

11

12

X 0e


439.3


400


360


320


280


240


200


160


120


80


40

0

Y 0e

11.4

33

52

68.5

84

96

107

114

120.5

125

128

129

Loc.

13

14

15

16

17

18

19

20

21

22

23

X 0e Y 0e

40

80

120

160

200

240

280

320

360

400

439.1

128

125

121

115

106.5

97

84

68.5

52

32.5

11.5

physical specimen are used as reference in geometry measurement under loaded condition. Projected part of the load connector is inserted into the hole of the master leaf and slotted discs (item 4) are placed on the load connector to impart static load on the leaf spring. A vertical guide rod is fitted to the load connector through the slots of loading discs and moves smoothly in vertical direction through a bush (item 2) fitted on guiding disc (item 3). Thus, the guide rod allows only vertical movement of the loading point and restricts horizontal movement of load application line. Hence, centre line of the applied static load always acts along the hole centre and produces eccentric loading on the master leaf spring with respect to mid-span vertical line. Due to the presence of ball bearings (item 13) in roller supports (refer Fig. 2(b)), application of static load causes frictionless horizontal movement of each half of the master leaf about vertical load line which experimentally simulates beam model with moving boundaries. During experimentation, static load is gradually increased through seven steps up to limit load which is calculated using analytical expression of Winkler-Bach beam theory. 75% of yield stress value of leaf spring material is considered as allowable load in this limit load calculation. First load step is coming from weight of vertical guide rod and load connector without application of any external dead weight and hence this step is termed as pre-load (8.5 N). The other six load steps are 84.1, 159.6, 235.2, 310.7, 385.2 and 470.6 (in N). Experimental measurement of deflection profiles of the master leaf under loaded conditions is not straight forward as in case of no-load profile and presented in the following sub-section. 2.1. Deflection profile Deflection profiles of the master leaf spring under loaded conditions are determined by using two different techniques – direct

measurement technique and image processing technique. Basic instrumentation required for direct geometry measurement under loaded condition includes steel rule (item C), plumb (item A) and a vernier height gauge (item B). A fixed plumb pointer (item 12) is attached to the pin of each roller supports by means of spacer and threaded fastening (refer Fig. 2(b)). To measure horizontal movements of these fixed plumb pointers, two steel rules are placed in grooves of the cast iron bed. These two steel rules form a Cartesian coordinate frame (X e ; Y e ), whose Y e axis is along vertical load line and X e axis is the virtual straight line joining two roller centres at mid-plane. At a particular load step, positions of the fixed plumb pointers in (X e ; Y e ) frame give span of the master leaf, whereas camber is measured using the vernier height gauge. Measurement of abscissa and ordinate of the marked points are made using a general purpose plumb and height gauge (Supplementary Table S1), which gives deflection profile of the master leaf under loaded condition. Deflection profiles of the master leaf coming from direct measurement under each of the seven load steps are shown in Fig. 3 by dots. Direct profile measurement for some intermediate points (1, 6–8, 15, 16 and 21) could not be made due to the presence of box structure. The figure also includes the corresponding deflection profiles obtained from photographs and the post processing of indirect measurement technique is discussed briefly in the following paragraph. During experimentation, photograph of the master leaf spring is captured at each load step, using a digital camera and stored (Supplementary Fig. S2). Each photograph is taken as background of AutoCADÒ, and the centre line of loaded master leaf spring is drawn by joining the prick punched marking points as a spline fitted polyline. As in case of direct measurement of deflection profile, origin of the Cartesian coordinate system (X e ; Y e ) is considered at intersection point of vertical load line and the horizontal line between roller centres. Now the AutoCADÒ drawing is scaled

S. Ghuku, K.N. Saha / Engineering Science and Technology, an International Journal 21 (2018) 408–420

Fig. 2. (a) Photograph and (b) schematic diagram of the experimental set-up.

Fig. 3. Deflection profile of the master leaf spring at the seven load steps.

411

412

S. Ghuku, K.N. Saha / Engineering Science and Technology, an International Journal 21 (2018) 408–420

considering origin of (X e ; Y e ) as base point in such a way that span and camber of loaded leaf spring in drawing match with actually measured span and camber at that particular load step. AutoCADÒ drawings of master leaf configurations corresponding to each of the seven load steps are presented in Supplementary Fig. S3. Coordinates of the centre line (X e ; Y e ) corresponding to twenty-three equidistant points along X e axis are exported (Supplementary Table S2) and plotted as shown in Fig. 3 through solid lines. It is obvious from the figure that deflection profiles coming from direct measurement and image processing technique through the software (AutoCADÒ) are matching quite well at each load step. Hence the inherent error of distance distortion of the photographs get corrected in the process of bidirectional scaling, within acceptable limits. The indirect measurement technique is advantageous as it gives continuous deflection result for the complete domain, whereas direct measurement give deflection measures at some limited points. 3. Mathematical formulation Any particular formulation method requires, (i) identification of problem parameters, (ii) identification of kinetic and kinematic constitutive relations and (iii) identification of boundary conditions; Or in other words, it is necessary to identify the relevant issues of the problem at hand, for which the mathematical model is being developed. To determine the load-deflection behaviour of the master leaf spring, we focus on free body diagrams of different components of the system. The free body diagrams are drawn in equilibrium position, attained after application of load in a particular load step. The configuration attained in this equilibrium position is a function of its initial geometry and initial locked-up stress field. However, at the beginning of loading, there is no locked-up stress in the beam. 3.1. Kinetic and kinematic constitutive relations Free body diagram of the loaded master leaf spring acted upon by external vertical and horizontal force components W and F is shown in Fig. 4 in OX 0 Y 0 Z0 , a local Cartesian coordinate system. A detail list of notation is provided in Appendix A. Vertical load W is apparent in all applications, but F is absent in many cases and even when it is there, its presence is not manifested directly. The reaction forces in left and right roller supports are RxL ; RyL and RxR ; RyR respectively. From force equilibrium conditions, W ¼ RyL þ RyR and neglecting the friction forces at the contact points F ¼ RxL
RxR . A global Cartesian coordinate system OXYZ , as used earlier in profile measurement, is also used and its origin is taken on the horizontal axis of OX 0 Y 0 Z0 , where geometric centre of curved beam drops a perpendicular. To simulate the actual working condition of leaf springs, line of action of load W is not assumed concentric with the origin of OXYZ , and eccentricity between them is denoted by eW (positive towards the positive X axis). The centre

distance between the two support rollers is L0 , commonly known as the ‘span’ of leaf spring. Another important specification of leaf spring is ‘camber’ (marked as point C 0 ), which is a measure of its height. For symmetric spring geometry and concentric loading, camber point will lie on mid-span vertical, but this is not the condition for present generalized curved beam model. Points A; B are roller centres, whereas points A0 ; B0 are contact points of the rollers with ground. At load step (i
1) we have the known configuration, on which an incremental load DW is applied. Hence at load step i, applied load W is known (¼ W i
1 þ DW), but the new configuration is unknown. To find out deformed configuration and L0 with known values of W and eW , moment of all the forces about the contact point of left roller (A0 ) is balanced to yield,

RyR L0 ¼ WL0W þ FH0W

ð1Þ

(L0W ; H0W )

In Eq. (1), is the coordinate of load application point in OX 0 Y 0 Z0 coordinate system. Hence, L0W ¼ L01
eW . Now to explore kinematic and kinetic details, schematic diagrams of the left and right roller supports are shown in Fig. 5(a) and free body diagram of the left roller is presented in Fig. 5(b). The figures show two concentric rollers, one virtual eye roller with radius eL and another physical support roller of radius RRoller . Point A1 is the virtual contact point between the roller and the leaf spring where the centre line of the curved beam intersects the pitch circle of the eye roller. It is assumed that, at this point the slope of the beam is tangential to pitch circle of the roller and hence the radial line AA1 is normal to the beam centre line. The force components acting at the contact point A0 and the components acting at point A1 are also shown separately. The direction of resultant reaction 1=2

force RL ½¼ ðR2xL þ R2yL Þ

, acting along line A0 A1 is indicated by angle

w0L . It should be noted that the direction of RL is not necessarily perpendicular to the deflection curve at point A1 . The normal AA1 to the curved beam profile at point A1 makes an angle wL with negative X 0 axis and it is termed as reaction force angle at left end. This reaction angle wL yields the slope of curved beam (hL ) at left end A1 and given by wL ¼ p=2
hL . The kinematic relationship between w0L and wL is established as cotðw0L Þ ¼ ðeL cos wL Þ=½eL ð1 þ sin wL Þ, which leads to the kinetic relation

RxL =RyL ¼ cos wL =ð1 þ sin wL Þ:

ð2aÞ

Similarly in case of right end of the curved beam, reaction force angle is given by wR ¼ p=2 þ hR at B1 . And hence similar kinematic relationship cotðw0R Þ ¼ ðeR cos wR Þ=½eR ð1 þ sin wR Þ yields

RxR =RyR ¼ cos wR =ð1 þ sin wR Þ:

ð2bÞ

Point A1 is considered as the origin of two different local coordinate systems: one Cartesian Oxyz and the other one is curvilinear Osnz as shown in Fig. 6. In the body fitted ‘s
n’ coordinate, reaction force RL has normal and tangential components indicated by N L and T L . In ‘x
y’ coordinate RL has horizontal and vertical components, indicated by PxL and PyL , where P yL ¼ RyL and

Fig. 4. Free body diagram of the leaf spring.

S. Ghuku, K.N. Saha / Engineering Science and Technology, an International Journal 21 (2018) 408–420

413

Fig. 5. (a) Schematic diagram of left and right roller and (b) free body diagram of left roller support.

model, the object of the present study, is shown in Fig. 7, where the horizontal length of problem domain gets modified from L0 to L. In Osnz system the problem domain is defined by arc length S of the curved beam A1 B1 . Relationship between spans of the domains in the two different coordinate systems is given by RS L ¼ L0 þ eL cos wL þ eR cos wR ¼ 0 ds cos h. Span L gets redefined at each load step and the corresponding stretching effect redefines arc length S as well. 3.2. Coordinate systems and their transformation

Fig. 6. Two different coordinate systems at point A1 .

P xL ¼ RL cos w0L . Free body diagram of the right roller is similar to that of the left roller and hence it is not shown, but locations of points B, B0 and B1 and corresponding angles wR and w0R are shown in Fig. 5(a). In case of the right roller, normal and tangential components of reaction force RR at B0 are indicated by N R and T R respectively. Similarly, horizontal and vertical components of RR are indicated by PxR and PyR , where PyR ¼ RyR and P xR ¼ RR cos w0R . The roller supports of the physical system are eliminated to obtain a simple curved beam model and the reaction forces at the roller-beam interface are imposed at the two end points of the curved beam model. Free body diagram of this curved beam

Profile of curved beam, at no load condition as well as at loaded condition is represented in global Eulerian coordinate system (GCS) OXYZ , as shown in Figs. 4 and 7. Major computations like force balance, energy balance, etc., of the curved beam is carried out in body fitted Lagrangian coordinate system Osnz , the origin of which is at the virtual point of tangency of left roller (A1 ). Another body fitted Cartesian Lagrangian frame Oxyz is also used, the origin of which coincides with the origin of Osnz . For transformation between global Eulerian and body fitted Lagrangian frames, another local Eulerian frame OX 0 Y 0 Z0 fitted at the centre of left roller (refer Fig. 6) is considered. Due to change in length of domain with load change, Lagrangian coordinates are calculated after each incremental load step thus rendering the solution process as an updated Lagrangian analysis. Transformation relation between local and global Eulerian coordinate systems are given by

X 0 ¼ X þ L01

and Y 0 ¼ Y

ð3Þ

The transformation relation between local Eulerian frame OX 0 Y 0 Z0 and Lagrangian frame Oxyz is given by

x ¼ X 0 þ eL cos wL

and y ¼ Y 0
eL sin wL

ð4Þ

414

S. Ghuku, K.N. Saha / Engineering Science and Technology, an International Journal 21 (2018) 408–420

Fig. 7. Free body diagram of the initially curved beam in updated Cartesian and body fitted curvilinear coordinate systems.

Using the above relations, transformation relation between OXYZ and Oxyz coordinates is established as

x ¼ X þ L01 þ eL cos wL ¼ X þ L1

and y ¼ Y
eL sin wL

ð5Þ

Once beam profile is transformed in Oxyz frame, from the GCS, the coordinate measurement in Lagrangian frame Osnz is given by

Z s¼

Z

s

ds ¼ 0

x

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 f1 þ ðdy=dxÞ g dx

ð6Þ

0

In Oxyz system, problem domain is normalized (0 6 nx 6 1) with beam span along x axis (0 6 x 6 L). On the other hand, in case of analysis in Osnz system, normalization (0 6 ns 6 1) is done with arc length (0 6 s 6 S). However, as span and arc length change with load change, the before said normalizations are carried out with updated span and arc length of beam. Since derivation of governing equation of the present beam problem and its solution are carried out in Osnz system, ns will be denoted by n from here onwards. 3.3. Boundary condition The boundary conditions of the curved beam (refer Fig. 7), prescribed at two ends, involve a displacement condition and a force condition. As the present theoretical analysis is carried out through displacement based method, displacement boundary conditions are identified and presented in this section. The displacement field in Osnz coordinate system is expressed as usn ¼ us ^ en . es þ ðwb þ un Þ ^ Here, wb ðsÞ is transverse body displacement and us ðsÞ is in-plane displacement of beam centre line due to bending and stretching respectively and these displacement components are only functions of s. On the other hand, un ðs; nÞ is deformation displacement field due to shear deformation of beam in (s; n) plane and this function is used to incorporate warping of cross-section. Boundary conditions of the displacement fields are presented below.

ðiÞ us ¼ unknown at s ¼ 0; S

ð7aÞ

ðiiÞ wb ¼ 0 at s ¼ 0; S

ð7bÞ

ðiiiÞ un ¼ 0 at s ¼ 0; S and at n ¼ h=2

ð7cÞ

The in-plane displacement field in both the boundaries is unknown but the transverse displacement fields coming from bending and shear deformations are prescribed. The displacement field causes a rigid body motion in deformed configuration of the curved beam at each stages of incremental loading, which is in addition to the deformation displacement. Addition of some rigid body displacements to the coordinate origins will take care of actual boundary condition. Thus this type of beam bending problem with moving boundary is analyzed through updated Lagrangian approach.

3.4. Governing equation Governing equation of the curved beam model, presented in Fig. 7, is derived in Lagrangian frame Osnz through minimization of total potential energy functional (p), which is implemented by

dp ¼ 0:

ð8Þ

Total potential energy (p) consists of two parts – strain energy stored in beam (U) and work potential of the external loads (V). Strain energy U is derived from volume integral of stress-strain R R product U ¼ 12 v rss ess dv þ 12 v rsn esn dv . Expressions of axial strain (ess ) and shear strain (esn ) in terms of displacement fields are obtained in body fitted curvilinear frame following geometrically exact beam theory [18]. The two base vectors defining undeformed beam geometry in body fitted curvilinear coordinate system, are taken along tangential and normal directions to beam centre line. The base vectors are allowed to deform with the body and deformation of such covariant base vectors gives strain–displacement n  2 relations as ess ¼ ð1
1jnÞ 12 @u@sn þ u0s þ 12 ðu0s Þ2 þ 12 ðw0b Þ2
nw00b þ  n  n g and esn ¼ ð1
1jnÞ @u . These two equations are valid for w0b @u @s @s the prevailing beam geometry corresponding to a particular load step, and hence when geometry is updated subsequent to a load increment, the equations remain valid as well. Now substituting the strain – displacement relations for axial strain (ess ) and shear strain (esn ), and using linear stress – strain relations, strain energy equation becomes

(  4 1 1 @un 1 1 2 4 4 þ ðu0s Þ þ ðu0s Þ þ ðw0b Þ ð1
4 4 4 j nÞ @s 0
h=2  2  2  2 3 @un 1 @un 2 2 @un 2 þ ðu0s Þ þ ðu0s Þ þ n2 ðw00b Þ þ ðw0b Þ 2 2 @s @s @s  2  3 @un @un 3 2
n w00b þ w0b þ ðu0s Þ þ ðu0s Þðw0b Þ
2nðu0s Þðw00b Þ @s @s     @un 1 @un 2 2 2 2 þ ðu0s Þ ðw0b Þ
nðu0s Þ w00b þ ðu0s Þ w0b þ 2ðu0s Þðw0b Þ 2 @s @s     0 2 00 0 3 @un 00 0 @un
2nwb wb dn ds
nðwb Þ wb þ ðwb Þ @s @s ( )   Z Z 2 Gb S h=2 1 @un þ dn ds ð9Þ @s 2 0
h=2 ð1
jnÞ

Eb U¼ 2

Z

S

Z

h=2

Expression of V in terms of external load W, and considering locked-up moment Ms , in-plane force T s and shear stress ssn generated in previous load step, is given by

V ¼
ðW cos hÞwb jsW
ðW cos hÞun jn¼0 s¼sW
ðW sin hÞus jsW  Z S Z S Z S Z h=2  @un þ Ms w00b ds þ T s u0s ds þ b ssn dn ds @s 0 0 0
h=2

ð10Þ

415

S. Ghuku, K.N. Saha / Engineering Science and Technology, an International Journal 21 (2018) 408–420

Energy equations are normalized by using normalization parameters n ¼ s=S and g ¼ ðn þ h=2Þ=h. In this normalized domain the unknown displacement fields are assumed as un ðn; gÞ ¼ Pns Pnw Pnn i¼1 c 1i a1i ðnÞa2i ðgÞ, us ðnÞ ¼ i¼1 c 2i bi ðnÞ and wb ðnÞ ¼ i¼1 c 3i ci ðnÞ. Now putting the normalized energy expressions in terms of the approximated displacement fields into system governing equation (Eq. (8)) and carrying out variation on the unknown coefficients, set of governing equation of the present curved beam bending problem is obtained as

½Kfcg ¼ ff g

ð11Þ

In the above equation ½K is stiffness matrix, fcg is set of unknown coefficients of displacement field and ff g is the load vector, having the following forms

2

k11

6 ½K ¼ 4 k21 k31

k12

k13

3

k22

7 k23 5;

k32

k33

8 9 8 9 > > < c1 > = = fcg ¼ c2 and ff g ¼ f 2 : > > : > ; : > ; c3 f3

The elements of stiffness matrix ½K and load vector ff g are given in Appendix B. 3.5. Geometry updation Due to application of incremental load DW on top of W i
1 at i-th load step, curved beam (X i
1 ; Y i
1 ) undergoes displacement fields i

i

uX and uY and occupies configuration (X ; Y ) with respect to global frame OXYZ as shown in Fig. 8(a). Hence X i ¼ X i
1 þ uX and

Y i ¼ Y i
1 þ uY . However, geometry calculation at loaded condition is not so straight forward as it requires post processing of displacement field usn , obtained through energy principle based iterative method in Lagrangian frame Osnz . Hence, geometry updation is first carried out in body fitted curvilinear coordinate system Osnz and then transformed into global frame. In this body fitted frame, geometry is defined by arc length (s) and slope (h) at each point along beam centre line and deformation is defined by change in slope and in-plane displacement of beam centre line with respect to the origin of Osnz . Change in slope of beam profile (si
1 ; hi
1 ), due to application of incremental load   DW, is given by D/ ¼ w0b þ @s@ui
1n  i
1 . On the other hand, in-plane n

¼0

displacement of beam centre line coming from stretching field (us ) is given in Osnz frame by us
us jsi
1 ¼0 . However, in addition to this stretching field, curved beam domain also changes due to shifts of tangency points, as shown in Fig. 8(c). This domain increment imposes a displacement field along centre line (um ) and calculated by um ¼ eL D/, when 0 6 si
1 6 si
1 D/¼0 and um ¼ eR D/, when i
1 6 Si
1 . Displacement of beam centre line coming from si
1 D/¼0 6 s um in Osnz frame is given by um
um jsi
1 ¼0 . Hence, deformed configuration of beam centre line at i-th load step in curvilinear frame is given by

si ¼ si
1 þ us
us jsi
1 ¼0 þ um
um jsi
1 ¼0  @un  hi ¼ hi
1 þ w0b þ i
1  @s i
1 n

and

¼0

Fig. 8. Deformed and undeformed configurations in (a) global and (b) local Cartesian frame, and (c) change in domain due to rotation of rollers.

ð12Þ

416

S. Ghuku, K.N. Saha / Engineering Science and Technology, an International Journal 21 (2018) 408–420

Now deformed geometry in curvilinear frame (si ; hi ) is transformed into local Cartesian frame Oxyz using the relations

Z xi ¼

si

i

cosðhi Þ ds

Z and yi ¼

0

si

i

sinðhi Þ ds

ð13Þ

0

In this local Cartesian frame Oxyz , deformed beam geometry will look like the deformed configuration shown in Fig. 8(b). Finally beam configuration in global frame OXYZ is obtained by using transformation relation between the frames as described in Eq. (5) and presented below

X i ¼ xi
ðL1 Þi

and Y i ¼ yi þ eL sin ðwL Þi

ð14Þ

4. Comparison and discussion The necessary starting functions for the unknown displacement fields are selected as a11 ¼ nð1
nÞ, a21 ¼ gð1
gÞ, b1 ¼ ðn
0:5Þ2 and c1 ¼ nð1
nÞ, satisfying the geometric boundary conditions as presented in Eq. (7). The higher order orthogonal functions are generated using Gram-Schmidt scheme. Number of coordinate functions are taken as nn ¼ ns ¼ nw ¼ 5, which provides satisfactory convergence to the computation scheme. As governing equation (Eq. (11)) involves strong nonlinearity, the solution process is implemented by successive relaxation scheme considering error limit for energy balance as l1 ¼ 1e
05. Once solution of unknown coefficients is converged, displacement fields become known and geometry is updated using the geometry updation scheme presented in previous section. With the updated beam geometry and known strain fields, incremental bending moment, in-plane force, shear force and shear stress are calculated from their fundamental definitions as given by

ðDM s Þi ¼
ðDQ s Þi ¼


Z Z

h=2


h=2 h=2


h=2

Ebess ndn;

ðDT s Þi ¼

Z

h=2
h=2

Ebess dn;

Gbesn dn and ðDssn Þi ¼ Gesn :

Now bending moment, in-plane force, shear force and shear strain fields at current load step are calculated by adding these incremental fields to the corresponding fields of previous load step: given by ðM s Þi ¼ ðM s Þi
1 þ ðDM s Þi , ðT s Þi ¼ ðT s Þi
1 þ ðDT s Þi , ðQ s Þi ¼ ðQ s Þi
1 þ ðDQ s Þi and ðssn Þi ¼ ðssn Þi
1 þ ðDssn Þi . These fields at current load step i will act as locked-up fields at the next load step i þ 1. If difference between the calculated average shear force and applied load, defined by l2 ¼ jðQ a Þi
W i j=W i , is greater than a predefined error limit, whole calculations of energy balance, geometry updation and force balance are repeated with DW i ¼ DW i =2. For the present computation, allowable limit of error l2 is considered as 7%. No-load profile of the master leaf spring presented in Table 1 is now transformed into global reference frame OXYZ and values of camber (C 0 ¼ 129:04 mm) and eccentricity with respect to geometric centre (eW ¼ 1:74 mm) are calculated. With the modified geometry of the experimental specimen, deflection profiles are generated at different load levels using the present numerical method and presented in Fig. 9(a). The contour plot in Fig. 9(a) indicates that specimen shape remain symmetric with increase in load. The computational three dimensional deflection profile surface is interpolated at the experimental load values, as marked on the load axis, and presented separately in Fig. 9(b). Fig. 9(b) also shows experimental deflection profiles and hence a comparison between numerical and experimental results are obtained. Variations of span (L) and arc length (S) at different experimental load values are presented in Fig. 10

Fig. 9. (a) Deflection profiles of the master leaf spring at different load obtained numerically and (b) comparison between numerically and experimentally obtained deflection profiles.

S. Ghuku, K.N. Saha / Engineering Science and Technology, an International Journal 21 (2018) 408–420

417

Fig. 10. Variations of numerical and experimental (a) span and (b) arc length with load.

(a, b). As mentioned earlier in connection with solution method, the present numerical approach is implemented through both energy and force balance at each incremental load step. Hence to assure satisfaction of these conditions, strain energy (U), potential energy (V) and average shear force (Q a ) are also calculated corresponding to each experimental load values. The numerically calculated energies (U and V) together with error in energy balance (¼ p
2U, which is V
U), are plotted against load and shown in Fig. 11(a). On the other hand, error in shear force l2 (in %), as observed at the experimental load values, is shown in Fig. 11(b). It is obvious from the comparisons, presented in Fig. 9(b) and Fig. 10(a, b), that numerical and experimental results are matching quite well at lower load values. However, discrepancy increases with increase in load and numerically obtained deflections are found less compared to experimental ones, indicating relatively more weakening (reduced stiffness) of physical structure with loading. Stiffness of physical system decreases with loading due to presence of hole in the master leaf spring, which is absent in numerical model. In addition, the hole causes local increment of stress and this stress concentration effect adds more weakness to the physical structure. Incorporation of such geometry variation in mathematical modeling may leads to more realistic simulation of the actual experimental system. However, implementation of such effects into the numerical scheme requires a different and broader research field and hence these effects are not included in the present model. However, within the framework of the present model, comparison of numerical results with experimental ones can be considered as satisfactory validation.

5. Conclusions Large deflection characteristics of initially curved beam with moving boundaries, subjected to combined bending-stretchingshear deformation, have been studied experimentally and theoretically using updated Lagrangian approach. A master leaf spring has been considered as equivalent physical system of the nonlinear curved beam problem and its deflection characteristics under static load are obtained experimentally using direct measurement and image processing technique in a specially designed three point bending testing rig. After identifying the relevant physical parameters and kinematic and kinetic constitutive relations, geometric nonlinear governing equation has been derived using minimization of total potential energy functional. The problem formulation is made generalized considering non-uniformly curved beam with asymmetric geometry and eccentric loading with respect to geometric centre. Iterative solution process is implemented in the computational environment of MATLABÒ. At each incremental load step during geometry updation, force balance condition is also satisfied in addition to energy balance. As the theoretical curved beam model includes several complicating effects in addition to nonlinear kinematics, it can be used in standalone mode for analysis of several complicated curved beam like structures. Moreover, the theoretical and experimental results can be of interest specifically in design applications of leaf spring system. Finally, several issues of the physical system have been addressed from the comparison between experimental and theoretical results. Some comments related to further improvement in simulation of large deflection behaviour of the slender physical structure have been reported.

Fig. 11. Variations of (a) energy together with error in energy balance and (b) error in shear force balance with load.

418

S. Ghuku, K.N. Saha / Engineering Science and Technology, an International Journal 21 (2018) 408–420

Appendix A. Notation A A; B A0 ; B0 A1 ; B1 b fcg C0 eL ; eR eW E ff g F G h HW H0W ½K L L0 LW L0W Ms ng nn; ns; nw NL ; T L ; NR ; T R Oxyz ; Ox0 y0 z0 ; Osnz OXYZ OX 0 Y 0 Z 0 PxL ; PyL ; P xR ; PyR Qa Qs RL ; RR RRoller RxL ; RyL ; RxR ; RyR sW S Ts un ; us usn U

v

V wb W W max Xe ; Y e

a; b; c

Cross sectional area of beam Roller centres Contact points of left and right rollers with ground Virtual contact points of left and right rollers with beam Width of beam cross section Set of unknown coefficients Camber of leaf spring in OX 0 Y 0 Z 0 or OXYZ Pitch circle radius of left and right rollers Eccentricity of load W measured in OXYZ Young’s modulus of beam material Load vector Internally generated friction force Modulus of rigidity of beam material Thickness of beam cross-section Ordinate of loading point in Oxyz Ordinate of loading point in OX 0 Y 0 Z 0 Stiffness matrix Span of problem domain in Oxyz Distance between roller centres of master leaf spring i.e., span Abscissa of loading point in Oxyz Abscissa of loading point in OX 0 Y 0 Z 0 Bending moment field Number of precision points used for computation Number of orthogonal functions for un ; us ; wb Normal and tangential components of reaction forces at points A1 and B1 Body fitted Lagrangian Cartesian and Curvilinear systems Global Eulerian coordinate system Local Eulerian coordinate system Horizontal and vertical components of reaction forces at points A1 and B1 Average shear force Shear force field Reaction forces at left and right rollers at points A0 ; A1 and B0 ; B1 Outer radius of left and right rollers Reaction forces at points A0 and B0 in left and right roller supports s coordinate of loading point in Osnz Arc length of curved beam (A1 B1 ) In-plane force field coming from stretching field Displacement fields due to shear and stretching in Osnz Displacement field in Osnz Strain energy Volume of beam Potential energy due to external loading Displacement field due to bending in Osnz Externally applied vertical load (at ith load step W ¼ W i ) Maximum value of externally applied load Cartesian frame used in experimental geometry measurement Set of orthogonal functions for shear, stretching and bending displacement fields

DW ess ; esn

Load increment Axial and shear strain Normalized n coordinate Slope of curved beam profile Curvature of curved beam centre line Error limit used in iterative energy balance Error in iterative shear force balance Normalized s coordinate Normalized s coordinate of loading point Total potential energy Radius of curvature of curved beam profile (¼ 1=j) Axial and shear stress Locked-up shear stress Reaction force angles at left and right roller centres (wL ¼ p=2
hL at s ¼ 0 and wR ¼ p=2 þ hR at s ¼ S) Reaction force angles at left and right roller support points

g h

j l1 l2 n nW

p q rss ; rsn ssn wL ; wR

w0L ; w0R

Appendix B Elements of stiffness matrix

k11 ¼

nn X nn Z X i¼1 j¼1

0

8 2 !2  < X nn 1 EA 0 4 4 c1k a1k a2k f1
jhðg
0:5Þg 2S3 : k¼1 2

Z

1

1

0

nw X þ3 c3k c0k

!2

ns X þS c2k b0k

k¼1

! þ

k¼1

!

nw X


hðg
0:5Þ

00 k

c3k c

þ2

!2 c2k b0k

k¼1 nw X

k¼1

#  # GA 0 0 a1i a2i a1j a2j dgdn þ S

k12

ns X

!

0 k

c 3k c

k¼1

nn X c1k a01k a2k

!)

k¼1

( ! nn X 1 0 ¼ 2 c1k a1k a2k 2S i¼1 j¼1 0 0 f1
jhðg
0:5Þg k¼1 ! ! !) # nw ns nw X X 2 X b0i a01j a2j dgdn c3k c0k þ c2k b0k c3k c0k þ 3S k¼1 k¼1 k¼1 ns X nn Z EA X

k13 ¼

nw X nn Z EA X

2S3

i¼1 j¼1

1

0

1

Z

Z

1 0

1

"

"

" ( ! nn X 1
hðg
0:5Þ c1k a01k a2k f1
jhðg
0:5Þg k¼1 !)

nw X þ hðg
0:5Þ c3k c0k

c00i a01j a2j

k¼1

8 !2 ! !2 nn ns ns < X X 1 X þ c1k a01k a2k þS c2k b0k þ c2k b0k : k¼1 3 k¼1 k¼1 3 !2 !9 nw nw = X X 0 00 0 0 þ c3k ck
hðg
0:5Þ c3k ck c a a 5dgdn ; i 1j 2j k¼1 k¼1

k21

( ! nn X 1 0 ¼ 2 c1k a1k a2k 2S i¼1 j¼1 0 0 f1
jhðg
0:5Þg k¼1 ! ! !) # nw ns nw X X 2 X 0 0 0 0 0 c3k ck þ c2 b c3k ck a1i a2i bj dgdn þ 3S k¼1 k k k¼1 k¼1 nn X ns Z EA X

1

Z

1

"

S. Ghuku, K.N. Saha / Engineering Science and Technology, an International Journal 21 (2018) 408–420 ns X ns Z EA X

k22 ¼

2S

3

Z

1 0

i¼1 j¼1

ns X 1 4 2S2 þ c2k b0k f1
jhðg
0:5Þg : k¼1

1

0

!2

nn X

c1k a01k a2k

þ

8 <

2

ns X

þ 3S

k¼1 nw X
hðg
0:5Þ c3k c00k

!

k¼1

k23 ¼

nw X ns Z EA X

Z

1

!

c2k b0k

nw X c3k c0k þ

k¼1 nw 2 X c 3 c0 þ 3 k¼1 k k

!

Appendix C. Supplementary data Supplementary data associated with this article can be found, in the online version, at https://doi.org/10.1016/j.jestch.2018.04.007.

!2

k¼1 nn X c1k a01k a2k

!)

# b0i b0j dgdn

References

!

[1] K.E. Bisshopp, D.C. Drucker, Large deflection of cantilever beams, Q. Appl. Math. 3 (3) (1945) 272–275. [2] M. Batista, Large deflections of a beam subject to three-point bending, Int. J. Nonlinear Mech. 69 (2015) 84–92. [3] C.M. Wang, K.Y. Lam, X.Q. He, S. Chucheepsakul, Large deflections of an end supported beam subjected to a point load, Int. J. Nonlinear Mech. 32 (1) (1997) 63–72. [4] L. Chen, An integral approach for large deflection cantilever beams, Int. J. Nonlinear Mech. 45 (3) (2010) 301–305. [5] M. Dado, S. Al-Sadder, A new technique for large deflection analysis of nonprismatic cantilever beams, Mech. Res. Commun. 32 (6) (2005) 692–703. [6] X.T. He, L. Cao, Z.Y. Li, X.J. Hu, J.Y. Sun, Nonlinear large deflection problems of beams with gradient: a biparametric perturbation method, Appl. Math. Comput. 219 (14) (2013) 7493–7513. [7] A. Banerjee, B. Bhattacharya, A.K. Mallik, Large deflection of cantilever beams with geometric non-linearity: analytical and numerical approaches, Int. J. Nonlinear Mech. 43 (5) (2008) 366–376. [8] H. Tari, On the parametric large deflection study of Euler-Bernoulli cantilever beams subjected to combined tip point loading, Int. J. Nonlinear Mech. 49 (2013) 90–99. [9] G.D. Angel, G. Haritos, A. Chrysanthou, V. Voloshin, Chord line force versus displacement for thin shallow arc pre-curved bimetallic strip, P. I. Mech. Eng. C-J. Mec. 229 (1) (2015) 116–124. [10] T. Beléndez, C. Neipp, A. Beléndez, Large and small deflections of a cantilever beam, Eur. J. Phys. 23 (3) (2002) 371. [11] S. Ghuku, K.N. Saha, A theoretical and experimental study on geometric nonlinearity of initially curved cantilever beams, Eng. Sci. Technol. Int. J. 19 (2016) 135–146. [12] M. Mutyalarao, D. Bharathi, B.N. Rao, Large deflections of a cantilever beam under an inclined end load, Appl. Math. Comput. 217 (7) (2010) 3607– 3613. [13] M. Mutyalarao, D. Bharathi, B.N. Rao, On the uniqueness of large deflections of a uniform cantilever beam under a tip-concentrated rotational load, Int. J. Nonlinear Mech. 45 (4) (2010) 433–441. [14] A.K. Nallathambi, C.L. Rao, S.M. Srinivasan, Large deflection of constant curvature cantilever beam under follower load, Int. J. Mech. Sci. 52 (3) (2010) 440–445. [15] S. Navaee, R.E. Elling, Equilibrium configurations of cantilever beams subjected to inclined end loads, J. Appl. Mech.-T. ASME 59 (3) (1992) 572–579. [16] K. Murty, Toward a consistent beam theory, AIAA J. 22 (6) (1984) 811–816. [17] M.A. Benatta, I. Mechab, A. Tounsi, E.A. Bedia, Static analysis of functionally graded short beams including warping and shear deformation effects, Comput. Mater. Sci. 44 (2) (2008) 765–773. [18] K. Washizu, Some considerations on a naturally curved and twisted slender beam, Stud. Appl. Math. 43 (1–4) (1964) 111–116. [19] M.A. Crisfield, G. Jelenic´, Objectivity of strain measures in the geometrically exact three-dimensional beam theory and its finite-element implementation, P. Roy. Soc. Lond. A Mat. 455 (1999) 1125–1147. [20] S. Ghuku, K.N. Saha, A review on stress and deformation analysis of curved beams under large deflection, Int. J. Eng. Technol. 11 (2017) 13–39. [21] A. Mitra, P. Sahoo, K. Saha, Large displacement of crossbeam structure through energy method, Int. J. Automot. Mech. Eng. 5 (2012) 520–544. [22] W. Lacarbonara, H. Yabuno, Refined models of elastic beams undergoing large in-plane motions: theory and experiment, Int. J. Solids Struct. 43 (17) (2006) 5066–5084. [23] M. Saje, Finite element formulation of finite planar deformation of curved elastic beams, Comput. Struct. 39 (3) (1991) 327–337. [24] A. Atig, R. Ben Sghaier, R. Seddik, R. Fathallah, A simple analytical bending stress model of parabolic leaf spring, P. I. Mech. Eng. C-J. Mec (2017). 10.1177/ 0954406217709302. [25] H.A.F.A. Santos, P.M. Pimenta, J.M. De Almeida, Hybrid and multi-field variational principles for geometrically exact three-dimensional beams, Int. J. Nonlinear Mech. 45 (8) (2010) 809–820. [26] P. Hansbo, M.G. Larson, K. Larsson, Variational formulation of curved beams in global coordinates, Comput. Mech. 53 (4) (2014) 611–623. [27] H. Sugiyama, A.A. Shabana, M.A. Omar, W.Y. Loh, Development of nonlinear elastic leaf spring model for multibody vehicle systems, Comput. Methods Appl. M. 195 (50) (2006) 6925–6941. [28] D.A. Miller, A.N. Palazotto, Nonlinear finite element analysis of composite beams and arches using a large rotation theory, Finite Elem. Anal. Des. 19 (3) (1995) 131–152. [29] G. Wu, X. He, P.F. Pai, Geometrically exact 3D beam element for arbitrary large rigid-elastic deformation analysis of aerospace structures, Finite Elem. Anal. Des. 47 (4) (2011) 402–412. [30] J.N. Reddy, I.R. Singh, Large deflections and large-amplitude free vibrations of straight and curved beams, Int. J. Numer. Meth. Eng. 17 (6) (1981) 829–852.

k¼1

"

1

!2

"(

nw X

1 c3k c0k S f1
jhðg
0:5Þg k¼1 ! ! !) nn ns nn X X 2 X þS c1k a01k a2k þ c2k b0k c1k a01k a2k c0i b0j 3 k¼1 k¼1 k¼1 ( !) # ns X 0 00 0
2Shðg
0:5Þ þ hðg
0:5Þ c 2k b k ci bj dgdn 2S3

i¼1 j¼1

0

0

Z

"

k¼1

k31 ¼

nn X nw Z EA X

2S3

1

0

i¼1 j¼1

1

0

" ( ! nn X 1 c1k a01k a2k
hðg
0:5Þ f1
jhðg
0:5Þg k¼1 !)

nw X þ hðg
0:5Þ c3k c0k

a01i a2i c00j

k¼1

8 !2 ! !2 nn ns ns < X X 1 X 0 0 0 þ c a a þS c2k bk þ c2 b : k¼1 1k 1k 2k 3 k¼1 k k k¼1 9 3 !2 ! nw nw = X X 0 00 0 05 c 3k c k
hðg
0:5Þ c 3k c k a a c dgdn þ ; 1i 2i j k¼1 k¼1

k32 ¼

"( ! nw X 1 0 c c S 3k k 2S3 i¼1 j¼1 0 0 f1
jhðg
0:5Þg k¼1 ! ! !) nn ns nn X X 2 X 0 0 0 b0i c0j c1k a1k a2k þ c2 b c1k a1k a2k þS 3 k¼1 k k k¼1 k¼1 ( !) # ns X 0 0 00
2Shðg
0:5Þ þ hðg
0:5Þ c 2k b k bi cj dgdn ns X nw Z EA X

Z

1

"

1

k¼1

k33 ¼

nw X nw Z EA X

2S3

i¼1 j¼1

0

1 0

nn X c1k a01k a2k

þ3

28 !2 nw < X 1 4 4 c3k c0k f1
jhðg
0:5Þg : k¼1 2

Z

1

!2

k¼1 nw X


2hðg
0:5Þ

ns X

þS ! c 3k c

00 k

(


hðg
0:5Þ

nw X c3k c0k

!2

k¼1

nw X

þ2

nn X c1k a01k a2k k¼1

ns X c2k b0k

þ

k¼1

k¼1


hðg
0:5Þ

! c2k b0k

! c 3k c

nn X c1k a01k a2k

0 k

k¼1

!

!)

c0i c0j

k¼1

c00i c0j þ 2h2 ðg
0:5Þ2 c00i c00j

!

nn X

þ hðg
0:5Þ

k¼1

!) c1k a01k a2k

#

c0i c00j dgdn

k¼1

Elements of load vector nn  X  f 1 ¼ ðW cos hÞ a1j a2j  j¼1

n¼nW ;g¼0:5

ns ns Z X X  f 2 ¼ ðW sin hÞ bj n
j¼1

 nw X  f 3 ¼ ðW cos hÞ cj  j¼1

W

nW

j¼1




1

nw 1X S j¼1

nn Z X
bh j¼1

1

0

T s ðnÞb0j dn

0

Z 0

1

419

M s ðnÞc00j dn

Z 0

1

ssn ðn; gÞa01j a2j dgdn

420

S. Ghuku, K.N. Saha / Engineering Science and Technology, an International Journal 21 (2018) 408–420

[31] P.F. Pai, T.J. Anderson, E.A. Wheater, Large-deformation tests and totalLagrangian finite-element analyses of flexible beams, Int. J. Solids Struct. 37 (21) (2000) 2951–2980. [32] L.N.B. Gummadi, A.N. Palazotto, Large strain analysis of beams and arches undergoing large rotations, Int. J. Nonlinear Mech. 33 (4) (1998) 615–645. [33] K.J. Bathe, S. Bolourchi, Large displacement analysis of three-dimensional beam structures, Int. J. Numer. Methods Eng. 14 (7) (1979) 961–986. [34] Z.Q. Chen, T.J.A. Agar, Geometric nonlinear analysis of flexible spatial beam structures, Comput. Struct. 49 (6) (1993) 1083–1094. [35] T. Belendez, C. Neipp, A. Beléndez, Numerical and experimental analysis of a cantilever beam: a laboratory project to introduce geometric nonlinearity in mechanics of materials, Int. J. Engng. Ed. 19 (6) (2003) 885–892. [36] A. Dorogoy, D. Rittel, Transverse impact of free–free square aluminum beams: An experimental–numerical investigation, Int. J. Impact Eng. 35 (6) (2008) 569–577.

[37] A.G. Rodríguez, J.M. Chacón, A. Donoso, A.G. Rodríguez, Design of an adjustable-stiffness spring: mathematical modeling and simulation, fabrication and experimental validation, Mech. Mach. Theory 46 (12) (2011) 1970–1979. [38] G.S. Shankar, S. Vijayarangan, Mono composite leaf spring for light weight vehicle–design, end joint analysis and testing, Mater. Sci. 12 (3) (2006) 220– 225. [39] A.B. Pippard, The elastic arch and its modes of instability, Eur. J. Phys. 11 (6) (1990) 359–365. [40] H. Ahuett-Garza, O. Chaides, P.N. Garcia, P. Urbina, Studies about the use of semicircular beams as hinges in large deflection planar compliant mechanisms, Precis. Eng. 38 (4) (2014) 711–727.