Line data Source code
1 : #include <cmath>
2 : #include <iostream>
3 : #include <numeric>
4 : #include <set>
5 : #include <tuple>
6 :
7 : #include "Alg/Flow/EdmondsKarp.hpp"
8 : #include "Alg/Graph.hpp"
9 : #include "Alg/ShortestPath/BFS.hpp"
10 : #include "Static/supply/BPRNetwork.hpp"
11 : #include "data/SUMO/Network.hpp"
12 : #include "data/SUMO/NetworkTAZ.hpp"
13 : #include "data/SUMO/SUMO.hpp"
14 :
15 : using namespace std;
16 : using namespace Static;
17 : using namespace utils;
18 :
19 : using Alg::Graph;
20 :
21 : typedef SUMO::Network::Edge::Lane Lane;
22 : typedef SUMO::Speed Speed;
23 : typedef SUMO::Length Length;
24 :
25 281263 : Time BPRNetwork::Loader<SUMO::NetworkTAZs>::calculateFreeFlowSpeed(const Time &maxSpeed) const {
26 281263 : return maxSpeed * 0.9;
27 : }
28 :
29 37856 : Time BPRNetwork::Loader<SUMO::NetworkTAZs>::calculateFreeFlowSpeed(const SUMO::Network::Edge &e) const {
30 37856 : return e.speed() * 0.9;
31 : }
32 :
33 254801 : Time BPRNetwork::Loader<SUMO::NetworkTAZs>::calculateFreeFlowSpeed(const SUMO::Network::Edge::Lane &l) const {
34 254801 : return l.speed * 0.9;
35 : }
36 :
37 7534 : Time BPRNetwork::Loader<SUMO::NetworkTAZs>::calculateFreeFlowTime(const SUMO::Network::Edge &e) const {
38 7534 : Length length = e.length();
39 7534 : Speed freeFlowSpeed = calculateFreeFlowSpeed(e);
40 7534 : Time freeFlowTime = length / freeFlowSpeed;
41 7534 : return freeFlowTime;
42 : }
43 :
44 0 : Time BPRNetwork::Loader<SUMO::NetworkTAZs>::calculateFreeFlowTime(const SUMO::Network::Edge::Lane &l) const {
45 0 : Length length = l.length;
46 0 : Speed freeFlowSpeed = calculateFreeFlowSpeed(l);
47 0 : Time freeFlowTime = length / freeFlowSpeed;
48 0 : return freeFlowTime;
49 : }
50 :
51 : const Flow SATURATION_FLOW = 1110.0; // vehicles per hour per lane
52 : const Flow SATURATION_FLOW_EXTRA_LANES = 800.0;
53 :
54 : /**
55 : * @brief Cost of having traffic stop once. This should be a time penalty that
56 : * averages all the negative effects of cars having to start after the light
57 : * turns green.
58 : */
59 : const Time STOP_PENALTY = 0.0;
60 :
61 7534 : Flow BPRNetwork::Loader<SUMO::NetworkTAZs>::calculateCapacity(const SUMO::Network::Edge &e) const {
62 7534 : Speed freeFlowSpeed = calculateFreeFlowSpeed(e);
63 7534 : Time adjSaturationFlow = (SATURATION_FLOW / 60.0 / 60.0) * (freeFlowSpeed / calculateFreeFlowSpeed(50.0 / 3.6));
64 7534 : Time adjSaturationFlowExtraLanes = (SATURATION_FLOW_EXTRA_LANES / 60.0 / 60.0) * (freeFlowSpeed / calculateFreeFlowSpeed(50.0 / 3.6));
65 7534 : Flow c = adjSaturationFlow + adjSaturationFlowExtraLanes * (Time)(e.lanes.size() - 1);
66 :
67 7534 : const vector<reference_wrapper<const SUMO::Network::Connection>> &connections = e.getOutgoingConnections();
68 7534 : if(!connections.empty()) {
69 14776 : vector<Flow> capacityPerLane(e.lanes.size(), 0.0);
70 21231 : for(const SUMO::Network::Connection &conn: connections) {
71 13843 : Flow cAdd = adjSaturationFlow;
72 13843 : if(conn.tl) {
73 2157 : SUMO::Time
74 2157 : g = conn.getGreenTime(),
75 2157 : C = conn.getCycleTime();
76 2157 : size_t n = conn.getNumberStops();
77 2157 : cAdd *= (g - STOP_PENALTY * (Time)n) / C;
78 : }
79 13843 : capacityPerLane.at(conn.fromLane().index) += cAdd;
80 : }
81 17338 : for(Flow &capPerLane: capacityPerLane)
82 13194 : capPerLane = min(capPerLane, adjSaturationFlow);
83 17338 : Flow cNew = accumulate(capacityPerLane.begin(), capacityPerLane.end(), 0.0);
84 :
85 7388 : #pragma GCC diagnostic push
86 7388 : #pragma GCC diagnostic ignored "-Wfloat-equal"
87 7388 : if(cNew != 0.0) {
88 7955 : c = min(c, cNew);
89 : }
90 7388 : #pragma GCC diagnostic pop
91 : }
92 :
93 14922 : return c;
94 : }
95 :
96 254801 : Flow BPRNetwork::Loader<SUMO::NetworkTAZs>::calculateCapacity(const SUMO::Network::Edge::Lane &lane) const {
97 254801 : Speed freeFlowSpeed = calculateFreeFlowSpeed(lane);
98 254801 : Time adjSaturationFlow = (SATURATION_FLOW / 60.0 / 60.0) * (freeFlowSpeed / calculateFreeFlowSpeed(50.0 / 3.6));
99 254801 : Flow c = adjSaturationFlow;
100 254801 : return c;
101 : }
102 :
103 0 : BPRNetwork *BPRNetwork::Loader<SUMO::NetworkTAZs>::load(const SUMO::NetworkTAZs &sumo) {
104 0 : clear();
105 :
106 0 : network = new BPRNetwork();
107 :
108 0 : addNormalEdges(sumo);
109 :
110 0 : addConnections(sumo);
111 :
112 0 : iterateCapacities(sumo);
113 :
114 0 : addTAZs(sumo);
115 :
116 0 : return network;
117 : }
118 :
119 3 : void BPRNetwork::Loader<SUMO::NetworkTAZs>::addNormalEdges(const SUMO::NetworkTAZs &sumo) {
120 6 : const vector<SUMO::Network::Edge> &sumoEdges = sumo.network.getEdges();
121 19357 : for(const SUMO::Network::Edge &edge: sumoEdges) {
122 19354 : if(edge.function == SUMO::Network::Edge::Function::INTERNAL) continue;
123 :
124 7534 : const auto &p = adapter.addSumoEdge(edge.id);
125 7534 : const NormalEdge::ID &eid = p.first;
126 7534 : Node u = p.second.first, v = p.second.second;
127 :
128 : // clang-format off
129 7534 : normalEdges[edge.id] = network->addNormalEdge(
130 : eid,
131 : u, v,
132 7534 : *network,
133 7534 : calculateFreeFlowTime(edge),
134 7534 : calculateCapacity(edge)
135 7534 : );
136 : // clang-format on
137 :
138 7534 : in[edge.to.value().get().id].push_back(edge);
139 7534 : out[edge.from.value().get().id].push_back(edge);
140 : }
141 3 : }
142 :
143 3 : void BPRNetwork::Loader<SUMO::NetworkTAZs>::addConnections(const SUMO::NetworkTAZs &sumo) {
144 3 : auto connections = sumo.network.getConnections();
145 :
146 19357 : for(const SUMO::Network::Edge &from: sumo.network.getEdges()) {
147 19354 : if(from.function == SUMO::Network::Edge::Function::INTERNAL) continue;
148 :
149 30202 : for(const SUMO::Network::Edge &to: from.getOutgoing()) {
150 15134 : if(to.function == SUMO::Network::Edge::Function::INTERNAL) continue;
151 :
152 15134 : addConnection(sumo, from, to);
153 : }
154 : }
155 3 : }
156 :
157 15134 : void BPRNetwork::Loader<SUMO::NetworkTAZs>::addConnection(const SUMO::NetworkTAZs &sumo, const SUMO::Network::Edge &from, const SUMO::Network::Edge &to) {
158 26528 : auto fromToConnections = sumo.network.getConnections(from, to);
159 :
160 15134 : if(fromToConnections.empty()) return;
161 :
162 11394 : Speed v = min(
163 11394 : calculateFreeFlowSpeed(from),
164 11394 : calculateFreeFlowSpeed(to)
165 11394 : );
166 11394 : Time adjSaturationFlow = (SATURATION_FLOW / 60.0 / 60.0) * (v / calculateFreeFlowSpeed(50.0 / 3.6));
167 :
168 22788 : vector<Time> capacityFromLanes(from.lanes.size(), 0.0);
169 22788 : vector<Time> capacityToLanes(to.lanes.size(), 0.0);
170 :
171 11394 : double t0 = 0;
172 :
173 25237 : for(const SUMO::Network::Connection &conn: fromToConnections) {
174 13843 : Time cAdd = adjSaturationFlow;
175 : /// Traffic lights
176 13843 : if(conn.tl) {
177 2157 : SUMO::Time g = conn.getGreenTime();
178 2157 : SUMO::Time C = conn.getCycleTime();
179 2157 : SUMO::Time r = C - g;
180 2157 : size_t n = conn.getNumberStops();
181 2157 : cAdd *= (g - STOP_PENALTY * (double)n) / C;
182 :
183 2157 : t0 += r * r / (2.0 * C);
184 : }
185 13843 : capacityFromLanes[conn.fromLane().index] += cAdd;
186 13843 : capacityToLanes[conn.toLane().index] += cAdd;
187 :
188 : /// Direction changes
189 13843 : #pragma GCC diagnostic push
190 13843 : #pragma GCC diagnostic ignored "-Wswitch-enum"
191 : // clang-format off
192 13843 : switch(conn.dir) {
193 289 : case SUMO::Network::Connection::Direction::PARTIALLY_RIGHT: t0 += 1.0; break;
194 2415 : case SUMO::Network::Connection::Direction::RIGHT : t0 += 2.0; break;
195 245 : case SUMO::Network::Connection::Direction::PARTIALLY_LEFT : t0 += 2.0; break;
196 2475 : case SUMO::Network::Connection::Direction::LEFT : t0 += 5.0; break;
197 7 : case SUMO::Network::Connection::Direction::TURN : t0 += 20.0; break;
198 : default:
199 : break;
200 : }
201 : // clang-format on
202 13843 : #pragma GCC diagnostic pop
203 : }
204 11394 : t0 /= (double)fromToConnections.size();
205 :
206 26862 : for(Time &c: capacityFromLanes) c = min(c, adjSaturationFlow);
207 26794 : for(Time &c: capacityToLanes) c = min(c, adjSaturationFlow);
208 26862 : double cFrom = accumulate(capacityFromLanes.begin(), capacityFromLanes.end(), 0.0);
209 11394 : double cTo = accumulate(capacityToLanes.begin(), capacityToLanes.end(), 0.0);
210 11394 : double c = min(cFrom, cTo);
211 :
212 11394 : ConnectionEdge *e = network->addConnectionEdge(
213 : adapter.addEdge(),
214 11394 : adapter.toNodes(from.id).second,
215 11394 : adapter.toNodes(to.id).first,
216 : *network,
217 : t0,
218 : c
219 11394 : );
220 :
221 22788 : connectionEdges[e->id] = make_tuple(e, from.id, to.id);
222 :
223 : // connectionMap[from.id][to.id] = e->id;
224 : }
225 :
226 3 : void BPRNetwork::Loader<SUMO::NetworkTAZs>::iterateCapacities(const SUMO::NetworkTAZs &sumo) {
227 3 : map<Edge::ID, Edge *> &netEdges = network->edges;
228 :
229 3 : const Flow EPSILON = 1.0 / 60.0 / 60.0 / 24.0;
230 3 : const size_t ITERATIONS = 100;
231 :
232 3 : bool changed = true;
233 :
234 11 : for(size_t i = 0; i < ITERATIONS && changed; ++i) {
235 8 : changed = false;
236 :
237 94499 : for(auto &[edgeID, edge]: netEdges) {
238 : // 1.1
239 94491 : {
240 94491 : const auto &nextEdges = network->adj.at(edge->v);
241 94491 : if(!nextEdges.empty()) {
242 : Flow c = 0.0;
243 207531 : for(const Edge *nextEdge: nextEdges)
244 113754 : c += nextEdge->c;
245 :
246 93777 : if(edge->c > c + EPSILON) {
247 : // cerr << " 1.1. | "
248 : // << "Capacity of edge " << edge->id
249 : // << " (SUMO edge " << (adapter.isSumoEdge(edge->id) ? adapter.toSumoEdge(edge->id) : "-") << ")"
250 : // << " was reduced from " << edge->c
251 : // << " to " << c
252 : // << " (delta=" << edge->c - c << ")"
253 : // << endl;
254 : // cerr << " Note: Adjacent edges' capacities are:" << endl;
255 : // for(const Edge *nextEdge: nextEdges)
256 : // cerr << " " << nextEdge->id << " (SUMO edge " << (adapter.isSumoEdge(nextEdge->id) ? adapter.toSumoEdge(nextEdge->id) : "-") << "): " << nextEdge->c << endl;
257 : // if(c / edge->c < 0.5)
258 : // cerr << " WARNING: Capacity reduced by more than 50%! ================================" << endl;
259 1803 : edge->c = c;
260 1803 : changed = true;
261 : }
262 : }
263 : }
264 : }
265 :
266 94499 : for(auto &[edgeID, edge]: netEdges) {
267 : // 1.2
268 94491 : if(adapter.isSumoEdge(edge->id)) {
269 75228 : const SUMO::Network::Edge::ID fromID = adapter.toSumoEdge(edge->id);
270 37614 : const SUMO::Network::Edge &from = sumo.network.getEdge(fromID);
271 :
272 75228 : const vector<reference_wrapper<const SUMO::Network::Connection>> connections = from.getOutgoingConnections();
273 :
274 37614 : if(!connections.empty()) {
275 73800 : map<SUMO::Network::Edge::ID, Graph::Node> sumoEdges2nodes;
276 36900 : map<SUMO::Network::Edge::Lane::ID, Graph::Node> sumoLanes2nodes;
277 :
278 73800 : Graph G;
279 36900 : Graph::Node incNode = 0;
280 36900 : Graph::Edge::ID incEdge = 0;
281 :
282 36900 : Graph::Node vSource = incNode++;
283 36900 : Graph::Node vSink = incNode++;
284 :
285 106002 : for(const SUMO::Network::Connection &conn: connections) {
286 69102 : const Edge &nextEdge = network->getEdge(adapter.toEdge(conn.to.id));
287 :
288 : // This from edge was never seen before
289 69102 : Graph::Node edgeSource;
290 69102 : if(!sumoEdges2nodes.count(conn.from.id)) {
291 36900 : edgeSource = incNode++;
292 :
293 36900 : sumoEdges2nodes[conn.from.id] = edgeSource;
294 :
295 36900 : G.addEdge(incEdge++, vSource, edgeSource, edge->c);
296 : } else {
297 32202 : edgeSource = sumoEdges2nodes.at(conn.from.id);
298 : }
299 :
300 : // This to edge was never seen before
301 69102 : Graph::Node edgeSink;
302 69102 : if(!sumoEdges2nodes.count(conn.to.id)) {
303 56877 : edgeSink = incNode++;
304 :
305 56877 : sumoEdges2nodes[conn.to.id] = edgeSink;
306 :
307 56877 : G.addEdge(incEdge++, edgeSink, vSink, nextEdge.c);
308 : } else {
309 12225 : edgeSink = sumoEdges2nodes.at(conn.to.id);
310 : }
311 :
312 : // This fromLane was never seen before
313 69102 : Graph::Node u;
314 69102 : if(!sumoLanes2nodes.count(conn.fromLane().id)) {
315 49055 : u = incNode++;
316 :
317 49055 : sumoLanes2nodes[conn.fromLane().id] = u;
318 :
319 49055 : G.addEdge(incEdge++, edgeSource, u, calculateCapacity(conn.fromLane()));
320 : } else {
321 20047 : u = sumoLanes2nodes.at(conn.fromLane().id);
322 : }
323 :
324 : // This toLane was never seen before
325 69102 : Graph::Node v;
326 69102 : if(!sumoLanes2nodes.count(conn.toLane().id)) {
327 69102 : v = incNode++;
328 :
329 69102 : sumoLanes2nodes[conn.toLane().id] = v;
330 :
331 69102 : G.addEdge(incEdge++, v, edgeSink, calculateCapacity(conn.toLane()));
332 : } else {
333 0 : v = sumoLanes2nodes.at(conn.toLane().id);
334 : }
335 :
336 69102 : G.addEdge(incEdge++, u, v, INFINITY);
337 : }
338 :
339 73800 : Alg::ShortestPath::BFS sp;
340 73800 : Alg::Flow::EdmondsKarp maxFlow(sp);
341 36900 : Alg::Graph::Edge::Weight c = maxFlow.solve(G, vSource, vSink);
342 :
343 36900 : if(edge->c > c + EPSILON) {
344 : // cerr << " 1.2. | "
345 : // << "Capacity of edge " << edge->id
346 : // << " (SUMO edge " << adapter.toSumoEdge(edge->id) << ")"
347 : // << " was reduced from " << edge->c
348 : // << " to " << c
349 : // << " (delta=" << edge->c - c << ")"
350 : // << endl;
351 375 : edge->c = c;
352 375 : changed = true;
353 : }
354 : }
355 : }
356 : }
357 :
358 : // 2.
359 94499 : for(auto &[edgeID, edge]: netEdges) {
360 : // 2.1
361 94491 : if(!adapter.isSumoEdge(edge->id)) {
362 56877 : const auto &[_, fromID, toID] = connectionEdges.at(edge->id);
363 :
364 56877 : const SUMO::Network::Edge &from = sumo.network.getEdge(fromID);
365 56877 : const SUMO::Network::Edge &to = sumo.network.getEdge(toID);
366 :
367 56877 : const Edge *prevEdge = netEdges.at(adapter.toEdge(from.id));
368 56877 : const Edge *nextEdge = netEdges.at(adapter.toEdge(to.id));
369 :
370 113754 : const auto &conns = sumo.network.getConnections(from, to);
371 :
372 56877 : if(!conns.empty()) {
373 113754 : map<SUMO::Network::Edge::ID, Graph::Node> sumoEdges2nodes;
374 56877 : map<SUMO::Network::Edge::Lane::ID, Graph::Node> sumoLanes2nodes;
375 :
376 113754 : Graph G;
377 56877 : Graph::Node incNode = 0;
378 56877 : Graph::Edge::ID incEdge = 0;
379 :
380 56877 : Graph::Node vSource = incNode++;
381 56877 : Graph::Node vSink = incNode++;
382 :
383 56877 : Graph::Node source = incNode++;
384 70899 : G.addEdge(incEdge++, vSource, source, min(min(prevEdge->c, nextEdge->c), edge->c));
385 :
386 125979 : for(const SUMO::Network::Connection &conn: conns) {
387 : // This fromLane was never seen before
388 69102 : Graph::Node u;
389 69102 : if(!sumoLanes2nodes.count(conn.fromLane().id)) {
390 67542 : u = incNode++;
391 :
392 67542 : sumoLanes2nodes[conn.fromLane().id] = u;
393 :
394 67542 : G.addEdge(incEdge++, source, u, calculateCapacity(conn.fromLane()));
395 : } else {
396 1560 : u = sumoLanes2nodes.at(conn.fromLane().id);
397 : }
398 :
399 : // This toLane was never seen before
400 69102 : Graph::Node v;
401 69102 : if(!sumoLanes2nodes.count(conn.toLane().id)) {
402 69102 : v = incNode++;
403 :
404 69102 : sumoLanes2nodes[conn.toLane().id] = v;
405 :
406 69102 : G.addEdge(incEdge++, v, vSink, calculateCapacity(conn.toLane()));
407 : } else {
408 0 : v = sumoLanes2nodes.at(conn.toLane().id);
409 : }
410 :
411 69102 : G.addEdge(incEdge++, u, v, INFINITY);
412 : }
413 :
414 113754 : Alg::ShortestPath::BFS sp;
415 113754 : Alg::Flow::EdmondsKarp maxFlow(sp);
416 56877 : Alg::Graph::Edge::Weight c = maxFlow.solve(G, vSource, vSink);
417 :
418 56877 : if(edge->c > c + EPSILON) {
419 : // cerr << " 2.1. | "
420 : // << "Capacity of edge " << edge->id
421 : // << " was reduced from " << edge->c
422 : // << " to " << c
423 : // << " (delta=" << edge->c - c << ")"
424 : // << endl;
425 433 : edge->c = c;
426 433 : changed = true;
427 : }
428 : }
429 : }
430 : }
431 : }
432 3 : }
433 :
434 3 : void BPRNetwork::Loader<SUMO::NetworkTAZs>::addTAZs(const SUMO::NetworkTAZs &sumo) {
435 122 : for(const auto &[id, taz]: sumo.tazs) {
436 119 : const auto &[source, sink] = adapter.addSumoTAZ(taz.id);
437 1205 : for(const SUMO::TAZ::Source &s: taz.sources) {
438 1086 : const Edge *e = network->edges.at(adapter.toEdge(s.id));
439 1086 : network->addNormalEdge(
440 : adapter.addEdge(),
441 : source,
442 1086 : e->u,
443 : *network,
444 : 0,
445 : 1e9
446 1086 : );
447 : }
448 1205 : for(const SUMO::TAZ::Sink &s: taz.sinks) {
449 1086 : const Edge *e = network->edges.at(adapter.toEdge(s.id));
450 1086 : network->addNormalEdge(
451 : adapter.addEdge(),
452 1086 : e->v,
453 : sink,
454 : *network,
455 : 0,
456 : 1e9
457 1086 : );
458 : }
459 : }
460 3 : }
461 :
462 3 : void BPRNetwork::Loader<SUMO::NetworkTAZs>::clear() {
463 3 : adapter.clear();
464 3 : in.clear();
465 3 : out.clear();
466 3 : normalEdges.clear();
467 3 : }
468 :
469 3 : BPRNetwork::Loader<SUMO::NetworkTAZs>::Loader::~Loader() {}
|
