*Corr. Author’s Address: University of Ljubljana, Faculty of Mechanical Engineering, Slovenia, domen.seruga@fs.uni-lj.si 3
Strojniški vestnik - Journal of Mechanical Engineering 68(2022)1, 3-13 Received for review: 2021-07-14
© 2022 The Authors. CC BY-NC 4.0 Int. Licensee: SV-JME Received revised form: 2021-09-29
DOI:10.5545/sv-jme.2021.7329 Original Scientific Paper Accepted for publication: 2021-12-16
Incorporation of a Simplified Mechanical Joint Model into 
Numerical Analysis
Glavan, M. — Klemenc, J. — Malnarič, V. — Šeruga, D.
Mitja Glavan1 – Jernej Klemenc2 – Vili Malnarič1 – Domen Šeruga2,*
1 TPV Automotive, Slovenia 
2 University of Ljubljana, Faculty of Mechanical Engineering, Slovenia
In the development of new joining technologies, incorporation of mechanical joints in computer analyses for the evaluation of structures can 
be carried out by a practical, simplified mechanical joint model. Here, two most frequently used joining technologies were analysed, a self-
piercing rivet joint and a clinch joint. Physical tests of static load capacity of the joints were performed and numerical models for simulations 
were set-up. An optimization method was designed for estimating the material parameters of the mechanical joint for the needs of numerical 
analyses. For optimization purposes during the plan of experiments, a range of possible parameter values was investigated using a response 
surface method, results of simulations, results of physical tests and a genetic algorithm. The results of simulations using the optimal values of 
the material parameters are comparable to the experimental observations for the both joints.
Keywords: mechanical joints, self-piercing riveting, clinch joint, FEM of mechanical joints, parameter optimization, simplified FE model, 
response surface, genetic algorithm
Highlights
•	 Numerical model for modelling a mechanical joint has been defined.
•	 Methodology of defining the properties of material parameters in numerical analysis has been developed.
•	 The material parameters of the mechanical joint for 2 mm thick S500MC material have been defined.
•	 The developed methodology is fast and easy to use and represents a universal tool for industrial needs.
0  INTRODUCTION
The current trend in the automotive, aerospace, nautical 
and other similar industries is primarily a decrease in 
material consumption and indirectly a reduction in the 
mass of these means of transportation. In this respect, 
today, the so-called lightweight materials such as 
high-strength steel, aluminium, composite materials 
and the use of magnesium are increasingly emerging 
as construction materials [1] to [3]. However, due 
to the different properties of these materials their 
joining with the traditional welding process is not 
possible. Today, a large number of alternative joining 
technologies is emerging, or is already in use, to meet 
the growing need in joining such materials [4] to [9]. A 
review of such technologies and their use can be seen 
in [10] to [12]. A wide range of technologies is thus 
in use, each with its advantages and disadvantages. 
The use of a particular technology hence depends 
on the strength requirements of the finished product, 
materials, the possibility of integration into production 
facilities and prices. Two joining technologies have 
been compared in our research, which possess the 
appropriate characteristics to represent an alternative 
to spot-welding. Self-piercing riveting is a technology 
that provides adequate static strength, i.e. load-
bearing capacity of the joint and a slightly higher 
dynamic strength in contrast to spot welding, while 
the advantage of the clinch joint is a significantly 
higher dynamic load-bearing capacity, with the loss 
of a significant proportion of static load capacity [13] 
to [16]. Table 1 presents the comparison of the above 
technologies [17]. However, for a proper consideration 
of these technologies in numerical analyses during 
the pre-development and development phases of 
a research and development (R&D) cycle, it is 
necessary to develop appropriate numerical models of 
these joints. Moreover, knowing the exact parameters 
of joints contributes to reducing development costs, 
shorter development time and optimization of 
geometry before the first prototype is manufactured 
[18].
As numerical structural analyses are carried out 
way before the production of prototypes and their 
physical testing, accurate knowledge of material 
parameters of the joint plays a key role in the quality 
of the simulations. For example, Bouchard et al. [19] 
used an approach with exact modelling of the actual 
rivet, which produced encouraging results even in 
the local area around the joint. Thus, following this 
approach, the numerical results are well comparable 
to physical testing, but such a model is not suitable 
for large structures and serves more for the purpose 
of understanding what happens in the joint and its 
surroundings. The disadvantages of such an approach 
are time-consuming simulations and possible 
Strojniški vestnik - Journal of Mechanical Engineering 68(2022)1, 3-13
4 Glavan,	M.	—	Klemenc,	J.	—	Malnarič,	V.	—	Šeruga,	D.
computational difficulties. For use in large structures, 
a reliable simplified numerical model is needed. 
Today, commercial software tools already offer 
various integrated numerical models that can be used 
for simulation of mechanical joints. Such an approach 
was used by e.g., Porcaro et al. [10], Yang et al. [20] 
or Kulawik and Wrobel [21]. Porcaro et al. [10] used a 
pre-set spot-weld model in LS-Dyna software, which 
joined the sheet metal with a simple line element that 
connected it on the sheet through contact conditions 
that ensured a transfer of forces and torques between 
the two sheets. The material parameters, which 
influenced the result of numerical analyses, were 
adjusted and modified iteratively by the researchers 
based on matching the results of simulations to the 
results of physical tests. In a similar way, Sommer and 
Maier [22] and Hanssen et al. [23], who were dealing 
with this issue using a more complex approach, used a 
calibration method in combination with physical tests 
to define the necessary joint parameters. Sommer and 
Maier [22] compared both several different geometric 
models of connecting elements and several different 
material models. They built various combinations of 
both geometric and material models and compared 
the results of numerical analyses with physical testing 
at chosen load cases. The material parameters were 
estimated on the basis of the calibration process that 
was based on the numerical force-displacement curve 
fitting to the curve obtained from the test. The results 
were promising and proved a good match with the 
experimental observations. Similarly, Hanssen et al. 
[23] focused on a numerical model and developed 
an algorithm to identify parameters of a numerical 
model, again based on physical testing. According to 
this methodology, the user defines the point between 
the sheets where the riveted joint is located and 
defines the cross-section of that joint. The algorithm 
searches for all the relevant nodes in the cross-section 
area of the rivet and connects the sheet metal. Through 
these nodes the loads are transferred and the joint 
response is calculated by the algorithm. However, 
for the operation of the model and the algorithm, it is 
necessary to define a wide range of parameters, which 
must be estimated through physical testing in various 
load cases. Nevertheless, the results obtained from 
such a model match well with the actual state.
In this study, the mechanical joint has also been 
modelled in a simplified way, whilst the number of 
the necessary material parameters has been minimised 
to the lowest possible amount which still ensured 
comparable results of the simulations and experimental 
observations. This way a high computational 
efficiency for analyses of large structures can be 
preserved. Mechanical joints have been modelled 
with a line element that is connected to the two sheet 
metals. A bilinear material model has been assigned 
to both the sheet metal and the joint, as the results of 
physical tests show an extreme bilinearity until the 
point of failure. Nevertheless, the material parameters 
of the sheet metal and the joint have been separately 
considered in the simulations. Namely, the physical 
features of the joint have been incorporated into its 
bilinear behaviour. Thus, both the set and the range of 
the material parameters for the joint have first been 
defined. Next, using an optimization method based on 
the full factorial test and application of the response 
surface method in combination with a real-valued 
genetic algorithm, an optimal combination of the 
values of the material parameters has been determined. 
Finally, the numerical model of the mechanical joint 
for a 2 mm thick sheet metal of S500MC material has 
been validated.
1  METHODOLOGY
1.1  Experimental Set-up
Two types of joints of a 2 mm thick sheet metal made 
of S500MC material were tested, a self-piercing rivet 
joint and a clinch joint. Self-piercing rivets with 
diameter of 5 mm were manufactured and mechanical 
clinching tool with 8 mm diameter was used for clinch 
joints. Static tests have been performed on the joints 
during experimental validation. 
Furthermore, two specimen types have been 
manufactured, the shapes of which were in compliance 
with the standard DVS/EFB 3480-1 [24], see Fig. 1. 
The width of the specimens was 45 mm.
Table 1.  An overview and comparison of self-piercing riveting (SPR) 
technology and clinch joint technology (adapted from [17])
Property SPR Clinching
Positioning on boreholes No No
Joining different materials Yes Yes
Joining three or more boards Good Poor
Sealing joint Yes Yes
Corrosion resistance Good Good
Visuality of the joint Good Medium
Access with the tool Bilateral Bilateral
Energy consumption Low Low
Environmental impact Small Small
Dynamic Strength Good Medium
Cycle time Short Short
Price High Low
Strojniški vestnik - Journal of Mechanical Engineering 68(2022)1, 3-13
5Incorporation of a Simplified Mechanical Joint Model into Numerical Analysis
Fig. 1.  The specimens for physical tests:  
a) overlapping joint and b) flanged joint for peel tests
Fig. 2.  Loading conditions: a) for shear strength test; b) for the 
peel test; c) for the combined test; clamps of the test rig and 
the compensation tabs are pointed out during application of the 
combined load
Specimen types corresponded to two load cases, 
which are also prescribed by the standard DVS/EFB 
3480-1 [24], and an additional combination of these 
two load cases, see Fig. 2:
• Static shear strength test for an overlapping joint.
• Static peel strength test for a flanged joint.
• Static combined tensile-shear test for a flanged 
joint, with a load at an angle of 45°. Additionally, 
compensation tabs to accommodate the 
misalignments of the specimens were adopted 
during the of tensile-shear test.
Fig. 3.  Layout of the test rig
Each loading condition was repeated five times, 
for both types of joints. In total, 15 specimens for each 
type of the joint were tested. The tests were carried 
out on a certified Zwick/Roell Z150 test rig, see Fig. 
3. The specimens were inserted in the clamps in the 
length of 40 mm. The load was applied by means 
of the clamp displacement at a rate of 10 mm/s 
with monitoring the force on the clamp at a given 
displacement. The results of each five repetitions were 
averaged, and the resulting force-displacement curve 
served as the basis for comparison with numerical 
simulations.
1.2  Numerical Model
Simulia Abaqus software was used to build the 
simplified finite-element model. The sheet metals 
were modelled with shell elements. The joint was 
modelled as a line connecting sheets by a kinematic 
or distributing connection as shown in Fig. 4. The 
size and the quality of the finite element mesh was 
generated according to the industrial requirements 
for crash simulations. Specifically, Timoshenko 
beams designated as B32 type in Abaqus were used 
to represent either the self-piercing rivet or the clinch 
joint. One beam connected the centres of the holes and 
several beams connected each centre of the hole with 
the circumferential nodes (Fig. 4).
A structured and denser mesh was created around 
the joint. Rectangular S4 type shell elements and 
triangular S3R type shell elements were applied. 
The size of the bore was equal to the diameter of the 
Strojniški vestnik - Journal of Mechanical Engineering 68(2022)1, 3-13
6 Glavan,	M.	—	Klemenc,	J.	—	Malnarič,	V.	—	Šeruga,	D.
connecting element representing the joint. A general 
contact was assigned between the sheets taking 
into account the friction between the sheets by a 
coefficient µ = 0.6 [25]. The boundary conditions 
represented the clamping of the sheets into the test 
rig. This involved fixed nodes on one sheet metal and 
application of prescribed displacement to the nodes on 
the other sheet metal to mimic the experimental setup 
as depicted in Fig. 2. The prescribed displacement was 
applied in a series of steps. The size of each step was 
limited upwards with 10 % of the final value of the 
displacement.
Fig. 4.  Model of the mechanical joint
Fig. 5.  Bilinear true stress-strain curve used to model  
S500 MC sheet metal behaviour
A bilinear true stress-strain curve for S500 MC 
is given in Fig. 5. The sheet metal was modelled 
using these properties. A bilinear material model 
and round cross-sections were also assigned to the 
line elements, which were used to model the joint. 
The diameter of the cross-section was equal to the 
cross-section of the connecting element. However, 
the parameters of the bilinear elastoplastic properties 
of the joint differed from the true material properties 
for S500 MC and rather represented its physical 
response under the loading. All the peculiarities of 
the joint were hence integrated and represented by 
the bilinear characteristic of the joint. The goal of the 
simulations was then to determine the optimal values 
of the parameters of the modelled joint. In doing so, 
a constant value of Poisson's ratio and a constant 
material density were assumed. Optimal values of the 
following parameters were searched for, which fully 
described the behaviour of the joint:
• elastic modulus, E [MPa],
• tangent modulus, ET [MPa], and 
• yield stress, Rp0.2 [MPa].
To define the optimal values of these parameters, 
a full factorial test was performed which used the 
response surface method in combination with genetic 
algorithm. In the first step, domain limits were defined 
within which the parameter values were searched 
for and these domains were then further divided into 
multiple values within the search area. In the next step, 
possible combinations p of the parameter values were 
defined. For each of the combinations, a numerical 
analysis of the tensile test was then performed 
including each of the three load cases. In total, 3p 
simulations for each type of joint were performed. 
The numerically calculated force-displacement curves 
for each parameter combination and load case were 
compared to the experimental results. The agreement 
between the measured and the modelled response of a 
joint was then estimated using the sum of the squared 
distances (SSQD) for the force differences at a certain 
displacement:
 SSQD F x F xj i i
