Abstract

Based on Timoshenko’s beam theory, this paper adopts segmented strategy in establishing the governing equations of a multibeam system subjected to various boundary conditions, in which free, clamped, hinged, and elastic constraints are considered. Meanwhile, Galerkin method is incorporated as a competitive alternative, in which a new set of unified, efficient, and reliable trial functions are proposed. A further optimization in regard to boundary distributions under forces is implemented and established on the least absorbed energy principle. High agreement is observed between the analytical results and the FEM results, verifying the correctness of the derivations. Complete comparisons between the analytical and the numerical results indicate the Galerkin method is beneficial when slender ratio is larger than 30, in which the continuity of the deformation is proved to be a crucial influencing factor. A modified numerical strategy about optimal boundary is employed and the remarks imply the algorithm can be availably used to reduce the energy absorption of the whole system.

1. Introduction

The static and dynamic properties of the beam structures in engineering possess great significance, in which the achievements can be directly implemented in wide fields including civil engineering, composite material manufacture, and wing design. In beam analysis, the neglect of the shear and the rotary inertia restricts the utilization of the Euler-Bernoulli theory in short beam cases [1], while the Timoshenko theory can solve these issues well and provide far more reliable solutions [24]. Papers reported in [5, 6] were such examples, in which the mode superposition method and the Laplace transform method, for the analytical characteristics, have become the leads. Nevertheless, in some cases that the analytical solution cannot be derived, for example, when system incorporates complicated boundaries, the typical numerical methods like the conventional residual method [7], the Rayleigh-Ritz method [810], and the Galerkin method [11, 12] can provide alternative ways. Further, on account of the great computing ability in computers, the finite element method (FEM) and the boundary element method (BEM) have been developed to calculate more sophisticated structural systems, for example, when system subjected to complicated loads and boundary conditions, which are also regarded as criterions to verify analytical and numerical methods.

The beam systems with various boundary conditions (abbreviated as BCs) and internal constraints (also called BCs), for example, clamped, free, hinged, and elastic BCs, are named as multibeam systems. Actually, such systems have been investigated constantly during recent years [1315]. Reference [5] has computed a beam system with hinge-joint constraint by employing BEM and FEM and compared the Euler-Bernoulli and the Timoshenko models, in which uniform descriptions for deflection, rotation, bending moment, and shear forces were derived at the first time. Reference [16] built up a FEM model in which the natural frequencies and the mode shapes of an Euler-Bernoulli pontoon system applied with arbitrary numbers of internal hinges were solved via the transfer matrix method. Through a limited superposed method in Euler-Bernoulli model, [17] studied the vibrational power transmission in a uniform and multisupported beam. Reference [18] calculated the eigenvalues by a proposed algorithm and the validity was verified through case studies. Reference [19] considered a multispan Euler-Bernoulli beam with a series of one-sided inner constraints to study the numerical simulations through assumed mode method. Further, [20] analyzed the multispan beams carrying multiple spring-mass systems with axial force based on Timoshenko theory. The free vibration of a multispan Timoshenko beam was studied by Rayleigh-Ritz method in [21], where the trial functions were composed of a series sine or cosine functions plus a set of polynomial functions. Furthermore, by using the transfer matrix method, [22] calculated the eigensolutions of a multispan beam with arbitrary numbers of flexible constraints and compared with the result of the Euler-Bernoulli beam. Reference [23] investigated the free vibrations of the Timoshenko beam with an internal hinge and computed the optimal location of this internal hinge maximized the fundamental frequency. Nevertheless, although most types of the constraints are mentioned in above literatures, few efforts has been made in counting more complex BCs. Additionally, a majority of references lack parallel comparisons with analytical method, the numerical approach, or the FEM. Hence, the primary goal of this article is to analyze a constructed multibeam system established in Timoshenko theory with arbitrary types and numbers of BCs, in which remarkable efforts have been made in derivations and case studies. Based on that, the optimization for BCs under given loads is put forward through the minimum energy principle.

