## Abstract

The member stiffness and pressure distribution in a bolted joint is significantly influenced by the contact area of the mechanical interface under a prescribed preload force. This research explores the influence of as-built surface profiles for nominally flat interfaces of a C-Beam assembly with two well-defined contact regions. A high-fidelity finite element model is created such that the model uncertainty is minimized by updating and calibrating the piece parts prior to the preload assembly procedure. The model is then assembled and preloaded to evaluate the contact stresses and contact area for both nominally flat and perturbed non-flat surfaces based on three-dimensional surface topography measurements. The predicted pressures are validated with digitized pressure-sensitive film measurements. The high-fidelity modeling reveals how the compliance and thickness of the pressure-sensitive film alter the measured pressures, leading to incorrect evaluations of the stresses and contact area in the joint. The resulting low-level dynamic behavior of the preloaded assembly is shown to be sensitive to the true contact area by linearizing the nonlinear finite element model about the preloaded equilibrium and performing a computational modal analysis. The resonant frequencies are validated with experimental measurements to demonstrate the effect of the contact area on the modal characteristics of the bolted assembly. Vibration modes and loading patterns exhibit varying levels of sensitivity to the contact area in the joint, leading to an improved physical understanding of the influence of contact mechanics on the low-level linear vibration modes of jointed assemblies.

## 1 Introduction

To aid in the manufacturing and assembly of typical engineering systems, such as aircraft, bridges, or automobiles, designers use shear joints to connect various components that comprise the full assembly. Shear joints are recognized as two or more components being connected by a compressive force, typically applied through a bolt or rivet [1]. The compressive contact force in the normal direction holds the components together while frictional forces arising from the compression restrict motion in the shearing direction. From a modeling and simulation perspective, the dynamic behavior of monolithic structures in the absence of joints can be accurately predicted using traditional analytical or numerical modeling techniques such as the finite element method [2]. The inclusion of joints in these mathematical representations introduces significant challenges to standard modal analysis practices and uncertainty on the dynamic behavior. For example, Beards [3] reports that joints can contribute up to 90% of the structural damping of an assembly due to microslip in interfaces. In another study by Brake et al. [4], the normal contact pressures in a bolted joint and the resonant frequencies of the global assembly are influenced by the surface geometry and finish, leading to significant uncertainty and variability in analytically characterizing the linear and nonlinear dynamics. The stress state in the joint not only influences the energy dissipation but also the contact area in the joint, which ultimately dictates the member stiffness of the joint and assembly.

Several studies provide evidence that the amount of bolt torque and preload in bolted joints influence the dynamic stiffness and resonant frequencies of a structural assembly [5–10]. This preload and pressure distribution also influence thermal conductivity across the interfaces since thermal conductivity is directly associated with the material contact area [11–13]. It is important to understand the transfer of heat and vibration energy through these connections to optimally design and accurately predict the performance of thermomechanical systems. The sensitivity of the mechanical stiffness and thermal conductance across a bolted joint is attributed to the true contact area of the mating surfaces under normal pressure. Classical solutions from nonadhesive contact mechanics, such as Hertzian contact [14] and receding contact [15], relate the contact area and traction forces of two contacting elastic bodies. In the case of bolted joints, the true contact area within the bolted joint must be predicted to obtain accurate predictions of the thermal and mechanical behavior.

Predicting the true contact area within a mechanical joint provides a significant modeling challenge since it is not a common practice to explicitly model the bolt preload prior to performing linear modal analysis using finite element methods. Several finite element modeling approaches have been proposed to capture the effect of the preloaded joint on the structural response. Kim et al. [16] compared several bolt modeling strategies of various fidelity and conclude that modeling the bolt with three-dimensional (3D) elements and the interfaces with surface-to-surface contact elements produces results that best agree with experiments. The stiffness of joints is often calibrated to experimental data using model updating techniques [17]. When calibrating linearized models of jointed structures, prior studies have optimized the values of springs, masses, and dampers at the model's interface to match modal characteristics measured from actual hardware [18–21]. Additional techniques include varying geometric parameters [22] as well as the material properties of a “doubly connective layer” [23]. The research by Ahmadian and Jalali [24] model the interface with linear and nonlinear springs and dashpots and calibrate the parameters to experimental data to capture the physics in the joint in the linear and nonlinear response regime. To the best of the authors’ knowledge, little research has been done to predict the structural dynamic characteristics of a bolted assembly that accounts for the true contact area within the preloaded assembly.

