Modal,static,harmonic,transient-dynamic and random vibration analyses of a two degree of freedom system by MatLAB and ANSYS

Abstract

Static, harmonic, transient-dynamic and random vibration analyses are widely used to design structures. Commercial finite element (FE) programs such as ANSYS are used to analyze structures with complex geometries. However, FE programs can produce erroneous results if the user makes even minor mistakes. Thus, it is important to understand the basic theory of FE solutions. The purpose of this work is to give the introductory theory of FE solutions by considering a two degree of mass-spring system. The MatLAB programs based on the analytical solutions are given. The MatLAB results are compared with the FE results by ANSYS. The system is modeled in SolidWorks and imported to ANSYS. Modal, static, harmonic frequency response, impulse response, and random vibration responses are considered. Modal and static analyses for unconstrained systems are also considered. The results obtained by MatLAB and ANSYS are in good agreement. The results of this study can be used by the users of other simulation programs for comparison.

Keywords

Modal analysis harmonic analysis transient-dynamic analysis random vibration analysis ANSYS

Introduction

Advanced mathematical modeling of mechanical systems requires application of the laws of mechanics and solve resulting linear or non-linear, algebraic or differential equations. Before the computer age, the solution capabilities of these equations by using slide rulers or graphical methods were very limited. Engineers started to write codes to solve the equations in the early stage of computers. Later, they started to use commercial simulation software such as ANSYS. Simulation software enables one to model the system only, and then the software generates the mathematical model and solves the equations. The users observe the results for design purposes. Developing mechanical systems followed produce, test, and revise procedures before the computer age changed to model, analyze, revise, and produce procedures today. However, it should be noted that simulation methods must be used with great care because a minor mistake in the process can give erroneous results.

Design of mechanical systems require solid modeling, assembly, rigid body dynamic, static, modal, harmonic, transient-dynamic and random vibration analyses. Analyses of unconstrained structures are used in the aerospace, automotive, and naval industries. The modal and static analyses of unconstrained systems can be performed in simulation software using the inertia relief option.

Akdağ et al.¹ presented an integrated approach to design a hexapod using computer aided engineering (CAE) simulation software, SolidWorks. Deepak and Shetty² studied static and free vibration of functionally graded rectangular plates using ANSYS. Imran et al.³ analyzed delaminated composites using ANSYS and ABAQUS to find natural frequencies. Pappalardo et al.⁴ modeled a single-link flexible robot arm in SolidWorks and then performed the modal analysis in ANSYS. Sahoo et al.⁵ investigated static, free vibration, and transient behavior of laminated composite curved shallow panels using ANSYS. Abramowicz et al.⁶ modeled laminated veneer lumber slabs of different thicknesses, width and types in ABAQUS and studied natural frequencies. Yavuz et al.⁷ used ANSYS to analyze the residual vibration of a curved non-uniform steel manipulator with a complex shape. Bachoo and Elishakoff⁸ analyzed random vibrations of laminated planar frames using ANSYS. Chiariello et al.⁹ used inertia relief option of NASTRAN to study the structural response of an aircraft wing. Hamedi et al.¹⁰ designed a solar car chassis using ANSYS. Fan et al.¹¹ studied dynamic characteristics of a rotating blade with a dovetail fixture in the ANSYS environment. Doori et al.¹² investigated static response of functionally graded porous circular plates using ANSYS. Yashar et al.¹³ studied natural frequencies and mode shapes of a rotating cracked beam using a simplified method and compared the results with ANSYS results. It should be noted that simulation programs are widely used in academia and industry and a limited number of example works in the literature is given in this introduction.

A two degree of freedom system is considered in this study. Modal, static, harmonic, transient-dynamic and random vibration analyses, which uses Laplace transform methods, are given in the text books.¹⁴ In this study, the solutions by the Laplace transform methods using MatLAB codes are given. The system is also modeled in SolidWorks and imported to ANSYS. The procedures to obtain solutions in ANSYS are described and the results obtained by MatLAB and ANSYS are compared. Users of simulation programs should confirm the validity of their procedures by comparing them with analytical solutions. The main purpose of this work is to introduce example analyses for comparison. The procedures for analyses in ANSYS are given. It is observed that the MatLAB and ANSYS results are in good agreement. The results given by the MatLAB solutions can be used by the users of other simulation programs for comparison.

Modeling and analyses of two degree of system

Creating model

The system considered is shown in Figure 1(a). Here m₁, m₂, and m_b are masses; k₁ and k₂ are spring constants; x_b, x₁ and x₂ are displacements; p₁ and p₂ are external forces. The following values are assigned for the study: m₁= 5 kg, m₂= 4 kg, m_b= 220.78 gr, k₁= 3000 N/m, k₂= 2200 N/m. The finite element (FE) analysis model created in ANSYS¹⁵ is shown in Figure 1(b). The masses are modeled in Solidworks¹⁶ as cubes having given mass values. The distances between the masses have no effect on the results. Modal Analysis is started in ANSYS Workbench. The SLDPRT model is saved in parasolid format in SolidWorks and imported to ANSYS Workbench as Geometry. The elasticity modulus, E, equals to 200 GPa, and the density, ρ, equals to 7850 kg/m³. ANSYS Mechanical is started in ANSYS Workbench, and the analyses are continued. The springs are assigned as Connections in Model. Face-1 of the mass, m_b, is fixed (constrained). See Appendix 1 for ANSYS modeling.

Modal analysis

Modal analysis in MatLAB

Modal analysis is performed to observe natural frequencies and mode shapes of a structure. This analysis is important to avoid resonances and to control vibration. The mathematical model of the system can be found either by applying the equations of motions and action-reaction law to the free body diagrams to each mass, or by using energy methods and applying the Lagrange equations.¹⁴ The equations can be transformed to matrix form as

