Modeling Human–Robot Impact Dynamics in Collaborative Applications

1. Introduction

Human–robot collaboration (HRC) has emerged as a transformative paradigm in manufacturing and industrial environments over the past decade. Driven by increasing demand for flexible, adaptive, and safe robotic systems capable of operating in close proximity to humans  [1], collaborative robots (cobots) have fundamentally altered field of industrial automation; in  [2], main standards were published to guarantee safe human-robot interactions, dividng human-robot co-tasks into 3 groups: Hand-guided control (HGC), speed and separation monitoring (SSM) and power-force limiting (PFL) . Unlike traditional industrial robots confined to protected industrial cells with perimeter fencing, collaborative robots are explicitly designed to share workspaces with human operators, enabling dynamic task allocation and real-time human–robot teaming. However, this physical proximity introduces unprecedented challenges in ensuring safety during both intentional and accidental contacts [3].

Especially in PFL approach, rather than eliminating human–robot contact, collaborative systems must render such contact inherently safe. This requires a rigorous understanding of contact mechanics, impact dynamics, and the mechanical response of both the robot and the human body. International standards such as ISO/TS 15066:2016 (Robots and Robotic Devices—Collaborative Robots) [4], and ISO 10218-2:2025 (Safety Requirements for Industrial Robots—Part 2: Robot Systems and Integration) [5] define stringent safety requirements for collaborative operations. These standards emphasize the critical importance of limiting both quasi-static contact forces and transient impact forces to prevent injury during accidental collisions.

Beyond the standards, other studies introduced different ways to determine biomechanical limits to avoid injuries. As example, in  [6] Behrens proposed a method to determine biomechanical limits based on the probability of injury, acquiring data from real collision experiments.

Impact scenarios are particularly challenging from a safety perspective, as their outcomes depend on multiple factors: robot control laws and trajectory planning, the effective mass distribution along the kinematic chain, relative velocity at contact, mechanical properties of the contact interfaces, and the anatomical and physiological characteristics of the human body segment involved [7]. A collision that may be safe from the perspective of force magnitude, may cause serious injury if the contact is concentrated over a small area or if the impact energy is dissipated in a way that amplifies local stress. Hence, collaborative systems must be designed in such a way that impact energy is absorbed and distributed across larger body regions  [8] (e.g., by employing round surfaces). Another way of increasing human safety is to employ algorithms and methods for collision detection and avoidance [9]. Therefore, understanding and modelling the dynamic response of the robot during collision events is essential for both validating safe physical interaction and designing control strategies that minimize injury risk. Recent findings also indicate that collision detection itself can be formulated as an optimization problem, enabling significantly faster and more scalable distance queries for real-time safety assessment [10]. Scoccia et al. [11] proposed a method to avoid obstacles during human-robot collaboration, using a depth camera to detect human presence and adapting robot trajectory in real time.

The scientific literature on human–robot impact dynamics reveals two complementary but distinct approaches, each with inherent strengths and limitations  [12]. In the first approach, experimental investigations provide empirical force–deformation data through carefully instrumented tests conducted with human surrogates (anthropomorphic dummies, tissue-equivalent materials) or compliant prototypes designed to replicate human body mechanics [13,14]. These experiments capture realistic, complex behaviour including nonlinear material response, friction effects, and anatomical detail. A representative example of such experimental work is provided by Fischer et al. [15], who performed extensive collision tests with a power-and-force measuring device to study how impact geometry, body-region stiffness, and damping affect measured forces and pressures relative to ISO/TS 15066 biomechanical limits. However, physical testing is often expensive, time-consuming, difficult to reproduce across different laboratories, and inherently limited to the specific test configurations actually performed. In the second approach, several analytical and numerical models have been developed to describe contact phenomena using mathematical formulations such as nonlinear stiffness–damping models, most notably the Hertzian elastic theory and the Hunt–Crossley viscoelastic model [16,17,18]. Analytical models are less time consuming, offer transparency and the ability to conduct parametric studies across a wide range of conditions. On the other hand, every model inherits some simplifications which yield to some errors in the prediction of the contact forces. Furthermore, model parameters are not easy to measure during online tasks, and the validation of the models with real data for collaborative robots remains limited in the literature.

A recent and promising research trend is the integration of physics-based formulations with data-driven and machine learning methods to predict impact outcomes under diverse and previously unseen conditions. These hybrid approaches aim to combine the interpretability and efficiency of analytical models with the fidelity and realism of experimental data.

In [19], Rezayati et al. trained a neural network to clusterize incidental and intentional collisions. Similarly,  [20] collected experimental data from the prototype presented in [14] and used them to train a neural network to predict the contact force during human–robot collisions, using contact velocity and distance from the robot base as input parameters.

