Abstract
The dynamic stiffness method (DSM) is proposed for vibration transmission from a machine with three degree of freedoms (DOFs) to its supporting beam structures. The machine is idealized as a rigid mass with three major degrees of freedom, namely, heave, roll, and pitch. It implies that the moments of inertia and torques due to unbalanced parts can be taken into consideration. The three types of vibration within beam structures, i.e., bending, longitudinal, and torsional motions, are formulated in terms of dynamic stiffness matrices. The finite element techniques can be applied similarly to assemble the developed dynamic stiffness matrices of both the machine and its supporting beam structures. A beam-like raft carrying a machine is designed to illustrate the accuracy of the proposed method in numerical simulation, where the differences brought by three-DOF modeling in vibration transmission analysis are discussed as well. This work would provide a novel and easy-to-use alternative to the existing mobility method and finite element method due to low discretization requirements, high efficiency, and high accuracy.
1. Introduction
For current passive isolation applications extensively found in a variety of engineering fields, three distinct configurations are employed, including a single-stage, a double-stage, or a floating raft isolation system [1, 2]. It is illustrated that two-stage isolation systems can achieve better vibration isolation performances than one-stage ones, especially at high frequencies [3, 4]. Distinguished from classical two-stage isolation, floating raft systems, which allow multiple hosts and auxiliary machines to be mounted on them, were developed and preferable for ships and submarines in recent years due to their effective isolation behaviors [5, 6]. Among different types of floating raft systems, beam-like rafts have received, increasing attentions currently. Compared with conventional plate-like rafts, the beam-like rafts have much more configuration flexibility and more space to accommodate machines and pipelines within them. Additionally, beam elements can be better maintained and readily replaced with dampers or active actuators [7].
Hence, it is of high necessity to investigate the isolation performance of a beam-like raft with operating machines. A typical floating raft system generally consists of three substructures, namely, a source (such as a machine), a raft, and a receiver (or a supporting structure), which are connected with each other by resilient mounts. A simple and intuitive solution is to idealize these three parts as a three-degree-of-freedom rigid body system [8]. However, with the development of raft systems towards large-scale, low-rigidity, and flexibility, larger errors from dynamic isolation models can be found if rigid approximation of the raft frame is still applied [9]. Besides, the existence of scattering of vibration energy within a raft implies that a small amount of damping material could potentially be very effective [10]. Based on the flexibility assumption, some researchers developed impedance, mobility, and transfer matrix methods to formulate the dynamical coupling interaction between substructures connected with each other in series with multiple inputs and outputs [5, 9, 11–25]. However, in impedance- or mobility-based method, the raft itself was commonly idealized as a simple structure, such as a free-free beam, to obtain the transfer matrix analytically between its upper input ends and lower output ones. Such simplification obviously cannot be applied in analyzing real engineering problems, as a raft can be complexly configured. Therefore, numerical methods are extensively utilized in such problems, as represented by the finite element method (FEM) [26]. Nevertheless, as the vibrations of these structures are mostly analyzed in a wide frequency range, it is always preferred to use an analytical method instead of cumbersome finite element modeling. To solve such a dilemma, a proper analytical method is strongly demanded.
Recently, the dynamic stiffness method (DSM) has attracted increasing attention worldwide as a powerful alternative to those methods mentioned above. DSM allows an assembly technique similar to that used in FEM, which enables free vibration analysis for either a single structural element or a combination of elements with different orientations. Besides, as an analytical method itself, the results derived from DSM are exact and independent of the number of elements used in the analysis [27]. There are abundant studies of DSM for different types of beams that can literally be found [28–34], yet seldom focused on the applications to floating raft systems. The Bernoulli–Euler beam carrying a two-degree-of-freedom spring-mass system was investigated to obtain the natural frequencies [27]. However, it paid attentions to only the influence of an attached mass on the vibration characteristics of a beam structure by replacing the mass with equivalent stiffness coefficients. An exact dynamic stiffness method was proposed for the free vibration analysis of multibody systems consisting of flexible structures and rigid bodies, where the attached rigid bodies are treated as specialized boundary conditions for beam or plate elements [35, 36]. Subsequently, another solution for beams carrying spring-mass systems was put forward and analyzed by DSM [37], where the assembling technique from FEM was simply but significantly utilized in the assembly procedure for dynamic stiffness matrices of a cantilever beam with a spring-mass system. Although only the heave DOF of a mass is taken into account, it lays the foundation of dynamic stiffness formulation for a beam structure carrying a machine. Inspired by Banerjee’s work, Li et al. developed a self-contained DSM for a machine supported with three-dimensional complex beam structures, where a machine was still simply idealized as a one-DOF spring-mass system [7]. However, one-DOF modeling apparently seems oversimplified for a real machine, as a machine normally idealized as a rigid body has three translation DOFs and three rotation DOFs and is supported on multiple mounts. Hence, there are practical demands to develop DSM for a floating raft carrying machines with multiple degree of freedoms. Meanwhile, questions may be raised and should be answered urgently, such as what are the differences in isolation performance prediction caused by modeling a machine as a mass with multiple DOFs instead of a single DOF?
From the view of assessment index for isolation performance, transmissibility, vibration level difference, and insertion loss are commonly used. Some researchers argued that the reaction force at each mount is typically different in terms of phase and direction, so one may find it difficult to directly use them to define an objective function for the purpose of design optimization [38]. Vibrational power analysis from a machine into its supporting structure seems attractive to researchers, which has been widely used to assess the effectiveness of passive and active vibration isolation systems [38]. Although the power flow is indeed a comprehensive index for isolation performance, the transmitted forces are underlying factors that trigger supporting structures such as ship hulls to radiate underwater noise.
Therefore, the analysis of vibration transmission from a machine through a beam-like floating raft to a receiver by the dynamic stiffness method is proposed in this paper. The first step is to build the dynamic stiffness formulation of beam structures carrying a machine with roll, pitch, and heave DOFs supported on four and more mounts. The roll, pitch, and heave DOFs are the main concerns instead of the other three DOFs, as they greatly affect the transmitted vertical forces to the wet surfaces, inducing sound radiation. The dynamic stiffness matrix for a given raft can be obtained by the same method described in the previous work [7]. In the subsequent numerical simulation, the capability and accuracy of the presented method are demonstrated by comparison with FEM. The differences between three-DOF modeling for a machine and one-DOF modeling are later illustrated in terms of transmitted forces, by taking an eccentric machine as an example. Finally, machines with six or eight mounts are investigated as well, proving the universality of the proposed method.
2. Modeling and Dynamic Stiffness Formulation for a Machine with Multiple Mounts
A machine can be idealized as a rigid cuboid body with mass m and heave, pitch, and roll DOFs, namely, illustrated in Figure 1. The center of mass is located at o, the distances from which to front, back, left, and right sides are l1, l2, l3, and l4, respectively. Any of the four resilient mounts can be idealized as a combination of a spring ki (i = 1, 2, 3, 4) and a damper ci (i = 1, 2, 3, 4). The upper ends of mounts connected with the machine have vertical displacements zsi (i = 1, 2, 3, 4), while zui (i = 1, 2, 3, 4) are the vertical displacements at the lower ends.