[M] {\ddot{x}} + [K] {x} = {p}

(1)

where

{x} = {\begin{matrix} x_{1} \\ x_{2} \end{matrix}}

[M] = [\begin{matrix} m_{1} & 0 \\ 0 & m_{2} \end{matrix}]

[K] = [\begin{matrix} k_{1} + k_{2} & - k_{2} \\ - k_{2} & k_{2} \end{matrix}]

{p} = {\begin{matrix} p_{1} \\ p_{2} \end{matrix}}

The base is fixed in this model. [M] is the mass matrix and [K] is the stiffness matrix. The external force term is zero for the modal analysis. A solution as {x} = {a}e^st is tried as the solution where t is the time, and s is eigenvalue, {a} = [a₁,a₂]^T. This solution is inserted in equation (1) and [s²M + K]{a} = {0} is obtained. This gives {a} = [s²M + K]⁻¹{0}. The trivial solution is {a} = 0, which means no vibration. The non-trivial solution is found when the determinant of [s²M + K] is zero. So, the eigenvalue equation is found as

| s^{2} M + K | = 0

(2)

The MatLAB program which gives the un-damped natural frequencies and mode shapes is given in the Appendix 2. The result of the program is given below:

The natural frequencies: f₁= 2.4934 Hz, f₂= 5.8359 Hz.

Modal matrix: $[U] = [\begin{matrix} 1 & 1 \\ 1.8058 & - 0.6922 \end{matrix}]$

a₂= 1.8058 if a₁= 1 for the first natural frequency, f₁ (First mode shape)

a₂= −0.6922 if a₁= 1 for the second natural frequency, f₂ (Second mode shape)

It should be noted that the two masses are modeled as cubes having the mass properties considered in the study. They are considered as elastic in ANSYS. Thus, there are more than two natural frequencies. 3rd and 4th natural frequencies are found as 7545.18 Hz and 13798 Hz, respectively. The natural frequencies beyond the 2nd natural frequency are much higher; thus the masses are approximately rigid, and only 2 modes are considered for ANSYS solutions.

Modal analysis in ANSYS workbench-mechanical

The users of simulation programs model systems and define constraints. The programs develop finite element model of the system after meshing. The degree of freedom of the system is very large for complex engineering structures. A limited number of modes are considered considering mode shapes which are expected to be important.

The details of the modal analysis process are given in Appendix 3. Face-1 in Figure 1(b) is fixed. The displacements in z and y directions on the front faces in Figure 1(b) are constrained as zero. Two modes are selected for the analysis. Modal analysis solution is obtained. Mode shapes can be animated and relative displacement values can be obtained by using Probe. The results are given in Figure 2.

Comparison of modal analyses by MatLAB and ANSYS

It is observed that the natural frequencies found by MatLAB and ANSYS are in good agreement. The mode shapes are un-definite and therefore a₂/a₁ can be compared. a₂/a₁= 0.42512/0.23541 = 1.8059 for Mode-1, and a₂/a₁= −0.2632/0.38024 = −0.6922 for Mode-2 for ANSYS results. The results can be compared quantitatively by percentage error, defined as

100 (MatLAB solution - ANSYS solution) / MatLAB solution

The percentage error for the first natural frequency is 100(2.4934−2.4934)/(2.4934) = 0. Similarly, the percentage errors are 0, −0.0055, and 0 for the second natural frequency, mode shape-1 and mode shape-2, respectively. It is observed that the MatLAB and ANSYS solutions are in good agreement.

Static analysis

Static analysis is conducted in many engineering applications where there is no change in the system and forces with time. The derivatives are zero for the static analysis. Therefore, equation (1) leads to [K]{x} = {p} for the static analysis. The static displacements are found by {x} = [K]⁻¹{p}. The following values are assigned for the analysis: p₁= 10 N, p₂= −5 N.

Static analysis in MatLAB

The MatLAB program which gives the static analysis is given in the Appendix 4. The result of the program is given as x₁= 1.6667 mm, x₂= −0.6061 mm.

Static analysis in ANSYS

The details of the analysis process are given in Appendix 5. Face-1 in Figure 1(b) is fixed. The displacements in z and y directions on the front faces in Figure 1(b) are constrained as zero. The results are given in Figure 2. The forces are defined on Face-4 as 10 N, and on Face-6 as −5 N in the x direction. The solution is obtained and the result is given in Figure 3.

The results of ANSYS are found as x₁= 1.6667 mm and x₂= −0.60606 mm.

The percentage errors are 0 and 0 for x₁ and x₂, respectively. It is observed that the results obtained by MatLAB and ANSYS are in good agreement.

Frequency response analysis

The harmonic input is defined as cos(2πf), where f is the frequency in Hz. There are applications for harmonic or combinations of harmonic excitations. It is important to evaluate frequency dependent behavior of a system to observe resonances.

Frequency response analysis in MatLAB

Proportional damping (Rayleigh damping) is considered in the frequency response analysis. Mass and stiffness matrices are found in FE analysis, but damping mechanisms are difficult to identify. The proportional damping is a common method of modeling viscous damping.¹⁷ The damping matrix is defined as [C] = α[M] + β[K] where α and β are the proportional damping constants. Modal damping ratios are given as ξ_i= α/(2ω_i) + βω_i/2, where i is the i-th mode, and ω_i is the i-th natural frequency. The damping constants can be determined for experimental data.¹⁷ The mathematical model of the vibratory system becomes

[M] {\ddot{x}} + [C] {\dot{x}} + [K] {x} = {p}

(3)

The following values are assigned as example values for the study: α = 0 and β = 0.005.

The force input is taken as {p} = {P}e^st where {P} = [P₁, P₂]^T. Then, the response can be written as {x} = {X}e^st where {X} = [X₁, X₂]^T, and {X} = [H]{P}. The following equations are obtained by inserting the assigned input and response: s²[M][H] + s[C][H] + [K][H] = 1,

