# Model-Based Optimization of Velocity Strategy for Lightweight Electric Racing Cars

Authors: Mirosław Targosz; Wojciech Skarka; Piotr Przystałka
Journal: Journal of Advanced Transportation
 (2018)

---

## Abstract

The article presents a method for optimizing driving strategies aimed at minimizing energy consumption while driving. The method was developed for the needs of an electric powered racing vehicle built for the purposes of the Shell Eco-marathon (SEM), the most famous and largest race of energy efficient vehicles. Model-based optimization was used to determine the driving strategy. The numerical model was elaborated inSimulinkenvironment, which includes both the electric vehicle model and the environment, i.e., the race track as well as the vehicle environment and the atmospheric conditions. The vehicle model itself includes vehicle dynamic model, numerical model describing issues concerning resistance of rolling tire, resistance of the propulsion system, aerodynamic phenomena, model of the electric motor, and control system. For the purpose of identifying design and functional features of individual subassemblies and components, numerical and stand tests were carried out. The model itself was tested on the research tracks to tune the model and determine the calculation parameters. The evolutionary algorithms, which are available in theMATLAB Global Optimization Toolbox, were used for optimization. In the race conditions, the model was verified during SEM races in Rotterdam where the race vehicle scored the result consistent with the results of simulation calculations. In the following years, the experience gathered by the team gave us the vice Championship in the SEM 2016 in London.


---

## Body

## 1. Introduction

The issue of energy saving in road transport is getting more and more important. It is especially significant in the context of electric drives in vehicles becoming widespread. Reduction of energy consumption together with the development of new energy sources of greater capacity is seen as a method of overcoming the main limitation of electric vehicles which is range.Shell Eco-marathon car race is the testing ground for new solutions in designing energy saving vehicles. A team of scientists and students from the Silesian University of Technology has been taking part in it since 2012, in energy saving vehicles with electric drives, which are designed and built by them in the following categories: Prototype, UrbanConcept with Battery Electric, and UrbanConcept with Hydrogen Fuel Cell Stack. Significant reduction in energy consumption is achieved by the usage of two methods [1]. The first one includes design changes which result in lower energy consumption whereas the latter one introduces new strategies of driving and drive steering, which allows minimalizing the energy consumption on a given route. By means of this method, it is possible to reduce considerably the energy consumption in vehicles and in particular in electric vehicles. Currently, electric vehicles (city cars: segment B) can reach the result of about 5km/kWh of energy, whereas respective vehicles of UrbanConcept Battery Electric category reach up to 200 km/kWh in simulated conditions of driving in a city. However, prototype battery electric vehicles which are different type as far as their structure is concerned reach the result of 1000 km/kWh. It shows the potential of the design solutions and the strategy of driving. In both categories numerical simulation vehicle models are used, which allow determining the direction of designing changes as well as planning the proper strategy of driving. These models are also used for optimization based on Model-Based Optimization methodology (MBO). This proprietary methodology used for vehicles development and planning the driving strategy has been described in this paper. Vehicles modelling and in particular electric vehicles modelling are widely known and become of an increased interest to research, in particular with hybrid drive or power source. The necessity to model phenomena and objects is commonly known and creates substantial branch of science and engineering activities. The research and analysis of mechanical systems, technological processes, and other phenomena in real world are possible thanks to the use of proper mathematical apparatus which is based on previously assumed mathematical models of these processes [2]. According to this, by means of mathematical model, one should understand mathematical description of an object, process, or phenomenon which it represents.Compared to cars with internal combustion engines, electric cars have a relatively simpler drive system; however, in order to achieve a significant reduction in energy consumption in the electric car's drive system, additional systems such as the energy recovery system, additional mechatronic subsystems, and advanced control systems are integrated, resulting in a continuous increase in the complexity of the electric drive systems currently used. The complexity of the technical means influences the complexity of the mathematical model of the vehicle. Among the models describing the dynamics of vehicles, which are used in the analysis and optimization of energy consumption, two approaches are generally observed due to the type of model used and the following models are distinguished:(i)
Analytical ones in which the model response to the signal is recorded as a motion equation(ii)
Models in the black box concept where the model is based on experimental dataThe created vehicle simulation model can be a combination of two model types. To some extent, the mathematical description of a given phenomenon is based on known dependencies, while another fragment of the model is created in the concept of a black box. Taking into account the direction of information processing in the model, the following models [3] are distinguished:(i)
Simple: where the calculations made by the model start with the engine, energy is transferred to the wheels of the vehicle, and vehicle behaviour can be analyzed(ii)
Reverse: in which the behaviour of the body of the vehicle is modelled, e.g., the speed, and on this basis the required torque and the speed of rotation of the motor shaft are determined.The choice of model type is mainly determined by the need for which the model is to be used. For modelling of dynamics of electric vehicles, many computer programs are used; among them noteworthy is theADVISOR (Advanced Vehicle Simulator) platform [4–6] written as a program in theMATLAB-Simulink package. The program in its libraries gives you the opportunity to analyze many of the solutions of the power system, drive system, or vehicle body. Another platform running in theMATLAB-Simulink environment isPSAT (The Powertrain System Analysis Toolkit), which simulates many predefined solutions for conventional, electric, and other vehicles [7]. The evolution ofPSAT isAutonomy [8] available onLMS Imagine.Lab. In addition to these simulation environments, where models are abstract and which take into account the mechatronic nature of the electric vehicle's propulsion system, there are many specialized dynamic analysis programs utilizing the multibody formalism, such as theLS-Dyna,Adams, orMotionSolve software. Simulation models designed to analyze and optimize energy consumption should allow for less time consuming calculations. This is especially important in optimization tasks, where simulation is performed many times (often over several tens of thousands of times). You can say that the faster the model is, the more abstract it is. In paper [9], this is demonstrated by the synchronous motor, written inMATLAB andSimplorer. In order to reduce the computational time, the model should be written in a low-level language, e.g.,C, as given in [10] and the simulation time can be reduced by up to 20 times compared to a program written inMATLAB's own language.


## 2. Research Methodology

In order to optimize the control strategy for a light, electric racing vehicle, the general methodology described below is presented, which explains the basic ideas and above all the definitions of the key concepts from energy consumption domain. The general functional structure of the drive system with particular attention to drive control and the optimization scheme are also presented below.Energy Consumption E of the Movement. Movement of the vehicle is a consequence of the longitudinal force, which overcoming the force of inertia and resistance against the movement performs work on a particular road section [12]. The energy that is associated with this work is called energy intensity and, in the case of wheeled vehicles, it can be expressed as (1)E=∫0DcFNsdswhereFNisdriving force andDC is travelled distance.Energy consumption of the movement determines the amount of energy supplied to the drive wheels anddoes not depend on the nature of the drive unit and the transmission.Total Energy Consumption E c. It determines total energy expenditure. In the case of an internal combustion engine unit expressed as the product of the fuel and its calorific value, or in the case of an electrical unit, the energy is taken from an electric source.Total energy consumption Ec is the sum of the energy consumptionE and the power losses in the motorΔEm and in the transmission systemΔETS [12]:(2)EC=E+ΔEm+ΔETSThe losses from (2) are usually taken into account inefficiencyη:(3)η=EECThe overall efficiency of the drive system is the product of the motor efficiencyηm and the efficiency of the transmission ηTS:(4)η=ηmηTSThe total energy consumptionEc, of the process, which is defined as driving of the given distanceDC, is influenced by many factors, the vehicle's structural parameters, the energy losses generated when generating and transmitting the driving torque, and the control method. A function describing the speed of a vehicle along a certain route is the so-called speed profile [12, 13]. The speed profile is usually determined by the driver of the vehicle, by means of appropriate accelerator pedal control and shift gear selection or by electronic control. A cruise control system controls the power delivered to the drive wheels in order to maintain a certain speed.If the vehicle is equipped with one drive system, it can be presented in the functional diagram as shown in Figure1.Figure 1
Functional diagram of the electric vehicle drive system.Specific blocks describe the performance of the individual vehicle systems:(i)
master control system: a speed profile generator, usually the driver, by adjusting the accelerator pedal, deciding on the profile, and depending on the road conditions;(ii)
motor control system: generator of appropriate motor control values (current, voltage);(iii)
motor: a machine that transforms the energy of a source into mechanical energy;(iv)
mechanical system of the vehicle: torque transmission from the motor shaft to the wheels of the vehicle (clutches, gears, etc.);(v)
measurement system: acting as feedback to provide information on the current position and speed of the vehicle.The optimization task can be reduced to the following:(i)
Optimization of design parameterspk and settingspn, e.g., power of the drive, TS ratio, i.e., optimizing the characteristics of the vehicle describing themechanical system andcurrent controller. The design of the system and its mathematical description in the form of the model are decisive for the optimization of these parameters.(ii)
Optimizing the control strategy, i.e., searching for the optimum speedVSET profile, or directing the control signal to the current controller USET.(iii)
Simultaneous optimization of control parameters and strategy.You can use the Model-Based Optimization method to solve the optimization task. The scheme of the optimization is shown in Figure2.Figure 2
Schematic diagram of the optimization using model [11].In this method, decision variables, whose value is determined by the optimization algorithm, are input variables for the simulation model. Numerical simulation of driving and obtained results during the simulation allows determining the result of the objective function. Repeating the simulation results for the subsequent driving parameters controlled by the optimization algorithm allows the calculation of successive results of the objective function. Frequently repeated results are aimed at finding the minimum of the target function and thus the optimal solution. Stopping the optimization algorithm loop depends on the type of the algorithm used. This usually occurs after a predetermined number of loops, after a lack of change in the value of the criterion function or other constraints.


## 3. Building the Model of the Electric Vehicle

To study the total energy consumptionEc of an electric vehicle, it is necessary to build a model to analyze the transfer of energy from the source of energy to the drive wheels of the vehicle. Figure 3 shows the functional diagram of the electric vehicle drive system. This diagram is the starting point for creating a simulation model.Figure 3
Functional diagram of an electric vehicle.According to the definition, the total energy consumption depends on the energy consumption of the movement E, that is, on the work generated by the driving force and the energy losses generated by the motorΔEm as well as on the energy loss of the transmission ΔETS. Accordingly, a mathematical model of the vehicle should be created describing the phenomenon of resistance and loss of energy in transmission system. The mathematical model of the electric vehicle is shown in this section. The model is based on equations of motion of a wheeled vehicle, biaxial with front axle steering, and rear wheel drive. The focus is on the representation of the driving forces and the resistance of movement acting in the direction of the longitudinal axis of the vehicle. The model described here is sufficient to simulate and optimize the parameters of an object that affect the total energy consumptionEc of an object. There was no need to determine the values of forces and displacements affecting individual vehicle components, as in the case of, for example, strength analyses. These issues are not the subject of this study, so the forces acting in the directions of the other axes of the vehicle are not dealt with, except for the forces which, as a result of their impact, affect the total energy consumption.


### 3.1. Mathematical Model of the Vehicle

