Line data Source code
1 : /****************************************************************************/
2 : // Eclipse SUMO, Simulation of Urban MObility; see https://eclipse.dev/sumo
3 : // Copyright (C) 2001-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 ROMAAssignments.cpp
15 : /// @author Yun-Pang Floetteroed
16 : /// @author Laura Bieker
17 : /// @author Michael Behrisch
18 : /// @date Feb 2013
19 : ///
20 : // Assignment methods
21 : /****************************************************************************/
22 : #include <config.h>
23 :
24 : #include <vector>
25 : #include <algorithm>
26 : #include <utils/common/SUMOTime.h>
27 : #include <utils/distribution/Distribution_Points.h>
28 : #include <utils/router/RouteCostCalculator.h>
29 : #include <utils/router/SUMOAbstractRouter.h>
30 : #include <router/ROEdge.h>
31 : #include <router/RONet.h>
32 : #include <router/RORoute.h>
33 : #include <od/ODMatrix.h>
34 : #include "ROMAEdge.h"
35 : #include "ROMAAssignments.h"
36 :
37 :
38 : // ===========================================================================
39 : // static member variables
40 : // ===========================================================================
41 : std::map<const ROEdge* const, double> ROMAAssignments::myPenalties;
42 :
43 :
44 : // ===========================================================================
45 : // method definitions
46 : // ===========================================================================
47 :
48 59 : ROMAAssignments::ROMAAssignments(const SUMOTime begin, const SUMOTime end, const bool additiveTraffic,
49 : const double adaptionFactor, const int maxAlternatives, const bool defaultCapacities,
50 : RONet& net, ODMatrix& matrix, SUMOAbstractRouter<ROEdge, ROVehicle>& router,
51 59 : OutputDevice* netloadOutput)
52 59 : : myBegin(begin), myEnd(end), myAdditiveTraffic(additiveTraffic), myAdaptionFactor(adaptionFactor),
53 59 : myMaxAlternatives(maxAlternatives), myUseDefaultCapacities(defaultCapacities),
54 59 : myNet(net), myMatrix(matrix), myRouter(router), myNetloadOutput(netloadOutput) {
55 59 : myDefaultVehicle = new ROVehicle(SUMOVehicleParameter(), nullptr, net.getVehicleTypeSecure(DEFAULT_VTYPE_ID), &net);
56 59 : }
57 :
58 :
59 59 : ROMAAssignments::~ROMAAssignments() {
60 59 : delete myDefaultVehicle;
61 59 : }
62 :
63 : // based on the definitions in PTV-Validate and in the VISUM-Cologne network
64 : double
65 161641 : ROMAAssignments::getCapacity(const ROEdge* edge) const {
66 161641 : if (edge->isTazConnector()) {
67 : return 0;
68 : }
69 161641 : const int roadClass = myUseDefaultCapacities ? -1 : -edge->getPriority();
70 : // TODO: differ road class 1 from the unknown road class 1!!!
71 161641 : if (edge->getNumLanes() == 0) {
72 : // TAZ have no cost
73 : return 0;
74 161641 : } else if (roadClass == 0 || roadClass == 1) {
75 159361 : return edge->getNumLanes() * 2000.; //CR13 in table.py
76 2280 : } else if (roadClass == 2 && edge->getSpeedLimit() <= 11.) {
77 0 : return edge->getNumLanes() * 1333.33; //CR5 in table.py
78 2280 : } else if (roadClass == 2 && edge->getSpeedLimit() > 11. && edge->getSpeedLimit() <= 16.) {
79 0 : return edge->getNumLanes() * 1500.; //CR3 in table.py
80 2280 : } else if (roadClass == 2 && edge->getSpeedLimit() > 16.) {
81 0 : return edge->getNumLanes() * 2000.; //CR13 in table.py
82 2280 : } else if (roadClass == 3 && edge->getSpeedLimit() <= 11.) {
83 0 : return edge->getNumLanes() * 800.; //CR5 in table.py
84 2280 : } else if (roadClass == 3 && edge->getSpeedLimit() > 11. && edge->getSpeedLimit() <= 13.) {
85 0 : return edge->getNumLanes() * 875.; //CR5 in table.py
86 2280 : } else if (roadClass == 3 && edge->getSpeedLimit() > 13. && edge->getSpeedLimit() <= 16.) {
87 0 : return edge->getNumLanes() * 1500.; //CR4 in table.py
88 2280 : } else if (roadClass == 3 && edge->getSpeedLimit() > 16.) {
89 0 : return edge->getNumLanes() * 1800.; //CR13 in table.py
90 2280 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() <= 5.) {
91 0 : return edge->getNumLanes() * 200.; //CR7 in table.py
92 2280 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 5. && edge->getSpeedLimit() <= 7.) {
93 0 : return edge->getNumLanes() * 412.5; //CR7 in table.py
94 2280 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 7. && edge->getSpeedLimit() <= 9.) {
95 0 : return edge->getNumLanes() * 600.; //CR6 in table.py
96 2280 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 9. && edge->getSpeedLimit() <= 11.) {
97 0 : return edge->getNumLanes() * 800.; //CR5 in table.py
98 2280 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 11. && edge->getSpeedLimit() <= 13.) {
99 0 : return edge->getNumLanes() * 1125.; //CR5 in table.py
100 2280 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 13. && edge->getSpeedLimit() <= 16.) {
101 2250 : return edge->getNumLanes() * 1583.; //CR4 in table.py
102 30 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 16. && edge->getSpeedLimit() <= 18.) {
103 0 : return edge->getNumLanes() * 1100.; //CR3 in table.py
104 30 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 18. && edge->getSpeedLimit() <= 22.) {
105 0 : return edge->getNumLanes() * 1200.; //CR3 in table.py
106 30 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 22. && edge->getSpeedLimit() <= 26.) {
107 0 : return edge->getNumLanes() * 1300.; //CR3 in table.py
108 30 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 26.) {
109 0 : return edge->getNumLanes() * 1400.; //CR3 in table.py
110 : }
111 30 : return edge->getNumLanes() * 800.; //CR5 in table.py
112 : }
113 :
114 :
115 : // based on the definitions in PTV-Validate and in the VISUM-Cologne network
116 : double
117 241444 : ROMAAssignments::capacityConstraintFunction(const ROEdge* edge, const double flow) const {
118 241444 : if (edge->isTazConnector()) {
119 : return 0;
120 : }
121 160422 : const int roadClass = myUseDefaultCapacities ? -1 : -edge->getPriority();
122 160422 : const double capacity = getCapacity(edge);
123 : // TODO: differ road class 1 from the unknown road class 1!!!
124 160422 : if (edge->getNumLanes() == 0) {
125 : // TAZ have no cost
126 : return 0;
127 160422 : } else if (roadClass == 0 || roadClass == 1) {
128 158160 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.*(flow / (capacity * 1.3)) * 2.); //CR13 in table.py
129 2262 : } else if (roadClass == 2 && edge->getSpeedLimit() <= 11.) {
130 0 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.*(flow / (capacity * 0.9)) * 3.); //CR5 in table.py
131 2262 : } else if (roadClass == 2 && edge->getSpeedLimit() > 11. && edge->getSpeedLimit() <= 16.) {
132 0 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.*(flow / (capacity * 1.)) * 2.); //CR3 in table.py
133 2262 : } else if (roadClass == 2 && edge->getSpeedLimit() > 16.) {
134 0 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.*(flow / (capacity * 1.3)) * 2.); //CR13 in table.py
135 2262 : } else if (roadClass == 3 && edge->getSpeedLimit() <= 11.) {
136 0 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.*(flow / (capacity * 0.9)) * 3.); //CR5 in table.py
137 2262 : } else if (roadClass == 3 && edge->getSpeedLimit() > 11. && edge->getSpeedLimit() <= 13.) {
138 0 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.*(flow / (capacity * 0.9)) * 3.); //CR5 in table.py
139 2262 : } else if (roadClass == 3 && edge->getSpeedLimit() > 13. && edge->getSpeedLimit() <= 16.) {
140 0 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.7 * (flow / (capacity * 1.)) * 2.); //CR4 in table.py
141 2262 : } else if (roadClass == 3 && edge->getSpeedLimit() > 16.) {
142 0 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.*(flow / (capacity * 1.3)) * 2.); //CR13 in table.py
143 2262 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() <= 5.) {
144 0 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.*(flow / (capacity * 0.5)) * 3.); //CR7 in table.py
145 2262 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 5. && edge->getSpeedLimit() <= 7.) {
146 0 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.*(flow / (capacity * 0.5)) * 3.); //CR7 in table.py
147 2262 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 7. && edge->getSpeedLimit() <= 9.) {
148 0 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.*(flow / (capacity * 0.8)) * 3.); //CR6 in table.py
149 2262 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 9. && edge->getSpeedLimit() <= 11.) {
150 0 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.*(flow / (capacity * 0.9)) * 3.); //CR5 in table.py
151 2262 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 11. && edge->getSpeedLimit() <= 13.) {
152 0 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.*(flow / (capacity * 0.9)) * 3.); //CR5 in table.py
153 2262 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 13. && edge->getSpeedLimit() <= 16.) {
154 2235 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.7 * (flow / (capacity * 1.)) * 2.); //CR4 in table.py
155 27 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 16. && edge->getSpeedLimit() <= 18.) {
156 0 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.*(flow / (capacity * 1.)) * 2.); //CR3 in table.py
157 27 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 18. && edge->getSpeedLimit() <= 22.) {
158 0 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.*(flow / (capacity * 1.)) * 2.); //CR3 in table.py
159 27 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 22. && edge->getSpeedLimit() <= 26.) {
160 0 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.*(flow / (capacity * 1.)) * 2.); //CR3 in table.py
161 27 : } else if ((roadClass >= 4 || roadClass == -1) && edge->getSpeedLimit() > 26.) {
162 0 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.*(flow / (capacity * 1.)) * 2.); //CR3 in table.py
163 : }
164 27 : return edge->getLength() / edge->getSpeedLimit() * (1. + 1.*(flow / (capacity * 0.9)) * 3.); //CR5 in table.py
165 : }
166 :
167 :
168 : bool
169 40576 : ROMAAssignments::addRoute(const ConstROEdgeVector& edges, std::vector<RORoute*>& paths, std::string routeId, double prob) {
170 : std::vector<RORoute*>::iterator p;
171 48083 : for (p = paths.begin(); p != paths.end(); p++) {
172 45457 : if (edges == (*p)->getEdgeVector()) {
173 : break;
174 : }
175 : }
176 40576 : if (p == paths.end()) {
177 2626 : paths.push_back(new RORoute(routeId, 0., prob, edges, nullptr, std::vector<SUMOVehicleParameter::Stop>()));
178 2626 : return true;
179 : }
180 37950 : (*p)->addProbability(prob);
181 : std::iter_swap(paths.end() - 1, p);
182 37950 : return false;
183 : }
184 :
185 :
186 : const ConstROEdgeVector
187 40595 : ROMAAssignments::computePath(ODCell* cell, const SUMOTime time, const double probability, SUMOAbstractRouter<ROEdge, ROVehicle>* router, bool setBulkMode) {
188 81110 : const ROEdge* const from = myNet.getEdge(cell->origin + (cell->originIsEdge ? "" : "-source"));
189 40595 : if (from == nullptr) {
190 0 : throw ProcessError(TLF("Unknown origin '%'.", cell->origin));
191 : }
192 81110 : const ROEdge* const to = myNet.getEdge(cell->destination + (cell->destinationIsEdge ? "" : "-sink"));
193 40595 : if (to == nullptr) {
194 0 : throw ProcessError(TLF("Unknown destination '%'.", cell->destination));
195 : }
196 : ConstROEdgeVector edges;
197 40595 : if (router == nullptr) {
198 15 : router = &myRouter;
199 : }
200 40595 : if (myMaxAlternatives > 0 && (int)cell->pathsVector.size() < myMaxAlternatives) {
201 40576 : router->compute(from, to, myDefaultVehicle, time, edges);
202 40576 : if (setBulkMode) {
203 40501 : router->setBulkMode(true);
204 : }
205 121728 : if (addRoute(edges, cell->pathsVector, cell->origin + cell->destination + toString(cell->pathsVector.size()), probability)) {
206 : return edges;
207 : }
208 : } else {
209 : double minCost = std::numeric_limits<double>::max();
210 : RORoute* minRoute = nullptr;
211 38 : for (RORoute* const p : cell->pathsVector) {
212 19 : const double cost = router->recomputeCosts(edges, myDefaultVehicle, time);
213 19 : if (cost < minCost) {
214 : minCost = cost;
215 : minRoute = p;
216 : }
217 : }
218 19 : minRoute->addProbability(probability);
219 : }
220 37969 : return ConstROEdgeVector();
221 40595 : }
222 :
223 :
224 : void
225 4 : ROMAAssignments::getKPaths(const int kPaths, const double penalty) {
226 8 : for (ODCell* const c : myMatrix.getCells()) {
227 : myPenalties.clear();
228 14 : for (int k = 0; k < kPaths; k++) {
229 67 : for (const ROEdge* const e : computePath(c)) {
230 57 : myPenalties[e] += penalty;
231 10 : }
232 : }
233 : }
234 : myPenalties.clear();
235 4 : }
236 :
237 :
238 : void
239 59 : ROMAAssignments::resetFlows() {
240 59 : const double begin = STEPS2TIME(MIN2(myBegin, myMatrix.getCells().front()->begin));
241 1999 : for (std::map<std::string, ROEdge*>::const_iterator i = myNet.getEdgeMap().begin(); i != myNet.getEdgeMap().end(); ++i) {
242 1940 : ROMAEdge* edge = static_cast<ROMAEdge*>(i->second);
243 1940 : edge->setFlow(begin, STEPS2TIME(myEnd), 0.);
244 1940 : edge->setHelpFlow(begin, STEPS2TIME(myEnd), 0.);
245 : }
246 59 : }
247 :
248 :
249 : void
250 200 : ROMAAssignments::writeInterval(const SUMOTime begin, const SUMOTime end) {
251 200 : if (myNetloadOutput != nullptr) {
252 96 : myNetloadOutput->openTag(SUMO_TAG_INTERVAL).writeAttr(SUMO_ATTR_BEGIN, time2string(begin)).writeAttr(SUMO_ATTR_END, time2string(end));
253 6491 : for (std::map<std::string, ROEdge*>::const_iterator i = myNet.getEdgeMap().begin(); i != myNet.getEdgeMap().end(); ++i) {
254 6459 : ROMAEdge* edge = static_cast<ROMAEdge*>(i->second);
255 6459 : if (edge->getFunction() == SumoXMLEdgeFunc::NORMAL) {
256 1219 : myNetloadOutput->openTag(SUMO_TAG_EDGE).writeAttr(SUMO_ATTR_ID, edge->getID());
257 1219 : const double traveltime = edge->getTravelTime(getDefaultVehicle(), STEPS2TIME(begin));
258 1219 : const double speed = edge->getLength() / traveltime;
259 1219 : const double flow = edge->getFlow(STEPS2TIME(begin));
260 1219 : const double timeGap = STEPS2TIME(end - begin) / flow;
261 1219 : const double spaceGap = timeGap * speed;
262 1219 : const double density = 1000.0 / spaceGap;
263 : const double laneDensity = density / edge->getNumLanes();
264 1219 : myNetloadOutput->writeAttr(SUMO_ATTR_TRAVELTIME, traveltime);
265 1219 : myNetloadOutput->writeAttr(SUMO_ATTR_SPEED, speed);
266 1219 : myNetloadOutput->writeAttr(SUMO_ATTR_SPEEDREL, speed / edge->getSpeedLimit());
267 1219 : myNetloadOutput->writeAttr(SUMO_ATTR_ENTERED, flow);
268 1219 : myNetloadOutput->writeAttr(SUMO_ATTR_DENSITY, density);
269 1219 : myNetloadOutput->writeAttr(SUMO_ATTR_LANEDENSITY, laneDensity);
270 1219 : myNetloadOutput->writeAttr("flowCapacityRatio", 100. * flow / getCapacity(edge));
271 2438 : myNetloadOutput->closeTag();
272 : }
273 : }
274 64 : myNetloadOutput->closeTag();
275 : }
276 200 : }
277 :
278 :
279 : void
280 55 : ROMAAssignments::incremental(const int numIter, const bool verbose) {
281 : SUMOTime lastBegin = -1;
282 : std::vector<int> intervals;
283 55 : int count = 0;
284 2081 : for (const ODCell* const c : myMatrix.getCells()) {
285 2026 : if (c->begin != lastBegin) {
286 196 : intervals.push_back(count);
287 196 : lastBegin = c->begin;
288 : }
289 2026 : count++;
290 : }
291 : lastBegin = -1;
292 251 : for (std::vector<int>::const_iterator offset = intervals.begin(); offset != intervals.end(); offset++) {
293 196 : std::vector<ODCell*>::const_iterator firstCell = myMatrix.getCells().begin() + (*offset);
294 : std::vector<ODCell*>::const_iterator lastCell = myMatrix.getCells().end();
295 196 : if (offset != intervals.end() - 1) {
296 141 : lastCell = myMatrix.getCells().begin() + (*(offset + 1));
297 : }
298 196 : const SUMOTime intervalStart = (*firstCell)->begin;
299 196 : if (verbose) {
300 0 : WRITE_MESSAGE(" starting interval " + time2string(intervalStart));
301 : }
302 : std::map<const ROMAEdge*, double> loadedTravelTimes;
303 26814 : for (std::map<std::string, ROEdge*>::const_iterator i = myNet.getEdgeMap().begin(); i != myNet.getEdgeMap().end(); ++i) {
304 26618 : ROMAEdge* edge = static_cast<ROMAEdge*>(i->second);
305 26618 : if (edge->hasLoadedTravelTime(STEPS2TIME(intervalStart))) {
306 2024 : loadedTravelTimes[edge] = edge->getTravelTime(myDefaultVehicle, STEPS2TIME(intervalStart));
307 : }
308 : }
309 196 : if (loadedTravelTimes.empty()) {
310 : // we don't have loaded travel times for this interval but we may have set ourselves some for the interval before but no need to warn
311 : ROEdge::disableTimelineWarning();
312 : }
313 4146 : for (int t = 0; t < numIter; t++) {
314 3950 : if (verbose) {
315 0 : WRITE_MESSAGE(" starting iteration " + toString(t));
316 : }
317 3950 : std::string lastOrigin = "";
318 : int workerIndex = 0;
319 44530 : for (std::vector<ODCell*>::const_iterator i = firstCell; i != lastCell; i++) {
320 40580 : ODCell* const c = *i;
321 40580 : const double linkFlow = c->vehicleNumber / numIter;
322 40580 : const SUMOTime begin = myAdditiveTraffic ? myBegin : c->begin;
323 : #ifdef HAVE_FOX
324 40580 : if (myNet.getThreadPool().size() > 0) {
325 60 : if (lastOrigin != c->origin) {
326 60 : workerIndex++;
327 60 : if (workerIndex == myNet.getThreadPool().size()) {
328 : workerIndex = 0;
329 : }
330 60 : myNet.getThreadPool().add(new RONet::BulkmodeTask(false), workerIndex);
331 : lastOrigin = c->origin;
332 60 : myNet.getThreadPool().add(new RoutingTask(*this, c, begin, linkFlow), workerIndex);
333 : } else {
334 0 : myNet.getThreadPool().add(new RoutingTask(*this, c, begin, linkFlow), workerIndex);
335 : }
336 60 : continue;
337 : }
338 : #endif
339 40520 : if (lastOrigin != c->origin) {
340 11680 : myRouter.setBulkMode(false);
341 : lastOrigin = c->origin;
342 : }
343 40520 : computePath(c, begin, linkFlow, &myRouter, true);
344 : }
345 : #ifdef HAVE_FOX
346 3950 : if (myNet.getThreadPool().size() > 0) {
347 40 : myNet.getThreadPool().waitAll();
348 : }
349 : #endif
350 44530 : for (std::vector<ODCell*>::const_iterator i = firstCell; i != lastCell; i++) {
351 40580 : ODCell* const c = *i;
352 40580 : const double linkFlow = c->vehicleNumber / numIter;
353 40580 : const SUMOTime begin = myAdditiveTraffic ? myBegin : c->begin;
354 40580 : const SUMOTime end = myAdditiveTraffic ? myEnd : c->end;
355 40580 : const double intervalLengthInHours = STEPS2TIME(end - begin) / 3600.;
356 40580 : const ConstROEdgeVector& edges = c->pathsVector.back()->getEdgeVector();
357 281588 : for (ConstROEdgeVector::const_iterator e = edges.begin(); e != edges.end(); e++) {
358 241008 : ROMAEdge* edge = static_cast<ROMAEdge*>(myNet.getEdge((*e)->getID()));
359 241008 : const double newFlow = edge->getFlow(STEPS2TIME(begin)) + linkFlow;
360 241008 : edge->setFlow(STEPS2TIME(begin), STEPS2TIME(end), newFlow);
361 241008 : double travelTime = capacityConstraintFunction(edge, newFlow / intervalLengthInHours);
362 241008 : if (lastBegin >= 0 && myAdaptionFactor > 0.) {
363 : if (loadedTravelTimes.count(edge) != 0) {
364 0 : travelTime = loadedTravelTimes[edge] * myAdaptionFactor + (1. - myAdaptionFactor) * travelTime;
365 : } else {
366 55200 : travelTime = edge->getTravelTime(myDefaultVehicle, STEPS2TIME(lastBegin)) * myAdaptionFactor + (1. - myAdaptionFactor) * travelTime;
367 : }
368 : }
369 241008 : edge->addTravelTime(travelTime, STEPS2TIME(begin), STEPS2TIME(end));
370 : }
371 : }
372 3950 : myRouter.reset(myDefaultVehicle);
373 : }
374 196 : writeInterval(intervalStart, (*(lastCell - 1))->end);
375 : lastBegin = intervalStart;
376 : }
377 55 : }
378 :
379 :
380 : void
381 4 : ROMAAssignments::sue(const int maxOuterIteration, const int maxInnerIteration, const int kPaths, const double penalty, const double tolerance, const std::string /* routeChoiceMethod */) {
382 4 : getKPaths(kPaths, penalty);
383 : std::map<const SUMOTime, SUMOTime> intervals;
384 4 : if (myAdditiveTraffic) {
385 0 : intervals[myBegin] = myEnd;
386 : } else {
387 8 : for (const ODCell* const c : myMatrix.getCells()) {
388 4 : intervals[c->begin] = c->end;
389 : }
390 : }
391 5 : for (int outer = 0; outer < maxOuterIteration; outer++) {
392 11 : for (int inner = 0; inner < maxInnerIteration; inner++) {
393 22 : for (const ODCell* const c : myMatrix.getCells()) {
394 11 : const SUMOTime begin = myAdditiveTraffic ? myBegin : c->begin;
395 11 : const SUMOTime end = myAdditiveTraffic ? myEnd : c->end;
396 : // update path cost
397 29 : for (std::vector<RORoute*>::const_iterator j = c->pathsVector.begin(); j != c->pathsVector.end(); ++j) {
398 18 : RORoute* r = *j;
399 18 : r->setCosts(myRouter.recomputeCosts(r->getEdgeVector(), myDefaultVehicle, 0));
400 : // std::cout << std::setprecision(20) << r->getID() << ":" << r->getCosts() << std::endl;
401 : }
402 : // calculate route utilities and probabilities
403 11 : RouteCostCalculator<RORoute, ROEdge, ROVehicle>::getCalculator().calculateProbabilities(c->pathsVector, myDefaultVehicle, 0);
404 : // calculate route flows
405 29 : for (std::vector<RORoute*>::const_iterator j = c->pathsVector.begin(); j != c->pathsVector.end(); ++j) {
406 18 : RORoute* r = *j;
407 18 : const double pathFlow = r->getProbability() * c->vehicleNumber;
408 : // assign edge flow deltas
409 201 : for (ConstROEdgeVector::const_iterator e = r->getEdgeVector().begin(); e != r->getEdgeVector().end(); e++) {
410 183 : ROMAEdge* edge = static_cast<ROMAEdge*>(myNet.getEdge((*e)->getID()));
411 183 : edge->setHelpFlow(STEPS2TIME(begin), STEPS2TIME(end), edge->getHelpFlow(STEPS2TIME(begin)) + pathFlow);
412 : }
413 : }
414 : }
415 : // calculate new edge flows and check for stability
416 : int unstableEdges = 0;
417 22 : for (const auto& it : intervals) {
418 11 : const double intervalLengthInHours = STEPS2TIME(it.second - it.first) / 3600.;
419 11 : const double intBegin = STEPS2TIME(it.first);
420 11 : const double intEnd = STEPS2TIME(it.second);
421 447 : for (std::map<std::string, ROEdge*>::const_iterator e = myNet.getEdgeMap().begin(); e != myNet.getEdgeMap().end(); ++e) {
422 436 : ROMAEdge* edge = static_cast<ROMAEdge*>(e->second);
423 : const double oldFlow = edge->getFlow(intBegin);
424 : double newFlow = oldFlow;
425 436 : if (inner == 0 && outer == 0) {
426 149 : newFlow += edge->getHelpFlow(intBegin);
427 : } else {
428 287 : newFlow += (edge->getHelpFlow(intBegin) - oldFlow) / (inner + 1);
429 : }
430 : // if not lohse:
431 436 : if (newFlow > 0.) {
432 137 : if (fabs(newFlow - oldFlow) / newFlow > tolerance) {
433 60 : unstableEdges++;
434 : }
435 299 : } else if (newFlow == 0.) {
436 299 : if (oldFlow != 0. && (fabs(newFlow - oldFlow) / oldFlow > tolerance)) {
437 0 : unstableEdges++;
438 : }
439 : } else { // newFlow < 0.
440 0 : unstableEdges++;
441 : newFlow = 0.;
442 : }
443 : edge->setFlow(intBegin, intEnd, newFlow);
444 436 : const double travelTime = capacityConstraintFunction(edge, newFlow / intervalLengthInHours);
445 436 : edge->addTravelTime(travelTime, intBegin, intEnd);
446 : edge->setHelpFlow(intBegin, intEnd, 0.);
447 : }
448 : }
449 : // if stable break
450 11 : if (unstableEdges == 0) {
451 : break;
452 : }
453 : // additional stability check from python script: if notstable < math.ceil(net.geteffEdgeCounts()*0.005) or notstable < 3: stable = True
454 : }
455 : // check for a new route, if none available, break
456 : // several modifications about when a route is new and when to break are in the original script
457 : bool newRoute = false;
458 10 : for (ODCell* const c : myMatrix.getCells()) {
459 5 : newRoute |= !computePath(c).empty();
460 : }
461 5 : if (!newRoute) {
462 : break;
463 : }
464 : }
465 : // final round of assignment
466 8 : for (const ODCell* const c : myMatrix.getCells()) {
467 : // update path cost
468 11 : for (std::vector<RORoute*>::const_iterator j = c->pathsVector.begin(); j != c->pathsVector.end(); ++j) {
469 7 : RORoute* r = *j;
470 7 : r->setCosts(myRouter.recomputeCosts(r->getEdgeVector(), myDefaultVehicle, 0));
471 : }
472 : // calculate route utilities and probabilities
473 4 : RouteCostCalculator<RORoute, ROEdge, ROVehicle>::getCalculator().calculateProbabilities(c->pathsVector, myDefaultVehicle, 0);
474 : // calculate route flows
475 11 : for (std::vector<RORoute*>::const_iterator j = c->pathsVector.begin(); j != c->pathsVector.end(); ++j) {
476 7 : RORoute* r = *j;
477 7 : r->setProbability(r->getProbability() * c->vehicleNumber);
478 : }
479 : }
480 8 : for (const auto& it : intervals) {
481 4 : writeInterval(it.first, it.second);
482 : }
483 4 : }
484 :
485 :
486 : double
487 0 : ROMAAssignments::getPenalizedEffort(const ROEdge* const e, const ROVehicle* const v, double t) {
488 : const std::map<const ROEdge* const, double>::const_iterator i = myPenalties.find(e);
489 0 : return i == myPenalties.end() ? e->getEffort(v, t) : e->getEffort(v, t) + i->second;
490 : }
491 :
492 :
493 : double
494 236 : ROMAAssignments::getPenalizedTT(const ROEdge* const e, const ROVehicle* const v, double t) {
495 : const std::map<const ROEdge* const, double>::const_iterator i = myPenalties.find(e);
496 236 : return i == myPenalties.end() ? e->getTravelTime(v, t) : e->getTravelTime(v, t) + i->second;
497 : }
498 :
499 :
500 : double
501 0 : ROMAAssignments::getTravelTime(const ROEdge* const e, const ROVehicle* const v, double t) {
502 0 : return e->getTravelTime(v, t);
503 : }
504 :
505 :
506 : #ifdef HAVE_FOX
507 : // ---------------------------------------------------------------------------
508 : // ROMAAssignments::RoutingTask-methods
509 : // ---------------------------------------------------------------------------
510 : void
511 60 : ROMAAssignments::RoutingTask::run(MFXWorkerThread* context) {
512 120 : myAssign.computePath(myCell, myBegin, myLinkFlow, &static_cast<RONet::WorkerThread*>(context)->getVehicleRouter(SVC_IGNORING), mySetBulkMode);
513 60 : }
514 : #endif
|