## Abstract

Riveted connections are widely used to join basic components, such as beams and panels, for engineering structures. However, accurately modeling joined structures with riveted connections can be a challenging task. In this work, an accurate linear finite element (FE) modeling method is proposed for joined structures with riveted connections to estimate modal parameters in a predictive manner. The proposed FE modeling method consists of two steps. The first step is to develop nonlinear FE models that simulate riveting processes of solid rivets. The second step is to develop a linear FE model of a joined structure with the riveted connections simulated in the first step. The riveted connections are modeled using solid cylinders with dimensions and material properties obtained from the nonlinear FE models in the first step. An experimental investigation was conducted to study accuracy of the proposed linear FE modeling method. A joined structure with six riveted connections was prepared and tested. A linearity investigation was conducted to validate that the test structure could be considered to be linear. A linear FE model of the test structure was constructed using the proposed method. Natural frequencies and corresponding mode shapes of the test structure were measured and compared with those from the linear FE model. The maximum difference of the natural frequencies was 1.63% for the first 23 out-of-plane elastic modes, and modal assurance criterion values for the corresponding mode shapes were all over 95%, which indicates high accuracy of the proposed linear FE modeling method.

## 1 Introduction

To build complex engineering structures, various types of fasteners are needed to join various components together. While the use of high-strength bolts becomes more and more popular, riveted joints are still one of the most commonly used fasteners nowadays. Specifically, riveted joints are mainly applied in industry such as aircrafts, bridges, and building frames [1,2]. Hence, assessing, maintaining, and retrofitting existing components with riveted joints are important for structural designers and engineers. They need to be knowledgeable about riveted joints [1]. A riveted joint can consist of one or multiple riveted connections. There are different types of rivets in daily applications, including solid, self-piercing, and blind rivets, and solid rivets are one of the most reliable types of rivets [1,3]. One typical application of solid rivets can be found in manufacturing structural components of aircraft by layering two or three sheets of materials and then joining them with solid rivets [4].

In recent studies, numerous finite element (FE) analysis results of structures with riveted connections have been presented. FE models of structures with riveted connections can be developed in two ways: three-dimensional and two-dimensional axisymmetric modeling. In three-dimensional modeling, detailed FE models with large numbers of elements are developed to simulate crashworthiness and joint failure tests of structures with riveted connections, where plastic deformations are considered. Langrand et al. [5] proposed a FE modeling method to investigate joint failure tests in airframe crashes using three-dimensional FE models. The FE modeling method can be used to optimize mechanical properties of joint components specifically for airframe crash simulations. Rans et al. [6] developed a FE model to study residual stresses beneath a rivet head. The developed FE model is three dimensional and quarter symmetric, where eight-node brick elements were used. By simulating a riveting process, residual stresses beneath a rivet head were studied to show that they could be influenced by compression of joined sheets. Sadowski et al. [7] investigated a FE model of a hybrid joint with adhesive and a rivet to study effects of adding a rivet to an adhesive joint. In a developed FE model, a large number of elements were used to accurately evaluate large plastic deformations and contacts during a riveting process. By using adhesive riveted joints, the energy absorption could be increased by about 35% compared with an adhesive joint without a rivet. Xiong and Bedair [8] proposed a FE modeling method for stress analysis of riveted lap joints. The method used cap and spring elements to model riveted connections, and it was capable of dealing with relatively complicated geometry and loading conditions. Bedair and Eastaugh [9] later extended the FE modeling method for structures with riveted splice joints, where effects of their secondary out-of-plane bending and interactions between plates and rivets were considered. Fung and Smart [10] presented numerical parametric studies on effects of clamping forces from riveted connections, contacts between rivets and clamped components, and frictions between contacted surfaces. However, riveting processes of the studied connections were not considered. Parameters in the numerical studies might not be realistic. One common issue of using three-dimensional FE models is that computation costs can be too high to afford due to large numbers of degrees of freedom in the models.

The use of two-dimensional axisymmetric FE models can reduce numbers of degrees-of-freedom and enable accurate studies of effects of dimensions of rivets and applied loads during a riveting process. Hence, they can be considered to be computationally efficient with high accuracy as three-dimensional FE models. Quality of riveted connections can be affected by many parameters such as the applied load during the riveting process, length of a rivet shank, diameter of the shank, and matching hole diameter tolerance. Cheraghi [11] investigated effects of various aforementioned parameters in a countersunk riveting process using two-dimensional axisymmetric FE models. By decreasing the countersunk depth to 8.13 × 10^{−4} m from a recommended depth of 1.07 × 10^{−3} m, higher allowable ranges of matching holes and rivet diameters could be obtained. The range of the applied load could also be increased without deteriorating quality requirements. Zhang et al. [12] proposed a mathematical model for a riveting process and verified it by using a two-dimensional axisymmetric FE model. The mathematical study of a riveting process could accurately estimate deformations of thin-walled sheet-metal parts. Blanchot and Daidie [13] conducted an intensive numerical parametric study on riveting processes and mechanical properties of joined structures. It was found that riveted connections could affect mechanical properties of a joined structure when it is subject to external loading.

