# Effect of Method Type on the Response of Continuum Vibro-Impact

Authors: Hongtao Wei; Gang Li; Pan Guo; Jun Zhao
Journal: Shock and Vibration
 (2019)

---

## Abstract

The force integration method (FIM) and the methods based on the mode transfer principle are the frequently used choices for solving the vibro-impact problem of continuum. Notably, there are different solving options such as numerical and semianalytical procedures. In this paper, a new modeling method based on the mode transfer principle, called the relative mode transfer method (RMTM), is proposed, and its semianalytical solution is obtained. A typical vibro-impact problem of continuum beam is studied. The time-history response under periodic excitation is obtained using the FIM, mode transfer method (MTM) through numerical procedure, and the RMTM through numerical procedure and semianalytical procedure. The effects of the method on the system steady-state response, amplitude-frequency response, and sticking motion are discussed. The impact type is divided into “hard impact” and “soft impact” when obtaining the results.


---

## Body

## 1. Introduction

The vibro-impact of continuum has always been a focus of nonlinear vibration research. Due to manufacturing errors and tolerances that occur when connecting through connectors, impact vibrations as a result of contact separation are present in a large number of practical mechanical systems, which leads to complex dynamics such as resonance peak shift, jump, sticking, and rising motion of the collision point. Many of these phenomena need to be avoided when the system is running. For example, in the ground vibration test of spacecraft, a resonance peak drift to low frequency caused by the gap in the connecting structure can lead to the coupling of the natural frequency and the rocket vibration frequency when the rocket is launching [1], which brings about a potential risk to the launch mission. In addition, the vibro-impact will lead the system to become unstable, produces additional wear and noise, and so on; thus, it can be said that the impact brings a lot of harm to the system. However, the vibro-impact can also be useful, especially in some emerging applications such as the tapping motion of microcantilever in an atomic force microscope [2] or in an energy sink device with impact vibration to improve collection efficiency [3]. Overall, the body of research on the impact vibration of continuum is growing.The key to this problem is the precise modeling and solution of the vibro-impact process. There are three main modeling methods that have been reported in the past literature, which are based on force integration [4–6], mode transfer [7–10], and velocity restitution coefficient [11–14], respectively. Furthermore, according to whether or not they ignore the impact process time, the above methods can be divided into the following two categories:(1)
Discontinuous analysis method: the method based on velocity restitution coefficient, which is called the coefficient of restitution (CoR), belongs to this category. The method considers the contact time of impact to be zero, which is suitable for large stiffness impact vibrations. It was first widely used for single-degree-of-freedom impact vibrations and then later extended to multidegree of freedom and continuum impact vibration research.(2)
Continuous analysis method: both the force integration method (FIM) and the mode transfer method (MTM) techniques belong to this category. Both methods consider the process of impact. FIM introduces the reaction force generated by impact into the equations without changing the initial mode of the continuum, while MTM regards impact as a change of boundary conditions, thus introducing a new mode to the governing equations. From this point of view, CoR and FIM adopt a unified mode to solve the equations in the whole process and so they are generally can be said of the same type.In general, numerical procedures are frequently used to solve the governing equations generated by the methods mentioned above, and some steps to improve the precision and efficiency of calculation have been proposed in the literature aiming at the nonlinearity of the vibro-impact problem [15–17]. There are few methods available to date that can solve the problem analytically, and FIM cannot obtain the analytical solution because it needs to use the displacement of the endpoint of the current step in the solution process.At present, much of the research on continuum vibro-impact is focused on the nonlinear dynamic response. Especially, the recent research on the “sticking motion” between beam and stop has become a new hotspot. In the impact vibration, the chattering phenomenon often occurs between the impact point and the stop, and the phenomenon can be divided into two types, namely, complete chattering and incomplete chattering. During one period, the complete chattering ends with an infinite amount of impacts in a finite time that eventually produces sticking motion. The incomplete chattering includes only finite number of impacts and does not produce sticking motion [18]. Toulemonde and Gontier were the first to study this phenomenon systematically, based on the P-C method; the characteristics of the sticking motion of the single-degree-of-freedom and multidegree-of-freedom systems were investigated [19]. Separately, Vyasarayani et al. used the CoR method to study the sticking motion in the vibro-impact of a beam and deduced the Lagrange multiplier when the system is in the sticking stage [14].In the past, the research on the influential factors of sticking motion has generally focused on the system parameters. For example, Li et al. used the CoR method to study the impact vibration of the driving foot of the linear ultrasonic motor, obtained the factors leading to sticking motion in the system, including the coefficient of restitution, pretightening force, and pretightening force spring stiffness, and offered the suggestion to reduce the sticking motion [20]. Wagg discussed the phenomenon of sticking motion in a series of studies on vibration systems with two degrees of freedom, including the influence of velocity recovery coefficient on the frequency range of sticking [21], rising phenomenon [22], and slip bifurcation phenomenon in sticking motion [23]. Nordmark and Piiroinen studied the chattering of the vibro-impact system, and the local discontinuous mapping is proposed [18]. Additionally, Luo et al. studied the relationship between dynamic response and system parameters for 2 DOF systems, focusing on the effect of excitation frequency and gap size on the response and discussed the “chattering-sticking” motion of the system [24]; he also studied the “chattering-sticking” motion of the traditional gear system and discussed the occurrence and disappearance of sticking motion, as well as the maintenance condition of the sticking motion [25]; meanwhile, Li et al. evaluated the control of sticking motion by means of the impulse control method, which shows that this method can effectively suppress sticking motion but also introduces sticking motion for systems with specific parameters [26]; Luo et al. studied the “chattering-sticking” motion of a vibro-impact system with large dissipation and discussed the relationship between the coefficient of restitution and the sticking motion characteristics [27].Although the solving procedures do affect the solution results, up until now, there has been no relevant research on the influence of the calculation method on the vibro-impact response of continuum. Janin and Lamarque contrasted the results of several numerical methods and theoretical solutions for single-degree-of-freedom vibro-impact systems and discussed the higher order convergence of these methods as well as the computational costs [28]. According to Tsai and Wu, the solution speed of the FIM is lower than that of the MTM, and, with an increase of the stiffness, the time for the solution of the FIM will increase as the transition from “soft impact” to “hard impact” is made. However, it is still considered that the FIM has the advantage of simple processing steps [7].In this paper, first, a semianalytic method based on the mode transfer principle is proposed. This method is then subsequently used to study the specific vibro-impact system of a cantilever beam, and the difference between the responses of the system obtained by numerical and then analytical procedures under the conditions of “soft impact” and “hard impact” are discussed using both the new method and the FIM. Finally, the influence of different methods on the specific nonlinear dynamic response of “sticking motion” is also investigated.


