Concept Modelling of Vehicle Joints and Beam-Like Structures through Dynamic FE-Based Methods

This paper presents dynamic methodologies able to obtain concept models of automotive beams and joints, which compare favourably with the existing literature methods, in terms of accuracy, easiness of implementation, and computational loads. For the concept beams, the proposed method is based on a dynamic finite element (FE) approach, which estimates the stiffness characteristics of equivalent 1D beam elements using the natural frequencies, computed by a modal analysis of the detailed 3D FE model of the structure. Concept beams are then connected to each other by a concept joint, which is obtained through a dynamic reduction technique that makes use of its vibration normal modes. The joint reduction is improved through the application of a new interface beam-to-joint element, able to interpolate accurately the nodal displacements of the outer contour of the section, to obtain displacements and rotations of the central connection node.The proposed approach is validated through an application case that is typical in vehicle body engineering: the analysis of a structure formed by three spot-welded thin-walled beams, connected by a joint.


Introduction
The virtual models obtained by computer-aided engineering (CAE) tools play a fundamental role in the development process of complex products, because they enable engineers to predict various performance attributes, avoiding the use of expensive physical prototypes and thus reducing the time of design process.Particularly in the field of automotive industry, the performances related to noise, vibration, and harshness (NVH) are hard to improve in the last steps of the development process without raising conflicts with others vehicle requirements.For this reason, many researchers have developed predictive concept modelling methodologies, which can be used to predict and improve the vehicle design from the concept phase onwards.
This paper focuses its attention on concept modelling techniques concerning the reduction in an equivalent simplified model of the detailed vehicle body in white (BIW) model, allowing to drastically reduce the required computational resources and the time needed for its modifications.
The reduction of detailed 3D FE models can rely on different commercial FE solvers available today (such as Nastran, Abaqus, and Ansys), which provide libraries of simplified elements (1D beam elements, superelements, etc.).In the literature, there are several approaches able to simplify the principal structural elements of a BIW, such as beam-like structures and joints.
Regarding element models having beam-like global behaviour, main methodologies to reduce a 3D element into a 1D element with equivalent characteristics can be grouped in three categories: geometric, static FE-based, and dynamic FEbased approaches.The geometric concept modelling methods rely on a geometric analysis of the beam cross-sections [1][2][3].The mass and stiffness properties of the equivalent 1D beam element are computed by analysing the mass distribution along the section and considering whether the section has a single or a multiconnected closed shape.
Instead, in the static FE-based methods [4,5], a set of static load cases is generated by applying bending, torsion, and axial loads at the end sections of each beam segment.For some beam cross-sections, central nodes are created and connected to the other nodes at the same cross-sections, by means of multipoint constraint (MPC) elements, which are of rigid type for the two end sections and of interpolation type for the intermediate sections.In this way, external and reaction loads are applied directly to the central node of the end sections, and the rigid elements transfer them to the rest of the structure, while the interpolating elements allow estimating the linear elastic deformation of the beam central line.Finally, the stiffness properties of the equivalent 1D beam are estimated by applying the linear elastic load-deformation relationships of the beam structure, starting from the static deformations predicted by analysing the detailed 3D model.A scheme illustrating the static FE-based method is shown in Figure 1.
This work discusses an original method pertaining to the last category, the dynamic FE-based approach [6], through which the stiffness characteristics (e.g., quadratic moments of inertia, torsional modulus, etc.) of equivalent 1D beam elements are estimated using the natural frequencies computed by a modal analysis of the detailed 3D FE model of the structure.Its main advantages came from the fact that any possible discontinuities and variations that may occur along a beam and that affect its stiffness are taken into account during the computation, so that the proposed method results in an accuracy benefit as compared to prior art methods.
Automotive joints are the second main part of a BIW to be conceptualized.Currently, the most usual techniques rely on the reduction of a joint in a superelement (SE), which is defined by reduced stiffness and mass matrices.In order to guarantee structural continuity between beams and joints in the whole concept structure, an interface between the concept 1D beam model and the detailed 3D joint model is created before applying joint reduction.The latter is then achieved by condensing the joint stiffness and mass properties to the nodes on the beam side of each beam/joint interface, as it will be explained in Section 2. Therefore, in the reduction process, the central nodes of the beam/joint interfaces represent the master degrees of freedom (DOFs) to be preserved, while the DOFs of the nodes belonging to the 3D FE structure of the joint are removed (slave DOFs).Figure 2 shows a typical BIW structure of a commercial vehicle and a typical joint and beam-like member used in automotive bodies.In the same figure, a graph showing the overall scheme of the technique proposed here for 1D modelling of beams and SE representation of joints is also illustrated.The reduction techniques on SE can be categorized into two types: static or dynamic.Guyan reduction [7] is the most common method for static condensation.It returns an exact reduced stiffness matrix and an approximated mass matrix, exploiting some static considerations between the master nodes of the joint.On the contrary, all the dynamic methods make use of the vibration normal modes of the structure, but they differ from each other for the applied boundary conditions and for the selection of enrichment vectors in addition to the normal modes.Two well-established examples are the Craig-Bampton fixed interface [8] and MacNeal's [9] approach.In the former approach, which is used in the research presented here, the normal modes are computed with the structure clamped at the connection interfaces, while the enrichment vectors are determined as constraint modes.The latter method uses the normal modes of the component in free-free conditions, while the enrichment vectors consist of residual flexibility modes.These dynamic reduction methods are reliable and easy to implement.However, in the specific application of thinwalled automotive joints, the accuracy of these methods is strongly dependent on the type of connection models that are used at each beam-joint interface.Rigid connection elements, such as Nastran RBE2 elements [10], can make the entire structure excessively stiff, while general-purpose interpolation elements, such as Nastran RBE3 elements [10], can lead to coarse inaccuracies especially with regard to the torsional stiffness.
For this reason, a new multipoint constraint (MPC) connection element is proposed and validated here, in order to overcome the limitations of standard connection elements and achieve more accurate concept models of automotive joints.The proposed MPC allows correlating displacements and rotations of the dependent node (i.e., the central node of the joint end section) with the displacements of the peripheral section nodes.These relationships are based on static considerations and return a transformation matrix, whose implementation process is explained in detail.To verify the improved accuracy of the proposed model, a comparative analysis has been carried out between a structure where concept joints are reduced by using rigid spiders as connection elements and another where the new MPC connection elements are employed.
The paper is organized as follows.Section 2 describes the mathematical definition of the new connection element.Section 3 describes the detailed 3D and the concept 1D FE models of the spot welded structure that are used for validation purposes.The dynamic validation results are reported in Section 4, demonstrating the improvements that can be obtained using the proposed connection elements.Section 5 concludes the paper, by reviewing the results achieved and providing an outlook on the foreseen next steps.