The theory of vehicle motion is widely published and has been presented many times in the scientific works [14, 15]. The following section provides formulas and symbols used to build the vehicle model. A detailed description of the presented issues should be sought, among others, in the above-mentioned bibliography. Figure 4 shows the forces acting on the vehicle moving along the path of the inclination angle α.Figure 4
The system of forces acting on the vehicle.The following designations were used:FN driving force, which is the quotient of the driving torque and the dynamic radius of the wheel; FT1 and FT2 of the rolling resistance of the wheels of front and rear axles, respectively;Qm vehicle gravity force;FP force of air resistance;FW = Qmsin⁡(α) uphill resistance;FB inertia force;Z1 andZ2 ground reactions to the front and rear axle wheels;mz =mδz reduced mass.The equation of motion of the vehicle in the driving direction can be written as the sum of the driving forces and drag forces parallel to the road surface as(5)mZx¨=FN-FT-FW-FPThe following simplifications were adopted. The drive torque applied to the rear wheel does not affect the pressure distributionZ1 andZ2. In fact, in vehicles with rear wheel drive, the rear is downforced whereas the front axle is relieved. The reaction valuesZ1 andZ2 were statistically determined. The coefficients of rolling resistance of all wheels are identical.Resistance forces depend on many factors and parameters; many of them are highly nonlinear. Modelling of wheel rolling resistance is the subject of separate research and the model used is very complex so in practice, a simplified model is often used.The rolling resistance of a single wheel is as follows:(6)FTi=crrZic rr is the coefficient of rolling resistance that depends on the type of road surface, the type of tire, or the pressure in the tire. This factor also depends on the speed of travel and increases with its growth; this is related to the increase in tire deformation. The value of the rolling resistance coefficient is usually given by the tire manufacturer, and exemplary values and empirical models on the dependence of this coefficient on the speed function can be found in the literature [16]. For low travel speeds (<40km/h), the value of the rolling resistance coefficient is assumed to be constant, independent of speed. In linear motion, the total rolling resistance is the sum of the rolling resistance of all vehicle wheels.If an identical rolling resistance coefficient for each wheel is assumed, which is often justified, given the similar construction of the front and rear axle wheels and identical tires and the same pressure operation, the sum of the rolling resistance shall be written as(7)FT=crrQmcos⁡αThe total wheel resistance should take into account the tire rolling resistance, bearing losses, and so forth.The essence of resistance in bearings is identical to the tire rolling resistance described previously. The bearing torque is proportional to the bearing load, friction coefficient, and bearing radius. Bearing friction is generally neglected, as the bearing resistance value is small compared to the tire rolling resistance [16].Similarly to the value of rolling resistance, the turning resistance can be expressed as the product of the weight of the car and the coefficient of turning, although, according to literature, at low speeds it plays negligibly smaller role than the rolling resistance [16].(8)Fs=mg⁡fzf s is the coefficient of turning and depends on the pressure force on the axles, tires, and so-calledintensity of turning, in other words, ratio of centrifugal and vertical acceleration v2/Rg.The value of air resistanceFP can be determined from(9)FP=γvsum22Acxwhere γ is the density of air.vsum [m/s] is the total air velocity resulting from the linear velocity of the vehiclevp and wind speedvw and is the vector sum,A is the vehicle front area [m2], andcx is the aerodynamic drag coefficient.In this study, a TS based on a constant gear ratio was used. PM BLDC (Permanent Magnet Brushless Direct Current) motors are very often used in electric power systems [17]. The mathematical description of all phenomena related to magnetism, electrical, thermal, or mechanical phenomena can lead to a considerable complexity of the model, so simplified models are often used. The electrical part of motor’s model can be found in work [9]. The electromagnetic torque generated by the motor is proportional to the flowing current and the so-calledthe motor torque constant Km and is given by [17](10)Me=2KmizIn steady state, the electromagnetic torqueMe is equal to the torque on the motor shaft, which is the sum of the load torque and the torque of the losses. The efficiency of an electric motor is defined as(11)ηm=PmechmPelectP m e c h m is the mechanical output power on the motor shaft andPelect is the electrical power input. Relations describing the phenomenon of loss can be found in professional literature [10].The motor torque is transmitted to the drive wheel via the TS transmission. The mechanical power on the drive wheel depends on the efficiency of the TS:(12)Pmechp=ηTSPmechm=ηTSηmPelectThe efficiencyηTS is a function of many factors; i.e., it depends on the powertrain design, the transmitted power, and quality of the system. For TS based on tension transmission, efficiency can be approx. 90%. The efficiency of the tie link transmission is a function of inter alia transmission ratioiTS, interaxial distance, belt angle, tension force, rotational speed, or torque [18–20]. In the case of lightweight vehicles in the TS design, it is often not necessary to use variable gear transmissions. Speed and torque control is achieved by the control system. On the one hand, it is an advantage because the TS is smaller in weight and less complicated, but on the other hand it forces the electric drive to operate in a wide range of torque and rotational speed, which affects the motor energy efficiency.


### 3.2. Building the Model of the Driving Route

The parameters of the path that the vehicle travels on have a significant effect on the resistance values shown in Section3.1. In the paper, the path of ride is determined by the exact path of the vehicle on a designated route. It is assumed that the best path is known and it is determined based on the team’s experience (it means that the quickest best path along the track is not calculated by the optimization method). The difference between the path and the track is as follows. A track is a system of physically existing roads between points A and B, and when the path does not exist physically, it is the trace that leaves the selected point of the vehicle (e.g., center of gravity) on the route surface. With this in mind, the journey along the route can take place on a number of paths. The route is only a limitation.Figure5 shows a schematic drawing of one path and its elevation profile. The description of the path is formed by its discretization, i.e., the division into finite numbers of elementary sections of Oti, i∈{1,2,...,n}.Each elementary segment is described by a parameter vector.Figure 5
Driving path and its inclination profile.O t i = [ L t i ,  R t i ,  S t i ,  E t i ], where Lti[m] is the length of the path, Rti[m] is the radius of curvature, Sti[%] is the inclination of the path, and Eti represents the type and condition of the path. The length of the path is the sum of the lengths of the Lti sections.These basic data are essential for the correct determination of motion resistance on each segment. It is assumed that any change in the inclination or radius of curvature makes it necessary to define further vectorOti.Bearing in mind that vehicles usually move in the open air, a suitable model should be created describing the external conditions. For this purpose, theOti vector describing the nature of the path should include additional columns such as velocity and direction of the wind vwi and vki, atmospheric pressure patmi, and air temperature tatmi. In this situation, the vector describing the path and external conditions takes the form: Oti=[Lti,Rti,Sti,Eti,vwi,vki,patmi,tatmi]. The speed and direction of the wind have a significant effect on the value of air resistance. Determining the wind speed is relatively difficult and not precise. Temporary winds may have different directions and a much higher amplitude than the average. Preliminary investigations have been carried out considering the uncertainty of weather conditions. During optimization, where optimal control strategies were sought, the wind speed was randomly changed. The concept of optimizing the control strategy from taking into account the changing weather conditions was to design a family of strategies. It would be possible to select one optimal strategy or part of it in the real world because of the existing and changing weather conditions.


### 3.3. Control System Design

Two modes of the control system are taken into consideration in this study. They reflect the physical implementations of the control system used in the built-in electric race vehicles of Prototype class [21] and UrbanConcept class [22]. Figures 6 and 7 represent block diagrams of the drive system of the vehicle for different control strategies. The control scheme A (Figure 6) is used in a case when a driver is able to control an acceleration of the vehicle by means of a binary signal. In this case, the goal of the optimization task is to find the control strategy taking into account two criteria, i.e., the determined travel time as well as the minimal power consumption of the vehicle. The control strategy is used to specify in which parts of the route the driver should press or release the acceleration treadle. When the setpoint signal USET∈{0,1} is assumed, then the optimization task is to search for logic values of the control signal as the function corresponding to the position of the vehiclex(t), USET =f(x(t)). The supervisory role of the driver is included in the functional scheme as a decision module. In this study, the behaviour of the driver is not modelled in the simulation. It means that he/she realises the strategy in the exact way. The input control signalUr is always different from the signal USET in real conditions. It means that the signal USET is mainly employed by the driver to know when the acceleration treadle should be used.Figure 6
Drive system diagram: control scheme A.Figure 7
Drive system diagram: control scheme B.In the case of the control scheme B (Figure7), the general goal of the optimization task is the same as in the previous case. However, the formal definition of the problem is quite different because in this case the optimal form of the velocity signal is searched according to the position of the vehicle on the routeVSET =f(x(t)). The velocity controller compares the setpoint signalVSET with the measured velocity of the vehicle and for certain controller setting as well as for control errorev it generates the control signal USET. As in the previous case, the supervisory role of the driver is included in the decision module. However, the main difference is that in the optimization process the setpoint signalVSET can be limited in the direct way.According to control schemes A and B, the input control signalUr is used by the motor control unit where the emended current regulator normally outputs the corresponding signals to the motor terminals. The driving torque of the motor shaft is transmitted through the transmission system to the wheels of the vehicle.


## 3.1. Mathematical Model of the Vehicle

The theory of vehicle motion is widely published and has been presented many times in the scientific works [14, 15]. The following section provides formulas and symbols used to build the vehicle model. A detailed description of the presented issues should be sought, among others, in the above-mentioned bibliography. Figure 4 shows the forces acting on the vehicle moving along the path of the inclination angle α.Figure 4
The system of forces acting on the vehicle.The following designations were used:FN driving force, which is the quotient of the driving torque and the dynamic radius of the wheel; FT1 and FT2 of the rolling resistance of the wheels of front and rear axles, respectively;Qm vehicle gravity force;FP force of air resistance;FW = Qmsin⁡(α) uphill resistance;FB inertia force;Z1 andZ2 ground reactions to the front and rear axle wheels;mz =mδz reduced mass.The equation of motion of the vehicle in the driving direction can be written as the sum of the driving forces and drag forces parallel to the road surface as(5)mZx¨=FN-FT-FW-FPThe following simplifications were adopted. The drive torque applied to the rear wheel does not affect the pressure distributionZ1 andZ2. In fact, in vehicles with rear wheel drive, the rear is downforced whereas the front axle is relieved. The reaction valuesZ1 andZ2 were statistically determined. The coefficients of rolling resistance of all wheels are identical.Resistance forces depend on many factors and parameters; many of them are highly nonlinear. Modelling of wheel rolling resistance is the subject of separate research and the model used is very complex so in practice, a simplified model is often used.The rolling resistance of a single wheel is as follows:(6)FTi=crrZic rr is the coefficient of rolling resistance that depends on the type of road surface, the type of tire, or the pressure in the tire. This factor also depends on the speed of travel and increases with its growth; this is related to the increase in tire deformation. The value of the rolling resistance coefficient is usually given by the tire manufacturer, and exemplary values and empirical models on the dependence of this coefficient on the speed function can be found in the literature [16]. For low travel speeds (<40km/h), the value of the rolling resistance coefficient is assumed to be constant, independent of speed. In linear motion, the total rolling resistance is the sum of the rolling resistance of all vehicle wheels.If an identical rolling resistance coefficient for each wheel is assumed, which is often justified, given the similar construction of the front and rear axle wheels and identical tires and the same pressure operation, the sum of the rolling resistance shall be written as(7)FT=crrQmcos⁡αThe total wheel resistance should take into account the tire rolling resistance, bearing losses, and so forth.The essence of resistance in bearings is identical to the tire rolling resistance described previously. The bearing torque is proportional to the bearing load, friction coefficient, and bearing radius. Bearing friction is generally neglected, as the bearing resistance value is small compared to the tire rolling resistance [16].Similarly to the value of rolling resistance, the turning resistance can be expressed as the product of the weight of the car and the coefficient of turning, although, according to literature, at low speeds it plays negligibly smaller role than the rolling resistance [16].(8)Fs=mg⁡fzf s is the coefficient of turning and depends on the pressure force on the axles, tires, and so-calledintensity of turning, in other words, ratio of centrifugal and vertical acceleration v2/Rg.The value of air resistanceFP can be determined from(9)FP=γvsum22Acxwhere γ is the density of air.vsum [m/s] is the total air velocity resulting from the linear velocity of the vehiclevp and wind speedvw and is the vector sum,A is the vehicle front area [m2], andcx is the aerodynamic drag coefficient.In this study, a TS based on a constant gear ratio was used. PM BLDC (Permanent Magnet Brushless Direct Current) motors are very often used in electric power systems [17]. The mathematical description of all phenomena related to magnetism, electrical, thermal, or mechanical phenomena can lead to a considerable complexity of the model, so simplified models are often used. The electrical part of motor’s model can be found in work [9]. The electromagnetic torque generated by the motor is proportional to the flowing current and the so-calledthe motor torque constant Km and is given by [17](10)Me=2KmizIn steady state, the electromagnetic torqueMe is equal to the torque on the motor shaft, which is the sum of the load torque and the torque of the losses. The efficiency of an electric motor is defined as(11)ηm=PmechmPelectP m e c h m is the mechanical output power on the motor shaft andPelect is the electrical power input. Relations describing the phenomenon of loss can be found in professional literature [10].The motor torque is transmitted to the drive wheel via the TS transmission. The mechanical power on the drive wheel depends on the efficiency of the TS:(12)Pmechp=ηTSPmechm=ηTSηmPelectThe efficiencyηTS is a function of many factors; i.e., it depends on the powertrain design, the transmitted power, and quality of the system. For TS based on tension transmission, efficiency can be approx. 90%. The efficiency of the tie link transmission is a function of inter alia transmission ratioiTS, interaxial distance, belt angle, tension force, rotational speed, or torque [18–20]. In the case of lightweight vehicles in the TS design, it is often not necessary to use variable gear transmissions. Speed and torque control is achieved by the control system. On the one hand, it is an advantage because the TS is smaller in weight and less complicated, but on the other hand it forces the electric drive to operate in a wide range of torque and rotational speed, which affects the motor energy efficiency.