[H] = [s^{2} [M] + s [C] + [K]]^{- 1}

(4)

{\begin{matrix} X_{1} \\ X_{2} \end{matrix}} = [\begin{matrix} H_{11} & H_{12} \\ H_{21} & H_{22} \end{matrix}] {\begin{matrix} P_{1} \\ P_{2} \end{matrix}}

(5)

Here, [H] is 2 × 2 transfer function matrix; X₁(s) and X₂(s) are the Laplace transforms of x₁(t) and x₂(t), respectively; P₁(s) and P₂(s) are the Laplace transforms of p₁(t) and p₂(t), respectively. Here, for example, H₂₁(s) is the transfer function from p₁ to x₂. The amplitude of x₂ for the harmonic input, p₁ is found as A₂₁= |H₂₁(iω)|. Here ω is the input frequency in rad/s. The input frequency in Hz is found from ω = 2πf where f is the input frequency in Hz. The frequency response is obtained by changing f from 0 to f_end, where f_end is the ending input frequency. f_end can be selected approximately 2f_max, where f_max is the highest natural frequency. The MatLAB program, freqra.m, which gives the plot of the frequency response, A₂₁(ω), by inverting the transfer function matrix is given in Appendix 6.

An alternative method to obtain the frequency response function, which is recommended for large number of degrees of systems, is the mode superposition method.¹⁴ First, the modal matrix [U] is found in this method. The transformation as {x} = [U]{q} is used to obtain uncoupled equations

m_{i} {\ddot{q}}_{i} + c_{i} {\dot{q}}_{i} + k_{i} q_{i} = Q_{i}

(6)

Here, i is the mode number (1 or 2 for 2 degrees of freedom system), q_i is the generalized coordinate for Mode-i. m_i, c_i, and k_i are the i-th row and i-th column of the diagonal matrices [U]^T[M][U], [U]^T[C][U], and [U]^T[K][U], respectively. [U]^T is the transpose of [U]. The generalized force,Q_i, is the i-th row of [U]^T{p}. The transfer function obtained from equation (6) is H_i(s) = 1/(m_is²+ c_is + k_i). In order to find A₂₁(ω) for the example studied, {p} is taken as [1,0]^T. Then, the generalized forces are obtained as Q₁= cos(ωt) and Q₂(t) = cos(ωt). q₁(t) and q₂(t) are found by taking the real parts of H₁(iω) and H₂(iω), respectively. x₁ (t) is obtained as the 1st row of [U]{q}. The amplitude of x₁(t) gives A₂₁(ω). The frequency response function plot is obtained by changing the frequency in Hz from 0 to f_end with a selected interval. The MatLAB program, freqrb.m, which gives the plot of the frequency response, A₂₁(ω), by the mode superposition method is given in Appendix 6.

Frequency response analysis in ANSYS

The details of the analysis process are given in Appendix 7. Face-1 in Figure 1(b) is fixed. The displacements in z and y directions on the front faces in Figure 1(b) are constrained as zero. Define Range Maximum as 10 Hz and Solution Interval as 100. Define Stiffness Coefficient as 0.005 for Damping Control. The amplitude (magnitude) of the input force is defined on Face-4 as 1 N in the x direction. The solution is obtained as the frequency response.

Transient response

Structures are excited by time dependent inputs in many engineering applications. These inputs may be impulse, step, or arbitrary time dependent inputs. The impulse response is considered as the study for the transient analysis. Impulse response is used for system identifications.

Transient response analysis in MatLAB

The inverse Laplace transform of H₂₁(s) is the impulse response of the system where p₁(t) is the input and x₂(t) is the output. It is required to define time step, Δt; and the steady state time, t_ss for the transient analysis. The time step can be defined as Δt = 1/fmax/20. The steady state time can be defined as t_ss= 0.6/f₁/ξ₁, where f₁ is the first natural frequency in Hz, ξ₁ is the modal damping for the first mode. The MatLAB program which gives the plot of the impulse response is given in the Appendix 8. The program gives Δt = 0.0086 s, and t_ss= 6.1439 s for the example values considered in this study.

Transient response analysis in ANSYS

The details of the analysis process are given in Appendix 9. Define Time Step as 0.0086 s, and Step End Time as 6.1439 s. Define Stiffness Coefficient as 0.005 for Damping Control. Define the impulse force on Face-4 (Figure 1(b)) using the numerical model given in⁷ as 1/Δt = 116.2791 N at t = Δt, and 0 as other sample values. The solution is obtained and the resulting Directional Deformation time history in the x direction on Face-6 (Figure 1(b)) is observed.

Comparison of transient analyses by MatLAB and by ANSYS

The MatLAB program which compares the MatLAB and ANSYS transient analysis results is given in Appendix 10. The impulse responses obtained by MatLAB and ANSYS are shown in Figure 4.

The number of peaks in Figure 4 are the same for MatLAB and ANSYS solutions. The percentage errors are 1.15 and −0.37 for the peak values and logarithmic decrements,¹⁴ respectively. It is observed that the results obtained by MatLAB and ANSYS are in good agreement.

Frequency response obtained from impulse response and comparisons of frequency responses

The magnitude of the FFT of the impulse response gives the frequency response. The MatLAB program which gives the frequency response obtained by taking the FFT of the impulse response obtained by ANSYS is given in Appendix 11. The program also compares the frequency responses obtained by MatLAB, ANSYS, and the frequency response obtained from the impulse response.

All the three frequency responses are compared in Figure 5(a).

Obtaining the frequency response from the impulse response by taking the FFT is practical in experimental applications. The natural frequencies and modal damping for the modes can be observed from the frequency response.¹⁴

