Comparison of 1D and 2D solutions for a beam under transverse impact

The problem of non-stationary vibration of an elastic beam caused by a transverse impact loading is studied in this work. In particular, two different approaches to the derivation of analytical solution of the problem are compared. The first one is based on the Timoshenko beam theory, the latter one follows the exact two-dimensional theory. Both mentioned methods are used for finding the response of an infinite homogeneous isotropic beam. The obtained analytical results are then compared and their agreement is discussed in relation to main factors, i.e. the beam geometry, the character of loading and times and points at which the beams responses are studied.


Introduction
Beams and beam-like structures are used in many engineering applications across the entire spectrum of industry. There exists a lot of simplified beam theories which have been developed in order to effectively investigate beams responses under different loading conditions. The basic assumptions under which these theories can be applied to static or dynamic beams problems are different. In the case of static problems, the ratio of beam length and its height plays the main role.
If we focus on dynamic problems of beams, there exist two main approaches: modal and wave. The modal method, making the use of eigen frequencies and eigen modes, is usually applied to vibration problems in which the maximal frequencies are of orders about ones or tens of kHz. This method is sometimes used also for problems of higher frequencies but it is usually associated with the problem of inaccurate determination of eigen values. On the other hand, the wave approach is applicable for nearly arbitrary range of frequencies but it becomes quite inefficient when used for lower frequency problems.
As can be seen from most of the published works dealing with problems of thin beams, the Bernoulli-Euler theory [1] is often used when the modal method is applied. In the case of wave approach, this elementary theory is inappropriate because it results in infinite phase velocity for increasing frequency. Another beam theory developed by Rayleigh in 1894 and taking into account the effect of rotary inertia [1] leads to finite phase velocity, but it is still limited to thin beams with negligible influence of shear forces. This effect has been added by Timoshenko in the early 20th century [2]. According to his theory, two characteristic wave velocities exist and they are very close e-mail: vadamek@kme.zcu.cz to those predicted by the exact two-dimensional formulation. Unfortunately, the Timoshenko theory still neglects the warping of the beam cross-section so it can be applied to thick beams problems only under certain assumptions. This restriction is removed by so called higher-order beam theories developed after the mid-20th century, e.g. [3]- [6]. Such refined theories enable to solve more complicated problems of beams, e.g. beams made of functionally graded materials, laminated beams, nanocomposite beams etc.
In this study, the application of the Timoshenko beam theory to impact problems of beams is discussed from the wave propagation point of view. The analytical solution for an elastic infinite beam under transverse impact obtained using the Timoshenko theory is compared to that derived using the exact two-dimensional theory. The measure of compliance of both solutions is discussed in relation to the beam geometry, to the Fourier spectra of applied load and to spatio-temporal parameters which the both solutions are compared for.

Formulation and assumptions of the problem
We will assume an infinite thin beam with uniform rectangular cross-section of dimensions b 0 × d 0 (width × height). The area of the beam cross-section will be denoted as A = b 0 d 0 in the following. The material of the beam will be isotropic and homogeneous and its elastic properties will be described by the modulus of elasticity E, shear modulus G and by the density ρ. The beam is in a rest at time t = 0 s. After that an external load perpendicular to the beam longitudinal axis x is applied to the upper edge of the beam. This load will be defined by the function q(x, t) MATEC Web of Conferences 148, 05005 (2018) https://doi.org/10.1051/matecconf/201814805005 ICoEV 2017 which will be considered as an even function of the coordinate x. With respect to this assumption, the problem will be solved as a symmetric one on the spatio-temporal domain 0, +∞) × 0, t max . Further, the bottom edge of the beam will be free of load and the beam deflection and the slope of the beam cross-section tend to zero for x → +∞. Finally, small displacements and rotations will be taken into account in the following sections.