## 3.2. Building the Model of the Driving Route

The parameters of the path that the vehicle travels on have a significant effect on the resistance values shown in Section3.1. In the paper, the path of ride is determined by the exact path of the vehicle on a designated route. It is assumed that the best path is known and it is determined based on the team’s experience (it means that the quickest best path along the track is not calculated by the optimization method). The difference between the path and the track is as follows. A track is a system of physically existing roads between points A and B, and when the path does not exist physically, it is the trace that leaves the selected point of the vehicle (e.g., center of gravity) on the route surface. With this in mind, the journey along the route can take place on a number of paths. The route is only a limitation.Figure5 shows a schematic drawing of one path and its elevation profile. The description of the path is formed by its discretization, i.e., the division into finite numbers of elementary sections of Oti, i∈{1,2,...,n}.Each elementary segment is described by a parameter vector.Figure 5
Driving path and its inclination profile.O t i = [ L t i ,  R t i ,  S t i ,  E t i ], where Lti[m] is the length of the path, Rti[m] is the radius of curvature, Sti[%] is the inclination of the path, and Eti represents the type and condition of the path. The length of the path is the sum of the lengths of the Lti sections.These basic data are essential for the correct determination of motion resistance on each segment. It is assumed that any change in the inclination or radius of curvature makes it necessary to define further vectorOti.Bearing in mind that vehicles usually move in the open air, a suitable model should be created describing the external conditions. For this purpose, theOti vector describing the nature of the path should include additional columns such as velocity and direction of the wind vwi and vki, atmospheric pressure patmi, and air temperature tatmi. In this situation, the vector describing the path and external conditions takes the form: Oti=[Lti,Rti,Sti,Eti,vwi,vki,patmi,tatmi]. The speed and direction of the wind have a significant effect on the value of air resistance. Determining the wind speed is relatively difficult and not precise. Temporary winds may have different directions and a much higher amplitude than the average. Preliminary investigations have been carried out considering the uncertainty of weather conditions. During optimization, where optimal control strategies were sought, the wind speed was randomly changed. The concept of optimizing the control strategy from taking into account the changing weather conditions was to design a family of strategies. It would be possible to select one optimal strategy or part of it in the real world because of the existing and changing weather conditions.


## 3.3. Control System Design

Two modes of the control system are taken into consideration in this study. They reflect the physical implementations of the control system used in the built-in electric race vehicles of Prototype class [21] and UrbanConcept class [22]. Figures 6 and 7 represent block diagrams of the drive system of the vehicle for different control strategies. The control scheme A (Figure 6) is used in a case when a driver is able to control an acceleration of the vehicle by means of a binary signal. In this case, the goal of the optimization task is to find the control strategy taking into account two criteria, i.e., the determined travel time as well as the minimal power consumption of the vehicle. The control strategy is used to specify in which parts of the route the driver should press or release the acceleration treadle. When the setpoint signal USET∈{0,1} is assumed, then the optimization task is to search for logic values of the control signal as the function corresponding to the position of the vehiclex(t), USET =f(x(t)). The supervisory role of the driver is included in the functional scheme as a decision module. In this study, the behaviour of the driver is not modelled in the simulation. It means that he/she realises the strategy in the exact way. The input control signalUr is always different from the signal USET in real conditions. It means that the signal USET is mainly employed by the driver to know when the acceleration treadle should be used.Figure 6
Drive system diagram: control scheme A.Figure 7
Drive system diagram: control scheme B.In the case of the control scheme B (Figure7), the general goal of the optimization task is the same as in the previous case. However, the formal definition of the problem is quite different because in this case the optimal form of the velocity signal is searched according to the position of the vehicle on the routeVSET =f(x(t)). The velocity controller compares the setpoint signalVSET with the measured velocity of the vehicle and for certain controller setting as well as for control errorev it generates the control signal USET. As in the previous case, the supervisory role of the driver is included in the decision module. However, the main difference is that in the optimization process the setpoint signalVSET can be limited in the direct way.According to control schemes A and B, the input control signalUr is used by the motor control unit where the emended current regulator normally outputs the corresponding signals to the motor terminals. The driving torque of the motor shaft is transmitted through the transmission system to the wheels of the vehicle.


## 4. Evolutionary Optimization of the Velocity Strategy

### 4.1. Formulation of the Optimization Problem

One of the ways to define and to solve the dynamic optimization task of the velocity strategy of the vehicle is to reduce such a problem to the static optimization task through the so-called discretization of the independent variable. If one assumes that the independent variable is the position of the vehicle (a travelled distance)x(t), then the velocity strategyst can be represented by a vector of instantaneous values USET orVSET for discrete parts of the route. Figure 8 illustrates the idea of dividing the route into segments.Figure 8
Concept of discretization of the route.The path is divided inton segmentsOi each of the lengthLi, where the sum of the lengths of the segments is equal to the length of the total path. For each segment, it is necessary to determine the USETi or VSETi setpoint. The task is to find the optimal velocity strategy for control scheme A in the following form:(13)stoptU=USET1USET2…USETnT,USETi∈0,1whereas in the second case for control scheme B(14)stoptV=VSET1VSET2…VSETnT,VSETi∈Vmin,Vmaxwhere Vmin and Vmax can be determined in different ways, e.g., based on the minimum/maximum vehicle speed or based on the simple calculation corresponding to guidelines of the race competition.Evolutionary algorithms can be used to solve this problem. In the case of the search for strategystoptU it is preferable to use a classical genetic algorithm. In this instance, the genes in a chromosome of the individual can be encoded directly as the binary values of the control signal USETi. Therefore, the length of the chromosomeNch is equal to the number of sections of the route.In accordance with the control scheme B, the velocity profile is determined as the vector including setpoint values of the velocitystoptV. The components of the stoptV vector are floating point values; hence it is possible to use an evolutionary strategy proposed in [23]. Each of the genes in the chromosome of the individual corresponds to the setpoint value of the velocity during a part of the route.


### 4.2. Chromosome Codding for Different Velocity Strategies

For a vehicle that is equipped with a binary control system of the motor, the chromosome can be represented as the optimal vector of the control signal in the following formstoptU = USET1USET2…USETnT,USETi∈{0,1} for each segmentOi of the route (for control scheme A) or the optimal vector of the control signal in the second form stoptV = VSET1VSET2…VSETnT,VSETi∈[Vmin⁡,Vmax⁡] for each segmentOi of the route (for control scheme B).In this optimization task, the total distance of the roadDC is composed of SEM race lapsNto∈{1,2,...,N}. The division of the whole route on a finite number of segmentsOi is not recommended, because it can lead to unnecessary increase in the length of the vectorst. If the route consists of a minimum of four laps, it is possible to limit the search of vectorst to only three types of laps:(i)
the initial lap marked byOis(ii)
the middle laps marked byOim(iii)
the last lap marked asOifSuch assumption can be proven for at least three reasons:( 1 ) In the practical implementation of velocity strategies, the driver should not implement different control rules for each lap. The master control unit informs the driver about the current velocity strategy, but despite all, the driver should remember this strategy to some extent.( 2 ) The limitation of the velocity strategy only to the three types of laps reduces the length of the vectorst. In this way, the number of elements of the vectorst corresponds to the number of decision variables of the optimization task.( 3 ) Such division is reasonable by the nature of the laps. In the initial lap, it is necessary to accelerate the vehicle. Within the middle laps, the velocity profile of the vehicle should be the same. During the last lap, the previously acquired kinetic energyEk is utilized.The task of finding the optimal vectorstopt can be solved using two types of evolutionary algorithms:(i) Classical genetic algorithm for control without a velocity regulator (according to control scheme A presented in Figure6)(ii) Evolutionary strategy for the control system of the vehicle with a velocity regulator (according to control scheme B shown in Figure7)The usage of two algorithms forces the appropriate coding of the velocity strategy, which is described in the next two paragraphs.


#### 4.2.1. Chromosome Coding For Control Scheme A

The phenotype of each individual determines the control vectorstu that means each gene of the individual of the population is the value of the control signal USETi for each section of the segmentOi of the route. Assuming that we are looking for a velocity strategy of only three types of laps, the phenotype can be written as (for short form U ≡ USET)(15)stU=U1sU2s…UnsU1mU2m…UnsU1fU2f…UnfT,where Uis, Uim, Uif are the control signals for the initial lap, middle laps, and the last lap.


#### 4.2.2. Chromosome Coding for Control Scheme B

In this case, the elements in the vectorVSETi change over a wide range and can assign real number values. In the evolutionary strategy, the phenotype of each individual takes the following form:(16)stV=V1sV2s…VnsV1mV2m…VnsV1fV2f…VnfT,V i s, Vim, Vif are the setting values of velocity for the initial lap, middle laps, and the last lap.In the control mode, according to scheme B, the control signal is generated by the velocity controller based on the control error:(17)eV=VSET-Vsim,where(18)USET=1foreV≥10otherwiseThe coding for control scheme B in comparison with the coding for control scheme A has undoubtedly one advantage, namely, the possibility of introducing a velocity limit during the simulation of the passage of a given section of the route. A velocity limitation is of great practical importance. In many situations, it is necessary to reduce the velocity of the vehicle during certain parts of the route, for example, during cornering. Accordingly, in case of control scheme B, two velocity vectors are designated as follows:(i)
stV15 for VSETi = [0,15] m/s: although it is formally limited to 15m/s, it can never be achieved because the maximum vehicle speed (based on maximum engine speed) is approximately 15m/s;(ii)
s t V 10 for VSETi = [0,10] m/s.


### 4.3. The Fitness Functions