The percentage errors are 0.00 and 3.13 for the first mode peak values for ANSYS and impulse response solutions, respectively. The percentage errors are 0.00 for the first mode resonance frequency values for ANSYS and impulse response solutions. The percentage errors are 0.00 and 5 for the first mode bandwidths¹⁴ for ANSYS and impulse response solutions, respectively.

The percentage errors are 0.00 and 0.18 for the second mode peak values for ANSYS and impulse response solutions, respectively. The percentage errors are 0.00 and 0.65 for the second mode resonance frequency values for ANSYS and impulse response solutions, respectively. The percentage errors are 0.00 and 1.94 for the second mode bandwidths for ANSYS and impulse response solutions, respectively.

The percentage errors are 0.00, 0.33 for the areas underneath of the frequency response curves for ANSYS and impulse response solutions, respectively.

It is observed from Figure 5(a) that the results obtained by MatLAB and ANSYS are in good agreement.

The response given in Figure 5(a) depends on the excitation and damping. The amplitudes p₁= 1, and p₂= 0; and the Rayleigh damping coefficients are α = 0, β = 0.005 for Figure 5(a). Resonances are expected close to natural frequencies, 2.4934 Hz and 5.8359 Hz. Excitation determines which modes are more effective in the response, and damping ratios increase for higher natural frequencies according to the Rayleigh damping model. The resonance is not observed for 5.8359 Hz as observed for 2.4934 Hz in Figure 5(a). The response obtained by ANSYS for the amplitudes p₁= 1, and p₂= −0.5; and the Rayleigh damping coefficients α = 0, β = 0.00001 is shown in Figure 5(b). It is observed that the resonance near 5.8359 Hz is observed with changing excitation and damping coefficients.

Random vibration analysis

The excitations are random in many engineering applications. The input is deterministic if the value of the input is known at a given time. The input is random if the value of the input cannot be predicted at a given time. Wind, road, earthquake excitations are examples of random inputs. Random vibration analysis is performed by using statistical quantities.¹⁴

Random vibration analysis in MatLAB

It is considered that Face-1 in Figure 1(b) is excited randomly. Then, it can be found that {p} = [k₁x_b,0]^T in equation (1). The following equation is used for the random vibration analysis¹⁴

S_{2} (ω) = | H_{r} (i ω) |^{2} k_{2} S_{b} (ω)

(7)

where S_b(ω) is the power spectral density of

{\ddot{x}}_{b}

and S₂(ω) is the power spectral density of (

\ddot{x}

₂ −

{\ddot{x}}_{b}

), H_r(iω) = k₁H₂₁(iω) − 1. Power spectral density of the input signal can be found by MatLAB commands such as pwelch. Let S_b(ω) is given as in Figure 6.

There are industry-based standards such as DO-160 and MIL-STD-810 to test the endurance of airborne equipment under severe levels of vibrations.¹⁸ The equipments are tested according to the standards with the random vibration test machines, which produce random vibrations having power spectral densities like in Figure 6. The MatLAB program which gives the plot of the power spectral density of the output ( $\ddot{x}$ ₂ − ${\ddot{x}}_{b}$ ) is given in Appendix 12. The area underneath of the power spectral density of the output gives the root mean square (rms) value of the response. The rms value is the expected value with one standard deviation (%68). The Matlab program gives 928.0075 mm/s² as the rms value.

Random vibration analysis in ANSYS

The details of the analysis process are given in Appendix 13. Define Stiffness Coefficient as 0.005 for Damping Control. Define The Power Spectral Density (PSD) of the input given in Figure 6 to the fixed support. The solution is obtained and the plot of the resulting PSD of the output, ( $\ddot{x}$ ₂ − ${\ddot{x}}_{b}$ ), is obtained as the Response PSD. The rms value (the expected value with one standard deviation) is read as 933 mm/s².

Comparison of random vibration analyses by MatLAB and ANSYS

PSD of the output obtained by MatLAB and ANSYS are compared in Figure 7. The MatLab program to obtain Figure 1 is given in Appendix 14.

The percentage error is −0.54 for the ANSYS solution in Figure 7. It is observed that MatLAB and ANSYS results are in good agreement.

Modal analysis of unconstrained system

Analysis of unconstrained structures is widely performed in the aerospace, automotive, and naval industries. There are no fixtures in the unconstrained systems. Consider the case where Face-1 in Figure 1(b) is unconstrained (free), so m_b is free to move.

Modal analysis of unconstrained system by MatLAB

To find the natural frequencies and mode shapes, the mathematical model of the unconstrained system can be derived and the mass and stiffness matrices in the MatLAB program given in Appendix 2 is changed as:

M = [m 1, 0, 0; 0, m 2, 0; 0, 0, mb]; K = [k 1 + k 2, - k 2, - k 1; - k 2, k 2, 0; - k 1, 0, k 1];

Modal analysis of unconstrained system by ANSYS

Fixed Support is suppressed and Max Modes to find is changed as 3 in Analysis Settings in Modal-ANSYS Mechanical to obtain modal analysis results.

Comparison of modal analyses of unconstrained system by MatLAB and ANSYS

The MatLAB and ANSYS results are given in Table 1.

Table 1.

Modal analysis results for unconstrained system.

MatLAB	ANSYS	Modal Matrix- MatLAB	Modal Matrix-ANSYS
f₁= 0 Hz	1.0802e-3 Hz	$[\begin{matrix} 1 & 1 & 1 \\ 1 & - 1.276 & - 0.0184 \\ 1 & 1.034 & - 49.0293 \end{matrix}]$	$[\begin{matrix} 0.33149 & - 0.29336 & - 6.3687 e - 2 \\ 0.33149 & 0.37432 & 1.1709 e - 3 \\ 0.33149 & - 0.30333 & 3.1226 \end{matrix}]$
f₂= 4.985 Hz	4.985 Hz
f₃= 27.7796 Hz	27.88 Hz