Governing equations
As mentioned in Section 1, the Timoshenko beam theory takes into account the effect of shear forces and rotary inertia to the beam deformation. The equations of motion in this case can be formulated as [1] where α(x, t) and v(x, t) denote the slope of the beam crosssection and the beam deflection, respectively, and J stands for the cross-sectional-area moment of inertia about the neutral axis. The functions Q(x, t) and M(x, t) denote the shear force and the bending moment present in the beam cross-section, respectively, and there hold The constant κ is the so-called Timoshenko's shear coefficient introduced to quantify the variable distribution of shear force in the cross-section of the beam. Its value is κ = 5/6 0.833 for the cross-section of a rectangular shape (see [1]).

Solution in Laplace domain
The solution of the system (1) can be found using the method of integral transforms, namely using the combination of Laplace transform in time and Fourier transform in space. The final formulas for the Laplace transforms of the unknown functions α(x, t) and v(x, t) can be derived from the formulas presented in [7]. In the mentioned work, the analytical solution for the problem of a simply supported viscoelastic beam under transverse dynamic loading is presented for the case of special orthotropy. If we do that we can write the solution of 1D problem as where the new variables p ∈ C and ω ∈ R have been introduced by the application of Laplace and Fourier transforms, respectively. The spectra of the Fourier integrals present in (3) can be expressed as complex rational functions in this case and there hold The complex polynomials P i (ω, p) in (4) have the form and the function Q(ω, p) represent the Fourier cosine integral of the Laplace transform of q(x, t).

Analytical solution -exact 2D theory 4.1 Governing equations
When the exact 2D theory of continuum is used for solving beam bending problems, the state of plain stress is assumed. In this case, the beam is considered as a thin 2D strip and the beam deflection and the slope of the beam cross-section can be derived from the components of displacement. Let the functions u(x, y, t) and v(x, y, t) denote the horizontal and vertical displacement component, respectively. Both are the functions of the horizontal coordinate x ∈ 0, +∞), the vertical coordinate y ∈ −d 0 /2, d 0 /2 and time t. Introducing the function of volume dilatation ∆(x, y, t) and rotation r(x, y, t) by relations the equations of motion for a 2D element of the strip (beam) can be written as an uncoupled system of PDEs for ∆(x, y, t) and r(x, y, t). This system of equations can be derived from the equations of motion presented in [8] for the case of an infinite viscoelastic strip. Doing so one obtains the following system: where the constants c 2 and c 3 represent the phase velocities of shear and dilatational waves, respectively, and there hold The symbol ν stands in (8) for the Poisson's ratio.

Solution in Laplace domain
The solution of the system (7) can be simply derived in Laplace domain using the solution given in [8]. Neglecting the material viscosity, the Laplace transform of displacement components can be expressed as Fourier integrals having the following form: where the positive sign of the vertical displacement transformv(x, y, p) corresponds to the positive value of the beam deflection transformv(x, p) from Section 3.2. The symbols sh and ch used in (9) represent for the function sinh and cosh, respectively. The unknown complex functions A 1 (ω, p), A 2 (ω, p), B 1 (ω, p) and B 2 (ω, p) follows from the application of boundary conditions of the problem. Similarly to [8] one can write The remaining complex functions λ i (ω, p) and k i (ω, p) present in (9) and (10) can be simply derived using the analogous functions given in [8].

Evaluation procedure
The evaluation process which leads to the results in spatiotemporal domain consists of two basic steps: (i) to evaluate the Fourier integrals contained in the derived formulas (3) and (9) and (ii) to perform the inverse Laplace transform (ILT). It is possible and advantageous to change the order of these steps and to perform the Laplace inversion of derived Fourier spectra first. The inversion can be done analytically or numerically. The first approach, which is usually based on the residue theorem, is feasible and quite straightforward in the case of 1D solution. For the 2D problem, this procedure is much more complicated and so the application of a suitable numerical algorithm is more effective. To guarantee the same evaluation conditions for both solutions compared, the numerical approach to ILT is used in this work. In particular, an algorithm based on the combination of FFT and Wynn's epsilon accelerator was used for this purpose. This algorithm is described in detail in [9]. It is precise and quite robust for solving different problems of elastodynamics, as shown in the work mentioned. The evaluation of the Fourier integrals has been done also numerically. The classical Simpson rule with constant integration step ∆ω has been used. The choice of ∆ω and of the upper limit of integration ω max depended on the value t max . If t max < 200 µs, the parameters ∆ω = 1 and ω max = 4000 were used. For greater values of t max , the choice was ∆ω = 0.25 and ω max = 8000 to guarantee high precision of obtained results.
In the following subsection, the results obtained for steel beams with material parameters E = 2.1 × 10 11 Pa, ν = 0.3 and ρ = 7800 kg/m 3 are presented. The effect of the beam height d 0 , the character of the external loading q(x, t) and the time and the position of observation are discussed in relation to the differences between both solutions.