In the case of evolutionary optimization, the fitness functionfval has to be designed. This section describes the formulas used to calculate the fitness function of an individual that is needed for solving a task of evolutionary optimization. Each chromosome is a potential solution. Each control strategyst is used during numerical simulation. The value of the fitness functionfval can be calculated based on the outcomes of the simulation such as(i)
the total energy consumptionrsim [km/kWh](ii)
the simulation timetsim(iii)
the distancexsimand limitations such as(i)
the maximum travel timeTmax(ii)
the total distance of the vehicleDCThe value of the function is also used to evaluate the velocity strategy of the vehicle. In the preliminary studies [24], the following fitness functionfval was used.Fitness Function 1 (19) f v a l 1 = w 1 f 1 + w 2 f 2 + w 3 f 3 ,where(20)f1=1+rsimλ1λ1,(21)f2=HDc-xsimDc-xsimDcλ2,(22)f3=Htsim-TmaxTmax-tsimTmaxλ2,where the weightsw1 = 0.50E05,w2 =w3 = 0.25 and coefficients λ1 =−2, λ2 = 1 were chosen experimentally. The penalty functionsf2 andf3 were introduced. The Heaviside step functionH(.) is used to penalize individuals representing unacceptable solutions [25, 26], e.g., those that exceed the maximum time of the simulationH(tsim− Tmax) or do not reach the specified distanceH(Dc− xsim).Fitness Function 2. In this study, the authors have proposed the fitness function which would lead to the convergence of the evolutionary algorithm without searching for optimal values of weights. Such a function is shown and described below: (23)fval2=1+f4+f5f1f2f3,(24)f1=rsim,(25)f2=cosDc-dsimDc,(26)f3=cosTmax-tsimTmax,(27)f4=χ·HDc-dsim,(28)f5=χ·HTmax-tsim,wheref1 is a measure related to the simulated energy consumption of the vehicle,f2 andf3 correspond to the limitations such as the maximum travel time and the total distance of the vehicle,f4 andf5 also correspond to the limitations but these play strong penalty role, and χ is arbitrarily large value, e.g., significantly higher than maximal valuesf2=f3 =1. In this variant of the fitness functiondeath penalty [26] is applied for unacceptable individuals in the form of factorsf4 andf5. Individuals who do not reach a specified distance or exceed the time limit are eliminated. These two criteria act as a constraint and correspond to the limitations formulated in the guidelines of the race competition. On the other hand, the criteriaf2 andf3 have continuous nature and play very important role in the fitness function. They are very needed to evaluate continuously the individuals that are very close to and very far from the extremum.The argument of cos(.) function is not related to an angle. However, the usage of the cos(.) function is needed to eliminate the possibility of dividing by zero. In addition, the form of the cos(.) function is preferred. At the initial stage of optimization, when the simulation timetsim and the distancexsim differ significantly from the maximum timeTmax and the specified distanceDc, a small change will have a greater effect on the value of the fitness function. At the final stage of optimization, the function fval2 is more dependent on the simulated energy consumption of the vehicle related to the factorf1, and thanks to the cos(.) function, it is less dependent to the factorsf2 andf3. The proposed function does not require the introduction of weights and therefore it is easier to use this function in other applications, e.g., to find optimal velocity strategies for different routes.It should be noted that fitness functionsfval have been written in such a way as to ensure that the optimization algorithm will look for a minimum of these functions.


## 4.1. Formulation of the Optimization Problem

One of the ways to define and to solve the dynamic optimization task of the velocity strategy of the vehicle is to reduce such a problem to the static optimization task through the so-called discretization of the independent variable. If one assumes that the independent variable is the position of the vehicle (a travelled distance)x(t), then the velocity strategyst can be represented by a vector of instantaneous values USET orVSET for discrete parts of the route. Figure 8 illustrates the idea of dividing the route into segments.Figure 8
Concept of discretization of the route.The path is divided inton segmentsOi each of the lengthLi, where the sum of the lengths of the segments is equal to the length of the total path. For each segment, it is necessary to determine the USETi or VSETi setpoint. The task is to find the optimal velocity strategy for control scheme A in the following form:(13)stoptU=USET1USET2…USETnT,USETi∈0,1whereas in the second case for control scheme B(14)stoptV=VSET1VSET2…VSETnT,VSETi∈Vmin,Vmaxwhere Vmin and Vmax can be determined in different ways, e.g., based on the minimum/maximum vehicle speed or based on the simple calculation corresponding to guidelines of the race competition.Evolutionary algorithms can be used to solve this problem. In the case of the search for strategystoptU it is preferable to use a classical genetic algorithm. In this instance, the genes in a chromosome of the individual can be encoded directly as the binary values of the control signal USETi. Therefore, the length of the chromosomeNch is equal to the number of sections of the route.In accordance with the control scheme B, the velocity profile is determined as the vector including setpoint values of the velocitystoptV. The components of the stoptV vector are floating point values; hence it is possible to use an evolutionary strategy proposed in [23]. Each of the genes in the chromosome of the individual corresponds to the setpoint value of the velocity during a part of the route.


## 4.2. Chromosome Codding for Different Velocity Strategies

For a vehicle that is equipped with a binary control system of the motor, the chromosome can be represented as the optimal vector of the control signal in the following formstoptU = USET1USET2…USETnT,USETi∈{0,1} for each segmentOi of the route (for control scheme A) or the optimal vector of the control signal in the second form stoptV = VSET1VSET2…VSETnT,VSETi∈[Vmin⁡,Vmax⁡] for each segmentOi of the route (for control scheme B).In this optimization task, the total distance of the roadDC is composed of SEM race lapsNto∈{1,2,...,N}. The division of the whole route on a finite number of segmentsOi is not recommended, because it can lead to unnecessary increase in the length of the vectorst. If the route consists of a minimum of four laps, it is possible to limit the search of vectorst to only three types of laps:(i)
the initial lap marked byOis(ii)
the middle laps marked byOim(iii)
the last lap marked asOifSuch assumption can be proven for at least three reasons:( 1 ) In the practical implementation of velocity strategies, the driver should not implement different control rules for each lap. The master control unit informs the driver about the current velocity strategy, but despite all, the driver should remember this strategy to some extent.( 2 ) The limitation of the velocity strategy only to the three types of laps reduces the length of the vectorst. In this way, the number of elements of the vectorst corresponds to the number of decision variables of the optimization task.( 3 ) Such division is reasonable by the nature of the laps. In the initial lap, it is necessary to accelerate the vehicle. Within the middle laps, the velocity profile of the vehicle should be the same. During the last lap, the previously acquired kinetic energyEk is utilized.The task of finding the optimal vectorstopt can be solved using two types of evolutionary algorithms:(i) Classical genetic algorithm for control without a velocity regulator (according to control scheme A presented in Figure6)(ii) Evolutionary strategy for the control system of the vehicle with a velocity regulator (according to control scheme B shown in Figure7)The usage of two algorithms forces the appropriate coding of the velocity strategy, which is described in the next two paragraphs.


### 4.2.1. Chromosome Coding For Control Scheme A

The phenotype of each individual determines the control vectorstu that means each gene of the individual of the population is the value of the control signal USETi for each section of the segmentOi of the route. Assuming that we are looking for a velocity strategy of only three types of laps, the phenotype can be written as (for short form U ≡ USET)(15)stU=U1sU2s…UnsU1mU2m…UnsU1fU2f…UnfT,where Uis, Uim, Uif are the control signals for the initial lap, middle laps, and the last lap.


### 4.2.2. Chromosome Coding for Control Scheme B

In this case, the elements in the vectorVSETi change over a wide range and can assign real number values. In the evolutionary strategy, the phenotype of each individual takes the following form:(16)stV=V1sV2s…VnsV1mV2m…VnsV1fV2f…VnfT,V i s, Vim, Vif are the setting values of velocity for the initial lap, middle laps, and the last lap.In the control mode, according to scheme B, the control signal is generated by the velocity controller based on the control error:(17)eV=VSET-Vsim,where(18)USET=1foreV≥10otherwiseThe coding for control scheme B in comparison with the coding for control scheme A has undoubtedly one advantage, namely, the possibility of introducing a velocity limit during the simulation of the passage of a given section of the route. A velocity limitation is of great practical importance. In many situations, it is necessary to reduce the velocity of the vehicle during certain parts of the route, for example, during cornering. Accordingly, in case of control scheme B, two velocity vectors are designated as follows:(i)
stV15 for VSETi = [0,15] m/s: although it is formally limited to 15m/s, it can never be achieved because the maximum vehicle speed (based on maximum engine speed) is approximately 15m/s;(ii)
s t V 10 for VSETi = [0,10] m/s.


## 4.2.1. Chromosome Coding For Control Scheme A

The phenotype of each individual determines the control vectorstu that means each gene of the individual of the population is the value of the control signal USETi for each section of the segmentOi of the route. Assuming that we are looking for a velocity strategy of only three types of laps, the phenotype can be written as (for short form U ≡ USET)(15)stU=U1sU2s…UnsU1mU2m…UnsU1fU2f…UnfT,where Uis, Uim, Uif are the control signals for the initial lap, middle laps, and the last lap.


## 4.2.2. Chromosome Coding for Control Scheme B

In this case, the elements in the vectorVSETi change over a wide range and can assign real number values. In the evolutionary strategy, the phenotype of each individual takes the following form:(16)stV=V1sV2s…VnsV1mV2m…VnsV1fV2f…VnfT,V i s, Vim, Vif are the setting values of velocity for the initial lap, middle laps, and the last lap.In the control mode, according to scheme B, the control signal is generated by the velocity controller based on the control error:(17)eV=VSET-Vsim,where(18)USET=1foreV≥10otherwiseThe coding for control scheme B in comparison with the coding for control scheme A has undoubtedly one advantage, namely, the possibility of introducing a velocity limit during the simulation of the passage of a given section of the route. A velocity limitation is of great practical importance. In many situations, it is necessary to reduce the velocity of the vehicle during certain parts of the route, for example, during cornering. Accordingly, in case of control scheme B, two velocity vectors are designated as follows:(i)
stV15 for VSETi = [0,15] m/s: although it is formally limited to 15m/s, it can never be achieved because the maximum vehicle speed (based on maximum engine speed) is approximately 15m/s;(ii)
s t V 10 for VSETi = [0,10] m/s.


## 4.3. The Fitness Functions

In the case of evolutionary optimization, the fitness functionfval has to be designed. This section describes the formulas used to calculate the fitness function of an individual that is needed for solving a task of evolutionary optimization. Each chromosome is a potential solution. Each control strategyst is used during numerical simulation. The value of the fitness functionfval can be calculated based on the outcomes of the simulation such as(i)
the total energy consumptionrsim [km/kWh](ii)
the simulation timetsim(iii)
the distancexsimand limitations such as(i)
the maximum travel timeTmax(ii)
the total distance of the vehicleDCThe value of the function is also used to evaluate the velocity strategy of the vehicle. In the preliminary studies [24], the following fitness functionfval was used.Fitness Function 1 (19) f v a l 1 = w 1 f 1 + w 2 f 2 + w 3 f 3 ,where(20)f1=1+rsimλ1λ1,(21)f2=HDc-xsimDc-xsimDcλ2,(22)f3=Htsim-TmaxTmax-tsimTmaxλ2,where the weightsw1 = 0.50E05,w2 =w3 = 0.25 and coefficients λ1 =−2, λ2 = 1 were chosen experimentally. The penalty functionsf2 andf3 were introduced. The Heaviside step functionH(.) is used to penalize individuals representing unacceptable solutions [25, 26], e.g., those that exceed the maximum time of the simulationH(tsim− Tmax) or do not reach the specified distanceH(Dc− xsim).Fitness Function 2. In this study, the authors have proposed the fitness function which would lead to the convergence of the evolutionary algorithm without searching for optimal values of weights. Such a function is shown and described below: (23)fval2=1+f4+f5f1f2f3,(24)f1=rsim,(25)f2=cosDc-dsimDc,(26)f3=cosTmax-tsimTmax,(27)f4=χ·HDc-dsim,(28)f5=χ·HTmax-tsim,wheref1 is a measure related to the simulated energy consumption of the vehicle,f2 andf3 correspond to the limitations such as the maximum travel time and the total distance of the vehicle,f4 andf5 also correspond to the limitations but these play strong penalty role, and χ is arbitrarily large value, e.g., significantly higher than maximal valuesf2=f3 =1. In this variant of the fitness functiondeath penalty [26] is applied for unacceptable individuals in the form of factorsf4 andf5. Individuals who do not reach a specified distance or exceed the time limit are eliminated. These two criteria act as a constraint and correspond to the limitations formulated in the guidelines of the race competition. On the other hand, the criteriaf2 andf3 have continuous nature and play very important role in the fitness function. They are very needed to evaluate continuously the individuals that are very close to and very far from the extremum.The argument of cos(.) function is not related to an angle. However, the usage of the cos(.) function is needed to eliminate the possibility of dividing by zero. In addition, the form of the cos(.) function is preferred. At the initial stage of optimization, when the simulation timetsim and the distancexsim differ significantly from the maximum timeTmax and the specified distanceDc, a small change will have a greater effect on the value of the fitness function. At the final stage of optimization, the function fval2 is more dependent on the simulated energy consumption of the vehicle related to the factorf1, and thanks to the cos(.) function, it is less dependent to the factorsf2 andf3. The proposed function does not require the introduction of weights and therefore it is easier to use this function in other applications, e.g., to find optimal velocity strategies for different routes.It should be noted that fitness functionsfval have been written in such a way as to ensure that the optimization algorithm will look for a minimum of these functions.