The focus of this study is to present a modeling approach that first predicts the true contact area within a bolted joint by performing a nonlinear static preload analysis, followed by the computation of the undamped (real) vibration modes of the assembly by linearizing the contact interface about the equilibrium state. This approach is shown to be sensitive to the as-built geometry of a structure, particularly to the surface profile of the mechanical interface, through validation with experimental measurements of pressure distributions in the joint (via pressure-sensitive film) and modal frequencies. The so-called C-Beam assembly, shown in Fig. 1, is used in this research to investigate the predictive modeling approach. The system consists of two individual beams, designated B12A and B12B, henceforth referred to as beam A and beam B, respectively, with raised pads on each end with a well-defined contact region for a bolted joint. Each pad is 50.8 mm long and 31.8 mm wide, with an 8.4 mm diameter thru-hole to accommodate the bolts. The beams have nominally the same dimensions and were designed such that a flatness tolerance of 25.4 *µ*m was imposed on the nominally flat surfaces of the contact pads. The beams were assembled using 5/16-14 UNF-2B bolts with corresponding nut and washers between the nut-to-beam and bolt head-to-beam interfaces, as highlighted by the inlet in Fig. 1. Each beam and bolt hardware was serialized such that the mechanical properties could be measured prior to assembly. In addition to characterizing the masses and linear elastic material properties of the individual parts, the surface profiles of the nominally flat contact pads were characterized with high-resolution 3D surface topography and updated in the finite element model. With the well-characterized, as-built piece parts, a nonlinear finite element model was created and used to predict the preloaded equilibrium state of the assembly and the linearized vibration modes.

This research builds upon previously presented work by several of the authors in Refs. [25,26]. The first study [25] evaluated a model updating procedure that changed the effective radius of a tied multipoint constraint (MPC) to best match the first six vibration modes of the C-Beam assembly. The latter study [26] used the contact surface profile data to update the nodal geometry in the joint region of the finite element mesh and evaluated the predicted vibration modes by linearizing about the preloaded equilibrium state. This paper builds on these approaches through critical evaluation of the contact pressure within the joint, investigation of the effect of thin pressure-sensitive films on contact pressure distribution, and evaluation of the influence on the dynamic stiffness of the low-level linear vibration modes. The remainder of the paper is organized as follows: Sec. 2 describes the testing and calibration efforts to characterize the individual piece parts of the C-Beam assembly by separately measuring their mass, vibration modes, and surface topography of the contact pads. In Sec. 3, the nonlinear finite element model is used to predict the stress state within the joint using the perturbed surface geometry of the mating contact interfaces, which are then validated with Fujifilm pressure-sensitive film measurements. The vibration mode frequencies about the linearized preload state are also predicted and validated against experimental measurements to understand the influence of the contact area on the stiffness of the joint. These results are benchmarked against predictions from nominally flat interfaces to reveal the influence of as-built geometry. Conclusions of the research are drawn in Sec. 4.

## 2 Individual Piece Part Testing and Calibration

The following subsections describe the efforts to calibrate the properties of the individual piece parts in the C-Beam assembly. One of the objectives of this research is to determine the influence of surface waviness and localized asperities on the ability to predict experimental modal test results using traditional finite element approaches. As such, it is important to capture the surface profile of the contacting interfaces, as well as the mass and linear elastic material properties of each individual part, as accurately as possible. Sections 2.1–2.3 discuss the measurements necessary for obtaining these results.

### 2.1 Single-Beam Experimental Modal Analysis.

Experimental modal tests were conducted on each of the individual beams to provide baseline data for model updating and correlation prior to assembly. Free–free boundary conditions were approximated for model comparisons in an effort to mitigate the effect of these uncertainties on the modal results. Noncontact response measurements were utilized to minimize the mass-loading effects from mounting sensors to the test articles. All data were collected under ambient laboratory conditions.

Images of the test setup are shown in Fig. 2. A Polytec MPV-800 Multipoint Laser Vibrometer system was used as the noncontact response measurement system, while the beams were suspended from soft bungee cords on both ends to approximate a free–free boundary condition. A foam block on a lift platform was raised underneath the beam to the point that contact was made at either end of the beam, leaving the mass of the beam primarily supported by the bungees. This was necessary to prevent large rigid body motions when the beam was impacted, which would invalidate the laser measurements. A Maul-Theet vImpact-60 automatic modal hammer, instrumented with a load cell to measure input force, was used to mechanically excite the structure. The auto hammer allowed for repeatable impacts of approximately 4.45 N-peak. The beam was excited in both the out-of-plane (Fig. 2(b)) and in-plane directions (not shown for brevity) at one end of the beam to capture the first five elastic modes. Each impact was repeated 25 times for averaging.