i
K j
( )
exp( ( ) ( )) , 

 sim 2
1
 (1)
where Fexp represents the measured force from the 
tensile test, Fsim is the calculated force as a result of 
numerical simulations, j is the associated load case, 
Kj is the number of repetitions of a load case, and x 
represents the displacement. Since the main objective 
was to determine the material parameters that would 
optimally fit all three load cases, a multi-criteria 
objective function was defined as:
 F w SSQDj
j
j
n
C
T
 

 ( ) ,
1
 (2)
where wj represents a weight by which the significance 
or impact of each load case is attributed and nT stands 
for the number of load cases. For each combination 
of parameters defined in Table 2, the value of the 
objective function was calculated as:
 
F F E E R
m k n l o
m n o m n o
C C T p
( , . ) ( ) ( )
.
( )( , , ;
,..., ; ,..., ; ,.
=
= = =
0 2
1 1 1 .., ,r  (3)
where m, n, and o are the running indices through 
the levels of the material parameters, which are 
listed in Table 2. The most suitable values of the 
Strojniški vestnik - Journal of Mechanical Engineering 68(2022)1, 3-13
7Incorporation of a Simplified Mechanical Joint Model into Numerical Analysis
material parameters are those that give the lowest 
value of the objective function FC. The known values 
of combinations of parameters, together with the 
associated value of the objective function, represent 
the points in a four-dimensional space. These points 
are distributed relatively evenly throughout the 
space. Since the optimal combination of parameters 
lies within this range, it is not necessary for the 
optimal combination to be one of the combinations 
that were defined initially; the optimal combination 
might lie between the chosen points. Therefore, a 
response surface method was used, which enables an 
approximation of the surface through the calculated 
points in the four-dimensional space. This step then 
makes possible to determine the value of the objective 
function for any combination of material parameters, 
even for that combination which was not used during 
the numerical analysis.
The response surface for the objective function 
FC was modelled as the sum of global trend FC,glob  and 
local specifics FC,loc of the objective function [10], [15], 
and [26]. The local deviations of the objective function 
from the global trend were modelled using a mixture 
of multivariate elliptical Gaussian functions whereas 
the whereas the global trend was modelled using a 
third-degree polynomial for the three independent 
variables (E, ET and Rp0.2) and one dependent variable 
(objective function FC):
F x x x a a x a x x
a x x a
r q q
r
qr
C,glob ( 1 2 3 0 2 4 1 2
1
3
1
3
2 5 1 3
, , ) , ,
,
  
 


2 6 2 3 3 4 1
2
2
3 5 1 2
2
3 6 1
2
3 3 7 1 3
2
3 8 2
2
3
, ,
, , ,
,
x x a x x
a x x a x x a x x
a x x

  
  a x x a x x x3 9 2 3
2
3 10 1 2 3, , ,  (4)
where:
 x
x
x
x
E
E
R






















1
2
3 0 2
T
p .
.  (5)
Parameters a0 and ar,q are the parameters of 
the polynomials determined by the polynomial 
approximation for the p points from Eq. (3). The 
coefficients of the polynomial were estimated from 
the experimental data using the pseudo-inverse matrix 
as follows:
 a = X X X y( ) ,-T T⋅ ⋅ ⋅1  (6)
where the individual vectors and the matrix are the 
following:
a y















a
a
a
F
F
F
C
C
C
k l r
0
1 1
3 10
1 1 1
1 1 2
,
,
( , , )
( , , )
( , , )
,
















,
( ) ( ) ( ) ( ) ( ) ( )
(
x
1
1
1 1 1 1 1 1
1
E E R E E R
E
T p0.2 T p0.2
) ( ) ( ) ( ) ( ) ( )
( ) ( ) (
E R E E R
E E Rk l r
T p0.2 T p0.2
T p0.2
1 2 1 1 2
1

       
) ( ) ( ) ( )
.
 E E Rk l rT p0.2














(7)
The displacement points, in which the forces 
were compared, were selected based on the physical 
test results. As different load cases withstood different 
breakage forces at different final displacements, 
there was no possibility to choose exactly the same 
displacement points for all load cases. Based on the 
displacement-force diagram from the physical tests, 
the breakage point was established for each load 
case. This displacement was then used in numerical 
simulation as a maximum displacement by which the 
samples were loaded. The whole displacement range 
was divided equidistantly into 10 points, the final 
point being the breakage point from physical test. 
Those points were the displacement points in which 
the measured and the calculated forces were compared 
and used in Eq. (1).
Having determined the parameters of the 
polynomials that describe the global trend of the 
response surface, the residual RC between the actual 
value of the objective function FC and the value of 
its global trend of the FC,glob could be determined for 
each point from Eq. (3) as:
 
R F F
m k n l o
m n o m n o m n o
C C C,glob
( , , ) ( , , ) ( , , ) ,
,..., , ,..., ,
 
  1 1 1,..., .q  (8)
Then, the RC values were integrated into 
a conditional-average estimator using local 
multidimensional Gaussian functions to obtain a 
model of local specifics of the objective function 
FC,loc [10], [15], and [25]: 
F
R W X X
W X X
m n o m n o m n o
o
q
n
l
m
k
C,loc
C

  



( , , ) ( , , ) ( , , );
111
( , , ) ( , , );
,
m n o m n o
o
q
n
l
m
k
  
 111
 (9)
where R m n oC
( , , )  stands for residuals from Eq. (8), x are 
vectors from Eq. (5), W represent multidimensional 
Gaussian functions with median values vector in 
points x(m, n, o) and Σ(m, n, o) are covariant matrices 
Strojniški vestnik - Journal of Mechanical Engineering 68(2022)1, 3-13
8 Glavan,	M.	—	Klemenc,	J.	—	Malnarič,	V.	—	Šeruga,	D.
of individual Gaussian functions. The individual 
multidimensional Gaussian functions W were hence 
defined as:
W m n o m n o
m n o
dx x  
 



 
 
   
  
, , , ,
, ,
;
det
exp


1
22
1
2

  

    

      x x x xm n o m n o m n o, , , , , , ,
T

1
 (10)
where d is the dimension of the vector of independent 
variables, in this case the number of load cases (d = 3). 
Covariance matrices are diagonal with their variances 
equal to the distances to the nearest adjacent point in 
the direction of the individual components of vector 
x. In order to facilitate the data processing, individual 
variables (independent and dependent) were converted 
to logarithmic values before the process of global-
local modelling of the objective function.
In the last step, the objective function surface was 
calculated and the response surface was defined:
 F X F X F XC C,glob C,loc( ) ( ) ( ).   (11)
This way it was possible to model a multi-
dimensional surface that arbitrarily approaches any 
actual value of the measured point or, in this case, the 
result of the numerical analysis for each combination 
of parameters, while maintaining the trend without 
local deviations even outside of the measurement 
range [18]. The response surface thus assists in finding 
the domain of the optimal parameter values, but it 
is still difficult to accurately determine the actual 
optimum. Here, a real-valued genetic algorithm was 
hence utilised, which enabled the determination of 
the optimal combination of material parameters of 
the joint, minimizing Eq. (11). The in-house genetic 
algorithm used in this study combined the classic 
single-point crossover and linear crossover. The 
mutation-with-momentum of the child chromosomes 
enabled a better fine-tuning of the later generations 
of chromosomes to the optimal solution [26] to 
[28]. A population size of m = 20 randomly selected 
starting points on the response surface was used 
with the probability of crossover pcr = 0.5, fraction 
of linear crossover flin = 0.6 and probability of 
mutation pmut = 0.05. The algorithm was stopped 
after 10000 generations which is usually sufficient 
for a satisfactory optimum [26]. The above-described 
procedure for determining the optimum values of the 
parameters E, ET and Rp0.2 was repeated for both joint 
types.
Fig. 6.  Broken-down specimens of: a) SPR joint, and b) clinch joint
2  RESULTS
During the test, the force-displacement dependence 
was monitored and recorded until failure. In Fig. 
6, the broken specimens are shown. The results of 
the measurements were filtered and represented in 
diagrams for all test repetitions of each load case, and 
then, based on the curves obtained, the average force-
displacement curve was defined. This curve later 
served as the comparative indicator in the objective 
function. The average curve was generated to the 
point where the joint started to loosen and the value of 
the force began to decline, which is also seen in Figs. 
7 and 8. The measurements showed good repeatability 
for both joint types and all three load cases. Only one 
of the results of the clinch joint measurements for the 
shear load case was later excluded as that specific 
measurement deviated incomparably from the other 
measurements. The final displacement of the average 
curve was also the displacement value up to which the 
numerical model of the joint was loaded.
Using the developed methodology, an optimal 
combination of the material parameters of the joint 
was generated, which on average can simulate the 
force-displacement responses best for all the three 
load cases, so they are hence suitable to use in a 
general model for any other load case. Additionally, 
the methodology was used to determine the material 
parameters considering only a single load case. This 
served for verification purposes and confirmation of 
Strojniški vestnik - Journal of Mechanical Engineering 68(2022)1, 3-13
9Incorporation of a Simplified Mechanical Joint Model into Numerical Analysis
the best-fit ability of the method itself. During the 
parameter-estimation process, a starting point in the 
set of parameters was defined and a range of values 
within the algorithm performed its operation.
The range of values of the material parameters 
was divided approximately equidistantly in the 
logarithmic scale within the defined range. Table 
2 shows the values and the range of parameters. In 
total, p = 175 combinations of material parameters 
from Table 2 are possible for a single load case, which 
results in a total of 3p = 525 numerical analyses.
a)   
b)   
c)   
Fig. 7.  Measurement results for SPR joint: a) shear load case,  
b) peel load case, and c) combined load case
a)   
b)   
c)   
Fig. 8.  Measurement results for clinch joint: a) shear load case,  
b) peel load case and c) combined load case
Table 2.  Values	of	parameters	E [MPa], ET [MPa] and Rp0.2 [MPa]
Parameter Search domain The values of the parameters taken into account
Elastic modulus - E 200-300000 200 500 1000 3000 10000 30000 300000
Tangent modulus - ET 0.1·E- 0.001·E 0.1·E 0.0286·E 0.01·E 0.00286·E 0.001·E
Flow limit - Rp0.2 50-1000 50 100 200 500 1000
Strojniški vestnik - Journal of Mechanical Engineering 68(2022)1, 3-13
10 Glavan,	M.	—	Klemenc,	J.	—	Malnarič,	V.	—	Šeruga,	D.
which serve as a comparative value for numerical 
simulations. Choosing a more complex material 
model might lead to a more accurate understanding of 
the joint behaviour, while on the other hand, results 
in a more complex definition of the objective function 
and extends the numerical calculation time, which is 
of crucial importance in practice. Likewise, it is not 
possible to accurately determine the values of the 
material parameters if the search area is too wide 
or too narrow. If the search domain differs from the 
final area presented in Table 2, it may happen that the 
positions of the optimal parameter values lie outside 
of the defined search range. An example of such an 
inadequately defined search range is presented for the 
SPR joint in Fig. 11. The comparison of the results 
between the inadequately and adequately defined 
3  DISCUSSION
The force-displacement results of the measurements 
indicated that the response of each joint had a 
pronounced bilinear characteristic. Accordingly, the 
stress-strain behaviour of the connecting element in 
the joint was hence modelled using a bilinear material 
model. The failure of both the SPR and the clinch joint 
occurred in the very connection or in the connective 
element (Fig. 6). In the case of the SPR joint, the rivet 
was pulled out, and in the case of the clinch joint, the 
mechanical connection between the two sheets was 
broken. From this point of view, attention could be 
paid only to the connective element and its modelling, 
whilst the base material was modelled with a simple 
elastic-plastic material model.
It can be noted that in the case when each load case 
is treated separately (Figs. 9 and 10 – specific result), 
it is possible to accurately predict the combination of 
parameters that produce a numerical response of the 
joint almost identical to the measured response. By 
changing the range of the input parameters or even 
the type of the material model, a better agreement 
might also be possible to achieve. However, the 
method is primarily intended to determine the material 
parameters of the joints for general use in large 
structures with a large number of such joints whilst 
they are subjected to a random load. This response 
should therefore be the best approximation of the 
observed response of the joint regardless of the size 
and the direction of the applied load.
The method thus enables consideration of the 
agreement between the numerical and experimental 
response of the joint for several load cases. 
Consequently, the resulting values of the material 
parameters are always a compromise for several 
load cases and the deviation between the numerical 
response and the experimental response. A better 
agreement is hence always observed when only a 
single load case is considered during the estimation of 
the optimal material parameters.
Furthermore, other types of joints will have 
other properties therefore it is necessary to perform 
some initial tests under shear, peel and combined 
tensile-shear load conditions to be able to search for 
the optimal values of the parameters in the numerical 
model. After these have been performed however, the 
numerical model of a complex structure with such a 
joint can be used under arbitrary load conditions.
For the correct performance of the method, it 
is essential to correctly define the material model, 
the range of possible values of associated material 
parameters and to perform accurate measurements 
a)   
b)   
c)   
Fig. 9.  Comparison of numerical and actual response of SPR joint: 
a) shear load case, b) peel load case, c) combined load case
Strojniški vestnik - Journal of Mechanical Engineering 68(2022)1, 3-13
11Incorporation of a Simplified Mechanical Joint Model into Numerical Analysis
search ranges reveals that the force-displacement 
description can drastically differ if the initial search 
domain of the parameters is too wide. 
The material model of the joint used in this 
study succeeds in describing the force-displacement 
response as long as the force increases with increasing 
displacement and shows a pronounced bilinearity. 
When the maximum force is achieved in practice, the 
joint begins to loosen, the force begins to decline, and 
the joint then gradually breaks down. In numerical 
analyses however, such behaviour of the joint cannot 
be modelled using a bilinear material model. Hence, 
only the “functional” part of the joint behaviour can 
be simulated, i.e., the increasing part of the force-
displacement characteristics until the maximum 
force. For the correct performance of the developed 
methodology this fact has to be taken into account 
in numerical simulations so that the final simulated 
displacement coincides with the experimental 
displacement at the maximum value of the force. If 
the displacements exceeding this limit are considered 
during the optimization of the material parameters, the 
error of the method will increase as a consequence of 
the diverging experimental and simulated curves.
The final values of the material parameters can 
also be influenced by weights assigned to the load 
cases in the objective function. The importance of a 
specific load case can hence be considered. This is 
most easily explained if the values of forces that occur 
in individual load case are observed. The force in the 
case of the shear test is approximately twice as high as 
in the other load cases. Consequently, at higher force 
values, the relative difference of the numerical and the 
experimental result will be comparable to other load 
cases, whilst the absolute value of the difference will 
be significantly higher than for the other load cases.
a)   
b)   
Fig. 11.  Force-displacement results for: a) inadequately defined 
search range, and b) adequately defined search range
Hence, the value of the objective function 
strongly depends on the tensile case and less on the 
other two cases. Thus, by using weights, individual 
contributions of various load cases to the objective 
function can be balanced and their importance for the 
final result can be evenly distributed. Nevertheless, 
the above-mentioned differences did not provide a 
significant impact on the final result during this study, 
a)   
b)   
c)   
Fig. 10.  Comparison of numerical and actual clinch joint 
response: a) shear load case, b) peel load case,  
c) combined load case
Strojniški vestnik - Journal of Mechanical Engineering 68(2022)1, 3-13
12 Glavan,	M.	—	Klemenc,	J.	—	Malnarič,	V.	—	Šeruga,	D.
therefore equal weights were attributed to the various 
load cases. The optimal values of the parameters to 
simulate the response of the self-piercing rivet joint 
were determined as E = 8980 MPa, ET = 260 MPa 
and Rp0.2 = 217 MPa. Similarly, the optimal values of 
the parameters for the clinch joint were determined as 
E = 1042 MPa, ET = 1.5 MPa and Rp0.2 = 195 MPa. 
Nevertheless, these values of the parameters represent 
synthetic properties of the joint and not the material 
properties of the assembly.
4  CONCLUSIONS
The force-displacement response of both the self-
piercing rivet joint and the clinch joint tested during 
this study exhibited distinctly bilinear characteristics 
up to the highest measured force before failure. The 
presented methodology could be applied to determine 
optimal values of material parameters for the 
bilinear material model. Numerical analyses of large 
structures with a number of such joints, each having 
the same values of the material parameters, can thus 
be performed to simulate their stress-strain or force-
displacement behaviour under variable loading. The 
optimal values of the material parameters represent a 
compromise between the quality of fit and the variety 
of load cases.
5  ACKNOWLEDGEMENTS
The authors acknowledge financial support from the 
Slovenian Research Agency (research core funding 
No. P2-0182 entitled Development evaluation).
6  REFERENCES
[1] Barimani-Varandi, A., Aghchai, A.J., Lambiase, F. (2021). 
Failure behavior in electrically-assisted mechanical clinching 
joints. Journal of Manufacturing Processes, vol. 68 p. 1683-
1693, DOI:10.1016/j.jmapro.2021.06.072.
[2] Mori, K., Abe, Y. (2018). A review on mechanical joining 
of aluminium and high strength steel sheets by plastic 
deformation. International Journal of Lightweight Materials 
and Manufacture, vol. 1, no. 1, p. 1-11, DOI:10.1016/j.
ijlmm.2018.02.002.
[3] Yang, Y., Xiong, X., Chen, J., Peng, X., Chen, D., Pan, F. (2021). 
Research advances in magnesium and magnesium alloys 
worldwide in 2020. Journal of Magnesium and Alloys, vol. 9, 
no. 3, p. 705-747, DOI:10.1016/j.jma.2021.04.001.
[4] Abe, Y., Kato, T., Mori, K. (2007). Joining of aluminium alloy 
and mild steel sheets using mechanical clinching. Material 
Science Forum, vol. 561 1043-1046, DOI:10.4028/www.
scientific.net/MSF.561-565.1043.
[5] He, X., Zhao, L., Yang, H., Xing, B., Wang, Y., Deng, C., Gu, 
F., Ball, A. (2014). Investigations of strength and energy 
absorption of clinched joint. Computational Materials Science, 
vol. 94, p. 58-65, DOI:10.1016/j.commatsci.2014.01.056.
[6] Maiwald, M., Thiem, J. (2012). Joining without a pilot hole 
using friction welding. ATZ Production Worldwide, vol. 5 10-17, 
DOI:10.1365/s38312-012-0045-0.
[7] Tozaki, Y., Uematsu, Y., Tokaji, K. (2010). A newly developed tool 
without probe for friction stir spot welding and its performance. 
Journal of Materials Processing Technology, vol. 210, no. 6-7, 
p. 844-851, DOI:10.1016/j.jmatprotec.2010.01.015.
[8] Özdemir, U., Sayer, S., Yeni, Ç. (2012). Effect of pin penetration 
depth on the mechanical properties of friction stir spot welded 
aluminium and copper. Material Testing, vol. 54, no. 4, p. 233-
239, DOI:10.3139/120.110322.
[9] Perrett, J.G., Martin, J., Threadgill, P.L., Ahmed, M.M. (2007). 
Recent Developments in Friction Stir Welding of Thick Section 
Aluminium Alloys. 6th World Congress, The aluminium Two 
Thousand. p. 1-11.
[10] Porcaro, R., Hanssen, A.G., Aalberg, A., Langseth, M. 
(2004). Joining of aluminium using self-piercing riveting: 
Testing, modeling and analysis. International Journal of 
Crashworthiness, vol. 9, no. 2, p. 141-154, DOI:10.1533/
ijcr.2004.0279.
[11] Abe, Y., Kato, T., Mori, K. (2006). Joinability of aluminium 
alloy and mild steel sheets by self piercing rivet. Journal of 
Materials Processing Technology, vol. 177, no. 1-3, p. 417-
421, DOI:10.1016/j.jmatprotec.2006.04.029.
[12] Li, D., Chrysanthou, A., Patel, I., Williams, G. (2017). Self-
piercing riveting-a review. The International Journal of 
Advanced Manufacturing Technology, vol. 92, p. 1777-1824, 
DOI:10.1007/s00170-017-0156-x.
[13] Di Franco, G., Fratini, L., Pasta, A., Ruisi, V.F. (2013). On the 
self piercing riveting of aluminium blanks and carbon fibre 
composite panels. International Journal of Materials Forming, 
vol. 6, p. 137-144, DOI:10.1007/s12289-011-1067-2.
[14] Su, Z.-M., Lin, P.-C., Ali, W.-J., Pan, J. (2015). Fatigue analyses 
of self-piercing rivets and clinch joints in lap-shear specimens 
of aluminum sheets. International Journal of Fatigue, vol. 72, 
p. 53-65, DOI:10.1016/j.ijfatigue.2014.09.022.
[15] Abe, J., Mori, K., Kato,T. (2012). Joining of high strength steel 
and aluminium alloy sheets by mechanical clinching with 
dies for control of metal flow. Journal of Materials Processing 
Technology, vol. 212, no. 4, p. 884-889, DOI:10.1016/j.
jmatprotec.2011.11.015.
[16] Martinsen, K., Hu, S.J., Carlson, B.E. (2015). Joining of 
dissimilar materials. CIRP Annals, vol. 64, no. 2, p. 679-699, 
DOI:10.1016/j.cirp.2015.05.006.
[17] Kaščak, L., Spišak, E. (2013). Mechanical joining methods in 
car body construction. Transfer inovacii, vol. 28, p. 168-171.
[18] Škrlec, A., Klemenc, J. (2017). Parameter identification for a 
Cowper-Symonds material model using a genetic algorithm 
combined with a response surface. Engineering Computations, 
vol. 34, no. 3, p. 921-940, DOI:10.1108/EC-03-2016-0099.
[19] Bouchard, P.O., Laurent, T., Tollier,L. (2008). Numerical 
modeling of self-pierce riveting-From riveting process modeling 
down to structural analysis. Journal of Materials Processing 
Technology, vol. 202, no. 1-3, p. 290-300, DOI:10.1016/j.
jmatprotec.2007.08.077.
Strojniški vestnik - Journal of Mechanical Engineering 68(2022)1, 3-13
13Incorporation of a Simplified Mechanical Joint Model into Numerical Analysis
[20] Yang, L., Yang, B., Yang, G. W., Xiao, S. N., Zhu, T., Wang, F. 
(2020). S-N Curve and Quantitative Relationship of Single-Spot 
and Multi-Spot Weldings. International Journal of Simulation 
Modelling, vol. 3, p. 482-493, DOI:10.2507/IJSIMM19-3-CO11.
[21] Kulawik, A., Wrobel, J. (2020). Numerical Study of Stress 
Analysis for the Different Widths of Padding Welds. Strojniški 
vestnik - Journal of Mechanical Engineering, vol. 66, p. 567-
580, DOI:10.5545/sv-jme.2020.6771.
[22] Sommer, S., Maier, J. (2011). Failure modeling of a self-
piercing riveted joint using LS-DYNA. 8th	 European	 LS-DYNA	
Conference, Strasbourg.
[23] Hanssen, A.G., Oolsson, L., Porcaro, R., Langseth, M. (2010). 
A large-scale finite element point-connector model for self-
piercing rivet connections. European Journal of Mechanics- 
A/Solids, vol. 29, no. 4, p. 484-495, DOI:10.1016/j.
euromechsol.2010.02.010.
[24] DVS/EFB 3480-1:2010 (2010). Testing of Properties of Joints 
- Testing of Properties of Mechanical and Hybrid Joints. DVS	
Media, Düsseldorf.
[25] Wittel, H., Muhs, D., Jannasch, D., Voßiek, J. (2013). 
Roloff/Matek Maschinenelemente, Springer Fachmedien, 
Wiesbaden, DOI:10.1007/978-3-658-02327-0.
[26] Škrlec, A., Klemenc, J. (2020). Estimating the strain-rate-
dependent parameters of the Johnson-Cook material model 
using optimisation algorithms combined with a response 
surface. Mathematics, vol. 8, no. 7, p. 1-18, DOI:10.3390/
math8071105.
[27] Wright, A.H. (1991). Genetic Algorithms for Real Parameter 
Optimization. Foundations of Genetic Algorithms, Rawlins, 
G.J.E. (ed.), vol. 1, Morgan Kaufman, San Mateo, DOI:10.1016/
B978-0-08-050684-5.50016-1.
[28] Temby, L., Vamplew, P., Berry,A. (2005). Accelerating real-
valued genetic algorithms using mutation-with-momentum. 
The 18th Australian Joint Conference on Artificial Intelligence, 
Sydney, p. 1108-1111, DOI:10.1007/11589990_149.