Alongside learning-based approaches, empirical formulations have been proposed to estimate impact forces using a limited set of physical parameters. In [21], Tae-Jung Kim proposed an empirical formulation to predict the maximum impact force using only effective mass and contact velocity as input parameters, consistently with the quantities defined in ISO/TS 15066:2016. However, this formulation was derived from only 50 collision experiments performed on a dummy representing solely the human face.

The quantitative assessment of safe collaborative operation requires accurate estimation of the robot’s effective mass along its kinematic chain. The effective mass, a fundamental concept introduced by Khatib in the context of operational space control [22,23,24], represents the apparent inertia of the robot end-effector or of a generic body point as experienced in a given direction and configuration. During impact, the effective mass directly determines the contact force for a given collision velocity and contact stiffness.

In [25], Kirschner proposed a method to determine the effective mass, and consequently in [26], Steinecker developed the concept of mean reflective mass (MRM). Furthermore, lumped-parameter models and impedance-based approaches have been recently refined to accurately characterize effective mass in realistic robotic systems, enabling precise yet computationally efficient prediction of impact dynamics [27], with improved consistency with ISO/TS 15066:2016.

Despite these advances, reliable and well-documented experimental data remain essential to calibrate advanced models and to ensure reproducible compliance with normative safety requirements.

This work aims to bridge the gap between empirical accuracy and analytical simplicity by developing a physics-guided modelling framework that combines the interpretability and efficiency of analytical contact models with the fidelity and realism of experimental data.

To this end, a large experimental dataset was collected by designing a low-cost and effective impact testing prototype compliant with PD ISO/PAS 5672:2023. An entire framework was then developed to compare classical contact models and machine learning approaches, in terms of their predictive accuracy for collaborative robot impacts. Then a predictive model was trained using joint angles and contact velocity as input parameters, enabling online execution directly on the robot and allowing the proposed framework to be employed as a real-time risk assessment tool.

In this work, the robot’s effective mass was estimated using an approximate yet computationally efficient method [27] based on the robot’s kinematic and dynamic parameters, enabling accurate computation of impact energy and safety metrics as defined in ISO/TS 15066:2016.

The primary contributions of this work are summarized as follows:

• The design, validation and implementation of a low-cost impact test prototype capable of reproducing controlled human–robot contact scenarios, enabling synchronized acquisition of force data during collision events.
• A quantitative comparative evaluation of classical contact models (Hertzian, Hunt–Crossley) and new machine learning approaches in terms of their predictive accuracy for robot impacts.
• An integrated framework combining experimental data and analytical modelling, leading to an approach that can be used for online risk assessment during PFL working mode, demonstrating compliance with international standards.

The paper is structured as follows. Section 2 reviews the current state of knowledge in contact mechanics and human–robot impact dynamics. Firstly, Section 3 provides a detailed description of the testing collision prototype and validation test protocol. After that, Section 4 describes the collaborative robot system employed, the test setup and the data collection process. Then, Section 5 reports quantitative comparisons between simulated predictions (using physical and machine learning models) and experimental measurements, discusses the implications for safety evaluation, and analyzes the robustness of model predictions. Finally, Section 6 summarizes the main findings, highlights the practical implications for collaborative robot design and deployment, and outlines promising directions for future research, including the integration of machine learning methods.

2. Human-Robot Impact Dynamics

This section establishes the safety framework for human-robot interaction in industrial environments, focusing on the Power and Force Limiting (PFL) mode regulated by the ISO TS 15066:2016 standard, which provides the theoretical context for the prototype developed. Human–robot impact dynamics is a multidisciplinary research area that integrates principles from robotics, biomechanics, materials science, and safety engineering. The primary objective is to understand and predict the mechanical interactions that occur during collisions between robots and humans, with the ultimate goal of ensuring safe and effective collaboration in shared workspaces; then regulate manipulator parameters to prevent or limit pain upon collision.

2.1. Contact forces

Safety protocols focus on limiting robot performance parameters such as forces, velocities, power, and momentum. Collision scenarios are primarily categorized as either quasi-static contacts or transient contacts. The first refers to situations where the robot applies a continuous force on the human body, typically during slow movements or when the human body part is spatially costrained, while the latter involves tipically brief, spatially uncostrained high-impact events where the robot collides with the human body at higher speeds. Tipically an impact event can be subdivided in a transient contact part (duration of 0.5 s according to standard), followed by a quasi-static contact part, when the robot remains in contact with the human body after the impact.

The normative framework models the human body by partitioning it into 12 regions [4]. To each region is assigned a specific maximum permissible force $F_{max}$ and is modeled, from a mechanical perspective, by an effective mass $m_{body}$ and an effective stiffness $k_{body}$. The force limit is critical as it relates to the thresholds for pain onset.

2.2. Energetic Approach

For defining transient contact, the standard model [4] applies a simplified physical-mathematical model of a completely inelastic collision between the robot and the effective mass representing the human body region. The maximum allowed energy transfer E at the human-robot interface is calculated using the biomechanical limits $F_{max}$ and the region’s effective stiffness $k_{body}$, equal to the elastic energy stored in the equivalent spring:

$$E={\displaystyle \frac{F_{max}^2}{2k_{body}}}$$

(1)

This energy is equated to the kinetic energy of the system, expressed via the reduced mass $\mu$ and the relative velocity $v_{rel}$:

$$E={\displaystyle \frac{F_{max}^2}{2k_{body}}}={\displaystyle \frac{1}{2}}\mu v_{rel}^2$$

(2)

The relative velocity $v_{rel}$ is normally calculated as the difference between the robot’s velocity (considered as the point of the robot that will impact the human) and the human body region velocity, both projected along the impact direction.

The ISO standard [4] calculates the reduced mass $\mu$ using the effective mass of the human region $m_{body}$ and the effective mass of the robot $m_{robot}$:

$$\mu ={\left({\displaystyle \frac{1}{m_{body}}}+{\displaystyle \frac{1}{m_{robot}}}\right)}^{-1}$$

(3)

,

which can be substituted into eq. 2, so the force at impact $F_{max}$ can be derived as:

$$F_{max}=\sqrt{2k_{body}\mu }{v}_{rel}$$

(4)

According to the ISO/PAS 5672:2023, when test collisions with rigidly fixed impact prototypes are performed, it is not necessary to adopt a mass $m_{body}$ that matches the body part under analysis, but instead a scaling factor is used to convert the measured force $F_{meas}$ to an equivalent value of force $F_{eq}$ similar to the one which would have occurred when impacting the real body part. The scaling factor is defined as:

$$S=\sqrt{{\displaystyle \frac{m_{body}}{m_{body}+m_{robot}}}}$$

(5)

The force measured during the test $F_{meas}$ is converted to the equivalent force $F_{eq}$, according to eq. 5, as:

$$F_{eq}=F_{meas}S=F_{meas}\sqrt{{\displaystyle \frac{m_{body}}{m_{body}+m_{robot}}}}$$

(6)

The ISO/TS 15066 defines the robot effective mass as $m_{robo{t}_{ISO}}=M/2+m_L$, where M is the total mass of the robot’s moving parts and $m_L$ is the load mass. Furthermore, the model assumes a perfectly flat contact surface between the robot and the human body region.

Accurate estimation of $m_{robot}$ is essential to predict the forces generated during a collision and assess potential injury risks. A more in-depth analysis can be performed by calculating the effective mass $m_{robot}$ of the robot at the contact point from the configuration-dependent Cartesian mass matrix $M_C\left(q\right)$:

$$M_C\left(q\right)={\left(J\left(q\right)M{\left(q\right)}^{-1}J{\left(q\right)}^T\right)}^{-1}$$

(7)

$$m_{robot}\left(q,\widehat{n}\right)={\displaystyle \frac{1}{{\widehat{n}}^T{M}_C{\left(q\right)}^{-1}\widehat{n}}}={\displaystyle \frac{1}{{\widehat{n}}^{T}J\left(q\right)M{\left(q\right)}^{-1}J{\left(q\right)}^T\widehat{n}}}$$

(8)

where $M\left(q\right)$ is the joint-space inertia matrix which can be computed via the Lagrange approach [28], $J\left(q\right)$ is the geometric Jacobian of the contact point, and $\widehat{n}$ is the unit vector along the impact/force direction.

Figure 1. Lumped masses configuration in KUKA IIWA 14 R820. The robot configuration is the first one experimented in the Section 4.

Figure 1. Lumped masses configuration in KUKA IIWA 14 R820. The robot configuration is the first one experimented in the Section 4.

This approach approximates the effective mass by considering only the major contributors along the kinematic chain, significantly reducing computational complexity while maintaining accuracy [27].

2.3. Classical Contact Models

Several classical contact models have been developed to describe the mechanical interactions during collisions. Most of them are derived from the Hunt–Crossley viscoelastic model [16].

This model is able to describe the contact force between two bodies during an impact event, incorporating both elastic and damping effects, accounted for energy dissipation during impact. The contact force is expressed as:

$$F(\delta ,\dot{\delta })=k_{HC}{\delta }^n+{\chi }_{HC}{\delta }^m\dot{\delta }$$

(9)

where $\delta$ is the penetration distance of the contacting bodies, $\dot{\delta }$ the penetration velocity, $k_{HC}$ the Hunt–Crossley stiffness coefficient, which accounts for stiffness of the colliding bodies, ${\chi }_{HC}$ the interaction damping coefficient of the colliding bodies, and m and n are nonlinear exponents. The exponents m and n are typically set to 1.5 for spherical contacts [18], but can be adjusted based on experimental observations to fit specific contact scenarios.

The coefficient ${\chi }_{HC}$ in (9) is defined as:

$${\chi }_{HC}={\displaystyle \frac{3}{2}}\alpha k_{HC}$$