## 2. Modeling

Figure1 shows the cantilever beam impact model; the clearance size between the beam end and the stop is Δ. The stop is approximated as linear spring and the base movement is w0t. First, the beam is modeled under the condition of noncontact vibration. Vibration equations exist as follows:(1)ρAw¨x,t+EIw″″x,t=−ρAw¨0t,wL,t<Δ,where · denotes derivation with respect to time variable t and ′ denotes derivation with respect to spatial variable x. wx,t is the transverse deformation of the beam relative to the clamped end, ρ is the mass density of the beam, A is the cross section area, E is Young’s modulus, and I is the area moment of inertia. At this point, the boundary conditions are(2)wx,tx=0=0,w′x,tx=0=0,w″x,tx=L=0,w‴x,tx=L=0.Figure 1
Cantilever vibro-impact model.The displacement at any point on the beam can be discrete in the form of modal superposition as follows:(3)wx,t=∑i=1∞aitφix.In the formula above,φix is the ith mode function, and ait represents the ith modal coordinates. Substitute equation (3) to boundary conditions in equation (2) and you get(4)φxx=0=0,φ′xx=0=0,φ″xx=L=0,φ‴xx=L=0,where the mode function can be written as(5)φx=Asinkx+Bcoskx+Csinhkx+Dcoshkx.If you substitute equation (4) into equation (5) and let B=−1, the ith mode function φix can be written as(6)φnx=coshknx−cosknx−λnsinhknx−sinknx,where(7)λn=cosknL+coshknLsinknL+sinhknL.If you substitute equation (3) into equation (1), multiply both sides of the equation with φjx, and then integrate on [0, L] because of the orthogonality of modal modes, you obtain(8)ρA∫0Lφixφjxdx=mj,i=j,0,i≠j,EI∫0Lφixφj″″xdx=ωj2mj,i=j,0,i≠j,where ωj is the jth order natural frequency, and the driving force is in form of w¨0t=W0sinfet. There are discrete partial differential equations as follows:(9)a¨it+ωi2ait=bisinfet,i=1,…,N,where bi=W0fe2∫0Lφixdx/∫0Lφi2xdx. A linear damping is introduced into equation (9) and then the system equation can be(10)a¨it+2κa˙it+ωi2ait=bisinfet,i=1,⋯,N,where κ=ξωi, and ξ is the damping ratio. The initial values of the equations above can be used as follows:(11)ai0=∫0Lw0φixdx,a˙i0=∫0Lw˙0φixdx.In this way, the exact analytical solution of the system equation without contact can be obtained.


### 2.1. Relative Mode Transfer Method Description