## 5. Case Study

As described in Sections3 and 4, process of modelling and optimization was applied and verified using a prototype electric vehicle on street race circuit in Rotterdam. The methodology was also carried out on the experimental rides track of the Fiat Chrysler Automobiles Poland factory, located in Tychy. Based on (5)-(12), the simulation model was created in the MATLAB/Simulink® environment. Optimization was also carried out in MATLAB using the evolutionary algorithms available in theGlobal Optimization Toolbox. The library of functions for simulating, optimizing, analyzing, and presenting results has been developed.


### 5.1. Research Object

Light electric vehicle was a research object. The vehicle shown in Figure9 is a three-wheeled self-supporting structure. The vehicle fuselage is made of composites based on epoxy resin and woven of carbon and aramid fibre [21].Figure 9
Research object: prototype electric vehicle “MuSHELLka” (photo by M. Wyleżoł, 2012).The dimensions of the vehicle and basic data are shown in Figure10 and in Table 1.Table 1
Research object main parameters.
Type of the motor PM BLDC Transmission unit Toothed synchronous belt The minimum turning radius 8 [m] The battery type and capacity Lithium-ion polymer 2100[mAh] Mass of the vehicle 42 [kg] Mass with the driver min. 92 [kg] Diameter of wheels 0.478 [m] Moment of inertia of the front wheelIkP 0.023 [kgm2] Moment of inertia of the rear wheelIkT 0.023[kgm2] Drag coefficientcx(see Table 2) 0,22 – 0,23 Frontal areaA[m2] (see Table 2) 0.275 -0.297 The total rolling coefficientcrr (see Table 4) 0,0035 – 0,0058,Table 2
The air drag coefficient and frontal area.
c x: numerical simulation c x: experimental research Frontal areaA[m2] 0.23 0.22 { 0.275,0.297 }Figure 10
Dimensions of the vehicle.Identification of parameters for the mathematical model of a research object was carried out by numerical simulations, tests on specially designed test bench [27] and during conducted experiments in real environment [28].The vehicle can be driven by two types of drive unit:(i)
DT1: drive unit based on PM BLDC motors and belt transmission. Drive wheel is equipped with a one-way clutch.(ii)
DT2: a high-torque PM BLDC motor mounted inside the drive wheel. The torque generated by the motor is directly transmitted to the rear drive wheel.The catalogue data can be used for creating the mathematical model of the DT1 motor, but for identifying the exact parameters, it was necessary to carry out experimental tests. The manufacturer of the DT2 motor does not provide much information, so it was necessary to test the drive unit as well. DT1 drivetrain efficiency was estimated using a BG75x50 motor with toothed synchronous belt transmission ( z1= 28 teeth, z2 = 200 teeth,iTS = 7.14). The map of the efficiency was shown on Figure 11(a). The same test was carried out for DT2 drive unit with a motor directly mounted on the drive wheel (Figure 11(b)). The obtained efficiency maps were used to build the simulation model of the tested drivetrain systems.Figure 11
Maps of the efficiency of tested drivetrain systems.
(a)
DT1 (b)
DT2


### 5.2. Determination of Air Drag Coefficient

The air drag coefficientcx was estimated as a result of CFD analysis using ANSYS software [24]; and it was also verified by the research in wind tunnel at the Institute of Aviation in Warsaw. Frontal area of vehicle A was determined from a three-dimensional model in CAD environment, Catia V5. The values of drag coefficientcx and frontal area A are shown in Table 2.


### 5.3. Determination of Rolling Resistance

According to formula (7), the rolling resistance depends on the weight of the vehicle (Table 1), the tires rolling resistance coefficient, and the resistance of the bearings. Depending on the drive unit (DT1, DT2) used in the vehicle, additional resistance has to be added. In the DT1 system, there is a one-way clutch which also generates some resistance torque. The DT2 system cannot be mechanically disconnected and therefore has a higher rolling resistance. As a result of the coast down test, the total rolling resistance of the vehicle has been identified for the four types of tires. The test was performed only with the DT1 system, i.e., the system with the clutch. Table 3 shows the basic parameters of the tires used in the experiments.Table 3
The types of tires used in the prototype vehicle.
Tire's marker Name of the tire Tire’s symbol Shape Type of tread Working pressure W 1 Michelin Shell Eco-marathon 45/75R16 Rectangular Slick max.500 [kPa] W 2 Schwalbe Ul-Termo 23-406 Round Slick 600-1100 [kPa] W 3 Continental Sport Contact 28-406 Round Slick 600 [kPa] W 4 Schwalbe Marathon 40-406 Round Symmetrical 400-700 [kPa]Table 4
The total rolling coefficient depending on the tires.
Tire's marker Name of the tire The total rolling coefficientcrr mMAPE [%] W 1 Michelin ShellEco-marathon 0.0035 7.84 W 2 Schwalbe Ultermo 0.0043 11.03 W 3 Continental SportContact 0.0058 12.52 W 4 Schwalbe Marathon 0.0053 4.45The prototype vehicle can be equipped with specially manufacturedMichelin tires with a low rolling resistant coefficient ofcrr = 0.0014.In addition to the specialMichelin tires, three other types of tires were tested. The identification of rolling resistance coefficients was also estimated using a coast down test. The measurements were made as follows:(i)
The driver accelerated the vehicle to a fixed speed.(ii)
At the moment of crossing the marked gate, the driver pressed the marker button to save the flag in the acquired data.(iii)
The run continued until free braking, after which the driver switched off the measuring system.(iv)
The measurements were made at both directions of the road, to minimize influence of the wind.(v)
Before the attempts, wind velocity, temperature, and atmospheric pressure were measured.(vi)
The distance of the road was measured using a measurement system mounted on the vehicle and it was also verified using distance measuring wheel.On the basis of the experimental data, the total rolling resistance coefficient of the vehicle was determined. Table4 shows the coefficients of rolling resistance of the vehicle depending on the tire used during the test. Total resistance includes the tire rolling coefficient, the bearings resistance, and the resistance of the clutch. Verification of the data derived from the simulation and the experiment was conducted usingmMAPE (modifiable Mean Absolute Percentage Error). For velocity simulationmMAPE can be determined by(29)mMAPE=100%1n∑i=1nVr-VsimVmax-Vmin.Figure12 shows the experimental data obtained in a coast down test compared with the results obtained during the computer simulation.Figure 12
Velocity of the vehicle during coast down test.


### 5.4. Model of the Race Circuit

In accordance with the methodology presented in Section3.2, three routes have been described:(i)
Route 1: a street race circuit at theShell Eco-marathon competition in Rotterdam.(ii)
Routes 2 and 3: on the experimental rides track of the Fiat Chrysler Automobiles Poland factory, located in Tychy, in the south of Poland.TheShell Eco-marathon competition is held on a street circuit in Rotterdam. One lap is 1630 meters long. Vehicles have to cover ten laps in no longer than 39 minutes so such requirements impose an average minimum speed of 25 km/h.The trajectory of the vehicle path was developed using satellite images. The vehicle trajectory was determined using the maximum curvature possible that could be obtained on the existing roads.Movement of the vehicle in a curved motion causes the centrifugal force. Centrifugal force creates additional resistance. This resistance is inversely proportional to the radius of the path (8), so it is justified to select the largest radius. In addition, if the centrifugal force exceeds the limit value, it could provide to the side slip of the vehicle. According to drawings inAutoCAD Civil 3D, the lengths of straight sections and the lengths and the radius of arcs were obtained. Figure 13 shows a plan view of the path of the drive.Figure 13
Plan view of the path of the drive on route 1 in Rotterdam.Driving the prototype vehicle on public roads is prohibited and can be dangerous. Vehicle testing was also carried out on the experimental rides track of the Fiat Chrysler Automobiles Poland factory.For the track in Tychy geodetic measurements were taken using theDGPS (Differential Global Positioning System) technology. The measuring points were collected along the axis of the road, with the distance between points amounting to 10 meters. If on the measured route, direction or elevation were changed, the frequency of measurement was increased. Geodetic coordinates were also imported intoAutoCAD Civil 3D and were used to prepare data concerning the route that was imported into the simulation environmentMATLAB-Simulink. Figures 14 and 15 show the path and elevation profile, respectively.Figure 14
Plan view of the path of the drive on route 3 in Tychy.Figure 15
Profile of the path on route 3 in Tychy.Selected routes 2 and 3 have significant differences compared to route 1, especially in the height profile. This set of routes was selected to verify the effect of the route profile on the achieved results and to verify the proposed optimization methods for various road conditions. The routes in Tychy are located in an open area; there are no large buildings, which guarantees stability of wind conditions. The technical condition of the surface on routes 2 and 3 is much worse than on route 1, which has a significant influence on the rolling resistance. Rolling resistance was also identified and introduced into the simulation model. For route 2, thecrr coefficient was 0.0056 and for route 3 was 0.0068.


### 5.5. Verification of the Simulation Model

In order to verify the compatibility of the identified simulation model, a series of numerical simulations and experimental tests were performed. The first experiment consisted of running 300 [m] with a control strategy of 50 [m] with the working motor, and after that the vehicle was moving freely for the next 50 [m]. The motor is turned on three times. Figures16(a) and 16(b) show the velocity of the vehicle and the power of the motor in function of the distance determined during the experiment (red curves) and during the numerical simulations (green curves).Figure 16
The velocity and the power determined with used DT1 drive unit.
(a)
DT1 (b)
DT1Verification of the data derived from the simulation and the experiment was conducted usingmMAPE (29). For motor power simulation,mMAPE can be determined by(30)mMAPEP=100%1n∑i=1nPr-PsimPmax-Pmin,where Pr is the motor power measured during experiment and Psim is the motor power conducted in simulation.As a result of the model validation tests, estimation error of the velocity was on level between 2.35% and 5.34%. However, the motor power error level was greater. It was between 8.02% and 9.98%.


### 5.6. Determination of the Optimal Control Strategy