In contact models for connections in most mechanical joints, slip behaviors need considering since they can introduce nonlinearities to structures with the joints. Segalman et al. [14] and Blake [15] investigated the mechanics of micro- and macro-slip behaviors and experimentally confirmed the existence of slip behaviors. Lei et al. [16] presented a numerical study for riveting processes of rivets that join components with countersunk holes. Nonlinear two-dimensional axisymmetric FE models were constructed to simulate the riveting processes, where surface contacts and plastic deformations were considered. Compared with three-dimensional FE models, two-dimensional axisymmetric FE models enable the use of finer meshes that can yield more detailed and accurate simulation results in an efficient manner. Currently, most studies that presented FE models of structures with riveted connections have focused on nonlinearity-related mechanical behaviors of structures with the connections, including crashworthiness of structures with riveted joint elements, large deformations of the structures, and effects of applied loads and joint dimensions during a riveting process. A few studies have focused on modal properties, such as natural frequencies and mode shapes, of structures with riveted connections. He et al. [17] presented a numerical investigation on natural frequencies of two beams joined by a self-pierce rivet and found that its natural frequencies increase with an increase of the elastic modulus of the beams, but they would slightly change when Poisson’s ratio of the beams was changed. Koruk and Sanliturk [18] and Altuntop et al. [19] developed and applied linear FE modeling methods for structures with riveted connections. They estimated dimensions and material properties of elements, which were used to model riveted connections, by model updating so that differences between modal parameters from FE models and experimental modal analysis were minimized. So far, a linear FE modeling method for structures with riveted connections, which can be used to accurately estimate their modal parameters in a predictive manner, is not available and needs developing. He and Zhu [20] developed a linear FE modeling method for joined structures with bolted connections, where solid cylinders were used to model bolted connections. In the developed method, nonlinear two-dimensional axisymmetric FE models were constructed to estimate dimensions and material properties of the solid cylinders, and a linear FE model was then developed with structural components joined by the cylinders. The developed linear FE model was capable of accurately estimating modal parameters of the structures with bolted connections in a predictive manner.

This work extends the linear FE modeling method developed by He and Zhu [20]. A linear FE modeling method is proposed for joined structures with riveted connections. A test structure, which consisted of two aluminum strips joined by six aluminum solid rivets, was prepared. A linearity investigation of the test structure was conducted to verify that it could be considered to be a linear structure. The proposed FE modeling method for a joined structure with riveted connections consists of two steps; solid cylinders are used to model the riveted connections. The first step is to develop nonlinear two-dimensional axisymmetric models to simulate riveting processes of the riveted connections and determine dimensions and material properties of the solid cylinders. In the nonlinear models, plastic deformations and surface contacts are considered. The second step is to develop a linear FE model of the joined structure, where the riveted connections are modeled using solid cylinders with dimensions and material properties obtained from the nonlinear models in the first step. Experimental investigation of the proposed linear FE modeling method was conducted, where natural frequencies and modes shapes of the test structure were estimated. Modal phase collinearity (MPC) values of the estimated mode shapes were calculated to validate that mode shapes of the test structure were real. The estimated natural frequencies and mode shapes were compared with those from a linear FE model of the test structure, which was developed using the proposed method.

The remaining part of this paper is outlined as follows. The preparation of the test structure with riveted connections is described in Sec. 2.1. The investigation of linearity of the test structure is presented in Sec. 2.2. Methodologies for developing nonlinear and linear FE models for a structure with riveted connections are described in Secs. 3.1 and 3.2, respectively. Experimental estimation of material properties of the joined strips is presented in Sec. 4.1. Developments of a nonlinear FE model and a linear FE model are presented in Secs. 4.2 and 4.3, respectively. Results from the developed linear FE model were validated by those from an experimental modal analysis in Sec. 5. Conclusions of this work are presented in Sec. 6.

## 2 Test Structure and Linearity Investigation

In this section, a joined structure with riveted connections was prepared and tested, and it had two identical aluminum strips to be joined by six identical aluminum solid rivets. An experimental investigation of linearity of the test structure was conducted.

### 2.1 Test Structure Preparation.

One strip to be joined was made of aluminum 6061-T6511. Dimensions and matching hole positions of the strip are shown in Fig. 1; the strip had a thickness of 4.75 × 10^{−3} m. The strip had six circular holes, and each of the holes had a diameter of 6.44 × 10^{−3} m. Dimensions of one solid rivet to join the strips are shown in Fig. 2. The rivet had a round-domed rivet head and a cylindrical rivet shank, and definitions of the rivet head and rivet shank are depicted in Fig. 2. Dimensions of the strips and rivets conformed to an engineering criterion for solid rivet installations: the diameter of a rivet shank needs to range from one to three times the thickness of a strip to be joined [15].

