**1. Introduction**

One of the crucial problems accompanying the automotive sector development are road accidents resulting in human and material loss. Such events involve material loss and a loss of life. In 2018, in Poland, as many as 31,674 road accidents were reported, in which 2862 people were killed and 37,359 were injured. However, in the corresponding period, the total of 436,414 road collisions resulted only in material loss. The events classified as collisions of vehicles in motion, which accounted for 53.8% of all the accidents, were found to be most common. As many as 44.3% of all the people were killed in car crashes, and 58.6% were injured. There were 3104 head-on collisions, resulting in 510 deaths and 4676 people su ffering injuries. In Poland, the death rate for every 100 accidents was found to be the highest in all the European Union member states in 2017 (8.6). The country which came right after Poland is Greece, with the rate of 6.7, whereas Germany had the lowest rate, at 1.1 [1], which suggests that Poland faces a serious road tra ffic safety problem to be urgently addressed.

Next to the human loss, one must also remove the material damage, repair the vehicles damaged in road accidents and collisions. Usually, only accidents involving material losses of road accidents and collisions are considered, whereas insurance companies are involved in compensation of both material and human loss related to the outcomes of accidents delayed in time. As such, they do not focus only on road accidents but also on road collisions and the subsequent costs of repairs. The costs of spare parts, painting, and man-hours for the post-accident vehicle repair are calculated from the car damage analysis.

In Europe, the USA, and Canada one must follow the requirements providing for the impacts of cars with a rigid barrier. Another study [2] discusses the vehicle superstructure designing method to develop the guidelines for the application of bumper beams and impact energy absorbers. The study also covered the head-on collision with rigid barriers. The low-speed collision requirements by the Insurance Institute for Highway Safety (IIHS) and Dunner tests, considered in this study, are demonstrated in Figure 1.

**Figure 1.** Requirements for low-speed collisions by the Insurance Institute for Highway Safety (IIHS) and Dunner tests.

For insurance companies, prior to the paymen<sup>t</sup> of damages for the car to be repaired, it is also very important to determine whether they are liable for damages or not. For that reason, in the process of insurance claim settlement, insurance companies usually use professional assessments provided by experts ordered by courts and public prosecutors to establish liability for the consequences of road accidents. Currently, the assessors use simulation programs to reconstruct road accidents, so manual calculations are hardly used. Adequate simulation results are of key importance as using incorrect simulation results leads to inadequate decisions in terms of liability for damages, which also triggers some legal effects for the parties to a car repair damages court case. The highest calculation accuracy is provided by programs with finite elements method modeling; the FEM convention, e.g., Abaqus, and LS-Dyna. Another study [3] uses LS-Dyna program for researching the vehicle structure stiffening elements during head-on collision. The test involved a Sport Utility Vehicle (SUV), very popular in China. The vehicle Finite Element (FE) model was initially validated against a full-frontal rigid barrier impact test by National Crash Analysis Center (NCAC). The study included also a head-on collision with a rigid barrier. A satisfactory compliance of the test and simulation deformations for each such collision was found (Figure 2).

**Figure 2.** Deformations of the car in test and simulation.

Other authors [4] report on using LS-Dyna program for the vehicle front protection system (VFPS). The study has been launched to address the problem of animals crashing into vehicles in Australia, causing damage to headlights, engine coolant radiators, and engines. The VFPS changes the vehicle crush characteristics during a head-on collision and affects the airbag activation characteristics. An inadequate VFPS structure results in a greater vehicle damage, hence resulting in higher repair costs and potential passenger injuries. Applying the FEM convention modeling program, a structure eliminating such defects has been developed (Figure 3).

**Figure 3.** Over the bumper type front protection system (FPS) mount and its mating parts.