The length of route 1 and the number of laps are determined by the organizers of the Shell Eco-marathon competition. The length of routes 2 and 3 is derived from the existing road plan, where it was possible to carry out trials in real conditions. The number of laps for routes 2 and 3 was assumed arbitrarily.The travel time for route 1 is based on the rules of the race, which determine the minimum average speed of the vehicle. Based on the average speed and the length of routes 2 and 3, the travel times for these routes were determined. The maximum travel times and the total length of the route, as well as the number and length of one lapD, are summarized in Table 5.Table 5
Travel times and total lengths of the routes.
Route Length of the lapD [m] The number of laps The total length of the routeDc [m] Travel timeTmax [s] Route 1 Rotterdam 1630 10 16 300 2340 Route 2 Tychy 640 5 3200 460 Route 3 Tychy 1000 5 5000 720Evolutionary algorithms which are available in theMATLAB Global Optimization Toolboxwere used. From the MATLAB workspace, the algorithm ran a simulation model in aSimulink environment. During optimization while using evolutionary algorithms, the maximum number of fitness function evaluations Nfval is the product of the population size Npop and number of generations Ngen.It was assumed that the total number of fitness function evaluations was equal to 104; therefore, should this value be exceeded, the calculation would be cancelled. Therefore, the total number of generations Ngen. is dependent on the value ofpopulation sizeNpop.As part of this work phase, the vectorsstU, stV15, and stV10 of the optimum control strategy for an electric vehicle for routes 1, 2, and 3 were determined. During evolutionary optimization, the heuristic method of crossing was used, and the probability of crossing was established on the basis of previous studies at the level of 0.8 [24, 29]. At the stage of reproduction, the elitist succession method was applied where the best adapted two individuals in the current population passed to the next population. An adaptive method of mutation was used. The remaining parameters of the algorithm remain the default for the MATLABGlobal Optimization Toolbox.The first case was the searchstU vector for route 1. Figures 17(a)–17(c) show the vehicle velocityVsim obtained by computer simulation. The graph of thestU control signal is graphically displayed. Figure 17 shows the results where the optimal control strategy stUopt was used, and population was Npop=100 and Ngen=100. The total energy consumption was estimated asrsim=527,633 [km / kWh]. Figure 17(a) shows the vehicle speed and control signalstU during the first lap while Figure 17(b) shows center laps and Figure 17(c) shows the last lap. Figure 17(d) shows the curves of the mean values of the fitness functionsfval, depending on the number of calculations of these functions, during the optimization of the control strategy. There are eight curves in the graph (Figure 17(d)) for each of the Npop parameters used in the algorithm.Figure 17
Vehicle speed and control signalstU for route 1.
(a)
Lap 1 (b)
Laps 2-9 (c)
Lap 10 (d)
Mean values of the fitness functionsWhen analyzing the vehicle velocity chart on the last lap of route 1 (Figure17(c)), it can be seen that, according to the determined control strategy, the propulsion system of the vehicle should no longer be used. The vehicle should use the stored kinetic energy.The analysis of vehicle velocityVsim during the last lap indicates that at the moment of reaching the given distance the speed was still about 2 [m/s]. This demonstrates that the simulated vehicle has kinetic energy that would have to be dispersed by the braking system or recovered by the electric motor. Having such a speed does not mean that there is still possibility for improvement. Lower speed at the end of the route has to cause higher speed of the earlier parts of the route. This speed is the optimum speed due to the chosen criterion. From the analysis of the shape of the fitness function (Figure 17(d)), it can be seen that, for route 1, it is preferable to use a smaller population of Npop=25 and increase the number of Ngen=400. The curve shape for this case is steep and it falls to the minimum value of the fitness function, which may indicate that the algorithm has reached a convergence.Results of the optimization where B type of coding was used were analyzed. The results obtained as the mean of thersim result together with the standard deviation from all the tests were summarized in Table 6.Table 6
The results for route 1.
The type of coding,A The type of coding,B s t U s t V 15 s t V 10 N p o p r s i m [km/kWh] std r s i m [km/kWh] std r s i m [km/kWh] std 25 524,65 1,459 428,45 17,41 515,69 5,4775 50 525,79 1,0053 444,31 22,486 514,92 3,4094 75 524,97 3,9913 438,42 28,173 512,56 8,0682 100 525,89 1,8128 422,32 17,97 511,66 4,9379 125 526,25 1,3862 427,33 39,616 513,41 7,0765 150 524,3 2,979 443,5 47,701 514,86 2,7086 175 525,13 2,511 414,28 27,299 518,38 1,9358 200 525,79 1,9666 423,37 25,776 512,04 7,8357By analyzing the results in Table6, it can be seen that for route 1 better results were obtained (reduced simulated total energy consumption) using A coding compared to B, especially for stV15, whose results are characterized by the highest standard deviation. For routes 2 and 3, the same as for route 1, the procedure was followed. Optimization strategies have been optimized for all variants of control vectorstU, stV15,st10.The results using the control strategystV10 are shown in Figures 18(a)–18(c) (Npop=100 wherersim=351,795 [km / kWh]). Figure 18(d) illustrates the graph of the value of the fitness function.Figure 18
Vehicle speed and control signalstV10 for route 2.
(a)
Lap 1 (b)
Laps 2-4 (c)
Lap 5 (d)
Mean values of the fitness functionsThe results obtained for route 2 are presented in Table7.Table 7
The results for route 2.
The type of coding,A The type of coding,B s t U s t V   15 s t V   10 N p o p r s i m [km/kWh] std r s i m [km/kWh] std r s i m [km/kWh] std 25 341,38 6,6734 342,94 2,3839 348,73 5,3532 50 345,73 4,6359 350,03 4,62 350,6 3,8862 75 345,12 5,5134 349,61 3,3476 354,19 0,67858 100 348,36 4,5332 352,53 2,6117 354,4 1,9295 125 346,32 5,2704 350,17 4,9039 354,83 1,4756 150 345,52 6,3714 353,26 1,4885 353,43 0,71015 175 350,63 1,8728 352,74 4,7141 355,64 0,38022 200 350,73 1,4496 354,0 1,3389 355,69 1,4942For route 2 compared to route 1, the best results were achieved using thestV10 vector. When using vectorsstU and stV15 the results were slightly worse. Having analyzed the size of the population used in the optimization process, it should be noted that for route 2 it is preferable to use higher value of Npop. For all the analyzed cases, the best average results for Npop=200 were obtained. Table 8 summarizes the results for route 3.Table 8
The results for route 3.
The type of coding,A The type of coding,B s t U s t V   15 s t V   10 N p o p r s i m [km/kWh] std r s i m [km/kWh] std r s i m [km/kWh] std 25 280,71 10,173 286,43 7,323 309,94 2,776 50 293,46 6,142 295,81 7,894 308,74 4,399 75 297,11 5,002 292,17 7,367 311,41 1,856 100 298,26 3,231 295,15 7,022 309,57 4,049 125 297,81 2,980 298,66 12,977 311,68 3,445 150 294,53 3,346 302,17 9,623 312,88 1,249 175 292,34 2,317 306,21 6,816 314,81 1,261 200 296,86 2,337 305,31 11,737 310,89 4,880For route 3, the average results were reported for a larger population ofNpop with B coding. In the case of using vector with A type coding, no characteristic change was noted. As expected, comparing the standard deviation between algorithms using B coding, the standard deviation is less for the less variance range.Evolutionary algorithms are not complete search algorithms; therefore, it is not possible to determine whether the solution found is best for the chosen criterion and the constraints assumed.


### 5.7. Applying the Optimal Control Strategy in Real Conditions

Electric vehicle energy consumption optimization was finally verified in real conditions at the European edition of theShell Eco-marathon held on 16-18 May 2014 in Rotterdam.Optimal control strategystU was displayed by a prototype visualization system. Visualization unit is installed in the steering wheel of the vehicle. The visualization system shows the basic data such as instantaneous speed, average speed, travel time, and distance. In addition, the visualization system provides information on which part of the track the driver should press the acceleration button to execute the optimal control strategy. Measurement of the position on the track is carried out by the odometrical measuring system. Due to the possibility of accumulating errors in the measurement of the road, the driver will reset the measuring system every time he/she exceeds the lap line.Table9 presents the results of all the measurements conducted during theShell Eco-marathon. The columns in the Table 9 contain(i)
rSEm [km/kWh]: the result recorded during the race by the measuring device provided by the organizer; the results are available on the Shell Eco-marathon official website;(ii)
rr [kW/kWh]: the result based on the indications of a prototype measuring system permanently installed in the vehicle;(iii)
rsim: the result estimated by computer simulation using recorded real control signalUr;(iv)
ropt: the result obtained by computer simulation using the optimal control strategyUopt.Table 9
Comparison of simulation results with experimental tests.
Number of attempts r S E – m r r r s i m r o p t e m e i e s t Attempt 1 365,4 360,8 - - 1,25 - - Attempt 2 481,3 508,2 514,8 526,5 5,58 1,29 3,6 Attempt 3 487,3 512,7 519,8 526,5 5,21 1,39 2,69 Attempt 4 Failed attempt, damaged tires during 1st lapTable9 also includes the following:(i)
Measurement error:(31)em=rSEm-rrrSEm100%.(ii)
Identification error:(32)ei=rr-rsimrr100%.(iii)
Strategy implementation error:(33)est=rr-roptrr100%.There were four measurement trials available. During run 1, a much weaker result was obtained than estimated by simulation. The driver performed run 1 without a defined control strategy. As part of the experiment, the driver and the support team were asked to select the optimal tires for the vehicle and to implement the proposed control strategy. During run 1, only the travel time of individual laps was checked. According to the choice of the team in run 1, tires W2 were selected (Table4). During this test the worst result was obtained about 24% and 25% worse than in the other trials. In subsequent runs 2 and 3, the vehicle was equipped with the best W1 tires (Table 4). The driver has implemented an algorithm-based optimal control strategy. Strategy information was provided by the visualization system. During run 3 the best result was recorded. This result was taken into account in the classification of competitions, where the team placed in 12 positions out of 28 classified teams. Trials 4, due to damage to the tires, were not successful.Measurement error was from 1.25% to 5.58%, depending on the sample. The measurement system installed permanently in the vehicle is a prototype system. This system has been calibrated using available laboratory equipment. Measurement errors indicate that other measuring devices of the higher accuracy class should be used.Based on the measurement data obtained during test drives, the model was tuned, thus reducing the identification error. The competition took place in stable weather; the divergence of results was small and was about 1.29% for sample 2 and 1.39% for sample 3.The difference between the results determined during numerical simulation and that obtained in the real run was 3.6% for sample 2 and 2.69% for sample 3. The difference at this level was considered to be small. Searching for optimal control strategy has not taken into account possible "errors" which could be made by the driver and which result from the current situation on the route, for example, selecting the path and using the route by other participants in the race. Including driver behaviour in the simulation model should be the subject of further work.Figure19 shows the vehicle velocity together with the control signals recorded during run 2. The black curves represent the real vehicle velocityVr and the control signal. Green curves are the optimum control strategystU and the optimal velocityVopt determined by the simulation. The graph shows the entire distance travelled (10 laps). The analysis of the graph shows that the driver did not perform the exact planned strategy. Errors in the implementation of the strategy are related to the disruptions. During the competition, the driver is not alone on the track; in certain situations the driver cannot overtake other participants and cannot take the optimal path. In addition, weather condition such as wind has influence on the motion resistance. In such cases, the driver has to modify the optimal control strategy.Figure 19
Comparison of speed and control signal during attempt 2.In Figure20, the velocity recorded during the raceVr (black) and the velocity obtained by the computer simulationVsim (red) were shown.Vsim was determined by computer simulation by using real control signal. The purpose of this simulation was to investigate the correctness of identification of the simulation model and to identify the identification error.Figure 20
Comparison of speed registered in real conditions with the speed determined by computer simulation.Compared to theShell Eco-marathon competitions held in 2012 and 2013, our performance has improved. The results in these years were equal to 425 [km / kWh] and 454 [km / kWh], respectively. In 2012 and 2013, the control strategy was chosen on the basis of the simulation model, but the model was not adequately identified, and the control strategy had to be modified. In previous editions of the race, the implementation of optimal control strategy was not supported by the visualization system, which would inform the driver of the need to take specific actions (pressing or releasing the accelerator button).


## 5.1. Research Object