The test structure with the six riveted connections was prepared using a rivet installation kit shown in Fig. 3, including a McMaster-Carr style-B air-power hammer, a rivet setter, and a bucking bar. The rivet setter was a bit to be inserted into the barrel of the air-power hammer. The rivet setter was selected based on the shape and size of the rivet head. The rivet setter had a convex-rounded tip, which could engage a rivet head in the riveting process, and continuous loads from the air-power hammer could be evenly applied to the rivet head in the riveting process. Otherwise, the rivet would be damaged due to unevenly applied loads. The bucking bar was a solid iron block to be placed at the end of the rivet shank to provide support for the rivet.

Prior to the riveting process, the two strips were clamped by two bench vices, and they were joined by six bolts with hand-tight nuts. As shown in Fig. 4, the riveting process of one rivet was started by replacing one of the bolts with a rivet, while keeping other bolts in place to prevent slips between the strips that could lead to misalignments of matching holes. The rivet setter inserted into the air-power hammer and the bucking bar were placed onto the rivet head and the other end of the rivet shank, respectively. By applying continuous loads in the form of impacts to the rivet head using the air-power hammer and fixing the bucking bar, the protruding part of the rivet shank was gradually formed into a driven head, which concluded the riveting process of one solid rivet, as shown in Fig. 5. Other five rivets were installed in the same way, and preparation of the test structure was finished, as shown in Fig. 6(a). The test structure had a length of 1.13 m, and an enlarged side view of the riveted connections is shown in Fig. 6(b). Thicknesses of the driven heads of the six installed rivets were measured, and their average was 2.60 × 10^{−3} m.

### 2.2 Linearity Investigation of Test Structure.

The linearity of the test structure in Fig. 6(a) was experimentally investigated by conducting a linearity check and a reciprocity check. Free boundary conditions of the test structure were simulated by two identical rubber bands that supported its two ends, as shown in Fig. 6(a). A grid with 3 × 13 measurement points was assigned to the test structure. Numbering and positions of the points are shown in Fig. 6(c). An impact hammer PCB 068C05 was used to excite the test structure in the form of single impacts. A scanning laser Doppler vibrometer Polytec PSV-500-HV was used to measure response of the test structure at the measurement points. The excitation and response measurement were in −*z* and +*z* directions, respectively, as indicated in Fig. 6(c). In the linearity check, the test structure was excited at point 4 with two different impact magnitudes, as shown in Fig. 7(a), and responses at point 20 were measured by the vibrometer. Two frequency response functions (FRFs) between the two single impacts and corresponding responses were estimated and shown in Fig. 7(b). It can be seen that magnitudes of the two FRFs almost overlap in the measured frequency ranging between 0 and 2300 Hz. The frequency response assurance criterion (FRAC) value between the two calculated FRFs was 96.5% [21], indicating that they were highly correlated. A reciprocity check for the test structure was also conducted by comparing reciprocal FRFs associated with points 4 and 17. As shown in Fig. 8, magnitudes of the reciprocal FRFs overlap well, and their associated FRAC value was 94.4%. Based on the linearity and reciprocity check results, the test structure could be considered to be linear within the measured frequency range. It was implied that existence of the riveted connections would not effect the linearity of the test structure, and modal parameters of the structure could be numerically and experimentally estimated from a linear FE model and modal analysis, respectively. This implication is similar to that for joined structures with bolted connections, which can be considered to be linear within a certain frequency range [20].

## 3 Finite Element Modeling Methodology

Even though linearity of the test structure had been experimentally validated, it could not be immediately modeled by a linear FE model with high accuracy due to the riveting process, where surface contacts and plastic deformations had occurred. To develop an accurate linear FE model of a joined structure with riveted connections like the test structure, a two-step FE modeling technique is proposed. The first step is to develop nonlinear FE models that simulate the riveting process of the riveted connections, where surface contacts and plastic deformations are considered. The second step is to develop a linear FE model of the structure, whose riveted connections are modeled using linear solid cylinders. Dimensions and material properties of the solid cylinders are determined using results from the nonlinear FE models. Joined components of the structure are modeled using linear elements that are tied to the solid cylinders. In this work, the commercial FE analysis software abaqus is used for FE modeling and analysis.

### 3.1 Nonlinear Finite Element Models of Riveted Connections.

The objective of developing nonlinear FE models that simulate riveting processes is to obtain dimensions and material properties of solid cylinders, which are used to model riveted connections in the linear FE model. To achieve the objective, the developed FE models need to have high accuracy. However, developing such FE models can be challenging for the following three reasons. The first reason is that nonlinear plastic deformations occur to a solid rivet during its riveting process, where the protruding part of its shank is formed into a driven head, as shown in Fig. 5. Due to the plastic deformations, the diameter of a rivet shank increases. As a result, the clearance, which exists between the rivet shank and matching holes of joined components before the riveting process, is filled up. Besides, the volume of the rivet is changed by applied loads, so is the density of the installed rivet. The second reason is that surface contacts occur between the rivet and joined components and between the joined components themselves in the riveting process. The third reason is that the required number of degrees of freedom of an accurate nonlinear three-dimensional FE model for the riveting process can be too high to be computationally efficient. To address the first and second reasons, plastic deformations and surface contacts are considered in nonlinear FE models to be developed. To address the third reason, the nonlinear FE models use axisymmetric two-dimensional elements.