Definition of New MPC Connection Element
In this paper, a new MPC connection element is defined, which enables to create an interface between the concept 1D beam model and the detailed 3D joint model, as shown in Figure 3.For such purpose, the kinematic relationships between the displacements and rotations of the beam node (dependent node of the connection element) and the displacements of the nodes of the detailed 3D FE model of the joint at the interface section (independent nodes of the connection element) are derived in the form of a transformation matrix [], by using a static approach based on equilibrium conditions [11].The basic idea is to obtain a second transformation matrix [] that defines the relationship between the total forces at the central node of the beam/joint interface and the nodal loads over the section, considering theoretical stress fields resulting from the Saint-Venant assumptions with respect to axial, bending, shear, and torsion load-cases for a beam-like structure [12].In linear elastic field, this loads correlation can be inverted, returning the searched kinematic relationship.

The transformation matrix [𝑆] is obtained by calculating the product of two submatrices:
(i) the stress recovery matrix [ 1 ], which correlates the sectional stresses applied to the nodes of shell elements and resultant loads applied on the central node of equivalent 1D beam element, by using loadstress relationships of Saint-Venant for beam-like structures; (ii) the form nodal load matrix [ 2 ], which links all the nodal forces and the nodal stresses of the section, by using the shape functions of linear finite elements.
In the next subsections each of these matrices is described and explained in detail.