(10)

where $\alpha$ is a parameter depending on geometry, material and kinematics of the two bodies. The coefficient of restitution $c_r$, representing the ratio between the relative velocity of the two bodies after ($v_r^{+}$) and before ($v_r^{-}$) the impact ($c_r={\displaystyle \frac{v_r^{-}}{v_r^{+}}}$), in the Hunt–Crossley formulation is defined as [18]:

$$c_r=1-\alpha |{\dot{\delta }}_0|$$

(11)

where $|{\dot{\delta }}_0|$ is the relative velocity just before impact.

Introducing (11) into eq. (9) gives:

$$F(\delta ,\dot{\delta })=k_{HC}{\delta }^n\left[1+{\displaystyle \frac{3(1-c_r)}{2}}{\displaystyle \frac{\dot{\delta }}{|{\dot{\delta }}_0|}}\right].$$

(12)

Is important to underline that, knowing the penetration distance and velocity parameters, $F(\delta ,\dot{\delta })$ depends on three parameters only: $k_{HC}$, $c_r$ and n.

When the damping effects are negligible respect to the elastic ones, hence $c_r$ tends to 1, the Hunt–Crossley model simplifies to the Hertz contact model [29], which describes purely elastic contact interactions. The Hertz contact model assumes that the contact area is small compared to the dimensions of the bodies and that the materials behave elastically. The contact force F is related to the penetration depth $\delta$ by a nonlinear relationship:

$$F\left(\delta \right)=k_H{\delta }^n$$

(13)

where $k_H$ is the Hertzian stiffness coefficient and n is the nonlinear exponent, which depends on the materials properties and geometries of the contact bodies. This model captures the nonlinear elastic behaviour observed in many contact scenarios, particularly for spherical or curved surfaces.

2.4. Data-Driven Contact Force Prediction

Physics-based models often struggle to provide accurate predictions in highly complex scenarios [19]. For this reason, data-driven approaches based on Machine Learning (ML) have been widely adopted in the literature [20]. These methods enable predictions to be obtained solely from experimental data collected in laboratory settings, without requiring an explicit physical model of the system.

A fundamental assumption is the existence of an unknown functional input-output relationship f between the space of input measurements $\mathcal{X}$ and the corresponding force values. The learning task therefore consists in approximating such a function from data. To assess the quality of a candidate function, a loss function$\ell :\mathbb{R}\times \mathbb{R}\to {\mathbb{R}}_{\ge 0}$ is introduced, quantifying the error between predicted and observed outputs.

Additionally, the input–output pairs $(x,y)\in \mathcal{X}\times {\mathbb{R}}^d$, where $\mathcal{X}$ is an arbitrary measurable input space and $d\in \mathbb{N}$ denotes the output dimension, are assumed to be drawn from an unknown joint probability distribution ${\mathcal{P}}_{x\times y}$, and the ideal objective is the minimization of the expected risk

$$L\left(f\right)={\mathbb{E}}_{(x,y)\sim \mu }\left[\ell \left(f\right(x),y)\right],$$

(14)

namely the expected value ${\mathbb{E}}_{(x,y)}$ of the loss function $\ell \left(f\right(x),y)$ with respect to the unknown distribution ${\mathcal{P}}_{x\times y}$. Hence, $L\left(f\right)$ defines the expected risk of the predictor f, namely the average loss incurred by f over the data-generating process. As a consequence, $L\left(f\right)$ evaluates the performance of f at the population level, rather than on individual samples, and represents the ideal learning objective. However, as the underlying distribution ${\mathcal{P}}_{x\times y}$ is not accessible, the expected risk cannot be minimized directly. In practical ML settings, the function $f:\mathcal{X}\to \mathbb{R}$ is instead estimated by restricting the search to a predefined hypothesis space $\mathcal{F}$ and by exploiting a finite dataset of labeled examples ${\mathcal{D}}_n=\{(x_1,y_1),\dots ,(x_n,y_n)\}$, where each pair $(x_i,y_i)$ is a sample from the input-output space. The expected risk is approximated by the empirical risk

$$\widehat{L}\left(f\right)={\displaystyle \frac{1}{n}}\sum _{i=1}^n\ell \left(f\left(x_i\right),y_i\right),$$

(15)

leading to the following optimization problem:

$$\begin{array}{c}\hfill \underset{f\in \mathcal{F}}{min}\widehat{L}\left(f\right).\end{array}$$

(16)

Within this framework, Support Vector Machines for Regression (SVR) [30] and Decision Trees (DT) [31] are employed to obtain accurate predictions of the force from experimental measurements. Specifically, the input space $\mathcal{X}$ is defined by a set of input within: the robot joint angles $q_n$ at impact , cartesian velocity at impact ${\dot{\delta }}_0$, $\delta$, the maximun penetration distance ${\delta }_{max}$, $m_{body}$, $k_{body}$ and $m_{robot}$. The output space is given by the measured contact force values $F_{eq}$, as in eq. 6.