Light electric vehicle was a research object. The vehicle shown in Figure9 is a three-wheeled self-supporting structure. The vehicle fuselage is made of composites based on epoxy resin and woven of carbon and aramid fibre [21].Figure 9
Research object: prototype electric vehicle “MuSHELLka” (photo by M. Wyleżoł, 2012).The dimensions of the vehicle and basic data are shown in Figure10 and in Table 1.Table 1
Research object main parameters.
Type of the motor PM BLDC Transmission unit Toothed synchronous belt The minimum turning radius 8 [m] The battery type and capacity Lithium-ion polymer 2100[mAh] Mass of the vehicle 42 [kg] Mass with the driver min. 92 [kg] Diameter of wheels 0.478 [m] Moment of inertia of the front wheelIkP 0.023 [kgm2] Moment of inertia of the rear wheelIkT 0.023[kgm2] Drag coefficientcx(see Table 2) 0,22 – 0,23 Frontal areaA[m2] (see Table 2) 0.275 -0.297 The total rolling coefficientcrr (see Table 4) 0,0035 – 0,0058,Table 2
The air drag coefficient and frontal area.
c x: numerical simulation c x: experimental research Frontal areaA[m2] 0.23 0.22 { 0.275,0.297 }Figure 10
Dimensions of the vehicle.Identification of parameters for the mathematical model of a research object was carried out by numerical simulations, tests on specially designed test bench [27] and during conducted experiments in real environment [28].The vehicle can be driven by two types of drive unit:(i)
DT1: drive unit based on PM BLDC motors and belt transmission. Drive wheel is equipped with a one-way clutch.(ii)
DT2: a high-torque PM BLDC motor mounted inside the drive wheel. The torque generated by the motor is directly transmitted to the rear drive wheel.The catalogue data can be used for creating the mathematical model of the DT1 motor, but for identifying the exact parameters, it was necessary to carry out experimental tests. The manufacturer of the DT2 motor does not provide much information, so it was necessary to test the drive unit as well. DT1 drivetrain efficiency was estimated using a BG75x50 motor with toothed synchronous belt transmission ( z1= 28 teeth, z2 = 200 teeth,iTS = 7.14). The map of the efficiency was shown on Figure 11(a). The same test was carried out for DT2 drive unit with a motor directly mounted on the drive wheel (Figure 11(b)). The obtained efficiency maps were used to build the simulation model of the tested drivetrain systems.Figure 11
Maps of the efficiency of tested drivetrain systems.
(a)
DT1 (b)
DT2


## 5.2. Determination of Air Drag Coefficient

The air drag coefficientcx was estimated as a result of CFD analysis using ANSYS software [24]; and it was also verified by the research in wind tunnel at the Institute of Aviation in Warsaw. Frontal area of vehicle A was determined from a three-dimensional model in CAD environment, Catia V5. The values of drag coefficientcx and frontal area A are shown in Table 2.


## 5.3. Determination of Rolling Resistance

According to formula (7), the rolling resistance depends on the weight of the vehicle (Table 1), the tires rolling resistance coefficient, and the resistance of the bearings. Depending on the drive unit (DT1, DT2) used in the vehicle, additional resistance has to be added. In the DT1 system, there is a one-way clutch which also generates some resistance torque. The DT2 system cannot be mechanically disconnected and therefore has a higher rolling resistance. As a result of the coast down test, the total rolling resistance of the vehicle has been identified for the four types of tires. The test was performed only with the DT1 system, i.e., the system with the clutch. Table 3 shows the basic parameters of the tires used in the experiments.Table 3
The types of tires used in the prototype vehicle.
Tire's marker Name of the tire Tire’s symbol Shape Type of tread Working pressure W 1 Michelin Shell Eco-marathon 45/75R16 Rectangular Slick max.500 [kPa] W 2 Schwalbe Ul-Termo 23-406 Round Slick 600-1100 [kPa] W 3 Continental Sport Contact 28-406 Round Slick 600 [kPa] W 4 Schwalbe Marathon 40-406 Round Symmetrical 400-700 [kPa]Table 4
The total rolling coefficient depending on the tires.
Tire's marker Name of the tire The total rolling coefficientcrr mMAPE [%] W 1 Michelin ShellEco-marathon 0.0035 7.84 W 2 Schwalbe Ultermo 0.0043 11.03 W 3 Continental SportContact 0.0058 12.52 W 4 Schwalbe Marathon 0.0053 4.45The prototype vehicle can be equipped with specially manufacturedMichelin tires with a low rolling resistant coefficient ofcrr = 0.0014.In addition to the specialMichelin tires, three other types of tires were tested. The identification of rolling resistance coefficients was also estimated using a coast down test. The measurements were made as follows:(i)
The driver accelerated the vehicle to a fixed speed.(ii)
At the moment of crossing the marked gate, the driver pressed the marker button to save the flag in the acquired data.(iii)
The run continued until free braking, after which the driver switched off the measuring system.(iv)
The measurements were made at both directions of the road, to minimize influence of the wind.(v)
Before the attempts, wind velocity, temperature, and atmospheric pressure were measured.(vi)
The distance of the road was measured using a measurement system mounted on the vehicle and it was also verified using distance measuring wheel.On the basis of the experimental data, the total rolling resistance coefficient of the vehicle was determined. Table4 shows the coefficients of rolling resistance of the vehicle depending on the tire used during the test. Total resistance includes the tire rolling coefficient, the bearings resistance, and the resistance of the clutch. Verification of the data derived from the simulation and the experiment was conducted usingmMAPE (modifiable Mean Absolute Percentage Error). For velocity simulationmMAPE can be determined by(29)mMAPE=100%1n∑i=1nVr-VsimVmax-Vmin.Figure12 shows the experimental data obtained in a coast down test compared with the results obtained during the computer simulation.Figure 12
Velocity of the vehicle during coast down test.


## 5.4. Model of the Race Circuit

In accordance with the methodology presented in Section3.2, three routes have been described:(i)
Route 1: a street race circuit at theShell Eco-marathon competition in Rotterdam.(ii)
Routes 2 and 3: on the experimental rides track of the Fiat Chrysler Automobiles Poland factory, located in Tychy, in the south of Poland.TheShell Eco-marathon competition is held on a street circuit in Rotterdam. One lap is 1630 meters long. Vehicles have to cover ten laps in no longer than 39 minutes so such requirements impose an average minimum speed of 25 km/h.The trajectory of the vehicle path was developed using satellite images. The vehicle trajectory was determined using the maximum curvature possible that could be obtained on the existing roads.Movement of the vehicle in a curved motion causes the centrifugal force. Centrifugal force creates additional resistance. This resistance is inversely proportional to the radius of the path (8), so it is justified to select the largest radius. In addition, if the centrifugal force exceeds the limit value, it could provide to the side slip of the vehicle. According to drawings inAutoCAD Civil 3D, the lengths of straight sections and the lengths and the radius of arcs were obtained. Figure 13 shows a plan view of the path of the drive.Figure 13
Plan view of the path of the drive on route 1 in Rotterdam.Driving the prototype vehicle on public roads is prohibited and can be dangerous. Vehicle testing was also carried out on the experimental rides track of the Fiat Chrysler Automobiles Poland factory.For the track in Tychy geodetic measurements were taken using theDGPS (Differential Global Positioning System) technology. The measuring points were collected along the axis of the road, with the distance between points amounting to 10 meters. If on the measured route, direction or elevation were changed, the frequency of measurement was increased. Geodetic coordinates were also imported intoAutoCAD Civil 3D and were used to prepare data concerning the route that was imported into the simulation environmentMATLAB-Simulink. Figures 14 and 15 show the path and elevation profile, respectively.Figure 14
Plan view of the path of the drive on route 3 in Tychy.Figure 15
Profile of the path on route 3 in Tychy.Selected routes 2 and 3 have significant differences compared to route 1, especially in the height profile. This set of routes was selected to verify the effect of the route profile on the achieved results and to verify the proposed optimization methods for various road conditions. The routes in Tychy are located in an open area; there are no large buildings, which guarantees stability of wind conditions. The technical condition of the surface on routes 2 and 3 is much worse than on route 1, which has a significant influence on the rolling resistance. Rolling resistance was also identified and introduced into the simulation model. For route 2, thecrr coefficient was 0.0056 and for route 3 was 0.0068.


## 5.5. Verification of the Simulation Model

In order to verify the compatibility of the identified simulation model, a series of numerical simulations and experimental tests were performed. The first experiment consisted of running 300 [m] with a control strategy of 50 [m] with the working motor, and after that the vehicle was moving freely for the next 50 [m]. The motor is turned on three times. Figures16(a) and 16(b) show the velocity of the vehicle and the power of the motor in function of the distance determined during the experiment (red curves) and during the numerical simulations (green curves).Figure 16
The velocity and the power determined with used DT1 drive unit.
(a)
DT1 (b)
DT1Verification of the data derived from the simulation and the experiment was conducted usingmMAPE (29). For motor power simulation,mMAPE can be determined by(30)mMAPEP=100%1n∑i=1nPr-PsimPmax-Pmin,where Pr is the motor power measured during experiment and Psim is the motor power conducted in simulation.As a result of the model validation tests, estimation error of the velocity was on level between 2.35% and 5.34%. However, the motor power error level was greater. It was between 8.02% and 9.98%.


## 5.6. Determination of the Optimal Control Strategy