In a vibration cycle, the system is assumed to start vibrating at rest, andt0 is the initial time. Transverse displacement at any point on the beam can be expressed as wx,t=w1x,t, as shown in Figure 2. Starting from t1 moment, the displacement of the endpoint is greater than the gap, and the system will be affected by the reaction of the endpoint stop. Eventually, the boundary conditions change so that the mode shape of the governing equations is changed too and the transverse displacement of any point on the beam can be regarded as the displacement superimposed on the displacement before the transition of the boundary condition, that is, wx,t=w1x,t1+w2x,t. By analogy to ti moment, ti is the time at which the boundary conditions of the system change and the transverse displacement of any point on the beam, which can be presented as follows:(12)wx,t=∑j=1i−1wjx,tj+wix,t,where(13)wix,t=∑n=1∞anitφnix,where φnix is the nth mode of the beam under this boundary condition in the ti−1∼ti period and anit is the corresponding modal coordinate. The form of φnix is consistent with equation (5), and method to compute its explicit expression can be found in reference [29]. While substituting equation (12) into the governing equations of motion, the equations of period ti−1∼ti become(14)ρAw¨ix,t+EIwi″″x,t=−ρAw¨0t−EI∑j=1i−1wj″″x,tj.Figure 2
System vibro-impact process sketch map.Upon substituting equation (13) into equation (14), by using the orthogonality of the mode and adding damping, the discrete N-order partial differential equation can be obtained as follows:(15)a¨jt+2κa˙jt+ωj2ajt=bjw¨0t+∑r=1i−1∑n=1Nhnjranrtr,j=1,…,N,where tr is the point at which the boundary conditions change and hnjr=−EI∫0Lφnr″″xφjxdx/ρA∫0Lφj2xdx, bj=∫0Lφjxdx/∫0Lφj2xdx. Following state transition, the initial value of the vibration equation is wix,ti+=0; thus, the initial value of the modal displacements of each order is also 0, and the modal velocity is w˙ix,ti+=w˙ix,ti−. Subsequently, the expression of the modal velocity of the nth order can be determined as follows:(16)a˙nix,ti+=∑P=1∞a˙Pi−1ti−∫0LφPi−1xφnixdx∫0Lφni2xdx,where φPi−1x denotes the Pth modal mode of ti− moment and a˙Pi−1ti− is the corresponding modal velocity.The biggest difference between the RMTM and the classical MTM [7] is that the expressions of displacement at arbitrary points of the beam before and after impact are different. The basic idea of the RMTM is to record the “form” of the current continuum system and introduce the concept of “relativity” each time mode transfer takes place. The vibration of the transformed system is regarded as an “incremental” movement relative to the form before the conversion, which is also the origin of the name of the RMTM. Using this idea, the expression of system displacement in equation (12) is obtained. However, the classical MTM only records the modal displacements before the conversion and makes use of the recorded values when the next conversion happens. Considering the system shown in Figure 2, at t1 moment system shifted from the cantilever beam state to the constrained beam state and at t2 moment system shifted from the constrained beam state to the cantilever beam state. When the RMTM is used, the displacement expressions before and after the state transition are expressed by formula (12). When the classical MTM is adopted, the modal coordinates ajt1− of each order at t1− moment are recorded at first and then system changed to a new vibration mode at moment t1+ and the initial values of the modal coordinates of each order are set to zero as follows:(17)ajt1+=0,j=1,…,N.Whent2 moment system reenters the cantilever beam state from the restricted motion state, the modal coordinate values recorded before the last state transition are back-substituted as follows:(18)ajt2+=ajt1−,j=1,…,N.The above steps are no longer applicable to more complicated boundary conditions such as hysteretic nonlinear boundary conditions and boundary conditions which can be discretized into piecewise linear ones, but the RMTM can still be adopted to solve these problems.


## 2.1. Relative Mode Transfer Method Description

In a vibration cycle, the system is assumed to start vibrating at rest, andt0 is the initial time. Transverse displacement at any point on the beam can be expressed as wx,t=w1x,t, as shown in Figure 2. Starting from t1 moment, the displacement of the endpoint is greater than the gap, and the system will be affected by the reaction of the endpoint stop. Eventually, the boundary conditions change so that the mode shape of the governing equations is changed too and the transverse displacement of any point on the beam can be regarded as the displacement superimposed on the displacement before the transition of the boundary condition, that is, wx,t=w1x,t1+w2x,t. By analogy to ti moment, ti is the time at which the boundary conditions of the system change and the transverse displacement of any point on the beam, which can be presented as follows:(12)wx,t=∑j=1i−1wjx,tj+wix,t,where(13)wix,t=∑n=1∞anitφnix,where φnix is the nth mode of the beam under this boundary condition in the ti−1∼ti period and anit is the corresponding modal coordinate. The form of φnix is consistent with equation (5), and method to compute its explicit expression can be found in reference [29]. While substituting equation (12) into the governing equations of motion, the equations of period ti−1∼ti become(14)ρAw¨ix,t+EIwi″″x,t=−ρAw¨0t−EI∑j=1i−1wj″″x,tj.Figure 2
System vibro-impact process sketch map.Upon substituting equation (13) into equation (14), by using the orthogonality of the mode and adding damping, the discrete N-order partial differential equation can be obtained as follows:(15)a¨jt+2κa˙jt+ωj2ajt=bjw¨0t+∑r=1i−1∑n=1Nhnjranrtr,j=1,…,N,where tr is the point at which the boundary conditions change and hnjr=−EI∫0Lφnr″″xφjxdx/ρA∫0Lφj2xdx, bj=∫0Lφjxdx/∫0Lφj2xdx. Following state transition, the initial value of the vibration equation is wix,ti+=0; thus, the initial value of the modal displacements of each order is also 0, and the modal velocity is w˙ix,ti+=w˙ix,ti−. Subsequently, the expression of the modal velocity of the nth order can be determined as follows:(16)a˙nix,ti+=∑P=1∞a˙Pi−1ti−∫0LφPi−1xφnixdx∫0Lφni2xdx,where φPi−1x denotes the Pth modal mode of ti− moment and a˙Pi−1ti− is the corresponding modal velocity.The biggest difference between the RMTM and the classical MTM [7] is that the expressions of displacement at arbitrary points of the beam before and after impact are different. The basic idea of the RMTM is to record the “form” of the current continuum system and introduce the concept of “relativity” each time mode transfer takes place. The vibration of the transformed system is regarded as an “incremental” movement relative to the form before the conversion, which is also the origin of the name of the RMTM. Using this idea, the expression of system displacement in equation (12) is obtained. However, the classical MTM only records the modal displacements before the conversion and makes use of the recorded values when the next conversion happens. Considering the system shown in Figure 2, at t1 moment system shifted from the cantilever beam state to the constrained beam state and at t2 moment system shifted from the constrained beam state to the cantilever beam state. When the RMTM is used, the displacement expressions before and after the state transition are expressed by formula (12). When the classical MTM is adopted, the modal coordinates ajt1− of each order at t1− moment are recorded at first and then system changed to a new vibration mode at moment t1+ and the initial values of the modal coordinates of each order are set to zero as follows:(17)ajt1+=0,j=1,…,N.Whent2 moment system reenters the cantilever beam state from the restricted motion state, the modal coordinate values recorded before the last state transition are back-substituted as follows:(18)ajt2+=ajt1−,j=1,…,N.The above steps are no longer applicable to more complicated boundary conditions such as hysteretic nonlinear boundary conditions and boundary conditions which can be discretized into piecewise linear ones, but the RMTM can still be adopted to solve these problems.