Apparently, the motions of , and zsi (i = 1, 2, 3, 4) are not mutually independent. In the case of small oscillation, the motional relationship between them should be
Meanwhile, the governing dynamic equations of the given machine in Figure 1 arewhere and are the corresponding moments of inertia about y and x, respectively. Besides, for any of the four massless mounts, the following equation is satisfied:
After substitution of Equations (1) into (2) and (3), one can yield the following equation in matrix form:where . M, C, and K are symmetric, which denote the mass, damping, and stiffness matrices, respectively. For simplicity, the upper triangular part of M, C, and K are shown in Equations. (5), (6), and (7). One can find that Equation (4) establishes the dynamic relations for the machine and the four lower mounts ends:
The elements in C and K are listed below:
Thus, for harmonic oscillation, the dynamic stiffness matrix for the above machine iswhere is circular frequency of oscillation .
The dynamic stiffness matrices for a machine with six or eight mounts can be similarly obtained, which are listed in Appendix A and B.
3. Dynamic Stiffness Formulation for Beam Vibrations
The dynamic stiffness matrix of a Bernoulli–Euler beam element can be readily obtained in [7, 39, 40]. Two beam elements are rigidly joined at their junctions, between which a lumped mass can possibly be deployed. Each element (shown in Figure 2) in a built-up beam structure can support at least four vibrational types: compressional, torsional, and two flexural vibrations [10].