In a nonlinear FE model that simulates a riveting process, a rivet, a rivet setter, a bucking bar, and components to be joined are parts that constitute a riveted connection, and they are modeled using axisymmetric two-dimensional elements assigned with respective material properties, including mass densities, elastic moduli, Poisson’s ratios, and plastic stress–strain relations. To simulate the process, surface contacts between the parts and their plastic deformations are enabled. The modeled parts are assembled to appropriate positions. Boundary conditions are defined: the bottom surface of the rivet setter is radially and vertically fixed, and the edge of the rivet along the axis of symmetry is radially fixed, as shown in Fig. 9, where *x*- and *y*-axes represent radial and vertical directions in the cylindrical coordinate system, respectively. The riveting process simulated in the nonlinear FE model can be divided into two stages: a load forming stage and a load removing stage. The load forming stage is achieved by slowly moving the bucking bar downward along the *y*-axis to a desired position. The movement of the bucking bar is so slow that the load forming stage can be considered to be quasi-static. When the bucking bar arrives at the desired position and the driven head is under no load from the bucking bar, the protruding part of the rivet shank is formed into a driven head. The load removing stage is achieved by moving the bucking bar upward along the *y*-axis to another desired position, which concludes the development of the nonlinear FE model. From the developed nonlinear FE model, one can quantitatively observe the aforementioned surface contacts and plastic deformations of the rivet and joined components, with which dimensions and material properties of solid cylinders for a linear FE model can be determined. A flowchart that describes the development procedure of the nonlinear model is shown in Fig. 10.

### 3.2 Linear Finite Element Model of Structure With Riveted Connections.

Based on results from the nonlinear FE models of the riveted connections, a linear FE model of the joined structure can be developed. In the linear FE model, the riveted connections are equivalently modeled using solid cylinders with solid elements. Dimensions of the solid cylinders are obtained from the nonlinear FE models. The radius of one solid cylinder is equal to the sum of the radius of the matching holes of the joined components and an effective dimension of the contact area associated with the joined components. Specifically, the effective dimension is defined as the difference between the outer and inner radii of the contact area. The height of the solid cylinder is equal to the sum of thicknesses of the joined components, and it is denoted by *L*_{s}.

*M*

_{r}is the mass of a portion of the installed rivet, which is embedded in the matching holes,

*M*

_{j}is the mass of portions of the joined components that are in contact, and

*V*

_{s}is the volume of the solid cylinder, which can be calculated using the dimensions of the solid cylinder. Meanwhile, mass densities of the formed driven head and round-domed rivet head are assumed to be unchanged and equal to the mass density of the rivet before the riveting process. The equivalent elastic modulus of the solid cylinder can be calculated using Hooke’s law [15]:

*A*

_{r}and

*A*

_{j}are the cross-sectional area of the installed rivet and the contact area, respectively, and

*k*

_{s}is the equivalent stiffness of the solid cylinder. The equivalent stiffness can be obtained as a sum of stiffnesses of two parallelly connected springs:

*k*

_{r}and

*k*

_{j}are stiffnesses of the portion of the installed rivet and that of the joined components that are in contact, respectively. The two stiffnesses can be calculated by

*E*

_{r}is the elastic modulus of the rivet, and

*E*

_{j}is the elastic modulus of the joined components. Meanwhile, elastic moduli of the formed driven head and round-domed rivet head are assumed to be unchanged and equal to those of the rivet before the riveting process. Portions of the joined components beyond the contact areas are tied to the solid cylinders. The developed linear FE model of the joined structure can be used to accurately estimate its natural frequencies and mode shapes.

## 4 Numerical and Experimental Investigations

### 4.1 Dimensions and Material Properties of Strips.

Dimensions and material properties of the two identical strips of the test structure, including the elastic modulus, Poisson’s ratio, and mass density, were experimentally studied before the riveting process. Dimensions of one of the strips were measured, as shown in Fig. 1; the thickness of the strips was measured to be 4.75 × 10^{−3} m. The mass density of the strip was measured to be 2.67 × 10^{3} kg/m^{3}, which was calculated by dividing its total mass by its volume that was calculated based on its measured dimensions. Poisson’s ratio of the strip was 0.33, and the elastic modulus was 69.0 GPa, which were nominal values of aluminum 6061-T6511. A FE model of the strip was developed with 43,400 quadratic triangular shell (STRI65) elements. Free boundary conditions were applied to the strip in the FE model. Natural frequencies of the strip were calculated and listed in Table 1.

Mode | Measured (Hz) | Before updating (Hz) | Error (%) | After updating (Hz) | Error (%) |
---|---|---|---|---|---|