The FEM convention modeling programs are most commonly used in research facilities. However, they have not been commonly applied for the reconstruction of road accidents as geometric and material data must be introduced, which is discussed also by other authors [5,6]. The authors of this study have performed numerical research in the LS-Dyna environment. The results of that research was experimentally validated with crash tests using a full-size road infrastructure support structure dimension. The structure satisfied the Polish PN standard - PN-EN 12767 requirements: "Passive safety of road infrastructure support structures—Study requirements and methods". The study has demonstrated as a good compliance of the simulation with the experiment.

In practice, however, programs for modeling impacts and post-impact movement of vehicles in multi-body system convention (MBS) are used for analysis of vehicle collisions. They are not so precise as FEM programs; however, they offer significantly shorter calculation time and easier operation without the need to enter geometric and material constants. In consequence, modeling with the use of these programs is simplified. As practice shows, the users of these programs, who are commissioned to provide the court with assessment of road accidents, are not aware of it and make mistakes accepting incorrect simulation results as final. Programs for modeling in MBS convention include, for example, V-SIM, PC Crash, and Virtual Crash. Another study [7] reports on the tests of head-on collision with a rigid barrier at 13–95 km/h and full overlap, as compared with the experimental results. V-SIM, however, has not been verified with a full-frontal rigid barrier impact test. In other studies [8,9], contact parameters and other data were identified from crash tests and used in simulation. The articles present both the results of computer simulation of a passenger cars crash with a non-deforming barrier and a pole with V-SIM and the experimental results published by Allgemeiner Deutscher Automobil-Club and AUTOBILD. The results were received for the same initial conditions, which has facilitated evaluating the credibility of computer simulation. The V-SIM program can also be applied to verify damage using the SDC (Static, Dynamic, Characteristic) analysis method to prevent the payments of undue damages for vehicle repairs. It uses the results of geometrical parameters and the appearance of damage, as well as the results of a dynamic collation of the collision. A practical application of this method is described in another study [10], presenting research methods which facilitate an efficient verification of the vehicle collision to determine whether the crash actually occurred or not. It also presents an IT tool automating the decision-making process and supporting the expert's operations, verifying the impact with the SDC method.

The application of this program in Poland is common. Since it is so common to use the V-SIM program for the reconstruction of accidents and for insurance claim settlements, the problem of reliability of the simulation and the effect of the input data on the simulation results is addressed.

### **2. Materials and Methods**

The V-SIM code is designed for the analysis of the impacts of mechanical vehicles (passenger cars, trucks, semi-trailer trucks, and semi-trailers and trailers). It also facilitates the analysis of the impact of vehicles with rigid barriers, such as a wall. The program models the movement of the vehicle, which is treated as a solid with 10 degrees of freedom, moving as a result of the external forces acting on it. The model involved the application of two systems of reference. The first is the global one, which describes the momentary positions of the simulation objects and the position of the movement environment elements designed, e.g., the terrain barriers. The coordinate axes of the reference system are marked *x*, *y*, *z*. The axis is oriented opposite to the force of gravity, and the beginning of the reference system is selected by the program operator. The second system of reference is related to the simulated vehicle, and its axes are marked *<sup>x</sup>'*, *y*', *<sup>z</sup>'*. In that system, the external forces acting on the vehicle are determined. Axis *x'* in that system is directed according to the movement of the vehicle, while axis *z'* is directed vertically upwards the unladen vehicle, and its direction and sense result from the orthogonality and clockwise nature of the system. The positioning of the center of the vehicle mass is determined by radius vector →*r*. (Figure 4).

**Figure 4.** Model of vehicle with co-ordinate systems [11].

Equations of motion for the vehicle in the reference system related to the center of mass of the vehicle have been provided with the following formulas [11]:

$$
\stackrel{\rightarrow}{F}' = \left( \stackrel{\cdots}{r}' + \stackrel{\cdots}{r}' \times \stackrel{\cdots}{\omega}' \right) \cdot m\_r \tag{1}
$$

$$
\stackrel{\rightharpoonup}{M}' = \stackrel{\rightharpoonup}{\Theta\_{\mathfrak{C}} \cdot \stackrel{\rightharpoonup}{\omega^{\prime}}} + \stackrel{\rightharpoonup}{\omega^{\prime}} \times \Theta\_{\mathfrak{C}} \cdot \stackrel{\rightharpoonup}{\omega^{\prime}},\tag{2}
$$