For beams undergoing transverse vibration in the z-direction, the following equation is satisfied:where is the transverse displacement in z-direction, E and represent Young’s modulus and density, respectively, denotes the flexural rigidity in y-direction, and is the mass per unit length. , wherein the damping of material is taken into consideration, with as the damping loss factor and as the damping ratio.
Without loss of generality, the solution can be readily expressed in the following form:where and can be written aswhere and .
For the sake of convenience, the time dependence of harmonic term is suppressed. Subsequently, one can derive the rotation and the force (moment) of flexural motion in z-direction:
Thus, the relationship between nodal forces and displacements can be given bywhere is the nodal displacement vector, is the nodal force vector, and denotes the dynamic stiffness matrix for the bending vibration in the x-z plane. Similarly, the dynamic stiffness matrix for transverse vibration in y-direction can be readily obtained.
For vibration in x-direction, the governing equation iswhere is the displacement in x-direction and is the longitudinal wave speed, equaling to .
As for beams with circular sections, the torsional vibration is normally uncoupled with bending vibration, which can simplify the derivation of the given governing equation:where is the torsional wave speed and is equal to . G is the shear modulus. By a similar mathematical manipulation, the dynamic stiffness matrices for longitudinal and torsional vibration in the x-direction mentioned above can be readily obtained.
By establishing the relationship between the nodal force vector and displacement vector for all three types of vibrations, the dynamic stiffness matrix for a beam element can be derived as
Withwhere is the vector for the general displacements at the two ends of a beam element, while is the corresponding general forces. The expression of the elemental dynamic stiffness matrix will not be illustrated in detail here for sake of simplicity, which can be still viewed in [7].
Suppose a lumped mass is attached at one end of a beam element; its dynamic stiffness matrix can be obtained by investigation of the general inertia forces [7]:where subscript lm denotes the lumped mass. The stiffness matrix in equation (20) of a lumped mass should be assembled into that of a beam as one lumped mass is deployed at the end of a beam.
In order to obtain the dynamic stiffness matrix for an entire beam structure, the local stiffness matrix shall be transformed into global coordinates and assembled into an overall stiffness matrix according to the layout of a given structure. The transformation and assembly procedures are similar to standard finite element techniques, which deduce the global dynamic stiffness matrix aswhere and are the global and elemental stiffness matrices, respectively, is the transformation matrix for from global coordinates to local ones, and M is the number of beam elements.
Meanwhile, the stiffness matrix assembly process for a machine and a raft should be the same as FEM, as long as one can discover that the lower mount ends are in fact overlapped with the corresponding joint nodes of the raft. It means the DOFs of are shared and coupled between the raft and the machine, while the DOFs of the machine are incorporated independently. The dynamic stiffness matrix assembly procedure is depicted in Figure 3, where the blue blocks represent a beam structure including its dynamic stiffness matrix and corresponding DOFs, while the green ones are for a given machine.

4. Numerical Simulation
4.1. A Machine with Four Mounts
A simple raft carrying a machine with four mounts is investigated, shown in Figure 4. The raft colored green is mounted on two supporting beams colored blue via four springs. Dampers are not taken into consideration. The machine, idealized as a cuboid is mounted on the raft. The geometrical and material parameters are listed in Table 1, where r is the radius of the circular section of any beam element. The nodes are numbered 1 to 24. Apart from the global coordinate O-X-Y-Z, a set of local coordinate o-x-y is located at the geometrical center of the machine. The translational and rotational displacements along or around X (x)-axis, Y (y)-axis, and Z-axis for the node i on the beams are denoted as . The vertical displacement at the mass center and the rotations about the X (x)-axis and Y (y)-axis for the machine are denoted as , , and , respectively.