1 | 65.91 | 68.04 | 3.14 | 66.56 | 0.99 |

2 | 180.44 | 187.14 | 3.58 | 183.00 | 1.39 |

3 | 353.15 | 366.38 | 3.61 | 357.90 | 1.32 |

4 | 471.94 | 467.91 | −0.86 | 462.94 | −1.95 |

5 | 583.41 | 605.14 | 3.59 | 590.57 | 1.21 |

6 | 690.25 | 710.37 | 2.83 | 685.53 | −0.69 |

7 | 871.72 | 903.67 | 3.54 | 880.68 | 1.02 |

8 | 947.93 | 938.73 | −0.98 | 940.77 | −0.76 |

9 | 1218.25 | 1262.40 | 3.50 | 1228.70 | 0.85 |

Mode | Measured (Hz) | Before updating (Hz) | Error (%) | After updating (Hz) | Error (%) |
---|---|---|---|---|---|

1 | 65.91 | 68.04 | 3.14 | 66.56 | 0.99 |

2 | 180.44 | 187.14 | 3.58 | 183.00 | 1.39 |

3 | 353.15 | 366.38 | 3.61 | 357.90 | 1.32 |

4 | 471.94 | 467.91 | −0.86 | 462.94 | −1.95 |

5 | 583.41 | 605.14 | 3.59 | 590.57 | 1.21 |

6 | 690.25 | 710.37 | 2.83 | 685.53 | −0.69 |

7 | 871.72 | 903.67 | 3.54 | 880.68 | 1.02 |

8 | 947.93 | 938.73 | −0.98 | 940.77 | −0.76 |

9 | 1218.25 | 1262.40 | 3.50 | 1228.70 | 0.85 |

Note: Percentage errors between the natural frequencies from the FE models and the measured ones are listed.

An impact test of one strip was conducted to estimate its natural frequencies. The strip was hung using two cotton ropes to simulate free boundary conditions. Impacts and response were generated and measured using the same impact hammer and vibrometer as in Sec. 2.1, respectively. One impact point and nine measurement points were randomly selected. Nine FRFs associated with the nine measurement points and the impact point were measured. Each FRF was obtained with five averages. The nine FRFs were then analyzed by a modal parameter estimation algorithm PolyMax in LMS Test.Lab 17 to estimate natural frequencies of the first nine out-of-plane elastic modes of the strip, which are listed in Table 1. The free boundary conditions of the strip were validated, since the highest measured natural frequency of the rigid-body modes was 1.89 Hz, which is less than 10% of the natural frequency of the first elastic mode (65.91 Hz) [22].

By comparing with the natural frequencies from the FE model and the impact test, the maximum percentage error was 3.61%. The reason was that the nominal values of the elastic modulus and Poisson’s ratio of the strip in the FE model were not accurate. To reduce the error, the elastic modulus and Poisson’s ratio were updated and applied to the FE model. After updating the elastic modulus to 66.0 GPa and Poisson’s ratio to 0.30, the maximum percentage error for the natural frequencies of the nine elastic modes decreased to 1.95%, as shown in Table 1.

### 4.2 Nonlinear Finite Element Model Development.

To develop a nonlinear FE model of the riveting process of a riveted connection of the test structure in Fig. 6(a), five structural parts were developed: a bucking bar, two strips to be joined, a solid rivet, and a rivet setter. The rivet is meshed using 850 axisymmetric four-node bilinear (CAX4R) elements. Elements close to edges of the meshed rivet were finer than those close to the axis of symmetry since the former would undergo surface contacts during a riveting process. Each of the two strips were meshed with 720 CAX4R elements. In the FE model, the solid rivet, which was made of aluminum 1100, had a nominal mass density of 2.70 × 10^{3} kg/m^{3}, Poisson’s ratio of 0.33, and a nominal elastic modulus of 69.0 GPa. The strips, based on results from Sec. 4.1, had a mass density of 2.68 × 10^{3} kg/m^{3}, Poisson’s ratio of 0.30, and an elastic modulus of 66.0 GPa. The bucking bar and the rivet setter, which were made of steel, had a mass density of 8 × 10^{3} kg/m^{3}, Poisson’s ratio of 0.25, and an elastic modulus of 200 GPa. In addition, the rivet and strips had a plastic stress–strain relation of aluminum alloy [23]. The meshed parts were assembled to appropriate positions for a riveting process, as shown in Fig. 11(a), and a clearance between the rivet and strips can be observed in Fig. 11(b). Surface contacts were defined to numbered surface pairs: (1)–(2), (1)–(3), (3)–(4), (3)–(5), (3)–(8), (6)–(7), (9)–(10), (9)–(12), and (11)–(13), as indicated in Fig. 12. Normal and frictional surface contact properties of the above contact surface pairs were defined as follows. For the normal surface contact property, any penetration between two surfaces in contact was not allowed, and separation of surfaces that had once been in contact was allowed. For the frictional surface contact property, friction coefficients were defined: 0.61 for contacts between aluminum and steel surfaces and 1.05 for contacts between the aluminum surfaces [24]. The boundary conditions were defined: radial displacements of nodes on the edge of the meshed rivet along the axis of symmetry were constrained and both of vertical and radial displacements of nodes on the bottom of the meshed rivet setter were constrained, as shown in Fig. 9. Since the nonlinear model was developed with axisymmetric elements, it had a cylindrical coordinate system, where *x*- and *y*-axes denote the radial and vertical directions of the nonlinear FE model in Fig. 9, respectively.