For comparison all the rows in the Modal Matrix-ANSYS are divided to the first raw. The 2^nd column becomes [1, 0.37432/(−0.29336), −0.30333/(−0.29336)]^T= [1, −1.276, 1.034]^T for example. It is observed that the MatLAB and ANSYS results are in good agreement. First mode is rigid body motion because the system is unconstrained.

The percentage errors are 0.00 for the natural frequencies, and mode shapes for ANSYS in Table 1. It is observed that MatLAB and ANSYS results are in good agreement.

Static analysis of unconstrained system

External forces applied to unconstrained systems impose rigid body motions and conventional static FE analyses cannot be performed. One way to overcome this problem to apply inertia relief method where distributed inertia forces due to the rigid body accelerations and dummy boundary conditions are used. Reaction forces at the dummy fixtures become zero because static equilibrium condition is satisfied without the need for the reaction forces.

Static analysis of unconstrained system in MatLAB

The following equation of motion can be written for the system in Figure 1 when the fixture at Face-1 is suppressed, and the masses are assumed to be connected rigidly.

p_{1} + p_{2} = (m_{b} + m_{1} + m_{2}) a

(8)

where a is the acceleration of the masses. Then the fixture is un-suppressed and the displacements are found by {x} = [K]⁻¹{p_sum} where [K] is the same as in equation (1), and

p_{sum} = [p_{1}, p_{2}]^{T} + [- m_{1} a, - m_{2} a]^{T}

(9)

The following values for the external forces are considered as in static analysis of constrained system: p₁= 10 N, p₂= −5 N. The following results are obtained for the displacements: x₁= 0.0399 mm, and x₂= −3.2187 mm.

Static analysis of unconstrained system in ANSYS

The inertia relief option in Solver Controls of Details Window of the Analysis Settings is changed as “On” in the Static Analysis of constrained system in ANSYS. ANSYS calculates the accelerations and assigns distributed inertia forces and gives the results after the solution. The following results are obtained for the displacements: x₁= 0.039907 mm, and x₂= −3.2187 mm.

Comparison of static analyses of unconstrained system by MatLAB and ANSYS

The MatLAB and ANSYS results are given in Table 2.

Table 2.

Static analysis results for unconstrained system.

[x₁, x₂] = [0.0399, −3.2187] by MatLAB

[x₁, x₂] = [−0.039907, −3.2187] by ANSYS

The percentage errors 0.02 for the ANSYS solution in Table 2. It is observed that the MatLAB and ANSYS results are in good agreement.

Results and discussions

A two degree of system for which analytical solutions can be found is considered in this study. The results for modal, static, harmonic, transient-dynamic, and random vibration analyses for a selected mass and spring values are obtained using analytical methods and implemented by MatLAB programs are given.

The authors used SolidWorks to model masses and their positions and transferred the geometry to ANSYS in parasolid format. The geometry modeling can be performed by ANSYS-SpaceClaim or by other simulation programs. The authors preferred SolidWorks with which they are familiar. Next, the joints (springs between the masses) are defined in ANSYS-Mechanical to do analyses. The results for modal, harmonic, static, transient-dynamic, and random vibration analyses for a selected mass and spring values are obtained using ANSYS are given. Also, modal and static analyses for unconstrained systems obtained using analytical methods implemented with MatLAB programs and ANSYS are presented. The MatLAB programs and the ANSYS procedures using Graphical User Interface (GUI) are given in the Appendices. The screenshot of the ANSYS Workbench Project Schematic is given in Figure 8.

MatLAB and ANSYS results are compared in terms of the percentage error, which is defined as 100(MatLAB solution-ANSYS solution)/MatLAB solution. It is observed that the results are in good agreement. Small differences result due to different solution methods.

It should be noted that commercial simulation programs are used in the analysis of complex engineering systems. Modeling of systems, constraints, analysis types, analysis parameters, inputs are defined to the programs by keyboard and GUI. Any minor or major error in the definitions may result in large erroneous solutions. For example, if the forces are defined as p₁ = 10 N, p₂ = 5 N instead of p₁ = 10 N, p₂ = −5 N in the static analysis of constrained system, the displacements are obtained as x₁ = 5 mm and x₂ = 7.2727 mm instead of x₁ = 1.6667 mm and x₂ = −0.6061 mm. Simulation programs allow users to analyze complex engineering systems, which have a large number of degrees of freedom, without the need to develop mathematical models and programs. However, they should be very careful about the definitions of models, constraints, analysis parameters, and inputs.

Conclusions

Mechanical systems with bodies with complex shapes can be modeled and analyzed in various simulation programs. Modeling of the geometry and analyses can be done with different programs and the geometry model is imported to the analysis program. As in this work, the geometry model can be created in SolidWorks and transferred to ANSYS. Then the joints are defined in the analysis program.

Simulation programs are widely used to perform modal, static, harmonic, transient-dynamic and random vibration analyses. Modal and static analyses of unconstrained systems by simulation programs are considered in aerospace automotive, and naval industries, also. The users of simulation programs must ensure about the procedures applied for the analyses. When first using simulation programs, it is advised that solutions of commercial simulation programs are compared with analytical solutions obtained for basic systems; and then continue to use simulation programs for analyzing complex engineering systems. A two degree of freedom vibratory system for which analytical solutions can be obtained is considered in this study. The analytical solutions with MatLAB and ANSYS simulations are given and the results are compared for modal, static, harmonic, transient-dynamic, and random analyses for the constrained system; and modal and static analyses for the un-constrained system. Mechanical systems with complex geometries cannot be analyzed analytically and simulation programs are used. The procedures given for simulation can be applied for mechanical systems with complex geometries. It is proposed that the users of simulation programs may study the two degree of systems in this work for comparison before studying mechanical systems with complex geometries.