As to examine the accuracy of the proposed method, FEM is utilized to model the structure and obtain the dynamic behaviors as well. In FEM modeling, the beam element type is B31 in ABAQUS, the mesh size of which is 0.02 m. It is worth mentioning that, in FEM modeling, the machine is modeled as a rigid cuboid mass as well.
4.1.1. A Machine with Uniform Mass Distribution
In this case, the center of mass of the machine is localized at the point o in the local coordinates, on which the external harmonic excitation is also applied. The geometrical and mass-related parameters are listed in Table 2, where MC is denoted as the mass center. In this case, MC is coincided with the geometrical center, denoted as GC.
It is shown in Figure 5 that the dynamic responses of the vertical displacement at MC, Node 4, and Node 22 when a concentrated unit force F are applied on MC. A perfect agreement can be found between DSM and FEM in the frequency range up to 1 kHz. It can be seen from Figure 5 that the response amplitudes are decreased with frequency from the view of overall trend. In addition, resonant peaks from an upper substructure can be transmitted to and found in the lower one, which can be observed from comparisons between the figures in Figure 5.

(a)

(b)

(c)
The machine can also be subjected to a unit torque My about Y-axis at MC to simulate the excitation due to an unbalanced component. The dynamic responses at MC and Node 9 obtained by DSM and FEM are illustrated in Figure 6 where the results from DSM are excellently overlapped with those from FEM. An overall decreasing trend can be still observed in the dynamic response curves. The vertical displacement at Node 9 shares the same resonant peaks with θx because is attributed to the rotation motion due to the torque excitation.

(a)

(b)

(c)
4.1.2. A Machine with Uneven Mass Distribution
(1) Dynamic Responses at the Machine and a Node. In this case, uneven mass distribution of a machine is taken into consideration, which greatly leads to the noncoincidence of applied forces and MC. The geometrical and physical parameters of a given machine are listed in Table 3, the size and mass of which are kept consistent with those in Table 2. Due to the uneven mass distribution, the moments of inertia and localization of MC have been changed accordingly. The concentrated force in Z-direction is still applied on (0,0), namely, GC. The force applied on the GC can be equivalent to a concentrated force and two torques about x- and y-axis on the MC. The dynamic responses at MC and Node 4 obtained by DSM and FEM are demonstrated in Figure 7, respectively. An excellent agreement can be found between results from DSM and FEM, especially below 600 Hz, which proves the accuracy of the proposed method in solving dynamic responses of a raft system carrying an eccentric machine.

(a)

(b)

(c)

(d)

(e)

(f)
(2) Spring Forces between the Raft and the Supporting Beams. As far as ship sound radiation is concerned, the forces applied on the ship-hull structures receive wide concern. As the ship hull is simplified as two supporting beams, the forces of the springs connecting the raft and the supporting beams are worth being investigated. As to answer the question in Section 1 about the influences on dynamics prediction brought by the differences of modeling for a machine, the results obtained from the one-DOF modeling are plotted as well as compared. As long as a concentrated force is applied on MC, only heave DOF is activated so that one-DOF modeling can be implemented.
The forces of the springs attached with Node 22, 18, 23, and 19 are denoted as Fd1, Fd2, Fd3, and Fd4. It can be observed from Figures 8(a) to 8(d) that each of the spring forces from three-DOF modeling are diverse from those from one-DOF modeling due to the additional torque excitation. However, the amplitude of resultant forces Fdr from three-DOF modeling agrees quite well with those from one-DOF modeling, which implies that the torque excitations do not make any difference in terms of resultant spring forces. It can be also proved from phase analysis in Figure 9, where the phase of the resultant force from three-DOF modeling is almost identical with that from one-DOF modeling. It means unbalanced torque excitations from three-DOF modeling affect the separate spring forces only and make little influence on the resultant force.

