SUMO - Simulation of Urban MObility
PHEMCEP.cpp
Go to the documentation of this file.
1 /****************************************************************************/
11 // Helper class for PHEM Light, holds a specific CEP for a PHEM emission class
12 /****************************************************************************/
13 // SUMO, Simulation of Urban MObility; see http://sumo-sim.org/
14 // Copyright (C) 2013-2014 DLR (http://www.dlr.de/) and contributors
15 /****************************************************************************/
16 //
17 // This file is part of SUMO.
18 // SUMO is free software: you can redistribute it and/or modify
19 // it under the terms of the GNU General Public License as published by
20 // the Free Software Foundation, either version 3 of the License, or
21 // (at your option) any later version.
22 //
23 /****************************************************************************/
24 
25 // ===========================================================================
26 // included modules
27 // ===========================================================================
28 #ifdef _MSC_VER
29 #include <windows_config.h>
30 #else
31 #include <config.h>
32 #endif
33 
34 #include <cmath>
35 #include <string>
36 #include <utils/common/StdDefs.h>
39 #include "PHEMCEP.h"
40 
41 #ifdef CHECK_MEMORY_LEAKS
42 #include <foreign/nvwa/debug_new.h>
43 #endif // CHECK_MEMORY_LEAKS
44 
45 
46 // ===========================================================================
47 // method definitions
48 // ===========================================================================
49 PHEMCEP::PHEMCEP(bool heavyVehicel, SUMOEmissionClass emissionClass, const std::string& emissionClassIdentifier,
50  double vehicleMass, double vehicleLoading, double vehicleMassRot,
51  double crossArea, double cWValue,
52  double f0, double f1, double f2, double f3, double f4,
53  double ratedPower, double pNormV0, double pNormP0, double pNormV1,
54  double pNormP1, std::string vehicelFuelType,
55  const std::vector< std::vector<double> >& matrixFC,
56  const std::vector<std::string>& headerLinePollutants,
57  const std::vector< std::vector<double> >& matrixPollutants,
58  const std::vector< std::vector<double> > matrixSpeedRotational) {
59  _emissionClass = emissionClass;
60  _resistanceF0 = f0;
61  _resistanceF1 = f1;
62  _resistanceF2 = f2;
63  _resistanceF3 = f3;
64  _resistanceF4 = f4;
65  _cwValue = cWValue;
66  _crossSectionalArea = crossArea;
67  _massVehicle = vehicleMass;
68  _vehicleLoading = vehicleLoading;
69  _massRot = vehicleMassRot;
70  _ratedPower = ratedPower;
71  _vehicleFuelType = vehicelFuelType;
72 
73  _pNormV0 = pNormV0 / 3.6;
74  _pNormP0 = pNormP0;
75  _pNormV1 = pNormV1 / 3.6;
76  _pNormP1 = pNormP1;
77 
78  std::vector<std::string> pollutantIdentifier;
79  std::vector< std::vector<double> > pollutantMeasures;
80 
81  // init pollutant identifiers
82  for (int i = 0; i < (int)headerLinePollutants.size(); i++) {
83  pollutantIdentifier.push_back(headerLinePollutants[i]);
84  } // end for
85 
86  // get size of powerPatterns
87  _sizeOfPatternFC = (int)matrixFC.size();
88  _sizeOfPatternPollutants = (int)matrixPollutants.size();
89 
90  // initialize measures
91  for (int i = 0; i < (int)headerLinePollutants.size(); i++) {
92  pollutantMeasures.push_back(std::vector<double>());
93  } // end for
94 
95  // looping through matrix and assigning values for speed rotational table
96  _speedCurveRotational.clear();
98  for (int i = 0; i < (int)matrixSpeedRotational.size(); i++) {
99  if (matrixSpeedRotational[i].size() != 2) {
100  throw InvalidArgument("Error loading vehicle file for: " + emissionClassIdentifier);
101  }
102 
103  _speedPatternRotational.push_back(matrixSpeedRotational[i][0] / 3.6);
104  _speedCurveRotational.push_back(matrixSpeedRotational[i][1]);
105 
106  } // end for
107 
108  // looping through matrix and assigning values for Fuel consumption
109  _cepCurveFC.clear();
110  for (int i = 0; i < (int)matrixFC.size(); i++) {
111  if (matrixFC[i].size() != 2) {
112  throw InvalidArgument("Error loading vehicle file for: " + emissionClassIdentifier);
113  }
114 
115  _powerPatternFC.push_back(matrixFC[i][0] * _ratedPower);
116  _cepCurveFC.push_back(matrixFC[i][1]);
117 
118  } // end for
119 
120 
121  // looping through matrix and assigning values for pollutants
122  double normalizingPower = 0;
123 
124  if (heavyVehicel) {
125  normalizingPower = _ratedPower;
126  } else {
128  } // end if
129 
130  const int headerCount = (int)headerLinePollutants.size();
131  for (int i = 0; i < (int)matrixPollutants.size(); i++) {
132  for (int j = 0; j < (int)matrixPollutants[i].size(); j++) {
133  if ((int)matrixPollutants[i].size() != headerCount + 1) {
134  return;
135  }
136 
137  if (j == 0) {
138  _powerPatternPollutants.push_back(matrixPollutants[i][j] * normalizingPower);
139  } else {
140  pollutantMeasures[j - 1].push_back(matrixPollutants[i][j]);
141  } // end if
142  } // end for
143  } // end for
144 
145  for (int i = 0; i < headerCount; i++) {
146  _cepCurvePollutants.insert(pollutantIdentifier[i], pollutantMeasures[i]);
147  } // end for
148 
149 } // end of Cep
150 
151 
153  // free power pattern
154  _powerPatternFC.clear();
155  _powerPatternPollutants.clear();
156  _cepCurveFC.clear();
157  _speedCurveRotational.clear();
158  _speedPatternRotational.clear();
159 } // end of ~Cep
160 
161 
162 double
163 PHEMCEP::CalcPower(double v, double a, double slope) const {
164  const double rotFactor = GetRotationalCoeffecient(v);
165  double power = (_massVehicle + _vehicleLoading) * GRAVITY_CONST * (_resistanceF0 + _resistanceF1 * v + _resistanceF4 * pow(v, 4)) * v;
166  power += (_crossSectionalArea * _cwValue * AIR_DENSITY_CONST / 2) * pow(v, 3);
167  power += (_massVehicle * rotFactor + _massRot + _vehicleLoading) * a * v;
168  power += (_massVehicle + _vehicleLoading) * slope * 0.01 * v;
169  return power / 950.;
170 }
171 
172 
173 double
174 PHEMCEP::GetMaxAccel(double v, double a, double gradient) const {
175  UNUSED_PARAMETER(a);
176  const double pMaxForAcc = GetPMaxNorm(v) * _ratedPower - PHEMCEP::CalcPower(v, 0, gradient);
177  return (pMaxForAcc * 1000) / ((_massVehicle * GetRotationalCoeffecient(v) + _massRot + _vehicleLoading) * v);
178 }
179 
180 
181 double
182 PHEMCEP::GetEmission(const std::string& pollutant, double power) const {
183  std::vector<double> emissionCurve;
184  std::vector<double> powerPattern;
185 
186  if (pollutant == "FC") {
187  emissionCurve = _cepCurveFC;
188  powerPattern = _powerPatternFC;
189  } else {
190  if (!_cepCurvePollutants.hasString(pollutant)) {
191  throw InvalidArgument("Emission pollutant " + pollutant + " not found!");
192  }
193 
194  emissionCurve = _cepCurvePollutants.get(pollutant);
195  powerPattern = _powerPatternPollutants;
196  } // end if
197 
198 
199 
200  if (emissionCurve.size() == 0) {
201  throw InvalidArgument("Empty emission curve for " + pollutant + " found!");
202  }
203 
204  if (emissionCurve.size() == 1) {
205  return emissionCurve[0];
206  }
207 
208  // in case that the demanded power is smaller than the first entry (smallest) in the power pattern the first two entries are extrapolated
209  if (power <= powerPattern.front()) {
210  double calcEmission = PHEMCEP::Interpolate(power, powerPattern[0], powerPattern[1], emissionCurve[0], emissionCurve[1]);
211 
212  if (calcEmission < 0) {
213  return 0;
214  } else {
215  return calcEmission;
216  }
217 
218  } // end if
219 
220  // if power bigger than all entries in power pattern the last two values are linearly extrapolated
221  if (power >= powerPattern.back()) {
222  return PHEMCEP::Interpolate(power, powerPattern[powerPattern.size() - 2], powerPattern.back(), emissionCurve[emissionCurve.size() - 2], emissionCurve.back());
223  } // end if
224 
225  // bisection search to find correct position in power pattern
226  int upperIndex;
227  int lowerIndex;
228 
229  PHEMCEP::FindLowerUpperInPattern(lowerIndex, upperIndex, powerPattern, power);
230 
231  return PHEMCEP::Interpolate(power, powerPattern[lowerIndex], powerPattern[upperIndex], emissionCurve[lowerIndex], emissionCurve[upperIndex]);
232 
233 } // end of GetEmission
234 
235 
236 double
237 PHEMCEP::Interpolate(double px, double p1, double p2, double e1, double e2) const {
238  if (p2 == p1) {
239  return e1;
240  }
241  return e1 + (px - p1) / (p2 - p1) * (e2 - e1);
242 } // end of Interpolate
243 
244 
245 double
247  int upperIndex;
248  int lowerIndex;
249 
250  PHEMCEP::FindLowerUpperInPattern(lowerIndex, upperIndex, _speedPatternRotational, speed);
251 
252  return PHEMCEP::Interpolate(speed,
253  _speedPatternRotational[lowerIndex],
254  _speedPatternRotational[upperIndex],
255  _speedCurveRotational[lowerIndex],
256  _speedCurveRotational[upperIndex]);
257 } // end of GetRotationalCoeffecient
258 
259 void
260 PHEMCEP::FindLowerUpperInPattern(int& lowerIndex, int& upperIndex, std::vector<double> pattern, double value) const {
261  if (value <= pattern.front()) {
262  lowerIndex = 0;
263  upperIndex = 0;
264  return;
265 
266  } // end if
267 
268  if (value >= pattern.back()) {
269  lowerIndex = (int)pattern.size() - 1;
270  upperIndex = (int)pattern.size() - 1;
271  return;
272  } // end if
273 
274  // bisection search to find correct position in power pattern
275  int middleIndex = ((int)pattern.size() - 1) / 2;
276  upperIndex = (int)pattern.size() - 1;
277  lowerIndex = 0;
278 
279  while (upperIndex - lowerIndex > 1) {
280  if (pattern[middleIndex] == value) {
281  lowerIndex = middleIndex;
282  upperIndex = middleIndex;
283  return;
284  } else if (pattern[middleIndex] < value) {
285  lowerIndex = middleIndex;
286  middleIndex = (upperIndex - lowerIndex) / 2 + lowerIndex;
287  } else {
288  upperIndex = middleIndex;
289  middleIndex = (upperIndex - lowerIndex) / 2 + lowerIndex;
290  } // end if
291  } // end while
292 
293  if (pattern[lowerIndex] <= value && value < pattern[upperIndex]) {
294  return;
295  } else {
296  throw ProcessError("Error during calculation of position in pattern!");
297  }
298 } // end of FindLowerUpperInPattern
299 
300 
301 double PHEMCEP::GetPMaxNorm(double speed) const {
302  // Linear function between v0 and v1, constant elsewhere
303  if (speed <= _pNormV0) {
304  return _pNormP0;
305  } else if (speed >= _pNormV1) {
306  return _pNormP1;
307  } else {
309  }
310 } // end of GetPMaxNorm
311 
312 /****************************************************************************/
std::vector< double > _powerPatternFC
Definition: PHEMCEP.h:288
bool hasString(const std::string &str) const
double GetEmission(const std::string &pollutantIdentifier, double power) const
Returns a emission measure for power[kW] level.
Definition: PHEMCEP.cpp:182
std::vector< double > _speedCurveRotational
Definition: PHEMCEP.h:298
double _massVehicle
vehicle mass
Definition: PHEMCEP.h:268
double _resistanceF1
Rolling resistance f1.
Definition: PHEMCEP.h:256
std::string _vehicleFuelType
Definition: PHEMCEP.h:300
const double NORMALIZING_SPEED
Definition: PHEMConstants.h:28
double _pNormP0
Step functions parameter for maximum rated power.
Definition: PHEMCEP.h:278
~PHEMCEP()
Destructor.
Definition: PHEMCEP.cpp:152
double _pNormV1
Step functions parameter for maximum rated power.
Definition: PHEMCEP.h:280
double _resistanceF3
Rolling resistance f3.
Definition: PHEMCEP.h:260
double _massRot
rotational mass of vehicle
Definition: PHEMCEP.h:272
#define UNUSED_PARAMETER(x)
Definition: StdDefs.h:39
double GetMaxAccel(double v, double a, double gradient) const
Returns the maximum accelaration for a vehicle at state v,a, slope and loading.
Definition: PHEMCEP.cpp:174
int _sizeOfPatternFC
Definition: PHEMCEP.h:284
double _pNormV0
Step functions parameter for maximum rated power.
Definition: PHEMCEP.h:276
const double NORMALIZING_ACCELARATION
Definition: PHEMConstants.h:29
double _crossSectionalArea
crosssectional area of vehicle
Definition: PHEMCEP.h:266
double _vehicleLoading
vehicle loading
Definition: PHEMCEP.h:270
void insert(const std::string str, const T key, bool checkDuplicates=true)
double _pNormP1
Step functions parameter for maximum rated power.
Definition: PHEMCEP.h:282
std::vector< double > _speedPatternRotational
Definition: PHEMCEP.h:296
const double GRAVITY_CONST
Definition: PHEMConstants.h:25
double GetRotationalCoeffecient(double speed) const
Calculates rotational index for speed.
Definition: PHEMCEP.cpp:246
double _resistanceF4
Rolling resistance f4.
Definition: PHEMCEP.h:262
PHEMCEP(bool heavyVehicel, SUMOEmissionClass emissionClass, const std::string &emissionClassIdentifier, double vehicleMass, double vehicleLoading, double vehicleMassRot, double crossArea, double cWValue, double f0, double f1, double f2, double f3, double f4, double ratedPower, double pNormV0, double pNormP0, double pNormV1, double pNormP1, std::string vehicelFuelType, const std::vector< std::vector< double > > &matrixFC, const std::vector< std::string > &headerLinePollutants, const std::vector< std::vector< double > > &matrixPollutants, const std::vector< std::vector< double > > matrixSpeedRotational)
Definition: PHEMCEP.cpp:49
double _resistanceF0
Rolling resistance f0.
Definition: PHEMCEP.h:254
std::vector< double > _cepCurveFC
Definition: PHEMCEP.h:292
double GetPMaxNorm(double speed) const
Calculates maximum available rated power for speed.
Definition: PHEMCEP.cpp:301
int _sizeOfPatternPollutants
Definition: PHEMCEP.h:286
std::vector< double > _powerPatternPollutants
Definition: PHEMCEP.h:290
const double AIR_DENSITY_CONST
Definition: PHEMConstants.h:26
double Interpolate(double px, double p1, double p2, double e1, double e2) const
Interpolates emission linearly between two known power-emission pairs.
Definition: PHEMCEP.cpp:237
double _ratedPower
rated power of vehicle
Definition: PHEMCEP.h:274
double _cwValue
Cw value.
Definition: PHEMCEP.h:264
double CalcPower(double v, double a, double slope) const
Returns the power of used for a vehicle at state v,a, slope and loading.
Definition: PHEMCEP.cpp:163
StringBijection< std::vector< double > > _cepCurvePollutants
Definition: PHEMCEP.h:294
T get(const std::string &str) const
void FindLowerUpperInPattern(int &lowerIndex, int &upperIndex, std::vector< double > pattern, double value) const
Finds bounding upper and lower index in pattern for value.
Definition: PHEMCEP.cpp:260
SUMOEmissionClass _emissionClass
PHEM emission class of vehicle.
Definition: PHEMCEP.h:252
double _resistanceF2
Rolling resistance f2.
Definition: PHEMCEP.h:258