The MPV-800 Multipoint Laser Vibrometer system was configured with 24 individual infrared laser heads, as shown in Fig. 2(a), which measured velocity at each location concurrently. The preliminary modal results from the uncorrelated finite element model were used to select measurement locations, such that the first five elastic modes could be uniquely extracted. A total of 18 measurement nodes were used; nine stations at 64 mm intervals along the back of the beam (side without interface pads) with nodes on the top and bottom of each surface. Retroreflective tape was applied at each node to provide better signal return to each of the laser heads. The top nodes at both ends and the center node were configured for triaxial measurements (three heads aimed at one node), whereas the remaining nodes were uniaxial measurements normal to the beam surface. The triaxial measurements were necessary to capture the first in-plane bending mode of the beam. The physical arrangement of the laser heads can be seen in Fig. 2(a) (two heads are not visible at the extreme right and left of the images).

Modes were fit to the experimental data using the Synthesize Modes And Correlate (SMAC) algorithm [27] and the resulting natural frequencies are reported in Table 2 in Sec. 2.2. Using the frame structure supplied with the MPV-800 to support the 24 laser heads, it was difficult to achieve a large enough angle between the triaxial head combinations to allow for a proper transformation of the velocity data into orthogonal three-dimensional measurements. Careful observation of the in-plane motion at the triaxial measurement locations (particularly at the ends) shows slight experimental errors in the shapes. These geometric issues did not affect the extracted natural frequencies that were used in the subsequent model updating efforts.

### 2.2 Model Calibration of Piece Parts.

Based on the data from the experimental modal tests, finite element models of the individual beams were calibrated to best match the extracted natural frequencies. The material densities of the beams were calculated directly from the measured mass and the volume of the discretized finite element model. The other properties were updated via a simple model updating procedure that sampled the elastic properties (i.e., Young's modulus and Poisson's ratio) over a grid of parameter space to minimize the *l*^{2}-norm of the frequency errors. The resulting calibrated, linear elastic material properties are provided in Table 1. The beams were manufactured from 4340 stainless steel and the calibrated properties for beam B was within the expected range for this specific alloy, but beam A was not. The material density for the latter beam is slightly lower than the minimum reported value of the material (7801 kg/m^{3} compared with 7850 kg/m^{3}), and the modulus is lower compared with the other beam (204.8 GPa compared with 209.6 GPa). Deviations in material properties between similarly manufactured components are often a result of using different material lots for each component. Although when determining material properties from the updating procedure used in this paper, other effects can also skew the properties such as inexact geometries.

Property | Beam A | Beam B |
---|---|---|

E, GPa | 204.8 | 209.6 |

ν, [–] | 0.285 | 0.285 |

ρ, kg/m^{3} | 7759 | 7812 |

Property | Beam A | Beam B |
---|---|---|

E, GPa | 204.8 | 209.6 |

ν, [–] | 0.285 | 0.285 |

ρ, kg/m^{3} | 7759 | 7812 |

Mode number | Beam A | Beam B | Mode Shape |
---|---|---|---|

1 | 185.1 Hz (0.60%) | 187.7 Hz (0.04%) | |

2 | 520.5 Hz (0.47%) | 525.7 Hz (0.31%) | |

3 | 609.4 Hz (−0.72%) | 611.4 Hz (−0.21%) | |

4 | 1034 Hz (0.30%) | 1044 Hz (0.23%) | |

5 | 1527 Hz (0.26%) | 1543 Hz (0.11%) |

Mode number | Beam A | Beam B | Mode Shape |
---|---|---|---|

1 | 185.1 Hz (0.60%) | 187.7 Hz (0.04%) | |

2 | 520.5 Hz (0.47%) | 525.7 Hz (0.31%) | |

3 | 609.4 Hz (−0.72%) | 611.4 Hz (−0.21%) | |

4 | 1034 Hz (0.30%) | 1044 Hz (0.23%) | |

5 | 1527 Hz (0.26%) | 1543 Hz (0.11%) |

The natural frequency comparison between the calibrated finite element model and the experimental modal data is provided in Table 2. The experimental natural frequencies are listed for each beam under consideration, and the frequency error of the model is listed in parentheses. The finite element model predictions used the material properties are listed in Table 1. Overall the agreement with experiment is within 0.3% frequency error for beam B; however, beam A has slightly higher error up to a maximum of 0.7%. This beam had an overall lower mass and calibrated Young's modulus compared with the expected range of the material. It should be noted that other beams not shown here (i.e., serial numbers B10A and B10B) were tested and the errors and material fits were consistent with those of beam B, suggesting that beam A has some slight anomaly associated with it. This evidence shows that there exists variation in the as-built response of even monolithic structures of the same nominal design.

In addition to the calibration of the beam material properties, the bolts, washers, and nuts used to preload the assembly were weighed to determine their mass. As with the beams, the material density was calculated based on these mass values and the volume from the discretized finite element model. The material model for the bolt hardware assumed linear elastic material properties with a density of 7492 kg/m^{3}, Young's modulus of 206.8 GPa, and Poisson's ratio of 0.3. It is assumed here that the mass of the bolt hardware is most critical to the vibration modes of the assembly, as these mass load the structure at the beam ends. The method of linearization used in this paper, in which the interfaces are directly connected with multipoint constraints, negates any significant contribution to the vibration modes that may arise from the elastic behavior of the bolts.