(a)

(b)

(c)

(d)

(e)

(a)

(b)
The phenomenon illustrated above is most likely caused by the symmetry of the whole structure. In other words, the machine is located on the center of the raft. If the machine is offset and hypothetically translated to Node 9, 13, 10, and 14 and the dynamic response is reinvestigated in the same loading case, the resultant force Fdr and its phase can still be analyzed, as shown in Figure 10. It can be found that the two curves are overlapped with each other only in the frequency range of 300 Hz to 450 Hz and are different from each other in the rest of the entire frequency range. This implies that, in many cases where the machine is not located at the center of the symmetric raft, overlooking the rotation DOFs of a machine motion will significantly lead to the inaccurate simulation of transmitted forces.

(a)

(b)
(3) Optimization of Isolators. The discrepancies of separate spring forces between three-DOF modeling and one-DOF modeling may bring different isolator optimization strategies. The case in which the machine is mounted on the center of the raft is taken as an example. It can be supposed that the raft is mounted on air springs, the stiffness of which can be changed via air inflation and deflation. The springs attached with Node 22, 18, 23, and 19 are designated to have stiffness of , respectively, where are the adaptive stiffness coefficients varying in the assumptive interval of [0.5, 1.2] with increment 0.1. As to control the peak value of radiated noise in far-field, the objective function of the optimization process is to minimize the maximum of Fdr in the analysis frequency range. By the exhaustion method, the optimal solution for in three-DOF modeling and one-DOF modeling can be obtained, as shown in Table 4, where Fdr-max denotes the maximum amplitude of Fdr. It can be discovered that the optimized results from three-DOF modeling are better than those from one-DOF modeling by naturally taking advantage of the phase differences between spring forces. In other words, even when the entire structure is geometrically symmetric so that Fdr is not affected by the additional torque excitation, the three-DOF modeling provides different and more precise optimization schemes and optimal outputs.
4.2. Two Operating Machines
Two identical machines listed in Table 3 are mounted on Node 9, 13, 10, and 14 and Node 3, 15, 12, and 16, respectively. The two machines are both operating, the excitation force of which can be simulated by applying a concentrated force on the GC of each machine, respectively. ξ changing from -180° to 180° is the phase difference between the two concentrated forces, which is hard to identify in real engineering problems. The resultant force transmitted from two machines to the supporting beams obviously varies with ξ. For the sake of conservative estimation of sound radiation, engineers always prefer to estimate the maximum values of the resultant force (denoted as Fdr as well). It is shown in Figure 11 that the dynamic response contours of Fdr with respect to ξ by three-DOF modeling and one-DOF modeling in DSM. The major contour is the result obtained from three-DOF modeling and the minor one on the upper-right corner is from one-DOF modeling. The black line on the major contour and the red line on the minor contour are the tracing lines of the maximum values of Fdr with respect to frequency. It can be observed that, in one-DOF modeling, the maximum values of Fdr are taken at ξ = 0 in a wide span of frequency, which means excitations with identical phases provide the maximum values of Fdr. On the contrary, in three-DOF modeling, the maximum values of Fdr are taken at varying ξ, instead of ξ = 0. It means the maximum values of Fdr are not always taken at ξ = 0 when the two rotation DOFs of machine motions are taken into consideration. The black dashed line and the red dotted line are the projections of black and red tracing lines and plotted together as a comparison. Two lines are partially overlapped with each other possibly owing to the symmetrical localization of the two identical machines. However, the differences can still be found in some frequency ranges, such as 700 Hz to 900 Hz. The red dotted line is generally below the black dashed line, which means the three-DOF modeling provides a more conservative estimation in case of applied resultant forces on the supporting structure.