This paper is organized as follows. In Section 2, the nondimensional governing equations of the multibeam system are derived through Hamilton’s principle and an elementary solution is obtained where several typical BCs are also listed. Section 3 performs the detailed derivations about the beam under arbitrary BCs through an analytical segmented method, in which the matrices are connected through transformation. As an alternative, Galerkin method is introduced in which a group of new trial functions with orthogonal properties are conducted. Section 4 puts forward the optimization process for optimal BCs based on Section 3, achieving the minimum energy absorption under specific loads. A total analysis is drawn in Section 5, in which validity, flexibility, and efficiency of the analytical and the Galerkin methods are measured through FEM in case studies. Based on that, a further optimization about BCs is conducted to verify the derivation in Section 4. Finally, the most relevant conclusions are summarized in Section 6.

2. Modeling Based on Timoshenko Theory

As depicted in Figure 1, an ordinary beam system is configured under external loads with various BCs (at , the coordinate is ), each BC herein can be free, simply supported, clamped, hinge-joint, or elastic constraint. As described, and denote the external force and moment respectively. is the longitudinal length of the entire beam. Suppose the numbers of elastic BCs are () and Hamilton’s principle is generally used to construct the governing equations of such multibeam system based on Timoshenko’s beam theory. Since this approach is widely used [2426], we refine the derivations as follows.

In consideration of the lateral rotation and the shearing effect, the kinetic energy and the internal energy can be expressed asHerein, , , and represent the bending deflection, bending rotation, and shear deformation respectively, and and are the elastic modulus and the shearing modulus. is the so-called shear coefficient [27, 28], giving the shear strain ratio on cross section. is the moment inertia in cross section, is the solid density, and are the elastic and torsional spring coefficients of the BC at , and is the Dirac delta function. In view of the inherent relationship within Timoshenko model [24], the shear in beam yields

For simplification, distributed loading functions in regard to the external forces are originally used and hold

In Hamilton’s manner, , and through the variation procedure herein and in conjunction with (2), the final governing differential equations are expressed as

For unification, we nondimensionalize (4a) and (4b) by introducing , , , , , , , , , and . Thus, (4a) and (4b) are rewritten aswhere symbols   ,, , and are employed to denote , , , and , respectively.

3. Derivations of Solution

3.1. Analytical Method (A-Method) through Segmented Strategy

Considering the big challenges in obtaining the analytical solution of this multiple systems with plenty different BCs, we observe, owning to the segmented geometries, one can solve the vibrational issues via (A.3) for each subdivided part. This inspiration seems falling in between the range of the matrix-transforming method and the FEM [29], but it is generally different from both of them. On the one hand, the arbitrary internal BCs will result in various transfer matrices. On the other hand, it is similar to the solving process of FEM but without any assumptive approximation. Therefore, the strategy to solve this multisystem with segments proposed here is that (the numbers of BCs are and we use to express it) we use the fundamental solution shown in Appendix A to solve first and connect them in sequence via physical boundary relations, by which we are able to get the nature frequency and its modal shape through eigenanalysis. Herein, a classified discussion about different BCs is listed as follows due to the arbitrariness of the BCs.

At first, a series of dimensionless transformations in localised coordinates are made here. For segment , we introduce , , , and to guarantee a uniform throughout the segments, and as a ripple effect, and in each fundamental solutions are changed into

Hence, the shape function of segment naturally holdsIn which , . To get the specific , each is to be determined. In view of the dissimilarity of BCs at , we generate the following: (1) When (BC at ), we enumerate the decompositions of the BCs as follows:(a)For free BC(b) For simply supported BC(c) For clamped BC(d) For elastic BC(2)When , the free constraint disappears, but the hinge-joint BC is added.(a) For simply supported BC(b) For hinge-joint BC(c) For clamped BC(d) For elastic BC(3)When , accordingly one has the following:(a) For free BC(b) For simply supported BC(c) For clamped BC(d) For elastic BC

In conclusion, all mentioned expressions have generated two linear equations at two edge nodes and four linear equations at inner nodes about (detailed information is shown in Appendix B), which constitute linear equations of nodes in total denoted by . That is,where the coefficient matrix just involves the independent unknown variable . We further define the characteristic function asBy which we can obtain each in sequence after solving this super sophisticated transcendental equation . Meanwhile, and are available to obtain through

Ultimately, for , through coordinate transformation, we generate

3.2. Numerical Solution Based on Galerkin Method (G-Method)