### 2.3 Interface Surface Measurements and Updating.

After calibrating the material properties of the individual beams, the surfaces of each contact pad were measured via high-resolution 3D surface topography, prior to the assembly process, to determine the true surface profile of the interfaces preloaded within the joint. These surfaces, when joined together via a tightened bolt connection, dictate the true contact area in the mechanical interface and hence the stiffness of the joint carried through the material in contact. A Zygo Nexview NX2 scanning white-light interferometer (SWLI) was used to make high-resolution 3D surface topography measurements using vertical scanning interferometry (VSI). VSI uses a white-light light-emitting diode (LED) source (coherence length typically less than 1 *µ*m), coupled to a beam splitter, reference mirror, digital camera, 3-axis motorized gantry, and Mirau interferometry objective. The setup of the measurement system is shown in Fig. 3. A digital camera records the vertical stage position that produces the largest contrast for each pixel, and commercial software determines the corresponding height of the pixel relative to surrounding ones. This process is repeated across the surface area of interest, and software stitching is used to generate a large-area topographical map. Typical vertical resolution in VSI mode is expected to be better than 10 nm, while lateral resolution was limited to 500 nm with the 20× objective used in this study [28,29].

The challenge with the surface topography measurements is measuring the correct datum of the contact pad. The alignment with the SWLI may induce unknown orientations with respect to yaw, pitch, or roll of the surface profile. A macro-scale topography map was taken for each beam using a 3D laser scanner attached to an eight-axis scan arm (Faro ScanArm and Prizm Laser 3D scanner attachment), allowing 3D real-time volume mapping. A coordinate measuring machine (CMM) is a commercial inspection device that uses a physical probe in contact with a surface to determine the surface's location in space, with respect to a reference or datum surface. In this study, all macro-topography measurements were conducted using a highly flat and smooth surface as the datum such that the contact pad surface orientation could be properly known. The CMM arm provides a reference location for the laser scanner attachment, allowing the software to relate the laser scan data to a datum with a high degree of precision [30]. A representative study of laser/CMM combinations indicates a systematic error of 10–15 *µ*m with other factors influencing measurement error [30] such as warm-up time, reflectivity, measurement distance, and speed of collection. Faro reports an overall error of 55 *µ*m for the system that was used to collect the data. The beam specimens were laid on the inspection surface with the contact pads facing up and all 3D measurements were taken with a spatial density high enough to enable overlaying of the SWLI scan maps onto the lower-resolution scans of the beams. The macro-topological maps were used to determine the orientation (yaw, pitch, and roll) of the pad surfaces with respect to the flat reference datum. The SWLI data was aligned to geometric features on the pad surfaces in both datasets (the corners and central hole). This allowed for the generation of a multiscale 3D model of each beam and the orientation of asperity level data on the pad surfaces where contact between beams was simulated/modeled. Figure 4 illustrates the methodology to collect these datasets and how they were combined to create the SWLI data with the correct orientation relative to a flat datum.

The SWLI was used in VSI mode to scan all contact pads of the individual beams and corrected to the correct datum to obtain the proper orientation. With these surface profile measurements, the data were decimated to a more workable size for postprocessing and mapping to the finite element mesh. With the in-plane camera resolution being about 6.1 *µ*m, every ten points are retained resulting in a new in-plane resolution of 61 *µ*m. This decimated dataset is further smoothed in matlab^{®} using a moving average filter with 5% of the data used in each window. Without this smoothing, the reduced dataset had peaks with a single node harsh rise. The decimated and smoothed point cloud is overlaid onto the finite element mesh of the interface surface. matlab's *k*-nearest neighbors search (knnsearch) finds the nearest point in the SWLI data to each node of the finite element mesh and projects the node to match the height of the SWLI data. The pre- and postsmoothed surfaces are shown in Fig. 5 along with the projection of the smoothed data onto the finite element mesh in the last row. This process is repeated for each surface being analyzed throughout this study.