The riveting process was simulated in two stages, including a load forming stage and a load removing stage, as described in Sec. 3.1. The load forming stage had a duration of 10 s, in which the bucking bar moved downward to a desired position, which was −7.22 × 10^{−3} m along the *y*-axis from its original position. In Figs. 13(a)–13(e), the load forming stage is visually depicted by the FE model at five different time instants. The load forming stage started at *t* = 0 s that corresponded to Fig. 13(a), and it ended at *t* = 10 s that corresponded to Fig. 13(e). Plastic deformations during the load forming stage could be observed at the protruding part of the rivet shank when the driven head was being squeezed by the bucking bar, as shown in Figs. 13(b)–13(d). Due to the plastic deformations, the clearance between the rivet and strips was filled up, and the radius of the rivet shank increased from 3.20 × 10^{−3} m to 3.22 × 10^{−3} m at *t* = 10 s. Contact pressure distributions between the two strips along the radial direction at the five time instants are shown in Fig. 14. Surface contacts were considered to occur within areas that had non-zero contact pressures. Note that the origin of the horizontal axis in Fig. 14 was at 3.42 × 10^{−3} m instead of 3.22 × 10^{−3} m, which was the radius of the hole on the joined strips at *t* = 0. The reason was that nodes of the two joined components at their holes were also in contact with the inserted part of the rivet shank, and one could not differentiate the contact pressure given by the rivet from that given by the strips. As shown in Fig. 14, the maximum contact pressure occurred at the nodes closest to the edge of the holes, and it gradually increased with time. At *t* = 10 s, the effective dimension of the contact area was measured to be 7 × 10^{−3} m, and the radius of the inserted part of the rivet shank was measured to be 3.328 × 10^{−3} m, which concluded the load forming stage.

The load removing stage had a duration of 10 s, in which the bucking bar moved upward to its original position in the first stage. As shown in Figs. 13(e) and 13(f), the load removing stage started at *t* = 10 s that corresponded to Fig. 13(e), and it ended at *t* = 20 s that corresponded to Fig. 13(f). Contact pressure distributions between the two strips along the radial direction at different time instants in the load removing stage are shown in Fig. 15. It could be observed that the contact pressures along the radial direction decreased with time, and the effective dimension changed from 7 × 10^{−3} m to 1.6 × 10^{−3} m at *t* = 20 s; the latter effective dimension was considered to be the final effective dimension of the contact area. The reason was that moving the bucking bar upward resulted in the removal of the applied loads and the formed rivet connection stretched to a new equilibrium position. The radius of the inserted part of the rivet shank decreased from 3.328 × 10^{−3} m to 3.315 × 10^{−3} m. Note that the origin of the horizontal axis in Fig. 15 was at 3.515 × 10^{−3} m instead of 3.315 × 10^{−3} m, which was the final radius of the inserted part of the rivet shank at *t* = 20 s in the load removing stage. The driven head had a thickness of 2.60 × 10^{−3} m at the end of the load removing stage (*t* = 20 s), as shown in Fig. 13(f), and the thickness was equal to the average measured thickness of the six driven heads of the test structure in Sec. 2.1. An extra stage of simulation was added here to study the effective dimension, and it had a duration of 10 s. In the extra stage, no parameters were changed since the load removing stage, and contact pressure distributions at *t* = 25 s and *t* = 30 s are shown in Fig. 15. The contact pressure distributions were observed to be unchanged throughout the stage.

A parametric study of the friction coefficient was conducted for the effective dimension of the contact area. The studied friction coefficient was the one associated with the contact between the aluminum parts, including the rivet and strips. The studied friction coefficient valued between 0.85 and 1.45 with an increment of 0.20. The contact pressure distributions between the two strips along the radial direction at *t* = 10 s and *t* = 20 s are shown in Figs. 16(a) and 16(b), respectively. In Figs. 16(a) and 16(b), the contact pressure distributions for four cases were compared well with each other along the radial direction, which indicates that the change of the friction coefficient in the range would not affect the contact pressure distributions. More importantly, the effective dimension of the contact area would not be changed by the change of the friction coefficient.

### 4.3 Linear Model Development.