3. Experimental Apparatus

The experimental setup is designed to systematically investigate human–robot impact dynamics using a collaborative redundant manipulator (KUKA LBR iiwa 14 R820) and a custom-designed compliant prototype for collisions. This section details the design of the impact test prototype, its validation, the experimental setup, instrumentation, test protocol, and data acquisition procedures used to collect impact measurements.

3.1. Prototype Design

An impact test prototype was developed to replicate the mechanical behaviour of various human body regions under controlled impact conditions. The prototype is designed to be compliant to PD ISO/PAS 5672:2023 [32] , allowing for realistic deformation and force response during collisions with the robot. The prototype consists of several key components (Figure 2):

• Compliant Damping Layer: A layer of compliant material placed on top of the impact piece to simulate the soft tissue response of human body regions. The material properties are selected to closely match the biomechanical characteristics of human tissues [13], ensuring realistic force-deformation behaviour during impacts.
• External Cover: A rigid outer shell that protects the internal components and permits a proper impact stimulation of the load cell. It works also as linear guide system that ensures the impact piece moves along a predefined path during impact, minimizing lateral movements and ensuring consistent contact conditions.
• Impact Piece: A rigid component that directly interfaces with the robot during impact.
• Spring Mechanism: A set of springs arranged to provide the desired stiffness, compliant to the ISO standard.
• Central Piece: A structural component that connects the impact piece to the load cell, ensuring proper alignment of springs and force transmission during impacts.
• Load Cell: A high-precision piezoelectric force sensor (PCB Piezotrnics, Model 208C03, capacity 2.224 kN) integrated into the prototype to measure the contact forces during impact events.
• Lower Piece: A structural component that supports the load cell and connects it to the base, ensuring stability during impact tests.
• Base: A stable platform that securely fix the prototype to the ground, ensuring consistent positioning during tests.
• Sticks: Proper sticks, dimensioned with finite element analysis approach, were positioned, in order to protect the prototype in case of exceeding loads or accidents. These pins are designed to fail in bending, triggering a controlled drop of the upper assembly when the total applied force exceeds a predefined threshold, avoiding damages to the internal components of the prototype.

Since the core function of the prototype is to replicate the mechanical response of human body regions, the design process focused on achieving biofidelity in terms of stiffness and damping characteristics. The prototype was engineered to simulate the effective stiffness values specified in ISO/TS 15066:2016 for different body regions, ranging from 10 N/mm to 150 N/mm. The spring constants can be adjusted inserting or removing springs, to reach a configuration of four springs arranged in parallel (three at the vertices of an equilateral triangle, one central), choosing within a set of spring constants of 2.5, 7.5, 12.5, and 37.5 N/mm, to match the effective stiffness $k_{body}$ of the target human body region. The prototype’s modular design allows for easy adjustment of spring configurations to simulate different body regions as needed. Regarding the damping contribute, a 7 mm thick layer of polyurethan with hardness shore A70 was used for the hardest body regions; a 14 mm thick layer of polyurethan with hardness shore A30 was used for the medium ones, and a 21 mm thick layer of silicone with hardness shore A10 was used for the softest body regions.

Given the stiffness contribution of soft tissues using standard compliant damping elements, an additional stiffness component ($k_{damp}$) in series with the effective body stiffness (k) must be considerd, leading to an equivalent stiffness ($k_{eff}$):

$$k_{eff}={\left({\displaystyle \frac{1}{k_{body}}}+{\displaystyle \frac{1}{k_{damp}}}\right)}^{-1}$$

(17)

Is relevant to notice that, using $k_{eff}$ in the energy equation 2, higher permissible relative velocities ($v_{rel}$) are permitted and still safe, potentially increasing productivity (e.g., up to $+63\%$ increase for the lower legs region).

The structural components were realized using fused deposition modeling (FDM). Tough PLA and PETG UHD materials were employed. To ensure strength, particularly in stress-prone parts like the cylindrical guide features, a 100% infill density was applied. The theoretical model assumes that impact piece is infinitely rigid, avoiding to add an additional stiffness element in series.

Cylindrical elements act as guides, ensuring the vertical movement of the impact piece during contact. Diverse lengths of those elements prevent the springs from reaching maximum compression (coil-to-coil contact) prematurely in less rigid configurations. This ensures springs to maintain linear behavior.

3.1.1. Frequency Analysis

A Fast Fourier Transform (FFT) analysis was carried out in design phase to identify the characteristic frequency range of impacts similar to the ones in literature. The importance of this analysis lies in the need to select a load cell with an appropriate frequency response, ensuring accurate force measurements during impact events.