The plots in Fig. 6 show the deviation of peaks and valleys on the surface about a plane passing through the average of all peaks and valleys for each raised C-Beam pad. For a C-Beam pad with the interface pointing upward, a peak is defined as any point in which both principle curvatures are negative and is greater than 6.1 *µ*m height. Valleys are defined in a similar fashion with the exception of the principle curvatures both being positive. With the manufacturing drawings calling out a flatness of 25 *µ*m on the interfaces, these plots indicate that beam A meets the specification, while beam B exceeds the allowable flatness in several locations. Figure 7 shows the gaps between relative points in the contact interface. The finite element mesh grid was constructed in such a way that grid points of opposing pads are exactly across from each other, which allows for a straight forward calculation of each grid points distance from it corresponding grid point on the opposite pad. With a maximum allowable deviation of ±25 *µ*m from the drawing, the maximum expected gap is 102 *µ*m. However, the measured surface data show that some of the deviations from the drawing observed from Fig. 6 for beam B are aligned in a nonfavorable way to allow a gap of greater than 102 *µ*m in some locations. This gap measurement can be deceiving, as interfaces are rarely perfectly aligned. As such, the interfaces were also shifted relative to each other by 254 *µ*m in each of the principal directions along the interface, creating four unique permutations of asperity engagement. Following the procedure for preload and modal analysis outlined in this paper, these shifts in asperity locations only produced a maximum 0.21% difference from the perfectly aligned interface's modal frequencies. With this negligible change, the authors proceeded with the perfectly aligned interfaces.

## 3 Assembled Beam Modeling and Validation

The calibration procedures now allow for the finite element model to be assembled such that the preloaded equilibrium state can be accurately computed using nonlinear static solvers. The nonlinear finite element model is then linearized about this state to compute the linear natural frequencies of the structure. In Secs. 3.1 and 3.2, the results of the updated model are presented and directly compared with measurements for validation.

### 3.1 Nonlinear Finite Element Model With Perturbed Surface Geometry.

A highly refined finite element mesh, shown in Fig. 8, is used to perform the numerical nonlinear preload analysis on the system. In this model, a subsection of the raised pads is discretized with linear eight-node hexahedral elements with a typical element characteristic length of 203 *µ*m, resulting in 250,128 elements. Each pad has 21,230 nodes and 20,844 element faces in the contact surface to resolve the normal contact forces within the bolted region. The remainder of the parts is more coarsely discretized with linear eight-node hexahedral elements with a typical element characteristic length of 762 *µ*m. To mimic the monolithic nature of the individual beams and also to reduce the overall computational burden, the finely meshed raised pads were connected to the more coarsely meshed portion of the beam with a tied MPC. Additionally, to focus the nonlinearity of the problem solely in the contact interface that occurs between the mating beam pads, all remaining contact interfaces were instead connected with tied MPCs. This includes the bolt head to top washer, top washer to top C-Beam top surface, bottom washer to bottom C-Beam bottom surface, washer to nut, and nut to bolt. Such constraints are reasonable as the pressures underneath the bolt heads, nuts, and washers are expected to be significantly large. The total element and node count for the entire finite element model is 1,759,376 and 1,947,296, respectively.

The perturbed surface mesh, as described in Sec. 2.3, and their connected beam structures were assembled together to determine their preloaded equilibrium state using nonlinear static solvers within SIERRA/SM [31], a highly parallelizable structural mechanics code developed by Sandia National Laboratories. The finite element model was preloaded by applying a negative artificial strain along the bolt shank to the preload section (i.e., small portion along the length of the shank). To experimentally determine the total applied bolt load, Fujifilm Prescale pressure-sensitive film^{2} was placed between the contact pads and the bolts were tightened to a prescribed level of 12.4 N-m using a calibrated torque wrench. The single-use film registered the pressure within the interface and was later digitized to show the pressure distribution. From these digitized distributions, each digitized pixel has an area and pressure associated with it, thus it can be used to approximate the total force acting on the interface. Using this technique, the total force was determined to be approximately 5.56 kN for each film measured under a nominal 12.4 N-m torque. As a result, the bolts in the finite element model were preloaded to 5.56 kN to best represent the true load in the experiment. The pressure-sensitive film measurements were also used to validate the distribution of the pressure in the contact pads predicted from the updated finite element model with a non-flat profile. In this study, the low (2.41 to 9.65 MPa) and medium (9.65 to 49.0 MPa) films were used to capture a range of distributions for a 12.4 N-m bolt torque, which were sufficient for capturing the pressure in the joint.

Initial results indicated the predicted pressure distribution in the preload analyses of the C-Beam interfaces did not match the contact shape nor pressure profile experimentally determined by the digitized pressure-sensitive film. It was hypothesized that the presence of the film directly in the contact interface was affecting the observed pressure distribution. To test this hypothesis, numerical preload analyses were conducted with a finite element model representation of the pressure-sensitive film inserted between the two contact interfaces. The films were finely meshed with four 8-node hexahedral elements through the thickness. As shown in Ref. [32], the approximate Young's modulus, Poisson’s ratio, and density of the pressure-sensitive film are 100.0 MPa, 0.45, and 1389 kg/m^{3}. A low and medium range films were analyzed. The low film has a measurement range of 2.41–9.65 MPa and a thickness of 191 *µ*m, while the medium film has a measurement range of 9.65–49.0 MPa and a thickness of 102 *µ*m. It is a well-known solution in contact mechanics [15] that two flat surfaces brought together by a centralized load tend to form a receding contact patch such that the contact area is independent of the load. As such, the C-Beam finite element model was analyzed with all interfaces being perfectly flat. As shown in Fig. 9, for an unimpeded flat-on-flat contact preload numerical simulation, a receding contact patch is calculated as shown by the contact status and normal contact pressure in the interface. However, by including the low-pressure-sensitive film, a nontypical receding contact patch is predicted. The area in contact is predicted to be larger, nearly consuming the entire patch, and the pressure distribution is more ovular in shape. This is a result of the additional thickness and compliance added by the film. When the preload occurs, the film fills in the regions that would otherwise recede, causing the contact pressure to redistribute.