## 3. Semianalytical Solution Procedure of Relative Mode Transfer Method

### 3.1. Semianalytical Solution and Parameters

As shown in Figure2, the system is in contact with or detached from the stop at tm moments. The term ∑r=1i−1∑n=1Nhnjranrtr in formula (15) changes with each contact and separation, but the term is a constant term, and therefore the system governing the equations of motion is always linear so that the system can obtain an analytical solution between two switching moments.Without losing generality, the analytical solution expression of the system equation (15) can be written for a complete contact or separation period tm−1∼tm as follows:(19)at=a1t+a2t+a3t,where a1t is the initial response solution of the following:(20)a¨jt+2κa˙jt+ωj2ajt=0,j=1,…,N,which is caused by the initial velocity and displacement of the equation introduced by the system after the contact or separation conversion. a2t is the zero initial value response of equation (21):(21)a¨jt+2κa˙jt+ωj2ajt=bjw¨0t,j=1,…,N.Meanwhile,a3t is the zero initial response of the following:(22)a¨jt+2κa˙jt+ωj2ajt=∑r=1i−1∑n=1Nhnjranrtr,j=1,…,N.In the process of solving, a single step-forward strategy is adopted, setting the step size toΔt. Subsequently, at tm−1 moment, the initial value of displacement of equation (20) is zero; the initial value of velocity is obtained by equation (16); and the initial values of equations (21) and (22) are all zero.For equation (20), when t∈tm−1tm, set the initial value of the order i equation to tm−1=t0, atm−1=x0, and a˙tm−1=x˙0. Unlike in Figure 2, since t0 represents the time value of the tm−1 moment, we have a general solution:(23)ait=e−ntC1coswdt+C2sinwdt.Substituting the initial value into formula (23), one can get C1 and C2 in order to solve the problem with time t increases without obtaining infinite numbers in the results, and the following treatment is imposed:(24)C1⋅e−nt0=x0−e−nt0⋅C2⋅sinwdt0coswdt0,C2⋅e−nt0=x˙0+x0wd⋅tgwdt0+n⋅coswdt0wd.Equation (23) can be written as(25)ait=e−nt−t0C1⋅e−nt0⋅coswdt+C2⋅e−nt0⋅sinwdt.In the same way, the solution of equation (21) can be found, that is,(26)ait=e−ntC1coswdt+C2sinwdt+A⋅sinwt−α.If one substitutes the initial value, we get(27)C1⋅e−nt0=x0−Asinwt0−α−e−nt0⋅C2⋅sinwdt0coswdt0,C2⋅e−nt0=x˙0−Awcoswt0−α+x0−Asinwt0−α⋅wd⋅tgwdt0+n⋅coswdt0wd.Equation (26) can be written as(28)ait=e−nt−t0C1⋅e−nt0⋅coswdt+C2⋅e−nt0⋅sinwdt+A⋅sinwt−α.The solution of equation (22) is(29)ait=−dωn2+e−nt−t0C1⋅e−nt0⋅coswdt+C2⋅e−nt0⋅sinwdt.After the initial values have been substituted, the following can be obtained:(30)C1⋅e−nt0=x0+d/ωn2−e−nt0⋅C2⋅sinwdt0coswdt0,C2⋅e−nt0=x˙0+x0+d1/ωn2⋅wd⋅tgwdt0+n⋅coswdt0wd.After the solutions of each independent equation are obtained, the total solution of the system can be secured by the combination of formula (19), and then the displacement of any point can be obtained by formula (12).


