Vibration Response Evaluation under Shock-Type Loading with Emphasis on Finite Element Model Updating

In some cases, impulseor shock-type excitations as the dynamic loading are inevitable, and obtaining proper response with the wellknown numerical methods is not easy. (is paper focuses on dynamic response estimation against short-time loading with an updated finite element model using frequency response functions (FRF) and particle swarm optimization (PSO) technique. Because there is not an analytical method for assessing the numerical responses under shock-type excitations, in this paper, experimental tests are designed on a laboratory scale to evaluate the numerical responses. (e vibration responses of the system against shock loading are compared with the Newmark average acceleration scheme and also with experimental data. (e results reveal that the unconditionally stable Newmark method against regular loads has an appropriate performance. Still, under short-time loading, it faces numerical damping error, and this method should not be blindly applied under shock-type loads.


Introduction
Direct numerical integration technique is commonly used for response estimation of any structural systems. e effective parameters of any dynamic problems are mass, damping, and stiffness matrixes and also applied load in the single-degree-of-freedom (SDF) or multi-degrees-of-freedom (MDF) systems. Usually, these matrixes are obtained from finite element modeling (FEM) of any structural systems, and modeling errors include idealization in FEM, size of the finite element (FE) mesh, material properties, and inaccuracies in boundary conditions. To reduce the mentioned errors and obtain positive responses from numerical modeling, updating the mass and stiffness matrixes is a crucial process because those are inherent properties of the system and they affect structural responses and behavior. e modal parameters, e.g., natural frequencies, mode shapes, and also modal damping, are required for updating the numerical models, and obtaining these parameters is called system identification. Various useful approaches have been presented to extract the modal parameters of a structure in the time or frequency domain [1]. In addition to the numerical errors, the measurement errors occur because of the environmental effects, sensor errors, random errors, and errors in signal processing. As the modeling errors and measurement errors can lead to poor structural parameter identification, it is required to validate the model updating experimentally. Model updating methods simultaneously utilize the structural response obtained by the finite element method (FEM) and the measured structural response to calibrate mathematical modeling. However, in this research, it is assumed that the experimental data are clean due to the appropriate controls during the experimental tests.
Model updating could be classified into modal-based and response-based methods. e modal-based model updating technique relies on the modal characteristics data obtained from an experimental modal analysis that is indirectly extracted from the measured FRF data. D'Ambrogio and Sestieri [2] emphasize that the numerical procedures used for modal identification using experimental vibration data can introduce errors exceeding the level of required accuracy to update FE models. In response-based finite element model (FEM) updating methods, the measured FRF data are directly utilized to identify the unknown structural parameters. Most of the FRF-based model updating techniques are used to minimize a residual error between the analytical and experimental input force and output response [3][4][5][6][7]. Many researchers utilized vibration data of structures for finite element (FE) model updating. Baghchi [8] used vibration data for updating the mathematical model. Lin and Zho [9] updated the FE model using vibration data under base excitation. Pradhan and Modak [10] modified the RFM and presented a normal FRF concept to update mass and stiffness matrices of structure. For the damped system model updating, several investigations are conducted. Lin and Zho [11], Arora et al. [12], and Yuan and Yu [13] applied the FRF approach to update the FE model in a damped system based on vibration information. Garcia and Santini [14] updated damped structural systems with a two-step model-updating algorithm.
Many kinds of researches are focused on model updating from a mathematical point of view. Sipple and Sanayei [15] used numerical sensitivity techniques for model updating using FRF. ey used mathematical-based sensitivity method for gradient calculations. Wei-Ming and Jia-Zhen [16] updated the FE model based on model reduction and iterative approach. Weng et al. [17] applied the inverse substructure method for model updating. In their research studies, a system was divided into several substructures, and updating is applied for each subsystem. Papadimitrious and Papdioti [18] used the component mode synthesis technique for the FE model updating. Wan and Ren [19] utilized a residual-based Gaussian process for model updating. Sarmadi et al. [20] implemented an iterative least-squares minimal residual technique for model updating, and Wei et al. [21] updated multi-degree-of-freedom systems based on the intrinsic chirp component decomposition method.
Also, model-updating techniques are usually involved with optimization methods or algorithms. In this regard, there are much related researches. Ntotsios and Papadimitrio [22] applied the multiobjective optimization technique to FE model updating based on eigenvalue strain energy residuals. Christodoulou et al. [23] used the Pareto optimization technique for structural model updating; Jung and Kim [24] used a hybrid genetic algorithm for finite element model updating. Shabbir and Omenzetter [25] utilized particle swarm optimization with a sequential niche technique for dynamic finite element model updating.
Various researchers studied vibration responses of the systems under impact loading. Tagarielli et al. [26] predicted the dynamic response of composite sandwich beams under shock loading. Wali et al. [27] obtained the vibration response of the sandwich plate under low-velocity impact loading. Wang et al. [28] estimated the dynamic responses of clamped sandwich beams subjected to impact loading with analytical, numerical, and also experimental studies. Mazurkiewicz et al. [29] compared the numerical responses with experimental testing data of individual structural elements under impulse loading.
To the knowledge of the authors, numerical studies were not conducted to calculate the vibration responses of MDF systems with the commonly used Newmark method against impact or impulse loads [30]. erefore, an experimental test is designed on a laboratory scale to evaluate and validate the numerical responses of the Euler-Bernoulli-type beam with experimental data. us, vibration responses of a cantilever aluminum beam were studied with analytical, updated numerical and also experimental methods. For this purpose, at first, the analytical formulation of the free vibration of the beam was derived. en the numerical model with the finite element technique was created, and modal data were extracted from the FRF response. e beam was excited with a hammer, and the responses were measured by using the acceleration sensors at four reference points. e updating method was based on the particle swarm optimization (PSO) technique. After updating the finite element model, the dynamic responses of the beam against impact loading were evaluated with the Newmark method and experimental data.

Analytical Free Vibration Formulation.
In this section, the analytical free vibration formulations for modal data extraction are presented based on the differential equation of vibration of the Euler-Bernoulli-type beams. For this purpose, according to Figure 1, a cantilever beam with length L, a rectangular cross section with width b and height h, a concentrated mass M at its free end, and stiffness equal to EI is considered.
In the analytical method, by solving the equilibrium equation of the beam for free vibration, the frequencies and the mode shapes are determined. e dynamic equation of the Euler-Bernoulli beams with a uniform cross section in the forced-vibration case is obtained as follows: where ρA is equal to the mass of unit length or m and p(x, t) is the applied force to the arbitrary point. e parameter w(x, t) also refers to the beam displacement response equation along the z direction at any point in the time domain. For free vibration analysis, the right side of equation (1) should be equal to zero. As a result, equation (1) is rewritten as follows: To solve equation (2), the method of separation variable is used as follows: is equation is transformed into two separate equations in the time domain (t) and the space domain (x) as follows: e new variables are defined as follows: By solving equation (4), the responses are obtained according to the following equation: where z(t) is the time-domain response and ψ(x) is the eigenvector of the beam. e boundary and initial conditions of the beam are required to determine the coefficients G i and A, B. According to Figure 1 and equation (7), in the left boundary of the beam, the displacement and slope are zero: Based on equation (7), the mode shape is rewritten as follows: e boundary conditions of the cantilever beam with a concentrated mass M at the free end are as follows: e result EIw ″ (l, t) � 0 and EIw ‴ (l, t) � M € w(l, t) could be rewritten as equations (10) and (11), respectively: By replacing the values ψ(l), ψ ″ (l) into equation (11), the following equation is obtained: By replacing λ � βl, α � ml/Mλ, and β 4 /ω 2 � EI/m in equation (10), equation (12) is expressed as follows: To obtain the eigenvalues or frequencies (λ), the determinant of equation (13) should be zero as follows:  (14) is expressed as equation (15). By solving this equation, the eigenvalues (λ n ) are available. is nonlinear equation requires special numerical techniques to obtain all roots of eigenvalues [31]: According to the following equation, the frequency (ω n ) in each mode is available: To obtain the eigenvectors of the beam, one should calculate the coefficients G 1 and G 2 . In this study, G 1 is set to the unit, and the general unknown factor G n is computed from the following equation: erefore, the eigenvectors or equivalently mode shapes of the studied beam in Figure 1 is available according to the following equation:

Numerical Modal
Parameters. e finite element method was used to obtain the modal parameters in the numerical method. Firstly, using the finite element method, the studied beam in Figure 1 was divided into twenty-five equal parts with 3 centimeters length. Each element has two nodes, and each node has a transitional degree of freedom (w i ) and a rotation degree of freedom (θ i ). Stiffness matrix K e and mass matrix M e of each element are determined as follows: where E is the modulus of elasticity, I is the moment of inertia of the cross section, and l is the length of each element (l � L/N elem ), where L is the length of the beam, N elem is the number of elements, ρ is the mass density of the material, and A refers to the cross section of the beam. It should be noted that the last element of the beam has a concentrated mass (M) at the second node. erefore, during the assembling process of the total mass matrix, this concentrated mass should be added to the mass matrix as follows: e beam's length was equal to 75 cm, and the number of elements was 25. To compare the results and update the finite element model, in the second stage, the number of elements reduced to four elements.

Experimental Modal Parameters.
e aluminum beam, as mentioned in Section 2.1, was also studied in an experimental study. As shown in Figure 2, the acceleration data were collected by installing unidirectional acceleration sensors at four points of the beam and applying the hammer force at the specified positions.
According to Figure 3, the input data (forces) and output data (accelerations) are recorded at positions P1, P2, and P3, and they are collected at 4 stations (CH1, CH2, CH3, and CH4). e entire data collection process was 30 seconds. Each test is repeated three times. Based on the data sampling frequency of the data acquisition system (400 Hz), for each test, 12,000 points are recorded.
Extraction of modal parameters was accomplished using frequency response function (FRF). For example, FRFs; in the case that the hammer exerted at the free end (Input 1), the acceleration data (Outputs 1, 2, 3, and 4) are presented in Figure 4. e required experimental modal data are extracted from FRFs and written as codes in MATLAB [32].

Updating the Finite Element Model.
e aim of updating the finite element model is correcting or improving the finite element model of a structure in terms of updating the global or local mass and stiffness matrices, boundary conditions, and geometrical and material properties. For this purpose, the updating process was accomplished with the modal assurance criterion (MAC) method [5] based on the recorded data at four stations in the experiment program to correcting the global mass and stiffness matrices. As well, the finite element model is reduced to coincide with the experimental setup with the SEREP algorithm [33]. Update has been conducted for mass and stiffness matrices using the particle swarm optimization (PSO) technique [34]. e objective function is defined based on the frequencies and the components of the mode shapes in the experimental and numerical data as follows: 4 Shock and Vibration where ω num

Frequencies and the Mode Shapes.
Modal parameters from analytical, numerical, and also experimental methods were extracted and presented in this section. Table 1 summarizes the geometrical and material properties of the studied aluminum beam.
In the numerical method, firstly, the beam is divided into 25 elements; after that, frequencies and mode shapes were extracted. en, to compare the numerical and experimental results due to the limitations of acceleration sensors, the numerical model was also reduced from 25 elements to four elements. e results of the frequencies are presented in Table 2. As well, in Figure 5, the mode shapes of the first four modes are illustrated with experimental [35], analytical, and also numerical methods.
As can be seen from Figure 5, the analytical and numerical mode shapes of the beam have good agreements. However, the experimental results in some modes differ from the other methods. e reason for the differences in some modes is related to the scaling of mode shapes. erefore, the numerical model has been updated to modify the mass and stiffness matrices of the finite element model [5], which is performed with the particle swarm algorithm based on the objective function defined in equation (21). e convergence of the objective is shown in Figure 6. As well, Table 3 presents the values of the numerical frequencies before and after the updating process. In this study, numerical model updating was conducted only with frequency information.
e updated values of mass and stiffness matrices are summarized in Table 4. ese values are applied to calculate the numerical responses of the beam under impact loading.
After calculating the updated mass and stiffness matrices, the damping matrix is obtained using the Rayleigh damping matrix as follows: e coefficients α and β are obtained, according to equation (23), where ω j and ω i are the first and last available vibration frequencies, respectively, and also, ξ is a critical damping ratio that was obtained equal to 0.4% from the experimental data: Loading). e unconditionally stable Newmark average acceleration method is recommended for the calculation of the response by modal superposition when filtering the high     frequencies is naturally carried out by using a truncated modal space. It is the most accurate method, and the absence of numerical damping is not a problem [36]. In this paper, the capabilities of the Newmark method with the direct integration method, instead of the modal superposition method, are investigated by studying an SDF system. For this purpose, consider the damped SDF model with the specification that is summarized in Table 5. e natural frequency of the SDF system is obtained equal to 12 rad/s. e structural response of the SDF system is presented in Figure 7 for harmonic loading with a frequency equal to 1.2 rad/s and damping ratio (ξ) equal to 0.4%. As seen in this figure, the numerical method has excellent capability to predict the response. In this case, the natural frequency of the SDF system is 12 rad/s, and from the structural dynamics point of view, an excitation with ω � 1.2 rad/s could be taken into account as a forced vibration-type loading (Figure 7). e structural response of the SDF system is presented in Figure 8 for harmonic loading with a frequency equal to 60 rad/s and damping ratio (ξ) equal to 0.5%. As seen in this figure, the numerical method has good capability to predict the responses. From structural dynamics points of view, an   Table 4: Mass and stiffness matrices of the beam before and after updating based on frequency data.

Dynamic Analysis Using the Newmark Method (Harmonic
Unit of mass is kg and stiffness is N/m   excitation with ω � 60 rad/s and ω � 12 rad/s could be taken into account as a free vibration. Figures 7 and 8 reveal that for regular loading, the Newmark average acceleration method has an appropriate performance. e free end displacement responses of the studied beam were investigated to show the abilities of the Newmark method under impact or very short-time loadings, and the results were compared with the experimental data.

Dynamic Analysis
Using the Newmark Method (Impact Loading). According to Figure 9, the studied beam is subjected to very short-time impact loading at its free end in the laboratory. As seen from Figure 10, the acceleration data at four stations (CH1, CH2, CH3, and CH4) are presented.
e Newmark-beta average acceleration technique [36] was used to determine the numerical response of the beam, and the results are compared with recorded experimental data at stations 1-4 (CH1 to CH4). As shown in Figure 11, the numerical and experimental values have no appropriate agreements even after updating the finite element model. e only possibility that comes to the authors' mind is that the Newmark is not suitable for shock-type loadings, and the dynamic responses of the low damping MDF structures should be obtained from other methods [37,38].

Conclusions
In this study, the modal parameters of the aluminum beam were determined by analytical, numerical, and also experimental methods. en, updating the finite element model is conducted based on the modal data and particle swarm optimization techniques. en, the updated FE model was used to estimate the response of the beam under a very shorttime impact load. Based on the results of this study, the following conclusions could be drawn: (i) In the numerical model (finite element), increasing the number of DOF leads to an accurate estimation of higher modes but does not improve the low mode values. As well, because of the technical limitations in measurement of the rotations, using model reduction techniques to eliminate the rotational DOF does not lead to a reduction of accuracy of the numerical model.
(ii) e results show that if the objective function is only defined based on frequency data, the accuracy of the model is slightly reduced. However, simultaneous use of mode shapes and frequencies data improve the accuracy, but extracting the mode shapes from the laboratory data is not an easy task. (iii) e results reveal that the Newmark average acceleration method has excellent capabilities to accurately predict the dynamic responses under regular type loading even for low damping systems without facing numerical damping error. (iv) e results show that the Newmark average acceleration method that inherently is unconditionally stable is not able to accurately estimate the dynamic responses under shock-or impact-type loading, and other numerical or analytical methods should be utilized.

Data Availability
No data were used to support the findings of the study.