The contents of this study are designed to be appropriate for undergraduate and graduate engineering students, as well as professionals participating in short-term continuing education programs. This material can be incorporated into courses such as Structural Dynamics, Finite Element Analysis, or Mechanical Vibrations. These courses are typically offered to third or fourth year undergraduate or first year graduate students. Additionally, the study is well-suited for 2–5 day continuing education courses aimed at engineers who wish to enhance their skills in finite element methods and vibration analysis. The study is relevant for both purely theoretical courses and those with practical components. For instance, the theoretical aspects of finite element solutions and mass-spring systems are suitable for formal lectures, while the MatLAB programs and ANSYS modeling are suitable for laboratory sessions or hands-on workshops.

By integrating this content into the curriculum or short-term courses, we aim to enhance participants’ understanding of the fundamental theory behind vibration analysis. The comparison between analytical and numerical results has potential to develop critical thinking skills and enable users to validate simulation outputs effectively. For professional engineers, this knowledge will help minimize errors in their use of commercial FE software, thereby improving the accuracy of their designs and analyses.

The results of this study can be used for various simulation programs. This work considers single-body analyses for simulation and can be extended for multi-body analyses where rotational and translational joints are available. Analysis of systems with nonlinearities or uncertainties are other areas for further exploration.

Footnotes

Acknowledgements

BLM Mekatronik Limited, who offers short courses on Computer Aided Engineering to companies which own licenses of simulation programs, is acknowledged for introducing the subject of this work.

Declaration of conflicting interests

The authors declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The authors received no financial support for the research, authorship, and/or publication of this article.

ORCID iDs

Hira Karagülle

İlker Bülbül

Murat Akdağ

Data availability statement

All data generated or analyzed during this study are included in this published article.

Modeling in ANSYS Workbench (WB)/Mechanical

Run ANSYS Workbench > Double-click Modal in Toolbox window and start Modal Analysis.

Right-click Geometry in Modal Analysis > Import Geometry > Browse > import parasolid file

Double-click Model in Modal Analysis and start Ansys Mechanical > continue in ANSYS Mechanical

Right-click Connection in Model > Insert > Spring > Assign k₁ value as Longitudinal Stiffness>…

Click Face-2 (Figure 1(b))> Apply Reference Scope in Details window (DW)

Click Face-3 (Figure 1(b))> Apply Mobile Scope in DW

Similarly, assign Spring, k₂, between Face-3 and Face-4

(Note: Quick Launch > Details > if DW does not appear)

Modal analysis in MatLAB

% modal.m

clc,clear,close all

%Define mass and stiffness matrices

m1 = 5;m2 = 4;mb = 0.22078;k1 = 3000;k2 = 2200;

M = [m1,0;0,m2];K = [k1 + k2,-k2;-k2,k2];

ndof = size(M);ndof = ndof(1,1);

% Obtain solutions using equation (2)

syms s;a = [s^2*M + K];d = det(a);p = double(solve(d));

p = abs(p);p = sort(p);

for k = 1:1:ndof;wn(k,1) = p(2*k−1);end

fn = wn/(2*pi);mod = [];s2 = -wn.^2;

for k = 1:1:ndof;

a = s2(k)*M + K;b = -a(2:ndof,1);a = a(2:ndof,2:ndof);

modk = inv(a)*b;modk = [1;modk];mod = [mod,modk];

end

fn,mod

Modal analysis in ANSYS mechanical

Analysis Settings in Modal > Assign 2 for Maximum Modes to Find

Right-click Modal > Insert > Fixed Support > Click Face-1 (Figure 1(b))> Geometry > Apply

Right-click Modal > Insert > Displacement > Press Ctrl and Click, and select all the front faces>… Geometry > Apply > Assign 0 for z and y components¹

Right-click Modal > Solve > Solution > Observe results in Tabular Data

Select all Frequencies in Tabular Data > Right Click > Create Mode Shape Results

Right-click Solution > Evaluate all Results > Play mode shapes in Graph

Right-click Solution > Insert > Deformation > Directional > X Axis for Orientation in DW>

Select 1 for Mode in DW > Right-click Solution > Evaluate all Results¹

Click Probe on Result Ribbon > Click on front faces in Figure 1(b) > Read mode shape values

Static analysis in MatLAB

% stc.m

clc,clear,close all

%Define mass and stiffness matrices

m1 = 5;m2 = 4;k1 = 3000;k2 = 2200;

M = [m1,0;0,m2];K = [k1 + k2,-k2;-k2,k2];

p = [10;−5];d = inv(K)*p*1000

Static analysis in ANSYS

Drag and drop Static Structural in Toolbox window in WB to Model in Modal Window

Select ANSYS Mechanical window

Drag and drop Fixed Support in Modal to Static Structural

Drag and drop Displacement in Modal to Static Structural

Right-click Static Structural > Insert > Force > Click on Face-4 (Figure 1(b))> Geometry in DW > Apply

Define Magnitude in DW as 10 N

Right-click Static Structural > Insert > Force > Click on Face-6 (Figure 1(b))> Geometry in DW > Apply

Define Magnitude in DW as −5 N

Right-click Static Structural > Solve > Solution

Right-click Solution > Insert > Deformation > Directional > X Axis for Orientation in DW>

Right-click Solution > Evaluate all Results

Click Probe on Result Ribbon > Click on front faces in Figure 1(b)> Read displacement values

Frequency response analysis in MatLAB

% freqra.m

% Obtaining frequency response function with inversion of transfer function matrix

clc,clear,close all

%Define mass, stiffness and damping matrices

m1 = 5;m2 = 4;k1 = 3000;k2 = 2200; alpha = 0;beta = 0.005;