The analysis was carried out on a force profile in [33] and revealed that the significant frequency components of the impact events fall between 2 Hz and 68 Hz (Figure 3). The frequency range identified through the analysis performed in MATLAB environment is consistent with literature [13], which are also referenced by commercial solutions [14]. In this context, it is indicated that signals arising from contact events, even when characterized by very rapid dynamics, typically exhibit frequency components extending up to approximately 100 Hz.

3.2. Pre-Deployment Test Protocol

To validate the prototype before using with a robotic manipulator, a mass drop test protocol was designed. The tests utilized three masses ($M=[21.59,28.30,41.24]g$) and two drop heights ($h_{fall}=[40,70]mm$).

Five force measurements were recorded for each combination of mass, height, and stiffness configuration. The impact velocity $v_M$ of the mass $M_i$ impacting on the prototype was calculated using the equation:

$$v_M=\sqrt{2g{h}_{fall}}$$

(18)

where g is the acceleration due to gravity. The theoretical impact energy $E_{fall}$ was determined by:

$$E_{fall}=Mg{h}_{fall}$$

(19)

The initial mathematical model for dynamic force ($F_{din}$) during contact after drops, assumes that the stiffness $k_{body}$ equals the effective elastic constant of the spring assembly. This model accounts for system masses (M for the falling mass, m for the impact piece) and the velocities post-impact ($V_M$, $V_m$) respectively of the falling mass $M_i$ and of the impact piece, which depend on the coefficient of restitution $c_r$. Infact, using the momentum conservation and the definition of $c_r$, the post-impact velocities can be expressed as:

$$V_m=v_m-M(1+c_r){\displaystyle \frac{v_m-v_M}{m+M}},V_M=v_M+m(1+c_r){\displaystyle \frac{v_m-v_M}{m+M}}$$

(20)

The dynamic force can then be derived using the principle of mechanical energy conservation. The maximum elastic deformation u of the prototype occurs when the sum of kinetic energy of the masses and the potential energy due to gravity is completely stored as elastic energy:

$${\displaystyle \frac{1}{2}}{k}_{body}{u}^2=(m+M)gu+{\displaystyle \frac{1}{2}}\left(m{V}_m^2+M{V}_M^2\right)$$

(21)

Solving this quadratic equation for the force gives:

$$F_{din}=k_{body}u=(m+M)g+\sqrt{{\left[(m+M)g\right]}^2+k_{body}(m{V}_m^2+M{V}_M^2)}$$

(22)

In Figure 4, the test for prototype validation is shown. The hypotheses of this model include assuming the impact piece infinitely rigid, a linear spring behavior, and neglecting air friction and static deflection due to the impact piece weight.

The results from the drop tests were compared to theoretical predictions based on the model in eq. 22. Model performance was quantified using the percentage deviation ($\Delta \%$) between the model value ($F_{din}$) and the corrected experimental value ($F_{exp}$):

$$\Delta \%={\displaystyle \frac{F_{din}-F_{exp}}{F_{exp}}}·100$$

(23)

The ISO/TS 15066:2016 model (with $c_r=0$, representing a completely inelastic collision) severely underestimated the measured forces. Negative deviations ranged from $41\%$ to $55\%$.

The model with $c_r=1$ (perfectly elastic collision) was selected as the most appropriate for estimating forces in HRC application, as it showed smaller percentage differences (always under 15%), especially for tests with larger masses. Based on the $c_r=1$ assumption and the conclusion that the prototype mass (m) is negligible for peak force magnitude, the simplified dynamic force formulation ($F_{din}^{*}$) for mass drop tests (M dropped from height h) was derived:

$$F_{din}^{*}=Mg\left(1+\sqrt{1+{\displaystyle \frac{2k_{body}{h}_{fall}}{M{g}^2}}}\right)$$

(24)

For most configurations using the larger test masses, the deviation ($\Delta \%$) between the simplified model and the corrected experimental data was deemed acceptable, generally falling within an interval of approximately $\pm 15\%$.

The maximum force $F_{max}$ recorded during each impact was analyzed to assess the prototype’s performance. The experimental data closely matched theoretical expectations, confirming the prototype’s suitability for simulating human–robot impacts.

4. Experimental setup and data collection

A comprehensive experimental setup was established to collect high-fidelity data on human–robot impact dynamics. The setup includes the collaborative robot, the impact test prototype, a data acquisition system, and control software to manage the robot’s trajectories and impact conditions. The KUKA LBR iiwa 14 R820, a 7-DOF lightweight manipulator, was selected for its advanced sensing capabilities and compliance features, making it suitable for HRC scenarios. The robot’s joint positions and torques were monitored using its internal sensors. The impact test prototype described in Section 3 was securely mounted to a reticular welded metal structure (compliant with ISO/PAS 5672:2023), ensuring fixed positioning during impact tests.

An additional aluminum disk was designed to be attached to the robot’s end-effector, providing a consistent and rigid impact surface. The disk had a diameter of 31.5 mm and a thickness of 8.35 mm, ensuring uniform contact during impact tests