4.3. A Machine with Six Mounts
In this section, a machine with six mounts is investigated to verify the feasibility and accuracy of the proposed modeling. The machine is mounted on Node 10, 14, 11, 15, 12, and 16, which is still idealized as a cuboid mass. The local coordinate is at the GC of the machine, the x- and y-axes of which are in the same directions with X- and Y-axes in global coordinate. The physical parameters of the machine are listed in Table 5. A concentrated force in the Z-direction is applied to the geometrical center of the machine, namely, (0, 0) in the local coordinate. Dynamic responses are demonstrated in Figure 12, including two DOFs at the MC of the machine and at Node 10. It can be observed that the results obtained from DSM agree well with those from FEM in a general trend, in spite of some discrepancies in the frequency range of 10 to 150 Hz.

(a)

(b)

(c)

(d)
4.4. A Machine with Eight Mounts
A machine with eight mounts is still idealized as a cuboid, which is mounted on Node 9, 13, 10, 14, 11, 15, 12, and 16. The local coordinate is at the geometrical center of the machine, the x- and y-axes of which are in the same directions with X- and Y- axes in the global coordinate. The physical parameters of the machine involved in the calculation are listed in Table 6. Similarly, a concentrated force in the Z-direction is applied on the geometrical center of the machine, namely, (0, 0) in the local coordinate. In Figure 13, the dynamic responses at the MC of the machine and Node 10 are illustrated. The perfect agreement can be found between results obtained from DSM and FEM, which proves the accuracy of the proposed method in modeling and analyzing the dynamics of the machine supported on multiple mounts.

(a)

(b)

(c)

(d)
5. Conclusion
The dynamic stiffness method is presented for vibration transmission from a 3-DOF machine to beam structures. Dynamic responses caused by complex loading cases including vertical forces and torques can be accurately predicted due to taking the rotation DOFs of a machine into consideration. The results in numerical simulation show the excellent agreement between DSM and FEM, as to prove the accuracy of the proposed method. In addition, it is illustrated that three-DOF modeling instead of one-DOF brings different but more precise prediction of dynamic responses including the resultant forces transmitted to the supporting structures. It is also shown that the optimization strategy on vibration control should make the corresponding changes, such as optimization for stiffness-adaptive isolators. It can be foreseen that such practical modeling for beam rafts carrying machines is welcomed in analyzing vibration transmission and isolation behaviors in real engineering problems.
Appendix
A. A Machine with Six Mounts
The six mounts of a machine are equidistantly distributed along the x-axis direction, as shown in Figure 14. The 5th and 6th mounts are in middle of 1st, 3rd and 2nd, 4th, respectively. The following equations can be derived based on linear interpolation:where the definitions of zs5, zs6, zu5, and zu6 are the same as the previous zsi and zui, (i = 1, 2, 3, 4). Similar to equation (2), the governing equations for dynamics of a given machine with six mounts are

Thus, M, C, K, and q in equation (4) are changed accordingly. Based on the observation on equations (6) to (9), great similarity can be found between the elements of C and K. Hence, C and its elements are not listed in detail anymore for sake of simplicity, which can be readily obtained by replacing ‘k’ in K with ‘c’:where diag is the diagonalization function. K and its elements are formulated in equations (A.1) and (A.2). The dynamic stiffness matrix can still be obtained according to the following equation:
B. A Machine with Eight Mounts
Based on the previous studies of a machine supported on four or six mounts, it is easily extended to the case in which the machine is supported on eight mounts, as shown in Figure 15. It is supposed that the isolators are equidistantly distributed along the x-axis, and the following equations can be obtained, including displacement relationship, stiffness, and mass matrices. It is still worth mentioning that, in cases of the mounts which are not uniformly distributed, the K, C, and M matrices can be deduced from linear interpolation of the displacements.

Data Availability
Some or all data, models, or code generated or used during the study are available in a repository or online in accordance with funder data retention policies (provide full citations that include URLs or DOIs.)
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
The authors wish to thank High-Tech Ship Fund from the Ministry of Industry and Information Technology (MIIT): Deepwater Semi-submersible Support Platform (Grant no. 2016 [546]).