The Stress Recovery Matrix
Starting from the sectional load resultant at the interface, the stress recovery matrix [ 1 ] is defined by the assumptions of linear elastic, homogeneous, and isotropic material.First, it is assumed that the thickness is sufficiently small with respect to the crosssection dimensions, which allows using the stress distribution properties of thin-walled beams.Furthermore, for a spot welded structure, which is open in some cross-sections but closed in those regions where spot welds are applied, the global behaviour can be assumed as similar to that of closed section; for this reason, the stress equations for closed thin-walled sections are used.For a 3D beam model, the stress-load relationship for each simple load case can be written at the interface section, with respect to a reference system placed on the centre of gravity, with axis directed along the longitudinal neutral axis and and -axes along the principal and secondary bending directions, respectively.The matrix relation can be given as follows: where {} is the vector of nodal stresses, [ 1 ] is the stress recovery matrix, and {} is the vector of total forces applied to the central node.In particular, (1) for a generic rectangular cross-section can be detailed as follows: (1)     (1)     (2)     (2)     (2)    . . .
The physical meaning of each parameter is given in Table 1.
The superscript () indicates the global node number.For the transition region between two different values of thickness, an average value is set at the node between two adjacent shell elements.Note that the effects of individual loads are considered uncoupled; in addition, the contribution of linear shear stresses, perpendicular to the shear direction, are not considered, because their global effect on displacements and rotations of the central node is zero, as well as the global effect of warping in torsion load cases [13].Figure 4 shows the stress distributions in these two cases.

The Form Nodal Load Matrix [S 2
].The form nodal load matrix [ 2 ] provides a relationship between nodal forces and nodal stresses of the 3D beam model.First, the nodal loads on one shell element  at the interface are examined.This is a bilinear 4-node shell element, then only two nodes ( and ) must be considered for the interface (Figure 5).
The stress distributions on the interface of  element are found by averaging the nodal stresses of nodes  and .Then, the average stresses of  element,    ,    , and    , in , , and  directions, respectively, it can be written as follows: Since the average stresses at the interface are constant, the nodal loads for  and  can be calculated by using shape functions, yielding where  () is the area of the element,  indicates the shape function, and  is the local coordinate, useful to measure the length of element.The same relations are valid for all the other elements along the interface and, then, an assembled matrix relationship can be obtained between nodal forces and nodal stresses: where {} is the vector of total nodal forces, [ 2 ] is the form nodal load matrix, and {} is the vector of total nodal stresses.
Note that, for open sections, the forces applied on internal nodes receive the contribution of average stresses relating to both adjacent elements; instead, external nodes belong to a single element and then their values are considerably lower.Therefore, by combining (1) and ( 5), it is possible to obtain the relation between the total forces at the central node, {}, and the nodal forces on the outer contour of the section {}: where  denotes the maximum number of nodes on the interface.
To obtain the [] matrix that relates 1D beam displacements and rotations and 3D displacement values at the interface, linear relations between forces and displacements are assumed.In this way, it is sufficient to transpose []: The searched kinematic relationship can be written as follows: where {} 6×1 is the nodal displacements and rotations vector of the beam node and {} 3×1 is the nodal displacements vector of the 3D finite element model.

Application Case
This section describes an application model, on which the proposed new connection element has been validated.The geometry of the reference 3D structure is described in Section 3.1.Thereafter, Section 3.2 describes the 1D concept model of the structure.