Each of the six riveted connections in the test structure was modeled using a solid cylinder with a round-domed rivet head and a driven head using 10,695 ten-node quadratic tetrahedral (C3D10) elements, as shown in Fig. 17. Dimensions of the solid cylinders were obtained from the riveted connection in the nonlinear FE model for the riveting process. Specifically, the height of one solid cylinder was equal to twice the thickness of the strips, and its radius is equal to the sum of the final radius of the inserted part of the rivet shank (3.315 × 10^{−3} m) and the final effective dimension of the contact area between the two strips (1.6 × 10^{−3} m); the sum was equal to 4.915 × 10^{−3} m. The mass density of the solid cylinder was calculated to be 2.722 × 10^{3} kg/m^{3} based on Eq. (1); the elastic modulus of the solid cylinder was calculated to be 67.34 GPa based on Eq. (2). The mass densities and elastic modulus of the formed driven head and domed rivet head remained the same as those before the riveting process, which were 2.79 × 10^{3} kg/m^{3} and 69 GPa, respectively. Shell elements were used to model the two strips. A total of 4146 linear triangular shell (S3) elements were used for one strip, which was also the case for the other strip. Material properties of the strips, including the mass density, elastic modulus, and Poisson’s ratio, in the linear FE model were those experimentally obtained in Sec. 4.1. The solid cylinders are tied to the joined strips based on two reasons. One reason was that non-zero contact pressure distributions along the contact surfaces between the rivet shank and two strips were observed, which could lead to permanent contacts between the rivet and joined strips. The other reason was that from the results of the linearity investigation of the test structure in Sec. 2.2, nonlinear macro-slip behavior did not occur after the riveting process and nonlinearity caused by micro-slip could be neglected [14,15]. Natural frequencies of the first 23 out-of-plane elastic modes of the test structure in the linear FE model are calculated and listed in Table 2. Mode shapes associated with the first four out-of-plane elastic modes are shown in Fig. 18.

Mode | Calculated (Hz) | Measured (Hz) | Error (%) | Mode | Calculated (Hz) | Measured (Hz) | Error (%) | Mode | Calculated (Hz) | Measured (Hz) | Error (%) |
---|---|---|---|---|---|---|---|---|---|---|---|

1 | 18.56 | 18.37 | 0.99 | 9 | 465.17 | 459.29 | 1.28 | 17 | 1299.71 | 1317.51 | −1.35 |

2 | 51.83 | 50.93 | 1.75 | 10 | 589.41 | 582.10 | 1.26 | 18 | 1436.11 | 1454.76 | −1.28 |

3 | 98.88 | 97.65 | 1.26 | 11 | 752.96 | 743.44 | 1.28 | 19 | 1541.01 | 1520.60 | 1.34 |

4 | 167.29 | 164.61 | 1.63 | 12 | 773.65 | 785.023 | −1.45 | 20 | 1684.51 | 1672.17 | 0.74 |

5 | 247.56 | 244.49 | 1.25 | 13 | 892.97 | 883.87 | 1.03 | 21 | 1840.34 | 1863.54 | −1.24 |

6 | 256.84 | 260.53 | −1.42 | 14 | 934.74 | 940.44 | −0.61 | 22 | 1969.48 | 1993.63 | −1.21 |

7 | 347.07 | 342.05 | 1.47 | 15 | 1111.53 | 1097.38 | 1.29 | 23 | 2040.90 | 2010.40 | 1.52 |

8 | 454.11 | 461.52 | −1.61 | 16 | 1257.84 | 1247.30 | 0.85 |

Mode | Calculated (Hz) | Measured (Hz) | Error (%) | Mode | Calculated (Hz) | Measured (Hz) | Error (%) | Mode | Calculated (Hz) | Measured (Hz) | Error (%) |
---|---|---|---|---|---|---|---|---|---|---|---|

1 | 18.56 | 18.37 | 0.99 | 9 | 465.17 | 459.29 | 1.28 | 17 | 1299.71 | 1317.51 | −1.35 |

2 | 51.83 | 50.93 | 1.75 | 10 | 589.41 | 582.10 | 1.26 | 18 | 1436.11 | 1454.76 | −1.28 |

3 | 98.88 | 97.65 | 1.26 | 11 | 752.96 | 743.44 | 1.28 | 19 | 1541.01 | 1520.60 | 1.34 |

4 | 167.29 | 164.61 | 1.63 | 12 | 773.65 | 785.023 | −1.45 | 20 | 1684.51 | 1672.17 | 0.74 |

5 | 247.56 | 244.49 | 1.25 | 13 | 892.97 | 883.87 | 1.03 | 21 | 1840.34 | 1863.54 | −1.24 |

6 | 256.84 | 260.53 | −1.42 | 14 | 934.74 | 940.44 | −0.61 | 22 | 1969.48 | 1993.63 | −1.21 |

7 | 347.07 | 342.05 | 1.47 | 15 | 1111.53 | 1097.38 | 1.29 | 23 | 2040.90 | 2010.40 | 1.52 |

8 | 454.11 | 461.52 | −1.61 | 16 | 1257.84 | 1247.30 | 0.85 |

Note: Percentage errors between the calculated and measured natural frequencies are listed.

## 5 Experimental Validation of Linear Finite Element Model