The length of route 1 and the number of laps are determined by the organizers of the Shell Eco-marathon competition. The length of routes 2 and 3 is derived from the existing road plan, where it was possible to carry out trials in real conditions. The number of laps for routes 2 and 3 was assumed arbitrarily.The travel time for route 1 is based on the rules of the race, which determine the minimum average speed of the vehicle. Based on the average speed and the length of routes 2 and 3, the travel times for these routes were determined. The maximum travel times and the total length of the route, as well as the number and length of one lapD, are summarized in Table 5.Table 5
Travel times and total lengths of the routes.
Route Length of the lapD [m] The number of laps The total length of the routeDc [m] Travel timeTmax [s] Route 1 Rotterdam 1630 10 16 300 2340 Route 2 Tychy 640 5 3200 460 Route 3 Tychy 1000 5 5000 720Evolutionary algorithms which are available in theMATLAB Global Optimization Toolboxwere used. From the MATLAB workspace, the algorithm ran a simulation model in aSimulink environment. During optimization while using evolutionary algorithms, the maximum number of fitness function evaluations Nfval is the product of the population size Npop and number of generations Ngen.It was assumed that the total number of fitness function evaluations was equal to 104; therefore, should this value be exceeded, the calculation would be cancelled. Therefore, the total number of generations Ngen. is dependent on the value ofpopulation sizeNpop.As part of this work phase, the vectorsstU, stV15, and stV10 of the optimum control strategy for an electric vehicle for routes 1, 2, and 3 were determined. During evolutionary optimization, the heuristic method of crossing was used, and the probability of crossing was established on the basis of previous studies at the level of 0.8 [24, 29]. At the stage of reproduction, the elitist succession method was applied where the best adapted two individuals in the current population passed to the next population. An adaptive method of mutation was used. The remaining parameters of the algorithm remain the default for the MATLABGlobal Optimization Toolbox.The first case was the searchstU vector for route 1. Figures 17(a)–17(c) show the vehicle velocityVsim obtained by computer simulation. The graph of thestU control signal is graphically displayed. Figure 17 shows the results where the optimal control strategy stUopt was used, and population was Npop=100 and Ngen=100. The total energy consumption was estimated asrsim=527,633 [km / kWh]. Figure 17(a) shows the vehicle speed and control signalstU during the first lap while Figure 17(b) shows center laps and Figure 17(c) shows the last lap. Figure 17(d) shows the curves of the mean values of the fitness functionsfval, depending on the number of calculations of these functions, during the optimization of the control strategy. There are eight curves in the graph (Figure 17(d)) for each of the Npop parameters used in the algorithm.Figure 17
Vehicle speed and control signalstU for route 1.
(a)
Lap 1 (b)
Laps 2-9 (c)
Lap 10 (d)
Mean values of the fitness functionsWhen analyzing the vehicle velocity chart on the last lap of route 1 (Figure17(c)), it can be seen that, according to the determined control strategy, the propulsion system of the vehicle should no longer be used. The vehicle should use the stored kinetic energy.The analysis of vehicle velocityVsim during the last lap indicates that at the moment of reaching the given distance the speed was still about 2 [m/s]. This demonstrates that the simulated vehicle has kinetic energy that would have to be dispersed by the braking system or recovered by the electric motor. Having such a speed does not mean that there is still possibility for improvement. Lower speed at the end of the route has to cause higher speed of the earlier parts of the route. This speed is the optimum speed due to the chosen criterion. From the analysis of the shape of the fitness function (Figure 17(d)), it can be seen that, for route 1, it is preferable to use a smaller population of Npop=25 and increase the number of Ngen=400. The curve shape for this case is steep and it falls to the minimum value of the fitness function, which may indicate that the algorithm has reached a convergence.Results of the optimization where B type of coding was used were analyzed. The results obtained as the mean of thersim result together with the standard deviation from all the tests were summarized in Table 6.Table 6
The results for route 1.
The type of coding,A The type of coding,B s t U s t V 15 s t V 10 N p o p r s i m [km/kWh] std r s i m [km/kWh] std r s i m [km/kWh] std 25 524,65 1,459 428,45 17,41 515,69 5,4775 50 525,79 1,0053 444,31 22,486 514,92 3,4094 75 524,97 3,9913 438,42 28,173 512,56 8,0682 100 525,89 1,8128 422,32 17,97 511,66 4,9379 125 526,25 1,3862 427,33 39,616 513,41 7,0765 150 524,3 2,979 443,5 47,701 514,86 2,7086 175 525,13 2,511 414,28 27,299 518,38 1,9358 200 525,79 1,9666 423,37 25,776 512,04 7,8357By analyzing the results in Table6, it can be seen that for route 1 better results were obtained (reduced simulated total energy consumption) using A coding compared to B, especially for stV15, whose results are characterized by the highest standard deviation. For routes 2 and 3, the same as for route 1, the procedure was followed. Optimization strategies have been optimized for all variants of control vectorstU, stV15,st10.The results using the control strategystV10 are shown in Figures 18(a)–18(c) (Npop=100 wherersim=351,795 [km / kWh]). Figure 18(d) illustrates the graph of the value of the fitness function.Figure 18
Vehicle speed and control signalstV10 for route 2.
(a)
Lap 1 (b)
Laps 2-4 (c)
Lap 5 (d)
Mean values of the fitness functionsThe results obtained for route 2 are presented in Table7.Table 7
The results for route 2.
The type of coding,A The type of coding,B s t U s t V   15 s t V   10 N p o p r s i m [km/kWh] std r s i m [km/kWh] std r s i m [km/kWh] std 25 341,38 6,6734 342,94 2,3839 348,73 5,3532 50 345,73 4,6359 350,03 4,62 350,6 3,8862 75 345,12 5,5134 349,61 3,3476 354,19 0,67858 100 348,36 4,5332 352,53 2,6117 354,4 1,9295 125 346,32 5,2704 350,17 4,9039 354,83 1,4756 150 345,52 6,3714 353,26 1,4885 353,43 0,71015 175 350,63 1,8728 352,74 4,7141 355,64 0,38022 200 350,73 1,4496 354,0 1,3389 355,69 1,4942For route 2 compared to route 1, the best results were achieved using thestV10 vector. When using vectorsstU and stV15 the results were slightly worse. Having analyzed the size of the population used in the optimization process, it should be noted that for route 2 it is preferable to use higher value of Npop. For all the analyzed cases, the best average results for Npop=200 were obtained. Table 8 summarizes the results for route 3.Table 8
The results for route 3.
The type of coding,A The type of coding,B s t U s t V   15 s t V   10 N p o p r s i m [km/kWh] std r s i m [km/kWh] std r s i m [km/kWh] std 25 280,71 10,173 286,43 7,323 309,94 2,776 50 293,46 6,142 295,81 7,894 308,74 4,399 75 297,11 5,002 292,17 7,367 311,41 1,856 100 298,26 3,231 295,15 7,022 309,57 4,049 125 297,81 2,980 298,66 12,977 311,68 3,445 150 294,53 3,346 302,17 9,623 312,88 1,249 175 292,34 2,317 306,21 6,816 314,81 1,261 200 296,86 2,337 305,31 11,737 310,89 4,880For route 3, the average results were reported for a larger population ofNpop with B coding. In the case of using vector with A type coding, no characteristic change was noted. As expected, comparing the standard deviation between algorithms using B coding, the standard deviation is less for the less variance range.Evolutionary algorithms are not complete search algorithms; therefore, it is not possible to determine whether the solution found is best for the chosen criterion and the constraints assumed.


## 5.7. Applying the Optimal Control Strategy in Real Conditions

Electric vehicle energy consumption optimization was finally verified in real conditions at the European edition of theShell Eco-marathon held on 16-18 May 2014 in Rotterdam.Optimal control strategystU was displayed by a prototype visualization system. Visualization unit is installed in the steering wheel of the vehicle. The visualization system shows the basic data such as instantaneous speed, average speed, travel time, and distance. In addition, the visualization system provides information on which part of the track the driver should press the acceleration button to execute the optimal control strategy. Measurement of the position on the track is carried out by the odometrical measuring system. Due to the possibility of accumulating errors in the measurement of the road, the driver will reset the measuring system every time he/she exceeds the lap line.Table9 presents the results of all the measurements conducted during theShell Eco-marathon. The columns in the Table 9 contain(i)
rSEm [km/kWh]: the result recorded during the race by the measuring device provided by the organizer; the results are available on the Shell Eco-marathon official website;(ii)
rr [kW/kWh]: the result based on the indications of a prototype measuring system permanently installed in the vehicle;(iii)
rsim: the result estimated by computer simulation using recorded real control signalUr;(iv)
ropt: the result obtained by computer simulation using the optimal control strategyUopt.Table 9
Comparison of simulation results with experimental tests.
Number of attempts r S E – m r r r s i m r o p t e m e i e s t Attempt 1 365,4 360,8 - - 1,25 - - Attempt 2 481,3 508,2 514,8 526,5 5,58 1,29 3,6 Attempt 3 487,3 512,7 519,8 526,5 5,21 1,39 2,69 Attempt 4 Failed attempt, damaged tires during 1st lapTable9 also includes the following:(i)
Measurement error:(31)em=rSEm-rrrSEm100%.(ii)
Identification error:(32)ei=rr-rsimrr100%.(iii)
Strategy implementation error:(33)est=rr-roptrr100%.There were four measurement trials available. During run 1, a much weaker result was obtained than estimated by simulation. The driver performed run 1 without a defined control strategy. As part of the experiment, the driver and the support team were asked to select the optimal tires for the vehicle and to implement the proposed control strategy. During run 1, only the travel time of individual laps was checked. According to the choice of the team in run 1, tires W2 were selected (Table4). During this test the worst result was obtained about 24% and 25% worse than in the other trials. In subsequent runs 2 and 3, the vehicle was equipped with the best W1 tires (Table 4). The driver has implemented an algorithm-based optimal control strategy. Strategy information was provided by the visualization system. During run 3 the best result was recorded. This result was taken into account in the classification of competitions, where the team placed in 12 positions out of 28 classified teams. Trials 4, due to damage to the tires, were not successful.Measurement error was from 1.25% to 5.58%, depending on the sample. The measurement system installed permanently in the vehicle is a prototype system. This system has been calibrated using available laboratory equipment. Measurement errors indicate that other measuring devices of the higher accuracy class should be used.Based on the measurement data obtained during test drives, the model was tuned, thus reducing the identification error. The competition took place in stable weather; the divergence of results was small and was about 1.29% for sample 2 and 1.39% for sample 3.The difference between the results determined during numerical simulation and that obtained in the real run was 3.6% for sample 2 and 2.69% for sample 3. The difference at this level was considered to be small. Searching for optimal control strategy has not taken into account possible "errors" which could be made by the driver and which result from the current situation on the route, for example, selecting the path and using the route by other participants in the race. Including driver behaviour in the simulation model should be the subject of further work.Figure19 shows the vehicle velocity together with the control signals recorded during run 2. The black curves represent the real vehicle velocityVr and the control signal. Green curves are the optimum control strategystU and the optimal velocityVopt determined by the simulation. The graph shows the entire distance travelled (10 laps). The analysis of the graph shows that the driver did not perform the exact planned strategy. Errors in the implementation of the strategy are related to the disruptions. During the competition, the driver is not alone on the track; in certain situations the driver cannot overtake other participants and cannot take the optimal path. In addition, weather condition such as wind has influence on the motion resistance. In such cases, the driver has to modify the optimal control strategy.Figure 19
Comparison of speed and control signal during attempt 2.In Figure20, the velocity recorded during the raceVr (black) and the velocity obtained by the computer simulationVsim (red) were shown.Vsim was determined by computer simulation by using real control signal. The purpose of this simulation was to investigate the correctness of identification of the simulation model and to identify the identification error.Figure 20
Comparison of speed registered in real conditions with the speed determined by computer simulation.Compared to theShell Eco-marathon competitions held in 2012 and 2013, our performance has improved. The results in these years were equal to 425 [km / kWh] and 454 [km / kWh], respectively. In 2012 and 2013, the control strategy was chosen on the basis of the simulation model, but the model was not adequately identified, and the control strategy had to be modified. In previous editions of the race, the implementation of optimal control strategy was not supported by the visualization system, which would inform the driver of the need to take specific actions (pressing or releasing the accelerator button).


## 6. Conclusion

The presented method was developed and tested for a prototype electric racing car that has been designed for theShell Eco-marathon race. As it was shown the results of the numerical simulations were very similar to the outcomes of the real races in the context of the final result of energy consumption and the profile of velocity of the car during the race. The method has been developed and applied since 2012 and led in 2014 to achieving the optimal result for this type of the vehicle.The optimization method was evaluated. The results achieved during the race attempts where the optimal control strategy was applied (481,3 and 487,3 km/kWh) were compared to the result where the strategy was proposed by the experts (365.4 km/kWh), who already had contact with the results of similar optimization models. Optimal attempts were better by over 30%, despite errors in the implementation of the strategy by the driver. In the case of experts who did not have such experiences, the differences were even greater.Further improvement of the control strategy is impossible for this type of a vehicle and the best way to reduce energy consumption is by improving the mechanical design. It was the reason to develop new racing cars for UrbanConcept Battery Electric and Hydrogen competitions. The developed research methodology for vehicle optimization and control strategies has been developed and implemented in the Shell Eco-marathon in 2015, which has improved the score to 504 km / kWh. However, the research has shown that the potential for improved performance of theMuSHELLka vehicle has been exhausted. No major structural changes in the vehicle can significantly improve the result, so in the following years, the race was abandoned and the focus was on the development of theHydroGENIUS, the electric vehicle of UrbanConcept class with fuel cell stack as power source, which eventually scored vice Championship inSEM 2016. In the further work, it is planned to improve the optimization results by using the presented methodology in the case of the models of newly constructed vehicles. The development of the optimization method should take into consideration the uncertainty resulting from the given velocity strategy considering the real conditions of the race. It is also possible to use this method for controlling vehicles in traffic. To do this, it is necessary to completely change the way of the route modelling and the driving conditions modelling, as well as to develop a new optimization method for obtaining results in real time.



---

*Source: 3614025-2018-06-07.xml*