A similar analysis was conducted for both contact patches in the C-Beam assembly with the perturbed non-flat surfaces, with results shown in Figs. 10 and 11. The top image in each figure shows the simulated normal contact pressure for direct metal-on-metal contact of the non-flat surface, and the remaining images show the results with either the low- or medium-pressure-sensitive film included. For each pressure-sensitive film simulation, the digitized experimental pressure-sensitive film is shown adjacently to highlight the agreement with measured pressure distributions. It is evident in these plots that the predicted contact area is greatly increased when including the pressure-sensitive film within the interface and the subsequent normal contact pressure is greatly reduced. As with the receding contact predictions shown in Fig. 9, the pressure-sensitive film tends to fill in the small gaps created by the surface waviness and asperities. When the gaps are filled with the soft film material, the film acts to redistribute the contact pressure and inherently changes the measured behavior, leading to incorrect measurements of the true pressure profile in the joint. To highlight this, the metal-on-metal contact analysis predicted peak contact pressure of 181.8 MPa, whereas the analyses with the model of the pressure-sensitive film predicted peak values of 21.2 MPa and 19.0 MPa for the medium and low-pressure-sensitive films, respectively. Including the pressure-sensitive film in the experiment and/or the analysis causes a reduction of 8.6× and 9.6× for the medium- and low-pressure-sensitive films, respectively, in the true peak pressure for the cases studied here.

The pressure-sensitive film digitization in Figs. 10 and 11 show localized regions of pressures such that vertical strips of contact area form due to the non-flat surfaces in contact. There is excellent agreement in the overall shape between the measured and predicted quantities, giving confidence that the preloading numerical techniques used are accurate in capturing the pressure distribution and contact area in the interfaces. For the left interface in the beam assembly (i.e., see Fig. 10), three distinct contacting regions carry the load across the joint, which persist across the low- and medium-pressure-sensitive film simulations as well as the metal-on-metal contact analysis. The non-flat surface pressure distribution is drastically different than the flat-on-flat receding contact prediction in Fig. 9. The localized pressure regions in the as-built interfaces will not only significantly influence the member stiffness in the joint, but this pressure pattern will also affect the wear patterns observed during operation. Such localized normal contact pressures will induce larger in-plane shear tractions under high dynamic forces, thus leading to greater potential for fretting and surface wear in the joints. These results demonstrate the design challenges associated with manufacturing a truly “flat” surface for contacting bodies and show how the true interface pressures may deviate from the quantities predicted by classical contact mechanics solutions. For example, the receding metal-on-metal simulation predicted peak stress of 13.6 MPa, whereas the perturbed surface predicted values up to 181.8 MPa, which is an amplification factor of 13.4×. Section 3.2 will evaluate the influence of the member stiffness between contacting surfaces in the bolted joint by investigating the resultant vibration mode frequencies.

### 3.2 Linearized Modal Analysis.

Modal testing was conducted on the beam assembly to validate the predicted vibration modes of the finite element model. The modal test setup is similar to the individual beam setup where free–free boundary conditions with soft bungee supports were utilized and noncontact response measurements were used to minimize the effects of mass-loading from mounting sensors to the test articles. All data were collected under ambient laboratory conditions.

The test setup for the assembled beams is shown in Fig. 12. The NV-Tech SAM-1 auto hammer, used in the assembled beam test, allowed for repeatable impacts of approximately 4.45 N-peak. Prior work confirmed that such input force level excites linear dynamics of the beam assembly [26]. Four impact degrees of freedom were used to excite the first seven elastic modes of interest and a total of ten measurement sets were taken for each excitation location for data averaging. For the assembled beam test, a Polytec PSV-500 3D scanning laser Doppler vibrometer system was used, which comprises three laser heads, each with mirrors that allow the laser spot to be moved around the surface of the test article. By measuring at a coincident location with all three lasers simultaneously, a triaxial velocity measurement is obtained at that location. Thus, measurements are made at each node one at a time and a final set of data are compiled at the end of the test. The same node locations used in the individual beam tests were used for the assembled test articles; the front and back of each assembly were scanned for a total of 36 nodes. Since all three translational directions are measured at each node, there were a total of 108 measured degrees of freedom. The physical arrangement of the laser system is shown in Fig. 12(a).