### 3.2. Solution Process and Zero-Crossing Detection

For a given beam vibro-impact system, the coefficients in formulas (10) and (15) are calculated, step Δt is set, and the displacement of any point on the beam is obtained by the solution process above and formula (12).In the process of solving, the displacement of the beam end of each step is always monitored, and state switching occurs once it crosses the zero point. In order to guarantee the efficiency of the solution and ensure that the difference between the displacement of the beam end and the zero is small enough at the zero-crossing point, the strategy of changing the calculation step around the zero-crossing point is adopted. Here, the zero-crossing moment is determined using a dichotomous step, as shown in Figure3.Figure 3
Zero-crossing detection using bisection method.In Figure3, it is assumed that, within the tl,tu step, the beam end displacement passes through zero, while the error is wtu,L−0>εy with εy being the error threshold. The calculation needs to be refined within this zero-crossing time step, with the time starting point still being tl, the last calculated step being halved, and the error estimated. If it is still greater than the threshold, then from tl start, the step size will be halved again, and the cycle will continue until the zero-crossing point tu is found to meet the error. Then, by changing the parameters in the analytical expression (15), the starting point of the calculation after zero-crossing is set to tu and the initial velocity value of each component equation is obtained from (16). The step is reset to the initial stride length to continue the calculation until the next zero-crossing and the beam endpoint displacement response in the specified time can be obtained through repeating the above steps. The semianalytical solution flow chart is shown in Figure 4.Figure 4
Flow chart of solving process.


## 3.1. Semianalytical Solution and Parameters

As shown in Figure2, the system is in contact with or detached from the stop at tm moments. The term ∑r=1i−1∑n=1Nhnjranrtr in formula (15) changes with each contact and separation, but the term is a constant term, and therefore the system governing the equations of motion is always linear so that the system can obtain an analytical solution between two switching moments.Without losing generality, the analytical solution expression of the system equation (15) can be written for a complete contact or separation period tm−1∼tm as follows:(19)at=a1t+a2t+a3t,where a1t is the initial response solution of the following:(20)a¨jt+2κa˙jt+ωj2ajt=0,j=1,…,N,which is caused by the initial velocity and displacement of the equation introduced by the system after the contact or separation conversion. a2t is the zero initial value response of equation (21):(21)a¨jt+2κa˙jt+ωj2ajt=bjw¨0t,j=1,…,N.Meanwhile,a3t is the zero initial response of the following:(22)a¨jt+2κa˙jt+ωj2ajt=∑r=1i−1∑n=1Nhnjranrtr,j=1,…,N.In the process of solving, a single step-forward strategy is adopted, setting the step size toΔt. Subsequently, at tm−1 moment, the initial value of displacement of equation (20) is zero; the initial value of velocity is obtained by equation (16); and the initial values of equations (21) and (22) are all zero.For equation (20), when t∈tm−1tm, set the initial value of the order i equation to tm−1=t0, atm−1=x0, and a˙tm−1=x˙0. Unlike in Figure 2, since t0 represents the time value of the tm−1 moment, we have a general solution:(23)ait=e−ntC1coswdt+C2sinwdt.Substituting the initial value into formula (23), one can get C1 and C2 in order to solve the problem with time t increases without obtaining infinite numbers in the results, and the following treatment is imposed:(24)C1⋅e−nt0=x0−e−nt0⋅C2⋅sinwdt0coswdt0,C2⋅e−nt0=x˙0+x0wd⋅tgwdt0+n⋅coswdt0wd.Equation (23) can be written as(25)ait=e−nt−t0C1⋅e−nt0⋅coswdt+C2⋅e−nt0⋅sinwdt.In the same way, the solution of equation (21) can be found, that is,(26)ait=e−ntC1coswdt+C2sinwdt+A⋅sinwt−α.If one substitutes the initial value, we get(27)C1⋅e−nt0=x0−Asinwt0−α−e−nt0⋅C2⋅sinwdt0coswdt0,C2⋅e−nt0=x˙0−Awcoswt0−α+x0−Asinwt0−α⋅wd⋅tgwdt0+n⋅coswdt0wd.Equation (26) can be written as(28)ait=e−nt−t0C1⋅e−nt0⋅coswdt+C2⋅e−nt0⋅sinwdt+A⋅sinwt−α.The solution of equation (22) is(29)ait=−dωn2+e−nt−t0C1⋅e−nt0⋅coswdt+C2⋅e−nt0⋅sinwdt.After the initial values have been substituted, the following can be obtained:(30)C1⋅e−nt0=x0+d/ωn2−e−nt0⋅C2⋅sinwdt0coswdt0,C2⋅e−nt0=x˙0+x0+d1/ωn2⋅wd⋅tgwdt0+n⋅coswdt0wd.After the solutions of each independent equation are obtained, the total solution of the system can be secured by the combination of formula (19), and then the displacement of any point can be obtained by formula (12).


