Dynamic Analysis of Encastre Beams by Modification of the System ’ s Stiffness Distribution

DOI: http://dx.doi.org/10.24018/ejers.2020.5.11.2220 Vol 5 | Issue 11 | November 2020 1 Abstract — Numerical and energy methods are used to dynamically analyze beams and complex structures. Hamilton’s principle gives exact results but cannot be easily applied in frames and complex structures. Lagrange’s equations can easily be applied in complex structures by lumping the continuous masses at selected nodes. However, this would alter the mass distribution of the system, thus introducing errors in the results of the dynamic analysis. This error can be corrected by making a corresponding modification in the systems’ stiffness matrix. This was achieved by simulating a beam with uniformly distributed mass with the force equilibrium equations. The lumped mass structures were simulated with the equations of motion. The continuous systems were analyzed using the Hamilton’s principle and the vector of nodal forces {P} causing vibration obtained. The nodal forces and displacements were then substituted into the equations of motion to obtain the modified stiffness values as functions of a set of stiffness modification factors. When the stiffness distribution of the system was modified by means of these stiffness modification factors, it was possible to predict accurately the natural fundamental frequencies of the lumped mass encastre beam irrespective of the position or number of lumped masses.


I. INTRODUCTION
Most real loads on structures are dynamic in nature [1], [2]. To make the analysis simpler the static equivalent of the dynamic load is used [3]- [5]. Dynamic loads cause structures that have mass and elasticity to vibrate [6]. When the exciting force is not continuous, the structure vibrates at one or more of its natural frequencies. The natural frequencies of a structure are perculiar to it. They are dependent on the mass distribution and the stiffness matrix of the structure [7].
The number of independent motions a structure can make is known as its degrees of freedom. Real structures have infinite degrees of freedom. It is however possible to analyse a structure with respect to only few degrees of freedom. This makes the analysis easier. A reduction in the structure's degrees of freedom can be achieved in a number of ways. One such way is by mass lumping. A more popular way of doing this is by breaking the structure into discrete elements and then using the finite element method. The finite element method has been largely used in structural analysis because of its ability to model complex structures. Its results are not however accurate as errors often arise from poorly known boundary conditions and oversimplification of very complex structural systems [8], [9]. While mass lumping and the finite element method both simplify the structure's mass distribution, mass lumping does so by ignoring the rotational inertia of the masses and treats them as point masses. The finite element method uses a consistency matrix to approximate the mass distribution. The consistency matrix is generated using a chosen shape function. The accuracy of the finite element method is improved upon by dividing the structure into smaller elements or by carefully selecting better shape functions [10], [11]. Another method that makes use of shape functions is the Rayleigh-ritz method. While the Rayleighritz method predicts more accurately the lower natural frequency of vibration, the finite element method gets a better prediction for higher frequencies of vibration [12]. The partial differential equations for the vibration of structures can be obtained using an energy method known as the Hamilton's principle. This is an extension of the principle of least work to accommodate dynamics. Its major drawback is that with the exception of very simple structures formulation of such equations often leads to very complicated mathematics. If, however the continuous mass of the system can be grouped as lumped/point masses, it is easy to analyze the structure using an energy method known as the Lagrange's equations. These equations enable the analysis of structures consisting of lumped masses connected by mass-less elements [13]. With a proper selection of representative masses, the results can be very close to the exact response. When the discretization of the continuous mass is increased by the use of a higher number of lumped masses, the accuracy of the dynamic response improves. Literature is replete with representations of continuous systems with lumped masses connected by massless elements [14]. Lin and Tsai [15] modeled a free vibration analysis of a uniform multi-span beam carrying multiple spring mass systems. A numerical model that simulates the vertical dynamic interaction between a lumped massed vehicle and a rail track has also been developed [16]. Its usefulness at offering a simplified solution was exemplified in the works of Katherine Reichi and Daniel Inman when they modeled a meta-structure as a onedimensional lumped mass [17]. The results from such analysis are usually approximate [18]. In the lumping of masses effort is made to preserve the total mass of the structure even though this alone cannot ensure the quality of the solution [19]. The use of lumped masses (diagonal inertia matrices) simplifies the program coding and results in significant reduction in computer memory and computational effort [20]. Mass lumping of finite elements @ @ @ has been applied in acoustics through the use of spectral elements [21]. It has been extended to the dynamic finite element analysis of plates and shells [22]. It is however a general belief that consistent matrix lead to a more accurate solution [23]. The diagonal mass matrix is still employed because of its lesser computational effort. The processes of mathematical modeling and numerical solutions have been transformed by the availability of modeling software [24]. Efforts have also been made at developing simplified models for predicting the natural frequencies of lumped mass structures. The concept of using shear wave in a solid prismatic bar was used by Sule [25] in studying and obtaining approximate natural frequency of vibration of beams under lateral vibration. Osadebe [26] also proposed a model for calculating the natural frequencies of some lumped mass beams after an attempt on finding effective masses for the system. The use of experimental methods through measured data to adjust the analytical mass and synthesize the mass matrix has also been explored [27], [28]. We observe that the approximate methods such as Lagrange's equations (by mass lumping), Rayleigh-ritz and the finite element methods all modify the mass distribution of the structure. While in the use of the Lagrange's equations on lumped masses, a simple redistribution on the basis of the centre of gravity can be done. The Rayleigh-ritz and finite element methods carry out a more complex mass formulation resulting in the consistency matrix. The effort of these methods is to obtain a mass distribution (inertia matrix) representative of the real structure and so produce exact results. In the use of the Lagrange's equations the masses are lumped. Efforts have been made at finding better patterns of lumping that will minimize errors due to the lumping of a continuous mass. Using evenly spaced lumpsand increasing the number of lumps is some of the results of such studies. We can see that all the approximate methods hover around the distribution of the mass of the structure; we therefore need to explore the possibility of modification in the structure's stiffness distribution.

