Line data Source code
1 : /****************************************************************************/
2 : // Eclipse SUMO, Simulation of Urban MObility; see https://eclipse.dev/sumo
3 : // Copyright (C) 2002-2024 German Aerospace Center (DLR) and others.
4 : // This program and the accompanying materials are made available under the
5 : // terms of the Eclipse Public License 2.0 which is available at
6 : // https://www.eclipse.org/legal/epl-2.0/
7 : // This Source Code may also be made available under the following Secondary
8 : // Licenses when the conditions for such availability set forth in the Eclipse
9 : // Public License 2.0 are satisfied: GNU General Public License, version 2
10 : // or later which is available at
11 : // https://www.gnu.org/licenses/old-licenses/gpl-2.0-standalone.html
12 : // SPDX-License-Identifier: EPL-2.0 OR GPL-2.0-or-later
13 : /****************************************************************************/
14 : /// @file HelpersMMPEVEM.cpp
15 : /// @author Kevin Badalian (badalian_k@mmp.rwth-aachen.de)
16 : /// @date 2021-10
17 : ///
18 : // The MMP's emission model for electric vehicles.
19 : // If you use this model for academic research, you are highly encouraged to
20 : // cite our paper "Accurate physics-based modeling of electric vehicle energy
21 : // consumption in the SUMO traffic microsimulator"
22 : // (DOI: 10.1109/ITSC48978.2021.9564463).
23 : // Teaching and Research Area Mechatronics in Mobile Propulsion (MMP), RWTH Aachen
24 : /****************************************************************************/
25 : #include <config.h>
26 :
27 : #include <utils/common/SUMOTime.h>
28 : #include <utils/common/StringUtils.h>
29 : #include <utils/geom/GeomHelper.h>
30 : #include <utils/emissions/CharacteristicMap.h>
31 : #include <utils/emissions/HelpersMMPEVEM.h>
32 :
33 :
34 : // ===========================================================================
35 : // method definitions
36 : // ===========================================================================
37 : /**
38 : * \brief Compute the power consumption of an EV powertrain in a certain state.
39 : *
40 : * The model assumes that the driver recuperates as much energy as possible,
41 : * meaning that he perfectly employs the mechanical brakes such that only what
42 : * lies beyond the motor's capabilities is dissipated. If the new state is
43 : * invalid, that is the demanded acceleration exceeds what the electric motor
44 : * can manage or the power loss map is not defined for a given point, the torque
45 : * and/or power are capped to the motor's limits or the power loss is set to
46 : * zero.
47 : *
48 : * \param[in] m Vehicle mass [kg]
49 : * \param[in] r_wheel Wheel radius [m]
50 : * \param[in] Theta Moment of inertia [kg*m^2]
51 : * \param[in] c_rr Rolling resistance coefficient [1]
52 : * \param[in] c_d Drag coefficient [1]
53 : * \param[in] A_front Cross-sectional area of the front of the car [m^2]
54 : * \param[in] i_gear Gear ratio (i_gear = n_motor/n_wheel) [1]
55 : * \param[in] eta_gear Gear efficiency [1]
56 : * \param[in] M_max Maximum torque of the EM [Nm]
57 : * \param[in] P_max Maximum power of the EM [W]
58 : * \param[in] M_recup_max Maximum recuperation torque [Nm]
59 : * \param[in] P_recup_max Maximum recuperation power [W]
60 : * \param[in] R_battery Internal battery resistance [Ohm]
61 : * \param[in] U_battery_0 Nominal battery voltage [V]
62 : * \param[in] P_const Constant power consumption of ancillary devices [W]
63 : * \param[in] ref_powerLossMap Power loss map of the EM + inverter
64 : * \param[in] dt Simulation timestep [s]
65 : * \param[in] v Vehicle speed at the end of a timestep [m/s]
66 : * \param[in] a Constant acceleration during a timestep [m/s^2]
67 : * \param[in] alpha Constant incline during a timestep [deg]
68 : * \param[out] ref_powerConsumption Power consumption during the last timestep
69 : * [W]
70 : * \returns true if the new state is valid, else false
71 : */
72 5280 : bool calcPowerConsumption(double m, double r_wheel, double Theta, double c_rr,
73 : double c_d, double A_front, double i_gear, double eta_gear, double M_max,
74 : double P_max, double M_recup_max, double P_recup_max, double R_battery,
75 : double U_battery_0, double P_const,
76 : const CharacteristicMap& ref_powerLossMap, double dt, double v, double a,
77 : double alpha, double& ref_powerConsumption) {
78 : const double EPS = 1e-6;
79 : const double RHO_AIR = 1.204; // Air density [kg/m^3]
80 : bool b_stateValid = true;
81 :
82 : // Mass factor [1]
83 5280 : const double e_i = 1.0 + Theta / (m * r_wheel * r_wheel);
84 : // Average speed during the previous timestep
85 5280 : const double v_mean = v - 0.5 * a * dt;
86 :
87 : // Force required for the desired acceleration [N]
88 5280 : double F_a = m * a * e_i;
89 : // Grade resistance [N]
90 5280 : double F_gr = m * GRAVITY * std::sin(DEG2RAD(alpha));
91 : // Rolling resistance [N]
92 5280 : double F_rr = m * GRAVITY * std::cos(DEG2RAD(alpha)) * c_rr;
93 5280 : if (std::abs(v_mean) <= EPS) {
94 : F_rr = 0;
95 : }
96 : // Drag [N]
97 5280 : double F_d = 0.5 * c_d * A_front * RHO_AIR * v_mean * v_mean;
98 : // Tractive force [N]
99 5280 : const double F_tractive = F_a + F_gr + F_rr + F_d;
100 :
101 : // Speed of the motor [rpm]
102 5280 : const double n_motor = v_mean / (2 * M_PI * r_wheel) * 60 * i_gear;
103 : // Angular velocity of the motor [1/s]
104 5280 : double omega_motor = 2 * M_PI * n_motor / 60;
105 5280 : if (omega_motor == 0) {
106 : omega_motor = EPS;
107 : }
108 :
109 : // Torque at the motor [Nm]. Please note that this model, like most real EVs,
110 : // utilizes the EM to hold the vehicle on an incline rather than the brakes.
111 5280 : double M_motor = F_tractive * r_wheel / i_gear;
112 5280 : if (F_tractive < 0) {
113 564 : M_motor *= eta_gear;
114 : } else {
115 4716 : M_motor /= eta_gear;
116 : }
117 : // Engine power [W]
118 5280 : double P_motor = M_motor * omega_motor;
119 :
120 : // Cap torque or power if necessary
121 : // Accelerating
122 5280 : if (M_motor >= 0) {
123 4716 : if (M_motor > M_max) {
124 : M_motor = M_max;
125 0 : P_motor = M_motor * omega_motor;
126 : b_stateValid = false;
127 : }
128 4716 : if (P_motor > P_max) {
129 : P_motor = P_max;
130 0 : M_motor = P_motor / omega_motor;
131 : b_stateValid = false;
132 : }
133 : }
134 : // Recuperating
135 : else {
136 : // Even when capping, the state is still valid here because the extra energy
137 : // is assumed to go to the brakes
138 564 : if (M_motor < -M_recup_max) {
139 : M_motor = -M_recup_max;
140 420 : P_motor = M_motor * omega_motor;
141 : }
142 564 : if (P_motor < -P_recup_max) {
143 : P_motor = -P_recup_max;
144 0 : M_motor = P_motor / omega_motor;
145 : }
146 : }
147 :
148 : // Power loss in the electric motor + inverter [W]
149 : double P_loss_motor =
150 5280 : ref_powerLossMap.eval(std::vector<double> {n_motor, M_motor})[0];
151 5280 : if (std::isnan(P_loss_motor)) {
152 : P_loss_motor = 0.0;
153 : b_stateValid = false;
154 : }
155 :
156 : // Power demand at the battery [W]
157 5280 : double P_battery = P_motor + P_loss_motor + P_const;
158 : // Total power demand (including the power loss in the battery) [W]
159 5280 : double P_total = (U_battery_0 * U_battery_0) / (2.0 * R_battery)
160 5280 : - U_battery_0 * std::sqrt((U_battery_0 * U_battery_0
161 5280 : - 4.0 * R_battery * P_battery) / (4.0 * R_battery * R_battery));
162 :
163 5280 : ref_powerConsumption = P_total;
164 5280 : return b_stateValid;
165 : }
166 :
167 :
168 :
169 :
170 : /**
171 : * \brief Constructor
172 : */
173 45820 : HelpersMMPEVEM::HelpersMMPEVEM()
174 45820 : : PollutantsInterface::Helper("MMPEVEM", MMPEVEM_BASE, MMPEVEM_BASE + 1)
175 45820 : { }
176 :
177 :
178 :
179 :
180 :
181 : /**
182 : * \brief Compute the amount of emitted pollutants for an emission class in a
183 : * given state.
184 : *
185 : * This method returns 0 for all emission types but electric power consumption.
186 : *
187 : * \param[in] c An emission class
188 : * \param[in] e An emission type
189 : * \param[in] v Current vehicle velocity [m/s]
190 : * \param[in] a Current acceleration of the vehicle [m/s^2]
191 : * \param[in] slope Slope of the road at the vehicle's current position [deg]
192 : * \param[in] ptr_energyParams Vehicle parameters
193 : *
194 : * \returns The electric power consumption [Wh/s] or 0 for all other emission
195 : * types
196 : */
197 27120 : double HelpersMMPEVEM::compute(const SUMOEmissionClass /* c */,
198 : const PollutantsInterface::EmissionType e, const double v, const double a,
199 : const double slope, const EnergyParams* ptr_energyParams) const {
200 27120 : if (e != PollutantsInterface::ELEC) {
201 : return 0.0;
202 : }
203 :
204 : // Extract all required parameters
205 : // Vehicle mass
206 5280 : const double m = ptr_energyParams->getDouble(SUMO_ATTR_MASS) + ptr_energyParams->getDoubleOptional(SUMO_ATTR_LOADING, 0.);
207 : // Wheel radius
208 5280 : const double r_wheel = ptr_energyParams->getDouble(SUMO_ATTR_WHEELRADIUS);
209 : // Internal moment of inertia
210 : const double Theta
211 5280 : = ptr_energyParams->getDouble(SUMO_ATTR_INTERNALMOMENTOFINERTIA);
212 : // Rolling resistance coefficient
213 : const double c_rr
214 5280 : = ptr_energyParams->getDouble(SUMO_ATTR_ROLLDRAGCOEFFICIENT);
215 : // Air drag coefficient
216 5280 : const double c_d = ptr_energyParams->getDouble(SUMO_ATTR_AIRDRAGCOEFFICIENT);
217 : // Cross-sectional area of the front of the car
218 : const double A_front
219 5280 : = ptr_energyParams->getDouble(SUMO_ATTR_FRONTSURFACEAREA);
220 : // Gear ratio
221 5280 : const double i_gear = ptr_energyParams->getDouble(SUMO_ATTR_GEARRATIO);
222 : // Gearbox efficiency
223 5280 : const double eta_gear = ptr_energyParams->getDouble(SUMO_ATTR_GEAREFFICIENCY);
224 : // Maximum torque
225 5280 : const double M_max = ptr_energyParams->getDouble(SUMO_ATTR_MAXIMUMTORQUE);
226 : // Maximum power
227 5280 : const double P_max = ptr_energyParams->getDouble(SUMO_ATTR_MAXIMUMPOWER);
228 : // Maximum recuperation torque
229 : const double M_recup_max
230 5280 : = ptr_energyParams->getDouble(SUMO_ATTR_MAXIMUMRECUPERATIONTORQUE);
231 : // Maximum recuperation power
232 : const double P_recup_max
233 5280 : = ptr_energyParams->getDouble(SUMO_ATTR_MAXIMUMRECUPERATIONPOWER);
234 : // Internal battery resistance
235 : const double R_battery
236 5280 : = ptr_energyParams->getDouble(SUMO_ATTR_INTERNALBATTERYRESISTANCE);
237 : // Nominal battery voltage
238 : const double U_battery_0
239 5280 : = ptr_energyParams->getDouble(SUMO_ATTR_NOMINALBATTERYVOLTAGE);
240 : // Constant power consumption
241 : const double P_const
242 5280 : = ptr_energyParams->getDouble(SUMO_ATTR_CONSTANTPOWERINTAKE);
243 : // Power loss map
244 : const CharacteristicMap& ref_powerLossMap
245 5280 : = ptr_energyParams->getCharacteristicMap(SUMO_ATTR_POWERLOSSMAP);
246 :
247 5280 : double P = 0.0; // [W]
248 5280 : bool b_stateValid = calcPowerConsumption(m, r_wheel, Theta, c_rr, c_d,
249 : A_front, i_gear, eta_gear, M_max, P_max, M_recup_max, P_recup_max,
250 5280 : R_battery, U_battery_0, P_const, ref_powerLossMap, TS, v, a, slope, P);
251 5280 : P /= 3600.0; // [Wh/s]
252 :
253 5280 : if (b_stateValid) {
254 : return P;
255 : } else {
256 0 : return std::nan("");
257 : }
258 : }
|