## 3.2. Solution Process and Zero-Crossing Detection

For a given beam vibro-impact system, the coefficients in formulas (10) and (15) are calculated, step Δt is set, and the displacement of any point on the beam is obtained by the solution process above and formula (12).In the process of solving, the displacement of the beam end of each step is always monitored, and state switching occurs once it crosses the zero point. In order to guarantee the efficiency of the solution and ensure that the difference between the displacement of the beam end and the zero is small enough at the zero-crossing point, the strategy of changing the calculation step around the zero-crossing point is adopted. Here, the zero-crossing moment is determined using a dichotomous step, as shown in Figure3.Figure 3
Zero-crossing detection using bisection method.In Figure3, it is assumed that, within the tl,tu step, the beam end displacement passes through zero, while the error is wtu,L−0>εy with εy being the error threshold. The calculation needs to be refined within this zero-crossing time step, with the time starting point still being tl, the last calculated step being halved, and the error estimated. If it is still greater than the threshold, then from tl start, the step size will be halved again, and the cycle will continue until the zero-crossing point tu is found to meet the error. Then, by changing the parameters in the analytical expression (15), the starting point of the calculation after zero-crossing is set to tu and the initial velocity value of each component equation is obtained from (16). The step is reset to the initial stride length to continue the calculation until the next zero-crossing and the beam endpoint displacement response in the specified time can be obtained through repeating the above steps. The semianalytical solution flow chart is shown in Figure 4.Figure 4
Flow chart of solving process.


## 4. Case Studies