A. Mathematical theory
The Lagrange's equations provide a way of analyzing multi-degree of freedom system, a similar approach for continuous structures is the Hamilton's principle.
From this principle the partial differential equation governing the free lateral vibration of a beam can be obtained as [29]: where EI is the flexural rigidity of the beam and µ is its mass per unit length. ̈ is the second time derivative of the lateral displacement u. is the fourth derivative of the lateral displacement u with respect to x (an axis parallel to the longitudinal axis of the beam).
By assuming a harmonic vibration and applying the boundary conditions of an encastre beam we obtain the equation of the mode shapes as: The values of βjL for the first seven modes of vibration are as presented in Table 1. By taking care of the initial conditions the general solution by mode superposition is therefore [30]: where the constants Aj and Bj can be determinedfrom the initial conditions using the expressions: The eigenfunctions also satisfy certain orthogonality relationships. Mj is the generalized mass for the j th mode of vibration [12].

II. METHODOLOGY
The parameters that determine the vibration of structural systems are the structure's mass distribution and the structure's stiffness [31], [32]. They can be seen in the structure's inertia matrix and stiffness matrix respectively. The roles they play can be seen in the equations of motion of a vibrating system or the structural dynamics' eigenvalue problem. If the j th mode shape and the j th circular frequency are kept constant, then any variation in mass distribution will have a corresponding change in the element flexural rigidity EI.
Two equations were compared. One is the force equilibrium equation written as: (where the force vector {P} acts at the element's nodes), {F} is the vector of fixed end forces generated when nodal displacements are restrained.
[k] is the element stiffness matrix and {D} a vector of nodal displacements [33]. The second is the equation of motion of a vibrating system written simply as: where the external force vector {P} acts at the element's nodes, [m] is the inertia matrix, [ ] is the element stiffness matrix and {u} a vector of nodal displacements u1 to u4. Equations (9) can be applied in structural dynamics if the equations for the vector of fixed end moments/forces {F} can be formulated. The continuous system was analyzed using the Hamilton's principle and the equations for the fixed end forces {F} and nodal displacements {D} formulated for any arbitrary segment of a vibrating beam at time t = 0. This was then substituted into equation (9) to get the vector of nodal force {P} that is causing the vibration.
The stiffness matrix [Kd] in equations (10) was taken as the stiffness matrix of the lumped-mass system. If a vibrating element of the continous system (beam with continuous mass) and that of a corresponding element of a lumped-mass system (beam with lumped masses) are to be equivalent then their deformation must be equal and the force acting on their nodes {P} will also be equal. Therefore: A segment of the encastre beam shownin Fig. 1 is being restrained by the fixed end forces F1 to F4. The reduced structure or basic system of the segment is shown in Figure  1b.
From equation (5) the acceleration at any point in the vibrating beam is given by mode superposition as: From the D'Alembert principle the force on the vibrating beam can be calculated from its inertia force. For an elementary part of the beam at a distance z from the origin this force is̈.
Using the principle of virtual work and the flexibility method we determine the fixed end forces F1 to F4 of the isolated element of the excited beam to be: To be able to evaluate these equations, there is need to derive an expression for Aj. and Mj.

A. Derivation of the expression for Aj for an encastre beam
Consider a uniform encastre beam under the action of its self weight. The equation of the bending moment at any distance x from the left support is: where µ is the mass per unit length of the beam and g is the acceleration due to gravity.
By substituting equation (22) into (23) and considering the boundary conditions we obtain the initial deflection of the beam under its self weight as: where a is dimensionless constant equal to 3 24 ⁄ . Equation (24) is substituted into Equation (6)and (8) to obtain an expression for Aj and Mj.
These equations enable us to calculate the fixed end forces F1 to F4 on a segment of a vibrating beam and consequently obtain the nodal forces {P}. An arbitrary segment of a vibrating element is identified by means of the normalized distances 1 and 2 of its nodes from an origin.
EI is the flexural rigidity of the beam. l is the length of the segment and is equal to x2 -x1 but since the length of the beam l has been normalized and is now equal to unity. L x u1, u2, u3 and u4 are the total nodal displacements at the coordinates of the forces F1, F2, F3 and F4 respectively (see Fig. 1). The total displacement is obtained by totaling the displacements due to all the modes of vibration.
If a segment of a vibrating beam is isolated, it will be in equilibrium with the application of the force vector {P}. The force vector {P} represents the effect of the removed adjourning elements on the isolated segment. Figure 2a shows a segment of the vibrating continuous beam. The nodal forces on the beam P1, P2, P3and P4are calculated from the equilibrium equations (30 -33). When the continuous beam is represented by a lumped mass beam (a beam that has its distributed masses lumped at selected nodes), the equivalent segment of the beam is shown in Figure 2b. Just like the real segment the equivalent segment is supported by the same nodal forces P1, P2, P3and P4and has the same nodal displacements as the continuous/real beam. This implies that for the lumped mass beam to be equivalent to the real or continuous beam they must have the same forces and displacements at the nodes. ( 1 , 0) ( 2 , 0) numerical demonstration of these steps is presented below (see Table 2). For clarity the calculations will be presented in a tabular form. Using the methods presented in Table 2 the values of stiffness modification factors at different values of 1 and 2 for the lateral vibration of the beam can be obtained.  The values in the tables were written at 7 significant figures except for the summary which are at 14 decimal points to capture the values that are near zero.
Using the methods presented in Table 2 the values of stiffness modification factors at different values of 1 and 2 for the lateral vibration of the beam can be obtained.

III. NUMERICAL APPLICATION AND DISCUSSION OF RESULTS
The lumped mass beam of Fig. 3a (if rigidly restrained at both ends) was analyzed at different values of L1 and L2 to obtain the values of the fundamental frequency. The analysis was done first with the Lagrange's equations, and then by using the method presented in Table 2, the stiffness factors were obtained and used with the Lagrange's equations to get an improved result. Finally, the beam was analysed using the finite element method. Their results are presented in Table 3. Table 3 shows that the values of the calculated frequency vary with the relative values of L1 and L2. When Lagrange's equations were used without stiffness modification factors the obtained values differed considerably. The difference between the highest and the lowest values was 48.121. This was obtained when L1 = 1/8, 7/8 and 1/2. The result from finite element analysis was better even though it exhibited the same trend. The difference between the highest and lowest obtained values was 9.298 and where obtained at same positions as that from the use of the finite element method. The application of stiffness modification factors gave approximately the exact result of 22.373.  The trend is better appreciated by taking a look at the plot of the fundamental frequency against L1/L in Fig. 4.
The plots in Fig. 4 are all symmetrical about L1/L = ½. This is expected because the lumping of the mass of the beam is symmetrical about the beam centre. The fixed supports are also of equal distance from the beam centre. From the plots of Fig. 4 we observe that it is possible to obtain a good approximate of the exact value of the fundamental frequency by a careful selection of the position of the lumped mass. The calculated values of the fundamental frequency were best near L1/L = 3/8 and 5/8 for the Lagrange's equations and L1/L = ½ for the finite element method. A single factor analysis of variance was used to compare the values obtained from the use of the Lagrange's equation (with stiffness modifications factors) and that from the finite element method. The results gave Fcrit = 4.747, F = 6.440 and P-value = 0.026. Since Fcrit< F and P-value < 0.05, we accept that there was significant improvement in the results obtained by the application of the modification factors.
By increasing the number of lumps to two as shown in Figure 3b, we now have two independent variables L1 and L2. The value of L3 can be obtained the moment the values of L1 and L2 are specified hence it is not independent. The results of the dynamic analysis of the beam of Fig. 3b is presented in Table 4.  Table 4 it would be seen that the application of the Lagrange's equations with stiffness modification factors gave the exact value of the fundamental frequency. It is also worthy of note that the value obtained was not dependent on the relative values of L1, L2 and L3. The values of the fundamental frequency ω obtained with the deployment of Lagrange's equations without the stiffness modification factors however varied with the relative values of L1 and L2. This same observation was also noticed in the results from a finite element analysis. The effect of the relative variation of L1 and L2 on the fundamental frequency was studied by carrying out a two-way analysis of variance ANOVA on the values in table 3 for L1/L and L2/L of 1/9 to 4/9. The result was presented in Table 5. In carrying out the ANOVA the rows represented L1/L and the columns the different values of L2/L.From Table 5, it would be seen that neither L1/L nor L2/L was found significant. Their p-values were 0.412 and 0.746 respectively and which are more than 0.05. This suggests that an interaction between them may be largely responsible for any observed variation in the values of the fundamental frequencies obtained. The interaction plots for L1/L and L2/L (see Fig. 5) however showed consideration interactions. The lines L1/L = 1/9, 2/9, 3/9 and 4/9 were not parallel and intercepted at some points. The graphs revealed that it is possible to obtain the exact value of the fundamental frequency by a careful selection of the values of L1 and L2. This  The lines in Fig. 6 are not parallel and intercept at some points. However the values of L2/L that gave the best approximation of the fundamental frequency was in the neighbourhood of 2/9 at L1/L = 3/9 and 4/9. At higher values of L2/L the results deviated more from the exact value.
The introduction of more lumps in Fig. 3c led to similar results. But unlike in Fig. 3b, the number of independent variables rose to three. They are L1/L, L2/L and L3/L. The value of L4/L is defined once the values of the others are known. A similar interaction plot of the type shown in Fig. 5 would result in bogus surface plots. Hence, we decided to compute the average value for fundamental frequency obtained at different values of L1/L, L2/L, L3/L and L4/L. The plots are presented in Fig. 7.
The graphs for span L1 and Span L4 are the same due to symmetry. The same applies to graphs for span L2 and span L3. From the graphs we observe that the results from finite element analysis gave average fundamental frequencies that are close to the exact frequency irrespective of the values of the span lengths. The values of the span L1 and L4 for best fit were such that made L1/L and L4/L approximately ¼. For the application of Lagrange's equations on lumped masses the best fit was obtained when L1/L and L4/L are approximately 5/16. The response varied more and so depends more on the relative values of the spans L1, L2, L3 and L4.

IV. CONCLUSION
It would be observed that the natural frequencies obtained from the use of Lagrange's equations on the continuous system that had been represented as lumped masses had some measure of errors as seen from its comparison with the exact results. The error was more when the spacing between the lumped masses was higher and when the number of lumps was less. However, when the stiffness of the system was modified using the stiffness modification factors, the use of Lagrange's equations were able to predict accurately the fundamental frequencies hence their percentage errors were zero. The rigors involved in the manual calculation of the stiffness modification factors can be substantially reduced by writing a simple computer program for its execution.
From this work we can infer that: 1) An accurate dynamic response from a lumped mass encastre beam can be obtained by a modification in the stiffness composition of the system.
2) The values of the calculated stiffness modification factors obtained for each segment are not equal as shown in Table 1. This implies that no linear modification of the stiffness distribution of lumped mass encastre beams under lateral vibration can cause them to be dynamically equivalent to the continuous beams.
3) It is possible to obtain an exact response from the lumped mass beam through a careful selection of the position and spacing between lumps.
The improvement in the accuracy of the estimated fundamental frequencies did not depend on the location of mass lumps or the number of lumped masses. The stiffness modification factors helped to correct the error introduced in the system by mass discretization or lumping. This is a work on the modification of a structure's stiffness in order to improve the results from a dynamic analysis a continuous beam represented as a lumped mass structure. Encastre beams (fixed at both ends) were used as a first demonstration but it can be easily implemented in beams of other restraints (fixed-pinned, fixed-free etc) and similar equations generated and with corresponding stiffness modification factors. These will reveal the relationship between the modification factors and the end restraints and since the stiffness modification factors are generated with respect to segments of these beams we will observe that they depend on the restraint provided by adjoining elements and hence on their stiffness. These will lay the foundation for the analysis of frames simply by lumping their continuous masses at their joints and then using modification factors evaluated as a function of the joint rigidities of member elements to carry out dynamic analysis of the frame with no compromise on the accuracy of results obtained.