Modes were fit to the experimental data using SMAC [27] and the averaged results from three separate disassemblies/assemblies are reported in Table 3. It was observed that the tests produced highly repeatable results, suggesting that this bolted assembly was relatively insensitive to assembly-to-assembly variability. The largest variation occurred in mode 6, where the three tests produced natural frequencies of 935.0 Hz, 935.6 Hz, and 936.2 Hz. The remaining modes were within ±0.1 Hz.

Elastic mode | Description | Frequency, Hz | Damping ratio, % | Mode shape |
---|---|---|---|---|

1 | 1st out-of-plane, out-of-phase bending | 282.7 | 0.12 | |

2 | 1st out-of-plane, in-phase bending | 352.4 | 0.10 | |

3 | 2nd out-of-plane, in-phase bending | 505.8 | 0.07 | |

4 | 1st in-plane, in-phase bending | 597.6 | 0.06 | |

5 | 2nd out-of-plane, out-of-phase bending | 774.7 | 0.06 | |

6 | 1st in-plane, out-of-phase bending | 935.6 | 0.05 | |

7 | Torsion | 1175.5 | 0.04 |

Elastic mode | Description | Frequency, Hz | Damping ratio, % | Mode shape |
---|---|---|---|---|

1 | 1st out-of-plane, out-of-phase bending | 282.7 | 0.12 | |

2 | 1st out-of-plane, in-phase bending | 352.4 | 0.10 | |

3 | 2nd out-of-plane, in-phase bending | 505.8 | 0.07 | |

4 | 1st in-plane, in-phase bending | 597.6 | 0.06 | |

5 | 2nd out-of-plane, out-of-phase bending | 774.7 | 0.06 | |

6 | 1st in-plane, out-of-phase bending | 935.6 | 0.05 | |

7 | Torsion | 1175.5 | 0.04 |

A linear modal analysis was conducted on the C-Beam assembly finite element model described in Sec. 3.1, which has been linearized about the preloaded states shown in Figs. 10 and 11. Note that only the metal-on-metal contact simulation is used to linearize the state of the joint since this scenario matches the configuration of the experimental modal analysis (i.e., no pressure-sensitive film was present in the interface during modal testing). The linearization is performed by (a) updating the finite element geometry based on the statically deformed preload shape, (b) calculating the geometric stiffness based on the preloaded stress state, and (c) determining which nodes in the surface are in contact (normal pressure > 0.000 Pa) and constraining them to the opposing surface via tied MPCs. With the linearized finite element model, a modal analysis is conducted using SIERRA/SD [33] to predict the real eigenmodes and eigenvalues of the system. This procedure, rather than a consistent Jacobian linearization, was necessary since the modal analysis is performed in a separate finite element code (i.e., Sierra/SM and Sierra/SD). This approach is only valid for small material and geometric nonlinearities, which is the case for the preloaded beam assembly. As with the experimental setup, the numerical analysis is conducted with free–free boundary conditions. In addition to using the MPCs from the as-calculated contact area, another set of MPCs are analyzed which reflect a slight increase in contact area over the as-calculated case. The two sets of MPCs are shown in Fig. 13.

Table 4 shows the results of the computational modal analysis, with comparisons between the experimentally measured natural frequencies and the simulated frequencies from various interface MPCs. The fully stuck interface refers to the case that ties the entire 50.8 mm by 31.8 mm pad without knowledge of the contact area from the preload, the most expedient treatment for contacting bodies. The receding flat interface MPC ties the nodes in contact per the results in Fig. 9 (without film). Both the calculated as-built non-flat surface MPCs are calculated based on the nonlinear preload simulations in Figs. 10 and 11. The increased contact area is evaluated by simply expanding the “calculated” MPCs to fill artificial rectangular boundaries bounding these areas of contact.

Nominal surfaces | As-built surface | ||||
---|---|---|---|---|---|

Mode number | Experimental (Hz) | Fully stuck (% error) | Flat-on-flat (% error) | Calculated (% error) | Increased (% error) |

1 | 282.7 | 3.82% | 1.66% | −0.68% | −0.13% |

2 | 352.4 | 0.92% | 0.37% | −0.58% | −0.30% |

3 | 505.8 | 0.35% | 0.25% | 0.21% | 0.22% |

4 | 597.6 | −0.48% | −0.58% | −0.58% | −0.58% |

5 | 774.8 | 4.06% | 1.69% | −0.81% | −0.22% |

6 | 935.6 | 4.63% | 1.87% | −2.23% | −0.99% |