Albeit the achievement in analytical method, the explicit is far more complex and costly in actual computation, especially when becomes larger. Actually, in allusion to the practicability, the static and dynamic properties of the multisystem, in many circumstances, are determined by the lower modes . Based on such consideration, the numerical solution, as an alternative way, may probably provide credible results with avoiding massive computing costs. Hence, the objective of this part is the derivations about numerical G-method.

Based on G-method, the residual functions of the entire system arewhere and are all differential operators about the physical residual and the natural residual. and separately, as the external forces and the physical constraints of BC, stand for the given restricting equations independent on the unknown function . Herein, is naturally satisfied. In general, the trail functions possess

Accordingly, the functions of residual errors are

To reduce the complexity, is usually satisfied in advance [11, 22, 25], which implies that all and should be constructed at first. According to Appendix C, the explicit and for free BCs are hard to work out because of the coupling in (C.1) about displacement and rotation from shear component. Nevertheless, if we reduce the demands, that is, to decouple (C.1) and replace it by a weaker Euler-Bernoulli condition in (27), not only will the difficulty be solved, but also the inner beam maintains Timoshenko precision

To guarantee an entire , the polynomials are employed to structure the fundamental trial functions and a series of and are proposed by expanding previous result respectively via polynomial, Legendre orthogonal polynomials, Chebyshenv orthogonal polynomials, and trigonometric functions, which is regarded as the explorations of preferable . Here, we have given the exact nature conditions in Appendix C that the trial functions should hold. Further, if the system has not clamped BCs in (actually, if clamped BCs exist in inner nodes, the total system can be divided into several individual beams), a set of explicit trial functions can be proposed aswhich include four kinds of : (1)Polynomial form (P-method)(2)Chebyshev form (C-method)(3)Legendre form (L-method)(4)Trigonometric form (T-method)

With the explicit , the final derivation is obtainable, and we conduct

and make

After substituting the specific expressions to (34) and a group of integrations [11, 12] (detailed in Appendix D), the special variables are vanished, leaving behind a time domain equation denoted as

The characteristic equation can be easily carried out through , by which it is capable of getting the approximate frequencies through . According to (here ), the corresponding eigenvector can be solved and consequently described as

4. The Optimization of BCs

Under external loads, the positions of BCs will determine the dynamic performance of the whole system, wherein a reasonable distribution about these nodes may improve the loading capacity and its stable margin [23]. Based on the work in Section 3, it is capable of acquiring the optimal BCs. For mode , through a series of modal transformations as reported in [30], we have

Under certain loads, we write , in modal superposition manner; (37a) and (37b) satisfy

Based on modal orthogonality, we obtain

To examine the quantitative performance, we further define the total absorbed energy of the whole system as

Hence, to achieve the minimum energy absorption, (40) must yield . In virtue of the derivations in Section 3, we get and , all of which are the functions of . Herein, we suppose the beam without initial disturbance and, after Duhamel integration, get the modal response as

After substituting it to (40), we obtain

Particularly for time-invariant loads, it is relatively easy to calculate the absorbed energy when system reaches into equilibrium state (), which holds

Further, to derive the optimal layout, (42) and (43) must satisfy (44), which comes down to a minimization issueBy which the optimal can be worked out.

5. Numerical Verification

To demonstrate the validity of the proposed solutions and make comparisons each other, we conduct 2 typical numerical examples illustrated in Figure 2 named as case 1 and case 2, respectively, in which different BCs are investigated and inherent variables are concerned, for example, the slender ratio and the BC positions . For the purpose of evaluation, the results of FEM are generated through Abaqus 6.13, whereby the verification is demonstrated and complete comparisons are made. Further, based on the effectiveness of Section 5.2, case 2 is chosen to carry out a further optimization study in Section 5.3.

5.1. Considerations during Computations

Before the numerations, several considerations are declared to obtain more reliable consequences. Strategy of algorithm: since the numerations in such issues involve massive matrix computations, Matlab programs are utilized. Whereas in special situations, for example, or , the bad condition-numbers of or may exist and distort the precision. Based on that, an examination and repairing step is incorporated during the calculation. Additionally, we reserve the real part of to solve some big eigenmatrices and acquire the eigenvalues and eigenvectors, holding which is proved with pretty effectiveness. “Mirror verification” (interpreted as algorithm flexibility analysis): for a set of symmetrical BCs, for example, , and , , the results of mode frequencies and mode shapes in each pair should be theoretically identical and symmetrical. This principle is regarded as a measurement about algorithm stability for different methods, whereby we call the “mirror verification.” Gridding scheme in FEM: basically, the accuracy of the FEM mostly depends on the mesh and the grid number, and hence the distributions of meshes as depicted in Figure 2 guarantee that no matter how long the proportion in each part, as sufficient as nodes are generated for each segment to improve the reliability of FEM.