Results comparison and discussion
The first comparison will be made for two beams with different heights d 0 = 2 mm and d 0 = 10 mm. Both beams are subjected to the load which is nonzero only for x ∈ 0, h and is defined as where σ a = 1 MPa, b 0 = 2 mm, h = 1 mm and t 0 = 2 µs. It means that the beam is excited by a pressure rectangular pulse lasting 2 µs and the amplitude of which is distributed along the x direction according to cosine function on a domain of total length 2 mm. If we compare the deflection histories v(t) at points x = {10, 50, 100} mm lying at beams longitudinal axes up to time 100 µs we find out that 1D and 2D approach lead for both values of d 0 to nearly identical results. To quantify the measure of this agreement, a relative error will be used. This error will be defined as The functions v 1D (t) and v 2D (t) in (12) represent the time histories of beam deflections determined by 1D and 2D theory, respectively, and the maximum in the denominator is looked for over the whole time interval 0, t max . Calculating δ(t) for all mentioned x one finds that it reaches its maximal value 5.8% in the case d 0 = 2 mm. When d 0 = 10 mm, the error is greater as expected. In this case the maximum error is about 19%. But it must be mentioned that this error is obtained only locally for time about 35 µs and for x = 100 mm. The comparison of vertical beam velocities dv(t)/dt is more important because it gives us deeper insight into the physics described by both approaches. This comparison is made in Fig. 1 and Fig. 2 for the same values of x as aforementioned. It is evident that the vertical beam velocity determined by the Timoshenko theory (thick solid curves) represents, let us say, a mean value around which the 2D solution oscillates (thin dashed curves). The amplitudes of oscillations for d 0 = 2 mm are quite small (see Fig. 1 the vertical velocity. On the other hand, the greater value of d 0 leads to considerable differences between both approaches, see Fig. 2. It could be expected because also in static beam problems the higher height of the beam means worse compliance of 1D and 2D solution. But in the elastodynamics of beams this fact is primarily connected with the wavelength of the waves propagated through the beam, as will be discussed later.   Fig. 2 correspond to real wave processes in the beam. For instance, the constant time-spacing between every two significant peaks corresponds to time during which the Swave propagates from the beam axis to the beam edge and back. To make a better picture of the wave propagation in the beam, the spatial distributions of total velocity in 100 mm long part of one half of the beam (strip) are used in Fig. 3 to show the positions of wave fronts at different times. Firstly, the P-wave and S-wave are prop-agating from the place of loading (denoted by the arrow) through the strip. At time t 1.6 µs, the front of P-wave is nearly reaching the lower strip boundary whereas the front of slower S-wave is located near the beam axis. At time t 2.7 µs, i.e. approximately 0.7 µs after the impact, the loading S-wave and reflected P-wave are interacting and the unloading P-and S-wave are propagated from the place of excitation. The process of waves interaction continues in time and the zone of transient state of stress moves farther from the vertical beam axis, see the state at t 10.6 µs. When t = 100 µs, the state in the observed part of the strip becomes quite homogeneous, the transient phenomena gradually disappear, so the amplitudes of oscillations of 2D solution around 1D results are decreasing, see Fig. 2. Another fact which we meet if we compare both theories is that the 1D solution derived using the Timoshenko theory is delayed in comparison with that obtained using the exact 2D approach. This time lag is obvious from Fig. 4 where the history of v(t) for x = 300 mm is depicted for the case d 0 = 2 mm. This is not an error but this is the consequence of different dispersion properties of both types of continuum being compared. In the case of Timoshenko beam, there exist only two dispersion curves [1] with the asymptotes corresponding to two characteristic wave velocities c 0 = E/ρ and c 2T = κG/ρ. On the other hand, the beam modelled as a 2D strip has an infinite number of dispersion curves and two different characteristic velocities c 3 and c 2 given by the relations (8). Calculating all of these velocities for parameters mentioned in Section 3.1 and 5.1, one obtains: c 0 = 5189 m/s vs. c 3 = 5439 m/s and c 2T = 2937 m/s vs. c 2 = 3218 m/s. The next factor which significantly influences the compliance of both types of results is the maximal frequency contained in the amplitude spectrum of the time part of q(x, t). From the elastodynamics point of view, a beam can be considered as thin only if the wavelength λ of the shortest wave present in the signal propagating through the beam is substantially larger than the beam height. For instance, when the Bernoulli-Euler beam theory is applied it must be true that k 0 /λ < 0.05, where k 0 is the cross-   . If we apply this criterion to the beam with rectangular crosssection b 0 × d 0 we approximately obtain the inequality λ < 5.77d 0 . It means that an Euler beam can be considered as thin only if its height is about six times smaller than the wavelength of the wave component with the highest frequency.
This fact will be demonstrated using two rectangular pulses of different duration. Let us assume that the external impact loading (11) takes 40 µs instead of 2 µs, as considered previously. It is clear that the high frequencies have much smaller contribution to the amplitude response than in the case of t 0 = 2 µs. If we compare the histories of the beam deflection v(t) determined for t 0 = 40 µs by 1D and 2D approach, we find out that the maximal value of δ(t) for all x = {10, 50, 100} mm and for d 0 = 10 mm is  only 0.7% contrary to 19% in the case t 0 = 2 µs, as mentioned above. Also the comparison of vertical velocities dv(t)/dt shows better agreement between 1D and 2D solutions, as clear from Fig. 5. The time-spacing of peaks in the 2D solution is naturally same as in the case of short pulse but the amplitudes of its oscillations around 1D results are significantly smaller then in Fig. 2.
All of the comparisons presented previously have been done for beams of quite small dimensions and for very short times. To present the behaviour of both solutions for a more realistic case, the transient response of a beam with cross-sectional parameters b 0 = 20 mm and d 0 = 200 mm to the pulse (11) with t 0 = 40 µs has been investigated as well. Fig. 6 and Fig. 7 show the comparison of 1D (thick solid curves) and 2D (thin dashed curves) solutions of v(t) and dv(t)/dt, respectively, for x = {1, 2, 3} m up to MATEC Web of Conferences 148, 05005 (2018) https://doi.org/10.1051/matecconf/201814805005 ICoEV 2017 time t max = 1 ms. It is clear from Fig. 6 that although the pulse takes 40 µs, the beam height is too large to satisfy the aforementioned criterion, so the agreement between the solutions v(t) is not so good as in the previous cases. However, the Timoshenko beam theory gives us deflections the character of which is very similar to those followed from the exact 2D theory. Unfortunately, the same cannot be said in the case of dv(t)/dt, see Fig. 7. It is clear from this figure that the height of the beam is too high compared to the wavelength of the shortest wave to use the Timoshenko theory for the investigation of vertical beam velocity.

Conclusions
This work presents the comparison of two types of analytical solutions for the problem of transverse impact on an infinite homogeneous elastic beam. Based on the analyses of the results obtained for different geometric and loading parameters, the main factors influencing the agreement between 1D Timoshenko theory and 2D exact theory are discussed from the wave propagation point of view. It is shown that the Timoshenko theory is applicable also to impact problems of beams but the credibility of its results significantly depends on the frequencies of waves propagated through the beam. This publication was supported by the project LO1506 of the Czech Ministry of Education, Youth and Sports.