3D Model Description.
The 3D application model consisted of three beams, connected with a joint.Each beam was one meter long and defined by two thin-walled sheets with Cshaped section, whose geometry and dimensions are shown in Figure 6(a).The cross-section had a vertical axis of symmetry and the centre of gravity has been considered coincident with the shear centre.This assumption is very important for two reasons.First, because it allows reducing detailed 3D beams with the dynamic method using uncoupled differential equations for flexural and torsional vibrations.Additionally, it also allows considering the effects due to the various stresses in the joint reduction as uncoupled.Since the sheet thickness was much smaller than the two transversal dimensions, it was possible to mesh each beam member by using 4-node finite shell elements.The three beams intersected in a joint, having the same cross-section and also modelled with shell elements.Therefore, this structure was composed of two parts, the upper and the lower parts, separated from each other.These parts have been connected by a set of equally spaced welding points along each of the longitudinal walls.The distance between spot welds has been chosen according to the typical layout in automotive beams, that is, equal to 100 mm.Each welding point has been created as a small Hexa solid element, connected to corner nodes of flanges by generic interpolation elements (Figures 6(b) and 6(c)).The Structures environment of LMS Virtual.Lab software [14] has been used to create the whole FE model.
The material assigned to the model, assumed homogeneous and isotropic, was typical steel with the following properties: (i) elasticity modulus:  = 210000 MPa; (ii) Poisson's ratio: ] = 0, 3; (iii) mass density:  = 7, 9 ⋅ 10 −9 ton/mm 3 . tot is the total area of the thin-walled cross-section Ω is the area encompassed by the middle perimeter line  () and  () are the position coordinates of the ith node  () is the thickness average value relating to the ith node   * () and   * () are the cross-section static moments about x-and y-axes for ith node   and   are the cross-section moments of inertia about xand y-axes p is a factor equal to 1 on the upper horizontal wall, −1 on the lower, 0.5 on the nodes at intersections, and 0 on the vertical walls.
q is a factor equal to 1 on the right vertical wall, −1 on the left, 0.5 on the nodes at intersections, and 0 on the horizontal walls.
In order to compare the detailed 3D FE model and the concept structure also in terms of modal shapes, a beam centre line has been created in the 3D structure by defining a proper number of central nodes, located at a distance of 100 mm from each other along the longitudinal direction of the beam.Each central node has been connected to the nodes of detailed 3D model in the same cross-section (boundary nodes) by an interpolation RBE3 element, so that the modal displacements of each centre node have been estimated by interpolation of the modal displacements of the boundary nodes.Figure 7 shows the complete detailed 3D model of the structure.

Concept Model Description.
To define the concept structure, simplified concept models of the joint and of the three beams were created.For the beams, a dynamic FEbased method for the estimation of equivalent cross-section properties was used [6].
The main advantage of the method is that it takes into account all possible discontinuities and variations (like spot welds) that may occur along a beam and that affect its stiffness, especially under torsional loads.This method consists of two principal steps: (1) firstly, the natural frequencies of a given beamstructure have been estimated by means of a modal analysis of the detailed 3D FE model, in free-free conditions; (2) secondly, the cross-sectional stiffness properties were obtained from the flexural and torsional frequencies, using the differential equations of beam vibrations [15].
In particular, an unconstrained nonlinear minimization algorithm (the Nelder-Mead simplex algorithm [16]) was used to estimate beam section properties.The implemented objective function consisted of minimizing the squared sum of differences between the frequencies reference vector, obtained from the detailed 3D FE model, and a frequencies vector, iteratively computed by applying non-linear equations derived from the beam modal model [6].For the spot-welded beam model under study, the equivalent stiffness parameters estimated for the 1D concept beams are listed in Table 2. Figure 8(a) shows the concept model, before the joint is reduced.
Instead, for the joint reduction, the Craig-Bampton fixed interface technique was applied as follows; in the final concept model the joint was represented by a SE, consisting of a stiffness and a mass matrix condensed to a set of master nodes, which includes one node on the beam side of each beam/joint interface.Therefore, the Craig-Bampton reduction was implemented with Nastran software [17], by keeping the DOFs of the central nodes at the three beam/joint interfaces as master and the DOFs of the nodes belonging to the detailed 3D FE model of the joint as slave.Figure 8(b) shows the concept structure, where an SE representation of the joint replaces the detailed 3D joint model while guaranteeing the structural continuity of the whole concept model.Two different concept models of the joint have been created: in one model, the Craig-Bampton dynamic reduction was applied to the detailed 3D FE model with rigid RBE2 elements at each beam/joint interface; in the other model, the proposed MPC elements were used to connect the central node to nodes placed on the cross-section at each interface.Note that the joint region has been defined in such a way that the distance of each interface from the joint centre is sufficiently large (100 mm in this case) to avoid any violation of the Saint-Venant beam assumptions.Note also that the detailed 3D model has over 60000 degrees of freedom (DoFs), while the concept model has 180 DoFs only, which allows a significant reduction of the computational time required by the FE simulations.