M = [m1,0;0,m2];K = [k1 + k2,-k2;-k2,k2]; C = alpha*M + beta*K;

fend = 10;ns = 100;

syms s;H = inv(s^2*M + s*C + K);Hr = H(2,1);

df = fend/ns;f = df:df:fend;v = 2*pi*f*i;

A = abs(double(subs(Hr,s,v)));plot (f,A)

v = [f’,A’];save frmla.txt v -ascii

‘——————————————————

% freqrb.m

% Obtaining frequency response function by modal analysis

% Run modal.m (Appendix 2) first to find modal matrices

modal; alpha = 0;beta = 0.005;C = alpha*M + beta*K;

% Frequency response function from p1 to x2

clc;pinp = [1;0];nout = 2;

% Obtain transformed mass, stiffness, damping, and force matrices

Mu = mod'*M*mod;Ku = mod'*K*mod; Cu = mod'*C*mod; pu = mod'*pinp;

% Assign frequency interval

fend = 10;ns = 200;f = 0:fend/(ns−1):fend;

for j = 1:1:length(f);

w = 2*pi*f(j);T = 1/f(j);

t = 0:T/64:T;ns = length(t);

% Find frequency response functions for transformed variables

for n = 1:1:ndof;

m = Mu(n,n);k = Ku(n,n);c = Cu(n,n);p = pu(n);

X = p/((k-m*w^2)+i*c*w);xu(n,1:ns) = real(X*exp(i*w*t));

end

% Calculate frequency response function

x = mod(nout,1)*xu(1,:)+mod(nout,2)*xu(2,:);A(j) = max(x);

end

plot (f,A);v = [f’,A’];save frmlb.txt v -ascii

Frequency response analysis in ANSYS

Drag and drop Harmonic Response in Toolbox window in WB to Solution in Modal window

Select ANSYS Mechanical window

Analysis Settings > Assign 10 Hz to Range Maximum > Assign 100 to Solution Intervals¹

Analysis Settings > Damping Controls > Assign 0.005 to Stiffness Coefficient

Right-click Harmonic Response > Insert > Force > Click Face-4 (Figure 1(b))> Geometry in DW > Apply

Define Magnitude in DW as 1 N²

Right-click Harmonic Response > Solve > Solution

Right-click Solution > Insert > Frequency Resp.>Deformation > Click Face-6 (Figure 1(b))> Geometry in DW > Apply

Frequency Response > Display in DW > Amplitude > Chart Viewing Style > Linear

Right-click on Tabular Data window > Export > Save given file name (frwb)

Open frwb.txt file with notepad editor > Delete the first line

Impulse response analysis in MatLAB

% impr.m

clc,clear,close all

%Define mass, stiffness and damping matrices

m1 = 5;m2 = 4;k1 = 3000;k2 = 2200;alpha = 0;beta = 0.005;

M = [m1,0;0,m2];K = [k1 + k2,-k2;-k2,k2];C = alpha*M + beta*K;

f1 = 2.4934; fmax = 5.8359;

dt = 1/fmax/20,w1 = 2*pi*f1;ksi1 = alpha/(2*w1)+beta*w1/2; tss = 0.6/f1/ksi1

dt = 0.0086;tss = 6.14;

% Define the transfer function and find impulse response

syms s;H = inv(s^2*M + s*C + K);Hr = H(2,1);

[nomHr,denHr] = numden(Hr);nHr = sym2poly(nomHr);dHr = sym2poly(denHr);

sysHr = tf(nHr,dHr);

t = 0:dt:tss;t = t’;v = impulse(sysHr,t);

vv = [t,v];save impml.txt vv -ascii;plot (t,v)

Transient response analysis in ANSYS

Drag and drop Transient Structural in Toolbox window in WB to Solution in Modal window

Select ANSYS Mechanical window

Analysis Settings > Assign 0.0086 to Time Step > Assign 6.1439 to Step End Time

Analysis Settings > Damping Controls > Assign 0.005 to Stiffness Coefficient

Right-click Transient > Insert > Force > Click Face-4 (Figure 1(b))> Geometry in DW > Apply

Magnitude > Tabular Data > Enter Force values 0, 116.2791, 0,0 for the Time values 0, 0.0086, 0.0172, 6.1439, respectively

Right-click Transient > Solve > Solution

Right-click Solution > Insert > Deformation > Click Face-6 (Figure 1(b))> Geometry in DW > Apply

Right-click on Tabular Data window > Export > Save given file name (imprwb)

Open frwb.txt file with notepad editor > Delete the first line

Comparison of transient analyses by MatLAB and ANSYS

% compimp.m

clc,clear,close all

load impml.txt;t1 = impml(:,1);v1 = impml(:,2);

load impwb.txt;t2 = impwb(:,2);v2 = impwb(:,3);t2 = [0;t2];v2 = [0;v2];

plot (t1,v1,t2,v2)

Frequency response from impulse response and comparisons

% imprfft.m; Freq.resp. from impulse resp. using FFT; Comparison of freuency responses

clc,clear,close all

load impwb.txt;t = impwb(:,2);v = impwb(:,3);t = [0;t];v = [0;v];

dt = t(2)-t(1);fcut = 10;nfft = 2048;

fs = 1/dt;df = fs/(nfft−1);

fv = fft(v,nfft);ampl = abs(fv)*dt;

ncut = round(fcut/fs*(nfft−1) + 1,0);

f = 0:df:(ncut + 1)*df;f = f’;

A = ampl(1:length(f));

%—- compare frequency responses

load frmla.txt;fml = frmla(:,1);Amla = frmla(:,2);

load frwb.txt;fwb = frwb(:,2);Awb = frwb(:,3);

plot (fml,Amla,fwb,Awb,f,A)