A National Instruments module NI9230 was connected to the load cell, integrated into the impact prototype, to record force measurements. A custom control application was developed using KUKA’s Sunrise.OS platform. The application allowed for a precise definition of the robot’s trajectories, enabling controlled impact tests at varying velocities and configurations.

Four different configurations were tested, varying the $\alpha$ angle of redundancy of the robot, as shown in Figure 5 and listed in Table 1.

For each configuration, a range of impact velocities was tested by adjusting the robot’s approach speed. The real impact velocity was calculated as the velocity of the end-effector at the instant of contact with the prototype. The set of experiments that results from this combination of configurations and velocities, was then repeated for 5 different body regions (Skull and forehead, Face, Neck, Chest, Back and shuoulders), changing the spring configuration of the prototype and the damping layer to match the effective stiffness $k_{body}$ of each region, as specified in ISO/TS 15066:2016.

The resulting multi-modal dataset was used to systematically identify contact parameters (stiffness and restitution coefficient, as seen in eq. 12), and to validate several data-driven approaches of force prediction. The performance of each model was assessed quantitatively in terms of its ability to accurately replicate measured maximum forces, force–time profiles, across the range of impact velocities and effective masses tested.

The set of experimental points collected for each body region, varying the contact velocity and the effective mass (configuration), is shown in Figure 6.

The effective mass used in this article is computed using the equation 8, considering lumped masses of the robot instead of link’s full inertia contributes and avoiding computation of rotational inertia [27]. In detail, three lumped masses are considered: one at the shoulder, one at the elbow and one at the wrist of the robot. An additional mass is positioned at the end-effector (Figure 1).

Each test, consisted in a controlled movement that can be divided in three phases:

1.

Approach phase: the robot moves forward mantaining the end-effector pose at a specified velocity

2.

Contact phase: in the last part of the approach phase, the end-effector contacts the impact prototype, generating a force peak that is recorded by the load cell.

3.

Retraction phase: after 0.5 seconds from the robot movement stop, still in contact, the robot retracts for 50 mm to the non-contact position.

A trapezoidal velocity profile was used to control the robot during the impact tests (Figure 7), ensuring consistent approach speeds and minimizing acceleration effects at the moment of contact. In this experimental scenario, for simplicity, the human body region is considered static, so the relative velocity $v_{rel}$ is equal to the robot’s velocity at the contact point.

5. Experimental results

5.1. Parameters Identification

Using the Hunt–Crossley formulation in eq. (12) and carrying out a costrained optimization algorithm in MATLAB, it is found out that the damping term of the Hunt–Crossley model is negligible for all body regions.

The optimization algorithm minimizes the following cost function:

$$\mathrm{fun}\left(\mathbf{w}\right)=\sum _{i=1}^{N_{\mathrm{test}}}\left[\sum _k{\left(F_{\mathrm{meas},i}\left(k\right)-w_1{\delta }_i{\left(k\right)}^{w_3}\left(1+{\displaystyle \frac{3(1-w_2){\dot{\delta }}_i\left(k\right)}{2|{\dot{\delta }}_0|}}\right)\right)}^2\right]$$

(25)

where $N_{test}$ is the number of tests performed for each body region, $F_{meas,i}\left(k\right)$ is the experimental force measured at time step k during test i, ${\delta }_i\left(k\right)$ is the penetration depth at time step k during test i, and ${\dot{\delta }}_i\left(k\right)$ is the penetration velocity at time step k during test i, and ${\dot{\delta }}_0$ is the initial penetration velocity k. The parameters to be identified are collected in the vector $\mathbf{w}=[w_1,w_2,w_3]=[k_{HC},c_r,n]$.

The estimated restitution coefficients $c_r$ are greater than $0.94$, leading to an almost elastic impact. Consequently, the model reduces to a Hertz-like formulation with a single dominant parameter $k_{\mathrm{HC}}$ and an exponent n slightly lower than the theoretical value $1.5$ ($0.98<n<1.13$). Since $c_r$ and n were very close to 1 for all body regions, the fallback Hertz model in eq. 13 was used to fit the experiment force profiles.Then, for each test setup in Figure 6, $k_{body}$ was evaluated (Figure 8).

5.2. Force Prediction

The identified parameters were used to predict the force–time profiles and maximum forces for each body region across the range of impact velocities and effective masses tested, using the Hunt–Crossley model.

In order to evaluate the performance of the Hunt–Crossley model against data-driven approaches, Support Vector Machine (SVM) and Decision Tree (DT) algorithms were implemented in MATLAB, using fitrsvm and fitrtree functions respectively. To randomize the training and test datasets, a 5-fold cross validation was performed for both algorithms, as described in [34].

A linear regression model was also implemented in MATLAB using the fitlm function, using as predictors ${\delta }_{max}$, $k_{body}$, $m_{robot}$ and $|{\dot{\delta }}_0|$, obtaining an $R^2=0.92$ and an $RMSE=36.1N$.