Dynamic Validation
A dynamic analysis has been carried out, to show the improvements on the predictive accuracy of the concept structure with the proposed MPC connection elements.A FE modal analysis in free-free conditions was performed using Nastran software as FE solver, in order to compare, in   terms of natural frequencies and mode shapes, the detailed 3D model to two different concept models.In both concept models, the properties of equivalent 1D beam elements have been estimated by using the dynamic FE-based technique described in Section 3.2.In the first concept structure, the joint reduction has been achieved by using Craig-Bampton technique and rigid RBE2 connection elements; in the second one, the proposed MPC connection elements have been used.Table 3 reports the natural frequencies of the first 12 global modes estimated for the detailed and for the two equivalent concept FE models, together with the percentage differences between each concept model and the reference structure.The diagonal values of the modal assurance criterion (MAC) matrix [18], obtained in Noise and Vibration environment of LMS Virtual.Lab, are reported as well.It can be observed that the concept model with rigid connection elements approximates very precisely the first vibration modes, in terms of frequencies and MAC values, concerning in particular axial and flexural vibrations; however, for flexural-torsional modes, significant differences between the concept and the detailed 3D model can be appreciated, with a maximum and an average value of 38.73% and 9.09%, respectively.This means that the stiffness of the joint is overestimated, especially when the interfaces undergo torsional deformation.Instead, the second concept structure, where the proposed MPC connection elements have been used, shows good accuracy for all modes, with a maximum and an average difference value of 4.53% and 1.91%, respectively.In particular, for modes involving torsional deformation at one or more joint/beam interfaces, the second concept model is up to 16 times more accurate than the first concept model.

Conclusions and Outlook
In this paper, new concept modelling methodologies, enabling an improved accuracy of beam and joint concept models, have been proposed and validated.For equivalent 1D beam elements, a dynamic FE-based method was applied, able to define equivalent characteristics of concept beam models, starting from the flexural and torsional natural frequencies of complex 3D beam models.Instead, for joint concept modelling, the Craig-Bampton dynamic reduction approach was used.This method is very fast and accurate, but it is strongly affected by the FE connection element between central and peripheral nodes of the joint/beam interface sections.For this reason, a new multipoint constraint connection element has been defined and implemented, in order to interpolate displacements of peripheral nodes and obtain displacements and rotations of the centre node with an increased accuracy, as compared to conventional connection elements, such as rigid spiders.
To assess the accuracy of the proposed method, an application case was analysed, consisting of a 3D structure where three spot welded beams are connected by a joint with the same cross-section geometry.By comparing two different concept models (one using rigid RBE2 and the second using the new MPC connection elements for the reduction of joint in SE) with the detailed 3D model in terms of natural frequencies and MAC values, it was proved that the proposed method has a good predictive accuracy.This confirmed that the proposed MPC connection elements permit to estimate the mass and stiffness characteristics of the reduced joint at each joint end section more accurately than with rigid spiders.
The proposed methodologies have been developed with the aim of enabling early predictions of static and dynamic behaviours in vehicle bodies already in the concept phase of the development cycle.However, it is worthy to notice that such concept modelling techniques can be exploited also in other application fields.
The next steps of this research will aim at extending the applicability of the proposed method to structures with general cross-section shape.If the assumption of doublesymmetry is removed, coupled effects of torsion and bending must be taken into account while calculating both the equivalent 1D beam properties and the coefficient of the MPC beam/joint connection elements.

Figure 1 :Figure 2 :
Figure 1: Application of static FE-based method: loads, constraint, and connection elements applied to the detailed 3D beam model.

Figure 3 :
Figure 3: The interface element (circled in red) between equivalent 1D beam model and detailed 3D joint model.The central node is in yellow (dependent node on the beam side) while the peripheral nodes are in red (independent nodes on the joint side).

Figure 4 :Figure 5 :
Figure 4: Stress distribution along rectangular cross-sections, due to vertical shear (  and   in (a)) and warping torsion (  in (b)): the global effect of linear parts in the central node is zero.

Figure 6 :
Figure 6: Application model: cross-section geometry of beams (a), mesh of detailed 3D beam model with welding zones (b), and spot weld model, with a central Hexa solid element connected to nodes of flanges by interpolation elements (c).

Figure 7 :Figure 8 :
Figure 7: Application case: detailed 3D FE model of the structure.

Table 2 :
Equivalent beam properties estimated by the dynamic FEbased method.

Table 3 :
Dynamic comparison between the original and the two-concept FE models, in terms of natural frequencies, % errors, and modal correlation factors.