To validate the linear FE model, an experimental modal analysis of the test structure was conducted with the same setup as that for its linearity investigation, which was described in Sec. 2.2.

### 5.1 Experimental Procedure.

The experimental modal analysis was conducted using the roving sensor technique [22], where one impact point and the 39 measurement points were selected. The location of the impact point was determined based on the first 23 out-of-plane elastic mode shapes from the linear FE model to avoid most nodal lines. Random errors, such as measurement noise and inconsistent impact locations, were minimized by spectrum averaging. Bias errors, such as leakage, were minimized by applying appropriate windows to excitation and response measurement data. At each measurement point, five sets of excitation and response measurement data were collected and used to calculate the associated averaged FRF with the H_{1} algorithm [21]. Each set of measurement data had a duration of 10.24 s and was collected with a sampling frequency of 5000 Hz; the frequency resolution of a FRF was 9.766 × 10^{−2} Hz. A force window and an exponential window were applied to the collected excitation and response measurement data, respectively.

### 5.2 Modal Parameter Estimation Results.

FRFs associated with the 39 measurement points were analyzed by the modal parameter estimation algorithm PolyMax in LMS Test.Lab 17. The sum of auto spectra of the 39 analyzed FRFs is shown in Fig. 19 to assist identification of the elastic modes of the test structure. The number of peaks in Fig. 19 could indicate that of the elastic modes within the Nyquist frequency associated with the sampling frequency (2500 Hz). Free boundary conditions had been validated in the linearity investigation in Sec. 2.2. Based on the measured FRFs, estimated natural frequencies and mode shapes associated with the first 23 out-of-plane elastic modes of the test structure were obtained. The estimated natural frequencies are listed in Table 2 and compared with those from the linear FE model. It can be seen that the maximum percentage error was 1.63% for the 23 out-of-plane elastic modes. The experimentally estimated mode shapes associated with the first four out-of-plane elastic modes are shown in Figs. 18(e)–18(h). Modal vector complexity associated with the 23 estimated mode shapes was quantified using MPC values, which served to check degrees of complexity of complex mode shapes [25]. The MPC values for the 23 mode shapes are calculated, and they were all over 99%. Such high MPC values indicated that the 23 out-of-plane elastic mode shapes were all real mode shapes, which was the case for the test structure in the linear FE model. A modal assurance criterion (MAC) matrix between the first 23 out-of-plane mode shapes from the linear FE model and the experimentally estimated ones was calculated, and all diagonal entries of the MAC matrix were over 95%, indicating that the mode shapes of the same modes were highly correlated. Based on the comparison of the natural frequencies and MAC matrix, accuracy of the linear FE model of the test structure was experimentally validated.

## 6 Conclusions

A linear FE modeling method for structures with riveted connections is proposed in this work. The proposed modeling method has two steps. The first step is to develop nonlinear FE models that accurately simulate riveting processes of riveted connections of a joined structure. In the nonlinear models, two-dimensional axisymmetric elements are used, and plastic deformations and surface contacts that occur during the riveting processes are considered. The second step is to develop a linear FE model of the joined structure. In the linear FE model, formed riveted connections were modeled using solid cylinders, driven heads, and rivet heads. Dimensions and material properties of the solid cylinders, driven heads, and rivet heads are estimated based on results from the nonlinear models. Joined components can be modeled using simple elements, such as plate/shell elements. The proposed method is novel in that solid cylinders are used to model riveted connections in a linear FE model of a structure with the connections, and modal parameters of the structure can be accurately estimated in a predictive manner. The proposed method can be used to estimate modal parameters of complicated structures with riveted connections, such as panels of ships and fuselages of airplanes.

To experimentally investigate the proposed method, a joined structure was prepared and tested, which had two identical aluminum strips. The two strips were joined by six identical solid rivets that were installed using a rivet installation kit. Linearity of the test structure was experimentally validated. Modal parameters, including natural frequencies and mode shapes, of the first 23 out-of-plane elastic modes of the test structure were estimated based on its measured FRFs. A linear FE model of the test structure was developed based on the proposed modeling method to numerically estimate its natural frequencies and mode shapes. A parametric study of the friction coefficient between aluminum parts in the test structure was conducted. It was found that changes of the friction coefficient would have negligible effects on contact pressure distributions and dimensions of a solid cylinder. The maximum error between the experimentally and numerically estimated natural frequencies was 1.63%, and all diagonal entries of a MAC matrix that corresponded to the experimentally and numerically estimated mode shapes were all over 95%. The relatively low percentage error and high MAC values indicate that the test structure could be accurately modeled to estimate its modal parameters. It was also experimentally shown that the estimated mode shapes had high MPC values, indicating that the estimated mode shapes were real mode shapes, which was the case for the test structure in the linear FE model.

## Acknowledgment

The authors are grateful for the financial support from the National Science Foundation (Grant Nos. CMMI-1762917 and CMMI-1763024). The first two authors are also grateful for the faculty startup support from the Department of Mechanical and Materials Engineering at the University of Cincinnati.