5.2. Cases Studies

As shown in Figure 2, a beam with two elastic edges and two internal hinged joints are considered in case 1, while in case 2 a beam with two free BCs and two simply supported internal constrains are investigated. Apparently, the rotational deformation in case 1 shows discontinuity, while case 2 is continuous as the opposite situation. Based on above specific configurations, it is able to construct the definite governing matrices according to (8a) to (19c) and (28a) and (28b). Herein, is employed in the study.

For A-method, the holds(1)case (2)case

And for G-method, it holds(1)case (2)case

During actual computations, (21) and (23) are applied to get modal frequencies and modal shapes. To solve the big eigenmatrix of (22), the proposed algorithm strategy is applied. Meanwhile, the gridding scheme is also applied to the FEM.

As presented in Figure 3, the frequency plots of FEM and A-method at and are depicted, respectively, in which three slender ratios are incorporated as the representatives of general short to long beams [1]. Results show that, no matter in case 1 or case 2, A-method shares same tendencies with FEM in which high coherence is observed in various . Figure 4 gives the corresponding modal shapes through normalization in each dimension, where only the natural displacement components are illustrated. It is observed that perfect coincidence appears except for the 1st and 3rd modes in case 1, demonstrating the overall correctness of the A-method.

Although the comparisons with FEM obtain satisfactory results, the high costs of A-method in computation, where the computation under even occupies tremendous costs, limits its application. Therefore, a total comparison is made between the G-method and the A-method to explore the potential substitutability. Figure 5 displays the outcome of the four G-methods based on polynomial, Chebyshev, Lagrange, and trigonometric functions, respectively, in which the computations consider a set of symmetric BC positions (, and , ) and various and are utilized. The comparisons within case 1, as exhibited in Figures 5(a) and 5(b), show there exist apparent errors between the A-method and the G-method. Wherein the 1st frequency generated by G-method lies between the 2nd and the 3rd modes calculated through A-method, which can be interpreted as, although the trial functions used in (13a), (13b), (13c), and (13d) already hold , the modal sectional angles at hinge-joint BCs ( and ) as demonstrated in Figure 6(b) inherit discontinuities, which exceeds the general ability of (13a), (13b), (13c), and (13d) and results in imprecise results. In case 2, the results turn better, as the displacement and the rotation at fix-hinge BCs transform smoothly (as shown in Figures 6(c) and 6(d)), and the roots of 1st frequency except for T-method fall in the adjacent range of exact value uniformly as pictured in Figures 5(c) and 5(d). The quantitative computing effect of G-method can be further measured through correlation coefficients relative to A-method [31].

Besides, it is undeniable that the slender ratio is regarded as an essential parameter in beam system [1, 32]. Hence, a large-scale variations about from 1 to 100 are carried out. Overall, we capture in Figure 7 that the curves of frequencies present irregular changes due to the discontinuity in case 1, while the trends are relatively uniform in case 2. Detailed information observed in Figures 7(a) and 7(b) shows that the numerical methods show surging changes in the plots within , demonstrating the numerical methods is not absolutely perfect for all slenderness and there still exist certain algorithm instabilities. For case 2, as illustrated in Figure 7(d), the C-method shows break-offs when which indicates the P-method and the L-method are more stable than the C-method and the T-method. Actually, the entire curve tendencies in Figure 7 demonstrate that the numerical precisions are significantly dependent on , where bad computing effect occurs within low , while in large the consequences of G-methods are fairly acceptable. Further, the curves inherit closing regions near to in case 1, after which the solving effect of G-methods slightly decreases. In case 2, although the G-method and the A-method share the consistent trends, the solving effect is indeed affected by . On the whole, we conclude the P-method and the L-method possess higher flexibility. It is also indicated that the G-methods can be deemed with satisfied precision compared with the A-method when .

5.3. BCs’ Optimization Based on Case 2