The predictive equation in this case is:

$$F_{\max }={\beta }_0+{\beta }_1{\delta }_{max}+{\beta }_2{k}_{\mathrm{body}}+{\beta }_3{m}_{\mathrm{robot}}+{\beta }_4\left|{\dot{\delta }}_0\right|+\epsilon$$

(26)

where ${\beta }_0$ is the intercept, ${\beta }_1,{\beta }_2,{\beta }_3,{\beta }_4$ are the coefficients for each predictor, and $\epsilon$ is the error term.

The magenta coloured line in Figure 9 and Figure 10 represents a +10% displacement from the 0-error line, used as a reference to evaluate the accuracy of the different models. The figures show that, with an overestimation of about 10%, the model can predict in advantage of security the force-time profile and the maximum force measured during the experiments.

The red coloured line represents the threshold of the force limit defined in ISO/TS 15066:2016 in the worst case scenario; in other words, values of force higher than 440 N, are never permitted.

For each type of algorithm and output ($F_{\mathrm{prof}}$ and $F_{\max }$), different combinations of input parameters were tested to identify the most relevant features for accurate force prediction. In the case of SVM, the best model for $F_{\mathrm{prof}}$ was built with the parameters $m_{robot}$, $k_{body}$, $m_{body}$, $\delta$, and shows $R^2=0.985$ and $RMSE=14.4N$, while the best model for $F_{\max }$ was built with the parameters $m_{robot}$, $k_{body}$, $m_{body}$, ${\delta }_{max}$, and shows $R^2=0.98$ and $RMSE=14.4N$.

Regarding Decision Tree, the best model for $F_{\mathrm{prof}}$ was built with the parameters $m_{robot}$, $k_{body}$, $m_{body}$, $\delta$, and shows $R^2=0.989$ and $RMSE=12.1N$, while the best model for $F_{\max }$ was built with the parameters $m_{robot}$, $k_{body}$, ${\delta }_{max}$, and shows $R^2=0.982$ and $RMSE=17.0N$, with only 3 parameters used, demonstrating the effectiveness of Decision Trees in capturing the relationships in the data with fewer features.

Using the decision tree algorithm, another model is built in order to use it during on-line tasks, to check simultaneously the force that could be exerted during an unexpected collision. Using ${\dot{\delta }}_0$ , $q_1$, $q_2$ and $q_4$, it is possible to estimate the maximum force that could be exerted during an impact, calculating the effective mass $m_{robot}$ from the joint positions and velocities. The result obtained is an $R^2$ of 0.965 and an $RMSE$ of 23.43 N, as shown in Figure 11.

6. Conclusions

This work presented an integrated experimental and modeling framework to analyze human-robot impact dynamics for collaborative applications. A dedicated impact prototype, compliant with ISO/TS 15066:2016 and PD ISO/PAS 5672:2023, was designed and validated to reproduce the mechanical behavior of different human body regions with tunable stiffness and damping. The resulting setup enabled the acquisition of a large and well-characterized dataset of impacts with a KUKA LBR iiwa 14 R820, spanning multiple effective masses, impact velocities and body regions.

The analysis of the experimental data showed that, under the tested conditions, collisions are almost elastic, with restitution coefficients close to one and a negligible contribution of the damping term in the Hunt-Crossley formulation. Consequently, a simplified Hertz-like model with a single dominant stiffness parameter was sufficient to capture the main characteristics of the force-penetration relationship, with identified stiffness values in good agreement with the normative ones. At the same time, the classical ISO/TS 15066 energetic model, when used with its standard effective mass definition, significantly overestimated peak forces, whereas the refined effective mass formulation adopted in this work led to more conservative and realistic predictions.

Data-driven models based on Support Vector Machine and Decision Trees further improved predictive performance, achieving high coefficients of determination and low RMSE on both instantaneous forces and peak values. In particular, models using a reduced set of physically meaningful inputs (effective mass, contact velocity and configuration-dependent parameters) proved suitable for on-line implementation, providing accurate estimates of impact forces that can support real-time risk assessment in Power and Force Limiting mode. This hybrid strategy, combining interpretable analytical formulations with flexible machine learning models trained on experimental data, offers a robust toolset for designing and validating safe collaborative robotic systems.

Future work will extend the proposed framework along several directions. On the experimental side, additional robot configurations, contact geometries and motion tasks will be investigated, including scenarios where the human body is not stationary and where multi-contact events may occur. On the modeling side, the integration of more advanced learning techniques, uncertainty quantification and formal safety margins will be explored, with the aim of embedding the proposed predictors directly into the robot controller for continuous monitoring of collision severity. Ultimately, the methodology introduced in this study can contribute to more informed standardization efforts and to the deployment of collaborative robots that are both productive and intrinsically safe.