Considering the system parameters in Table1, the FIM, MTM, and RMTM are adopted and both the numerical solution and semianalytical solution processing steps are used to obtain results. In the previous work, it was found that the results obtained by different methods are sensitive to the type of impact. The type of impact is defined by the parameter [13] k∗=k/3EI/L3, where k is the stop stiffness. When k∗=5 or 5000, the system is in the state of “soft impact” or “hard impact,” respectively. Parametric studies are carried out here to compare the differences in the results obtained by different solutions steps, with different methods for different impact types.Table 1
System parameters.
Symbol
Quality
Value
ρ
A
kg
/
m
Density
0.4649
L (m)
Length
0.2580
A (m2)
Area of the cross section
5.960 × 10−5
I
kg
⋅
m
2
Area moment of inertia
1.9867 × 10−11
E (N/m2)
Modulus
2.1 × 1011
Δ
m
Gap size
0
ξ
Damping ratio
0.02
W
0
Amplitude of excitation
5.05 × 10−4Firstly, the results of the endpoint response are compared between MTM numerical solution and RMTM numerical solution, and the following Figure5 is a comparison diagram of the results of 5-order and 20-order discrete with typical exciting frequency and different impact types. Line y=0 is on behalf of the impact line.Figure 5
Time-history response of soft impact and hard impact using MTM and RMTM. (a) Five-mode discretization: damping ratioξ=0.02, ω=60Hz, soft impact k∗=13k=1×104N/m, zero-crossing threshold εy<1×10−14, and step Δt=1×10−5. (b) Twenty-mode discretization: damping ratio ξ=0.02, ω=30Hz, hard impact k∗=13000k=1×107N/m, zero-crossing threshold εy<1×10−14, and step Δt=1×10−5.
(a)
(b)MTM numerical represents the numerical solution of the MTM, and RMTM numerical represents the numerical solution of the RMTM. Under the condition of different impact types, MTM and RMTM have good consistency when the solution is periodic both during noncontact and contact phase. It is important to note that when it is hard impact, the convergence of RMTM numerical is better than MTM numerical, as shown in Table2.Table 2
Mode number needed for convergence using different methods.
Exciting frequencyω (Hz)
FIM
RMTM analytical
RMTM numerical
MTM numerical
30
13
5
5
14
80
13
9
9
13
120
11
8
8
16Figure6 is the time-history response under the conditions of soft impact and hard impact with different modes and different methods. FIM represents force integral method and RMTM analytical represents the semianalytical solution of the RMTM. It can be seen from Figure 6(a) to 6(d) that the FIM and the RMTM have slightly different steady-state amplitudes under the conditions of different impact types and different number of modes. Using the results of the FIM as a reference, the difference of the steady-state amplitude obtained by the two methods is 8.23% in the case of the five modes discretized and soft impact type, while there are two peaks in one period with the minimum difference being 12.41% in the case of hard impact. The difference of the steady-state amplitude of the soft impact is 8.41% in the case of the 20-mode discretization, and the maximum difference between the two peaks of hard stiffness is 3.8%. Morphologically, the time-history responses obtained from the five-order and 20-order discretization are almost identical in the case of soft impact. When a hard impact occurs, the responses are consistent before the contact, while the responses after the contact show a greater difference. When five-order discretization is used, as shown in Figure 6(b), the response of FIM after contact is completely different from that of RMTM, while, when 20-order discretization is used, the response of FIM and RMTM after contact tends to be consistent. It can be seen from the above discussion that when the impact type is soft impact, both methods can be solved with less mode discretization; meanwhile, when the impact type is hard impact, more modes, such as 20 modes, must be used for discrete processing.Figure 6
Time-history response of soft impact and hard impact using FIM and RMTM. (a) Five-mode discretization: damping ratioξ=0.02, ω=60Hz; soft impact k∗=1.3k=1×103N/m; zero-crossing threshold εy<1×10−14; and step Δt=1×10−5. (b) Five-mode discretization: damping ratio ξ=0.02; ω=100Hz; hard impact k∗=13000k=1×107N/m; zero-crossing threshold εy<1×10−14; and step Δt=1×10−5. (c) Twenty-mode discretization: damping ratio ξ=0.02; ω=60Hz; soft impact k∗=1.3k=1×103N/m; zero-crossing threshold εy<1×10−14; and step Δt=1×10−5. (d) Twenty-mode discretization: damping ratio ξ=0.02, ω=70Hz; hard impact k∗=13000k=1×107N/m; FIM, RMTM analytical, and RMTM numerical I: zero-crossing threshold εy<1×10−14 and step Δt=1×10−5. RMTM numerical II: zero-crossing threshold εy<1×10−16 and step Δt=1×10−6.
(a)
(b)
(c)
(d)It should also be noted that, as shown in Figure6(d), RMTM numerical I uses a zero-crossing threshold εy<1×10−14 and step Δt=1×10−5, and the contacting response curve is significantly different from the contacting response using the analytical solving process. When a smaller zero-crossing threshold εy<1×10−16 and smaller step Δt=1×10−6 are used, the response curve of RMTM numerical II from the numerical solution almost coincides with the semianalytical solution curve. This shows that the semianalytical solution has the advantage of precision and speed in comparison with the numerical solution. At the same time, compared with the results of FIM and RMTM, MTM has higher consistency with RMTM, which is essentially because MTM and RMTM are both based on the principle of mode switching.Figure7 is the frequency spectrum obtained by different methods. The frequency spectrum is obtained from time-history response data through fast Fourier transform (FFT). The steady-state response data are taken when applying the FFT. It can be seen that under the same parameter settings, the frequency responses of different methods are only slightly different in terms of amplitude. Moreover, the frequency response of the semianalytical solution and the numerical solution of the RMTM are almost identical. Whether with hard or soft impact, the effect of mode number on response in the frequency domain can be almost ignored. As compared with the soft impact response in the frequency domain, the hard impact has a more abundant super-harmonic response, that is, there are some smaller wave peaks at the 60Hz octave frequency band. In addition, in the case of hard impact, the amplitude of the fundamental wave is larger, while in the case of soft impact, the energy is concentrated on the response of subharmonic 30Hz.Figure 7
Frequency spectrum of soft impact and hard impact using different methods. (a) Twenty-mode discretization: damping ratioξ=0.02, ω=60Hz, and soft impact k∗=1.3k=1×103N/m. (b) Twenty-mode discretization: damping ratio ξ=0.02, ω=90Hz, and hard impact k∗=13000k=1×107N/m.
(a)
(b)Figure8 shows the amplitude-frequency response of the system when different methods are used under soft impact condition. Using 20-mode discretization, the peak value of the system after steady state is taken as the amplitude. As can be seen in the figure, the amplitude-frequency response obtained by different methods is almost identical.Figure 8
Amplitude-frequency response of beam end. Twenty-mode discretization: damping ratioξ=0.02 and soft impact k∗=1.3k=1×103N/m.In soft impact, the time-history response of the system does not reflect multiperiodic motion; meanwhile, in hard impact, the system exhibits multiperiodic and quasiperiodic motion. Let the sweep frequency step be1Hz, and the local maximum value in the positive value of the time-history response after steady state is taken as the ordinate value. The bifurcation diagram of the results obtained by different methods is drawn as Figure 9.Figure 9
Bifurcation diagram with twenty-mode discretization: damping ratioξ=0.02 and hard impact k∗=13000k=1×107N/m. (a) Bifurcation diagram using FIM. (b) Bifurcation diagram using RMTM.
(a)
(b)The system bifurcations in Figures9(a) and 9(b) show complex quasiperiodic motions in the 101Hz−106Hz, 156Hz−171Hz, and 180Hz−200Hz intervals. In this state, the time-history responses of the FIM and the RMTM are completely different. It can also be seen that the bifurcation graphs obtained by the two methods are identical in shape and differ only in amplitude. The results of the numerical solution of the MTM and RMTM are consistent with the analytical solution when the step size and zero-crossing error are smaller. Therefore, the results of MTM are not discussed here.The convergence characteristics of different methods are studied. Considering the hard impact condition, the effect of mode number on convergence is discussed. Take RMTM analytical as example. The time history response using different discrete modes with exciting frequencyω=80Hz are shown in the following Figure 10. It should be noted that the solving process encountered problem in finding proper zero-crossing point when adopting 6-order discretization, and the simulation is not sustainable.Figure 10
RMTM analytical time-history response:ω=80Hz, ξ=0.02, and step=1×10−5. (a) Time-history. (b) Local amplification of the rectangle in (a).
(a)
(b)As shown in Figure10(b), when the discrete mode number increases, the curves converge gradually. Let PN stand for peak value of the steady-state response using N-order discretization. If PN−PN−1/PN−1×100%<1% and PN+1−PN/PN×100%<1% meet at the same time, then it is concluded the solution converges using N−1 modes. Using the step above and substituting the results in Figure 10, the solution converges when 9-order discretization is adopted.The results of the other methods adopting the same approach are listed in Table2. Specific frequencies when system experiences periodic motion are selected. As can be seen in Table 2, RMTM has the advantage of convergence on FIM and MTM.The computational efficiency of different methods is studied, through comparing the time required for different methods to simulate 0.5 s when the exciting frequencies are different under the condition of hard impact. At the same time, step length is1×105, and the 5 modes and 20 modes discretization are both adopted. It is important to note that different computer configurations and different memory states will affect the simulation time, and the time given in Table 3 may differ in different simulation environments.Table 3
Time consumption using different methods.
Exciting frequencyω (Hz)
Time (s), FIM-20 modes
Time (s), RMTM analytical-20 modes
Time (s), RMTM numerical-20 modes
Time (s), RMTM analytical-5 modes
Time (s), MTM numerical-20 modes
30
103.031468
136.740955
123.876582
37.812161
175.688769
60
108.737112
161.979584
148.081940
51.771183
248.205443
90
103.618053
147.316289
137.270971
40.824425
204.288178
120
105.332696
190.255675
177.958377
48.344781
203.100090
150
100.305788
130.061885
111.004888
44.388146
164.340308
180
101.289414
98.305552
75.175407
50.640532
104.304119As can be seen from the Table3, when using 20-order discrete, FIM has the shortest simulation time; in this point of view, the simulation efficiency of FIM is higher than RMTM and MTM. However, considering that RMTM analytical has better convergence and accuracy, RMTM can obtain the same result as the 20-order mode discretization of FIM by using 5-order discretization.The sticking motion is a hot topic in the discussion of the recent vibro-impact issue. The effects of two treatment methods and mode number on the sticking motion of the system are discussed as shown in Figure11.Figure 11
Effect of mode number and methods on sticking motion:ω=95Hz, damping ratio ξ=0.02, and hard impact k∗=13000k=1×107N/m. (a) Time-history. (b) Contacting vibration, local amplification of the rectangle in (a).
(a)
(b)As can be seen in Figure11, when the stiffness, damping, and other parameters are fixed, using five-mode discretization, the results from FIM just represent finite chattering after the endpoint contacts the stop, while the results obtained using RMTM show a direct sticking motion coupled with no chattering, as shown in Figure 10(b); when the discrete mode number increased to 20, the results from FIM demonstrate sticking motion, and the motion pattern is consistent with that of RMTM, which is different from the sticking motion pattern depicted in reference [18]. In fact, the “sticking-chattering” pattern obtained here is just the inverse type of the “chattering-sticking” in [18]. A possible reason may be that the study model in this work is continuum beam which is different from the model in [18].