7 | 1175.5 | 0.67% | 0.22% | −0.02% | 0.10% |

Mean % error | 2.13% | 0.95% | 0.73% | 0.36% |

Nominal surfaces | As-built surface | ||||
---|---|---|---|---|---|

Mode number | Experimental (Hz) | Fully stuck (% error) | Flat-on-flat (% error) | Calculated (% error) | Increased (% error) |

1 | 282.7 | 3.82% | 1.66% | −0.68% | −0.13% |

2 | 352.4 | 0.92% | 0.37% | −0.58% | −0.30% |

3 | 505.8 | 0.35% | 0.25% | 0.21% | 0.22% |

4 | 597.6 | −0.48% | −0.58% | −0.58% | −0.58% |

5 | 774.8 | 4.06% | 1.69% | −0.81% | −0.22% |

6 | 935.6 | 4.63% | 1.87% | −2.23% | −0.99% |

7 | 1175.5 | 0.67% | 0.22% | −0.02% | 0.10% |

Mean % error | 2.13% | 0.95% | 0.73% | 0.36% |

The fully stuck interface is a bounding case that produces the stiffest joint since the entire contact surface is tied. Reducing the contact area results in lower mode frequencies for the assembly. The results in Table 4 show that different mode frequencies have differing sensitivities to the contact area in the interface. The results show that modes 1, 5, and 6 are the most sensitive based on the change in frequency from the “fully stuck” model to the “as-built, calculated” model. The mode shapes for these are shown in Table 3, which reveal that these deformations act on the joint in a prying manner such that the joint opens/closes over an oscillation. Modes 2 and 3 impart a shearing load on the joint, and mode 4 does not appear to activate any loading on the joint (which is further supported by evidence that this mode was the least sensitive to the MPC). In general, the mean of the total frequency errors decreases with increased fidelity of the joint, i.e., going from fully stuck, to flat-on-flat, and to the as-built approach. The increased, as-built contact area gives the best prediction for the most sensitive modes (1, 5, and 6) with a percent error of −0.13%, −0.22%, and −0.99%, respectively. One interesting observation is that for mode 6, the frequency predicted from the calculated area of the as-built surfaces has an error of −2.23% while the slight increase in contact area reduces the magnitude of this error to −0.99%. The contact area in the latter case has not changed significantly, but clearly shows improvement on the results, suggesting this mode is the most sensitive and difficult to accurately predict with finite element modeling.

## 4 Conclusions

This study evaluates the influence of as-built properties of a bolted C-Beam assembly. A high-fidelity, nonlinear finite element model was developed such that each piece part was individually calibrated using the measured mass and vibration modes of the individual beams. The surface profiles of the contact pads were measured using 3D surface topography that was rotated to correctly identify the orientation with respect to a flat datum. The updated finite element models were then preloaded to evaluate the normal contact pressures that developed within the jointed region of the interface. These results were validated with pressure sensitive films and revealed two significant findings. First, it was discovered that the inclusion of pressure-sensitive film significantly influenced the pressure distribution in the joint. This was confirmed by explicitly modeling the pressure-sensitive film in the finite element model and comparing results with and without the film. The simulation results with the film agreed well with the measured data, validating the accuracy of the modeling approach. Second, it was confirmed that the true as-built surface profile had a significant influence on the stresses and the contact area within a nominally flat surface in a joint. The surfaces were designed to a low flatness tolerance and produced results different from expected classical mechanics solutions. The true as-built surfaces and resultant contact area are important to understand when modeling the contact area is important, such as cases for vibration or thermal analysis.

The influence of the contact area and member stiffness was further evaluated by predicting the vibration modes of the assembly linearized about the preloaded state. Several levels of fidelity for the multipoint constraint definition within the joint were evaluated. The results revealed that connecting the entire joint interface through a multipoint constraint produced the stiffest modal frequencies, whereas connecting the surfaces based on the as-built preload model, also through a multipoint constraint, produced the most accurate results. The modal investigation revealed that the frequency sensitivity of each mode to the contact area in the joint is dependent on its deformation shape at the joint. The results from the modal and stress analysis demonstrate the challenges associated with as-built assemblies based on designs with flat interfaces between bolted components.

## Footnote

## Acknowledgment

The views expressed in the article do not necessarily represent the views of the U.S. Department of Energy or the United States Government. Sandia National Laboratories is a multimission laboratory managed and operated by National Technology & Engineering Solutions of Sandia, LLC, a wholly owned subsidiary of Honeywell International Inc., for the U.S. Department of Energy's National Nuclear Security Administration under contract DE-NA0003525.

## Funding Data

This work was sponsored by Sandia National Laboratories’ Advanced Simulation and Computing (ASC) Physics and Engineering Models (PEM) subprogram and Delivery Environments (DE) Non-Contact Diagnostics subprogram.