Since the numerical method can predict the general frequency tendencies with reasonable and low-consumed properties, a further optimization based on case 2 is carried out in this part. Actually, the objective of the optimization conducted here is to obtain the optimal BCs proposed in Section 4 to achieve the minimum energy reservation in actual structural design. The A-method and the G-method are also applied together as a contrast. It is undeniable that the explicit expression of (46) is hard to get even for case 2 with only 2 fixed-hinge BCs. On the one hand, although the step performed in (21) can obtain the characteristic equations with variables and , to solve the frequency is almost impossible with two unknown parameters. On the other hand, the derivations of (46) about and have introduced a set of implicit characteristic equations, which are considered to be super complicated. Thus, a numerical strategy is put forward here on the base of Section 4.

Algorithm scheme: as shown in Figure 8, we traverse and in as small as possible step simultaneously to get reliable outcome during numerations and form the initial database. On such basis, a fitting procedure is drawn through Matlab’s griddata function (herein, means 4 grid points spline functions interpolation). In this process, mesh grids (1 million combinations in total) are generated for and , after which the minimum energy interpolating point is sought out. Sequentially, the dualistic algebraic interpolations are expanded at this point within the chosen grid area denoted as [33]. Afterwards, a procedure of derivations in (46) about and is made and the two-dimension equations are conducted to get the final optimal BC values. Actually, excessive integrals must be solved at each and combination during the traversal process to get the specific energy, which is considered to be a computationally expensive procedure. The balance of the traversal step must be considered due to the massive computing time, especially for A-method.

To conduct the numeration, we set and . To evaluate the numerical methods, the reservation of the first 20 modal energies in A-method is regarded as the accurate solution. According to the computing performance, we choose 0.05 in A-method, 0.01 in P-method, 0.05 in C-method, 0.02 in L-method, and 0.02 in T-method as the inherent steps for or , respectively. Figures 9 and 10 are the frequency spectrums and the energy maps of the A-method and the G-methods, in which the frequencies and the energy clouds are in keeping with each other. It is found that the BCs at higher frequencies tend to absorb more external energies generally in which, although the peaks and troughs are irregular, the G-methods perform good accordance to a certain degree. On the whole, we observe from Figures 9 and 10 that the P-method and the L-method present good “Mirror” properties, and the C-method is slightly inferior, while the T-method performs worse. The final optimal BCs of the two methods are listed in Table 1 in which the errors indicate the G-methods can effectively predict the optimal positions, demonstrating the substitutability of the numerical method.

6. Conclusion

Based on Timoshenko’s beam theory, this paper studies the vibrational performance of a continuous beam subjected to arbitrary types and numbers of constrains. Based on that, the optimization about boundary conditions, established on the minimum energy absorbed principle under loading, is proposed. During the derivations, analytical method is employed at first based on segmented strategy and the numerical Galerkin method is also adopted for comparison and as an alternative approach, in which a new type of trial functions is constructed. Further, the analytical equations are derived for boundary optimization considering the external loads. In case studies, two typical numerical examples representing the discontinuous and the continuous situations, respectively, are introduced. The results indicate that the analytical method is valid compared with the FEM but costly indeed. The comparisons between the analytical method and the Galerkin methods demonstrate the precisions of the numerical method are higher in continuous case, since the discontinuity in mode shapes hinders the entire simulation accuracy. It is also verified the Galerkin method is beneficial in the case of . Furthermore, a practical optimization is conducted based on case 2, in which a suitable computing strategy is put forward and the effectiveness of the G-methods are compared with A-method, and the optimization is finally demonstrated to be meaningful.

Appendix

A. Fundamental Solution

Herein, we consider the general solution of the beam with constraints on both sides, the typical equations areFor free harmonic vibration of mode , the modal response can write asWe directly give the fundamental solution of (A.1a) and (A.1b) under modeHere

B. Explicit Matrixes of Boundary Conditions

(1)(2), , , , (3)

C. Boundary Conditions in Trial Functions

The nature conditions of trial functions in Galerkin method are as follows. For free boundaries ( or ), the trail functions must yieldFor simply supported boundaries ( or )For clamped boundaries ( or )For simply supported constraints ()For hinge-joint constraints For clamped constraints

D. Explicit Matrixes of Time Domain Equation

Here, , , , and , in which

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

This work was supported in part by the National Key R&D Program of China under Grant no. 2016YFB1200100.