## 5. Conclusions

In dealing with the vibration problem of continuum with piecewise linear boundary conditions, the method of force integration is simple and widely used, but it can only be solved by numerical procedure. In this paper, by introducing the concept of relative displacement, a method based on mode transfer is proposed and a semianalytical solution step is obtained.The stop stiffness and the number of discrete modes are important parameters of the impact system and the solving procedure. Different methods and numerical and semianalytical solution steps are used to obtain and compare the time response, frequency response, and bifurcation diagram of specific parameters. The following can be concluded:(i)
In soft impact, the response of the system is relatively simple, and there is no complicated multiperiodic or quasiperiodic motion. No matter whether the method of force integration or the method based on mode transfer is adopted, a consistent result can be obtained using a small number of modes; there is only a small difference in amplitude. Therefore, one can choose the method which is most convenient to deal with when soft impact is present and so use fewer modes to solve the problem.(ii)
In contrast, when the system is hard impact, it is normally considered that such can satisfy the requirement of solving precision using five-mode discretization, but actually it can be seen from comparing the response of contact phase above, which is much different from that seen with more discrete modes such as the 20-mode version. Therefore, the conclusion is that five-mode discretization is insufficient for hard impact.(iii)
By comparing the results of different methods, the findings obtained by the force integration method and the method based on the mode transfer principle are almost identical under the conditions that enough modes and high precision are used. However, the difference of response amplitude and convergence rate between the two methods still exists under certain parameters. The RMTM has the advantage of convergence rate compared with FIM and MTM.(iv)
In this paper, a new “sticking-chattering” contact vibration pattern is obtained by using several different methods for the vibro-impact model of continuous beam. The number of discretization modes and the solution method have an effect on the sticking motion, and since sticking motion is an important response characteristic of hard impact, more modes need to be adopted in the study to make the results obtained convincing.



---

*Source: 2718502-2019-04-22.xml*