where:

→

→ *r* —radius vector,

ω—rotational speed of the vehicle in the system connected to the vehicle, and

Θ*c*—tensor of the mass moment of inertia.

For the rotational movement, there has been assumed the tensor of mass moment of inertia with zero mass moments of inertia, besides the main moments, which has resulted in the system of scalar equations describing the rotational movement, presented below [11]:

$$
\Theta\_{\mathcal{E}} = \begin{bmatrix} I\_{\mathbf{x'}} & 0 & 0 \\ 0 & I\_{\mathbf{y'}} & 0 \\ 0 & 0 & I\_{\mathbf{z'}} \end{bmatrix} \tag{3}
$$

$$I\_{\mathbf{x'}} \cdot \dot{\boldsymbol{\omega}}\_{\mathbf{x'}} = M\_{\mathbf{x'}} - I\_{\mathbf{z'}} \cdot \boldsymbol{\omega}\_{\mathbf{y'}} \cdot \boldsymbol{\omega}\_{\mathbf{z'}} + I\_{\mathbf{y'}} \cdot \boldsymbol{\omega}\_{\mathbf{y'}} \cdot \boldsymbol{\omega}\_{\mathbf{z'}},\tag{4}$$

$$I\_{\mathbf{y'}} \cdot \dot{\boldsymbol{\omega}}\_{\mathbf{y'}} = M\_{\mathbf{y'}} + I\_{\mathbf{z'}} \cdot \boldsymbol{\omega}\_{\mathbf{x'}} \cdot \boldsymbol{\omega}\_{\mathbf{z'}} - I\_{\mathbf{x'}} \cdot \boldsymbol{\omega}\_{\mathbf{x'}} \cdot \boldsymbol{\omega}\_{\mathbf{z'}} \tag{5}$$

$$I\_{\overline{z}'} \cdot \dot{\omega}\_{\overline{z}'} = M\_{\overline{z}'} - I\_{\overline{y}'} \cdot \omega\_{\overline{x}'} \cdot \omega\_{\overline{z}'} + I\_{\overline{x}'} \cdot \omega\_{\overline{x}'} \cdot \omega\_{\overline{y}'} \tag{6}$$

where:

*I x yz*—the main moments of inertia of the vehicle,

*M x yz*—components of the moment of external forces acting on the vehicle in the system connected to the vehicle, and

ω *x yz*—components of the rotational speed velocity of the vehicle against the axes selected.

The total force acting on the vehicle is the force of gravity, the forces of aerodynamics resistances, and the forces brought by the suspension of the vehicle wheels, which has been described with the following formula [11]:

$$
\overrightarrow{\overline{F}}' = \overrightarrow{\overline{F}}'\_{\mathcal{S}} + \sum \overrightarrow{F}'\_{i} + \overrightarrow{\overline{F}}'\_{ax} + \overrightarrow{\overline{F}}'\_{ay}.\tag{7}
$$

The gravitational force has been determined by transforming from the global system of reference according to the following formula [11]:

$$
\overrightarrow{F}\_{\mathcal{S}} = \begin{bmatrix} 0 \\ 0 \\ -m \cdot \mathcal{g} \end{bmatrix} \xrightarrow{\overrightarrow{F'}\_{\mathcal{S}}} \mathbf{x'}^{\overrightarrow{F}\_{\mathcal{S}}},\tag{8}
$$

where:

m—mass of the vehicle, and g—standardaccelerationduetogravity.

The forces of aerodynamic resistance are determined irrespective of the frontal and side area of the vehicle from the speed of the vehicle, while considering wind speed →*w* according to the formulas provided below [11]:

$$
\stackrel{\rightarrow}{F}'\_{\text{ax}} = \begin{bmatrix}
\pm \rho/2 \cdot \mathbb{C}\_x \cdot A\_x \cdot \left(\dot{r}\_{\text{x}'} - w\_{\text{x}'}\right)^2 \\
0 \\
0
\end{bmatrix} \tag{9}
$$

$$
\overrightarrow{F}\_{ay}^{\prime} = \begin{bmatrix} 0 \\ \pm \rho/2 \cdot \mathbb{C}\_{\mathcal{Y}} \cdot A\_{\mathcal{Y}} \cdot \left( \dot{r}\_{\mathcal{Y}} - w\_{\mathcal{Y}} \right)^2 \\ 0 \end{bmatrix} \tag{10}
$$

where:

ρ—air density,

Cx, y—coefficients of the frontal and side air resistance,

Ax, y—frontal and side areas of the vehicle,

. *rxy*—linear velocities of the vehicle in the system connected with the vehicle, and

wx, y—componen<sup>t</sup> wind speeds in the system connected with vehicle.

The vehicle suspension model also considers the independent vertical movement for each wheel of the vehicle. From the positioning of the vehicle and the local configuration of the road surface, deflection of the suspension Si is determined, whereas the normal force of the response of the suspension considering the elasticity with progressive characteristics and the damping values is calculated from the momentary deflection and the velocity of its changes according to the formula [11]:

$$F'\_{zi}(\mathbb{S}\_i, \dot{\mathbb{S}}\_i) = \max\{\max\{\mathbf{C}\_{3i} \cdot \mathbf{S}\_i^3 + \mathbf{C}\_{1i} \cdot \mathbf{S}\_i + F'\_{0i}, 0\} + \begin{cases} D\_{ii} \cdot \dot{\mathbb{S}}\_i^{\dot{\mathbb{S}}\_i \ge 0} \\ D\_{ii} \cdot \dot{\mathbb{S}}\_i^{\dot{\mathbb{S}}\_i \ge 0} \end{cases},\tag{11}$$

where:

*<sup>F</sup>zi*—vertical component of the response of the vehicle suspension, *Si*—reduced value of deflection of the vehicle wheel suspension, *C*3*i*—reduced coefficient of progression of stiffness of the vehicle wheel suspension,*C*1*i*—reduced coefficient of stiffness of the vehicle wheel suspension, *<sup>F</sup>*0*i*—normal force of the response of the vehicle wheel suspension, *Dci*—reduced damping coefficient for the compression of the vehicle wheel suspension,*Dri*—reduceddampingcoefficientforthestretchingofthevehiclewheelsuspension.

The program for determining the interaction of the tires with the road surface used a nonlinear model developed at the Highway Safety Research Institute (HSRI), University of Michigan, by the Dugoff team and TM-Easy. The forces of the tire reactions are calculated in the system connected with the tangent point of the tire with the road surface, and the normal force of the tire is assumed as the normal force of the reaction of the suspension. The normal force of the tire, the rotational speed of the wheel, and the parameters of the road surface are used to determine the forces of the reaction of the tire in the other axes of the system. For each wheel of the vehicle, additional degree of freedom has also been considered, and rotational movement has been considered. The model considers the driving torque and the movement resistances according to the formula [11]:

$$I\_i \cdot \dot{\omega}\_i = M\_{\rm ii} - \prescript{\prime}{F\_{\rm xi}^{\prime\prime}}{}{R\_{\rm di}} \pm (M\_{\rm hi} + M\_{\rm ci} + M\_{\rm ti} + M\_{\rm ri})},\tag{12}$$

 and

where:

*Ii*—mass moment of inertia of the vehicle wheels,

. <sup>ω</sup>*i*—rotational speed of the vehicle wheels,

*Mni*—reduced driving torque of the vehicle wheels,

*<sup>F</sup>xi*—componen<sup>t</sup> of the force of the response of the tire in the system connected with the tangent point of the vehicle tires with the road surface,

*Rdi*—rolling radius of the vehicle wheels,

*Mhi*—break torque of the vehicle wheel brake,

*Mei*—reduced engine brake torque which occurs on the vehicle wheels,

*Mti*—road surface rolling resistance torque which occurs on the vehicle wheels, and

*Mri*—vehicle-wheels own resistances torque.

The model of the brake system of the vehicle considers the corrector of the braking force of the rear axle for each load of the vehicle, and, optionally, it can consider also the effect of the anti-blocking system (ABS).

However, the driving torque of the engine is determined from its external characteristics from the following formula [11]:

$$M\_c(\omega\_c) \;=\; M\_m - \frac{M\_m - \frac{N\_n}{\omega\_n}}{\left(\omega\_n - \omega\_m\right)^2} \cdot \left(\omega\_c - \omega\_m\right)^2\tag{13}$$

where: *Me*—driving torque of the vehicle engine,

<sup>ω</sup>*e*—rotational speed of the vehicle engine,

*Mm*—maximum torque of the vehicle engine,

*Nn*—maximum power of the vehicle engine,

<sup>ω</sup>*n*—speed for the maximum power of the vehicle engine, and

<sup>ω</sup>*m*—rotational speed for the maximum torque of the vehicle engine.

The operation of the steering system of the vehicle was analyzed using the kinematic model following the Ackermann steering principles, considering the susceptibility of the real-life system and the correction of shearing forces of the reaction of the tires of the wheels of the steered axle. The program can also model additional tasks for the vehicle and for the driver; the vehicle acceleration and a turn with the steering wheel, braking with the foundation and secondary brakes, using the clutch and blocking the wheel, a pressure drop in the tire, and the deformation of suspension.

A force-based model has been developed for an impact with the impact forces developing in a constant manner while two objects are in contact [12]. The impact of vehicles is detected with a 2D or 3D model developed for collision detection, depending on the choice of the operator. This choice, however, should be informed and adequate to the collision analyzed. It mostly applies to the impacts of vehicles the geometric compatibility of which is inconsistent, which has also been discussed in other studies [13]. Solving problems of collision modeling is also discussed in Reference [14–16].

An impact force model of the V-SIM program suggests a development of force *F* in time during the impact compression, and then its collapsing, according to a certain function in the restitution phase. Such an approach assumes that all the external and mass forces, provided by the vehicle model and by the model of a wheel, occur during an impact [17,18]. Detailed methods for a vehicle impact modeling with a rigid barrier are used in another study and presented below.

Such modeling considers changes of force which appear during an impact and its influence on the vehicle movement, as well as the occurrence of two impact phases: the first is a compression phase (when the value of momentary *F* force increases), whereas the second one is a restitution phase (when the value of momentary *F* force decreases), and according to Newton's hypothesis, the above phases are associated with restitution coe fficient *k*.

On the basis of the concept of force impulse for the compression and restitution phases, formulas for both phases are used [19]:

$$\stackrel{\rightarrow}{S}\_k = \int\_0^{t\_k} \stackrel{\rightarrow}{F}(t)dt,\tag{14}$$

where:

→ *Sk*—represents impulse of force for the compression phase, ·*t*0—represents the compression phase beginning time, ·*tk*—represents the compression phase end time, and →

*<sup>F</sup>*(*t*)—represents a momentary force of impact.

$$\stackrel{\rightarrow}{S}\_{\mathcal{F}} = \int\_{t\_k}^{t\_r} \stackrel{\rightarrow}{F}(t)dt,\tag{15}$$

where:

→

 *Sr*—stands for impulse of the restitution phase, and

*tr*—stands for restitution phase end time.

The total impulse of the impact force is presented as the sum of force impulses for the phases of compression and restitution, expressed using the formula [19]:

$$
\overrightarrow{S} = \overrightarrow{S}\_k + \overrightarrow{S}\_r.\tag{16}
$$

In turn, joining the values of both impulses of the force by applying the restitution coefficient is described by the formula [19]:

$$S\_r = k \cdot S\_k.\tag{17}$$

The applied impact model assumes that, during the impact, the value of temporary impact force *F* is proportional to the volume of overlapping contours of the colliding simulation objects (Figure 5a) by the below formula [19]:

$$F(t) := \mathcal{c}f(t),\tag{18}$$

where:

*c*—stands for stiffness coefficient of the impact area [N/m3], and *f*(*t*)—stands for the impact area volume [m3].

After linearization of the course of action of forces that occur during the two earlier indicated phases, the impact force is presented with a diagram in Figure 5b. Modeling also includes the mean value of the force in the states of compression and restitution, expressing the coefficient of restitution by the following formula [19]:

$$k \cdot F\_{k \text{mid } t\_k} = F\_{\text{void}} \cdot t\_{\sigma} \tag{19}$$

where:

*Fkmid*—stands for the mean value of force *F* in the compression phase,

*Fomid*—stands for the mean value of impact force *F* in the restitution phase,

*to*—stands for the restitution phase duration time, and

*tk*—stands for the compression phase duration time.

**Figure 5.** Impact area of vehicles (**<sup>a</sup>**,**b**) impact force *F* in the compression and restitution phase: *Sk*—deformation depth in the compression phase, *So*—depth of the elastic disappearing deformation in the restitution phase, and *Fm*—maximum force value.

For the uniformly variable motion, the deformation depth has been described with the following formula [19]:

$$S\_{\mathbf{x}} = \frac{a \cdot t\_{\mathbf{x}}}{2},\tag{20}$$

where:

*Sx*—deformation depth, and

*a*—average deceleration. The course of action of impact force *F* during restitution is written by formula [19]:

$$F(t) := F\_m \left( 1 - \frac{1}{k^2} \left[ 1 - \frac{s(t)}{s\_k} \right] \right). \tag{21}$$

For k = 0 → F(t) = 0, where:

*Fm*—stands for maximal value of the impact force, and

*sk*—stands for depth of deformation at the end of the compression phase.

It has been assumed that having been hit by a vehicle, the wall does not undergo deformation. Following that approach, the vector of momentary impact force → *F*(*t*)is on the margin of a non-deforming barrier, in point *C*, found in the middle of the barrier length, in the area (Figure 6a), whereas, the position of this point for the vehicle is defined by radius vector.

Momentary impact force *F* is represented by the following components: tangent and normal to the surface of the non-deforming area. In its application point, the direction of normal component is perpendicular to the surface on the non-deformable area. The unit vector which defines the direction of this component is marked (Figure 6b) → *e n*1 , and it is defined as follows [19]:

$$
\overrightarrow{F}\_{m\_1}(t) = F\_n(t)\overrightarrow{\mathcal{e}}\_{m\_1}.\tag{22}
$$

where:

> *Fn*(*t*) stands for the normal component value.

**Figure 6.** Diagram depicting the calculation method for an impact with a rigid barrier (**a**) and momentary force of impact with a rigid barrier (**b**).

In this model, friction force occurring in the place of impact is represented by tangential component → *FT*(*t*). Its value must be such that the direction of the speed di fference is consistent with the direction of impact force resultant → *<sup>F</sup>*1(*t*). In Figure 6b, the angle between the direction of impact force component and the speed di fference direction is marked as δ. It was also accepted that the value of tangent component cannot exceed the product of normal component and friction coe fficient μ.

The value of tangent component is expressed with the formulas [19]:

$$F\_{T\_1} = F\_{n\_1} \tan \delta \leftrightarrow F\_{n\_1} \tan \delta \le F\_{n\_1} \mu\_\prime \tag{23}$$

$$F\_{T\_1} = F\_{n\_1} \mu \leftrightarrow F\_{n\_1} \tan \delta > F\_{n\_1} \mu. \tag{24}$$

The compression phase finishes when the impact area *f (t)* ceases to increase and the restitution phase begins. The value of normal component decreases from the maximal one down to zero, with the assumption of the above linear change in momentary impact force *F* to be a function of time and restitution coe fficient *k*.