Random vibration analysis in MatLAB

% rvr.m

clc,clear,close all

%Define mass, stiffness and damping matrices

m1 = 5;m2 = 4;k1 = 3000;k2 = 2200;alpha = 0;beta = 0.005;

M = [m1,0;0,m2];K = [k1 + k2,-k2;-k2,k2];C = alpha*M + beta*K;

% Define the transfer function

df = 0.01;fa = 0.4;fb = 1.4;fc = 8.4;fd = 9.4;

va = 5e−3;vb = 1e−2;vc = vb;vd = va;

syms s;H = inv(s^2*M + s*C + K);Hr = k1*H(2,1)−1;

f = fa:df:fd;f = f’;

n = length(f);v = 2*pi*f*i;aHr = abs(double(subs(Hr,s,v)));

% Define power spectral density of input as in Figure 6

for k = 1:1:n;

if f(k)<=fb;psdinp(k) = (f(k)-fa)/(fb-fa)*(vb-va)+va;

elseif f(k)>fa& f(k)<=fc;psdinp(k) = vb;

else f(k)>fc; psdinp(k) = (f(k)-fc)/(fc-fd)*(vc-vd)+vc;

end

k = k + 1;

end

% Find power spectral density of output, evaluate and save

psdout = psdinp’.*aHr.^2;

loglog (f,psdout);rmsv = sqrt(trapz(f,psdout))*1000

v = [f,psdout];save rvrml.txt v -ascii

Random vibration analysis in ANSYS

Drag and drop Random Vibration in Toolbox window in WB to Solution in Modal window

Select ANSYS Mechanical window

Analysis Settings > Damping Controls > Constant Damping > Manual > Assign Damping Ratio as 0

Assign Stiffness Coefficient as 0.005.

Right-click Random Vibration > Insert > PSD Acceleration

Boundary Condition in DW > Fixed Support

Tabular Data > Enter Acceleration values 0.005, 0.01, 0.01, 0.005 for the frequency values 0.4, 1.4, 8.4, 9.4, respectively¹

Direction > X axis

Right-click Random Vibration > Solve

Right-click Solution > Insert > Response PSD Tool

Right-click Response PSD Tool > Insert > Response PSD

Click Response PSD > Select a corner point on Face-6 (in Figure 1(b))>Geometry in DW > Apply

Response PSD > Result Type in DW > Acceleration > Reference > Relative to base motion

Right-click Solution > Evaluate all Results

Response PSD

Right-click on Tabular Data window > Export > Save given file name (rvrwb)

Open rvrwb.txt file with notepad editor > Delete the first line

Comparison of random vibration analyses by MatLAB and ANSYS

%comrvr.m

clc,clear,close all

load rvrml.txt;f1 = rvrml(:,1);a1 = rvrml(:,2);

load rvrwb.txt;f2 = rvrwb(:,2);a2 = rvrwb(:,3);

loglog (f1,a1,f2,a2)

References

Akdağ

Karagülle

Malgaca

. An integrated approach for simulation of mechatronic systems applied to a hexapod robot. Math Comput Simul 2012; 82: 818–835.

Deepak

Shetty

. Static and free vibration analysis of functionally graded rectangular plates using ANSYS. Mater Today Proc 2021; 45: 415–419.

Imran

Khan

Badshah

. Experimental, analytical, and finite element vibration analyses of delaminated composite plates. Sci Iran 2021; 28: 231–240.

Pappalardo

Lök

Şİ

Malgaca

, et al. Experimental modal analysis of a single-link flexible robotic manipulator with curved geometry using applied system identification methods. Mech Syst Signal Process 2023; 200: 110629.

Sahoo

Panda

Mahapatra

. Static, free vibration and transient response of laminated composite curved shallow panel–an experimental approach. Eur J Mech-A/Solids 2016; 59: 95–113.

Abramowicz

Chybiński

Polus

, et al. Free vibrations of sustainable laminated veneer lumber slabs. Sustainability 2023; 16: 166.

Yavuz

Akdağ

Karagülle

. A fast processing method to perform transient analysis for vibration control. Simul Model Pract Theory 2020; 104: 102152.

Bachoo

Elishakoff

. Random vibrations of laminated planar frames. J Sound Vib 2023; 562: 117860.

Chiariello

Vitale

Belardo

, et al. Non-linear analysis in post-buckling regime of a tilt rotor composite wing structure using detailed model and robust loading approach. J Compos Sci 2023; 8: 4.

10.

Hamedi

Marques

Vargo

, et al. Eco-friendly transportation: Solar car design for competition-part I, structural chassis. Open Access Repository 2023; 10: 16–35.

11.

Fan

Liang

Zeng

, et al. Dynamic characteristics of a rotating blade with a dovetail fixture. Machines 2023; 11: 643.

12.

Doori

SGM

Noori

Etemadi

. Static response of functionally graded porous circular plates via finite element method. Arab J Sci Eng 2024; 49: 141167–141181.

13.

Yashar

Ferguson

Ghandchi-Tehrani

. Simplified modelling and analysis of a rotating Euler-Bernoulli beam with a single cracked edge. J Sound Vib 2018; 420: 346–356.

14.

Rao

. Mechanical vibrations in SI units. 6th ed. Harlow, UK: Pearson Education Limited, 2018; p.1147.

15.

https://www.ansys.com/ .

16.

https://www.solidworks.com/tr .

17.

Schwarz

Richardson

. Proportional damping from experimental data. In: Topics in modal analysis, volume 7: proceedings of the 31st IMAC, a conference on structural dynamics. New York, USA: Springer, 2014, pp.179–186.

18.

Hacıömeroğlu

Uyulur

Güzel

. Automated vibration level assessment for airborne equipment. SAE Technical Paper (No. 2023-01-1143), 2023.