Line data    Source code

       1             : #include "Static/algos/ConjugateFrankWolfe.hpp"
       2             : 
       3             : #include <chrono>
       4             : #include <cmath>
       5             : #include <iomanip>
       6             : #include <iostream>
       7             : #include <memory>
       8             : #include <utility>
       9             : 
      10             : #include "Log/ProgressLogger.hpp"
      11             : #include "Log/ProgressLoggerTableOStream.hpp"
      12             : #include "Opt/QuadraticSolver.hpp"
      13             : #include "Opt/UnivariateSolver.hpp"
      14             : #include "Static/Solution.hpp"
      15             : #include "Static/algos/AllOrNothing.hpp"
      16             : #include "Static/algos/DijkstraAoN.hpp"
      17             : #include "Static/algos/FrankWolfe.hpp"
      18             : 
      19             : using namespace std;
      20             : using namespace Static;
      21             : 
      22             : typedef NetworkDifferentiable::Node Node;
      23             : typedef NetworkDifferentiable::Edge Edge;
      24             : 
      25             : typedef chrono::high_resolution_clock hrc;
      26             : 
      27             : const double EPSILON = 0.1;
      28             : 
      29           0 : ConjugateFrankWolfe::ConjugateFrankWolfe(
      30             :     AllOrNothing          &aon_,
      31             :     Opt::UnivariateSolver &solver_,
      32             :     Log::ProgressLogger   &logger_
      33           0 : ):
      34           0 :     FrankWolfe(aon_, solver_, logger_) {}
      35             : 
      36           0 : Solution ConjugateFrankWolfe::solve(
      37             :     const NetworkDifferentiable &network,
      38             :     const Demand                &dem,
      39             :     const Solution              &startingSolution
      40             : ) {
      41           0 :     supplyDifferentiable = &network;
      42             : 
      43           0 :     return FrankWolfe::solve(network, dem, startingSolution);
      44             : }
      45             : 
      46           0 : double ConjugateFrankWolfe::getExpectedIterations() {
      47           0 :     double linearWithIterations = pow(-log10(epsilon / zn), 12);  // This variable has a linear relation with number of iterations
      48           0 :     double expectedIterations   = linearWithIterations / 272014.433647;
      49           0 :     return min(expectedIterations, (double)iterations);
      50             : }
      51             : 
      52           0 : Solution ConjugateFrankWolfe::step1() {
      53           0 :     SolutionBase xAoN = aon.solve(*supply, *demand, xn);
      54             : 
      55           0 :     unordered_set<Edge::ID>        edgeIDs;
      56           0 :     const unordered_set<Edge::ID> &xnEdges   = xn.getEdges();
      57           0 :     const unordered_set<Edge::ID> &xAoNEdges = xAoN.getEdges();
      58           0 :     edgeIDs.insert(xnEdges.begin(), xnEdges.end());
      59           0 :     edgeIDs.insert(xAoNEdges.begin(), xAoNEdges.end());
      60             : 
      61             :     // Update lower bound
      62           0 :     Time zApprox = zn;
      63           0 :     for(const Edge::ID &eid: edgeIDs) {
      64           0 :         Network::Edge &e = supply->getEdge(eid);
      65             : 
      66           0 :         Flow xna   = xn.getFlowInEdge(eid);
      67           0 :         Flow xAoNa = xAoN.getFlowInEdge(eid);
      68             : 
      69           0 :         zApprox += e.calculateCost(xn) * (xAoNa - xna);
      70             :     }
      71           0 :     lowerBound = max(lowerBound, zApprox);
      72             : 
      73             :     // Get xStar
      74             : 
      75             :     // Store previous xStar value in xStarStar
      76           0 :     const Solution &xStarStar = xStar;
      77             : 
      78           0 :     double top = 0.0, bot = 0.0;
      79           0 :     for(const Edge::ID &eid: edgeIDs) {
      80           0 :         Edge &e = supplyDifferentiable->getEdge(eid);
      81             : 
      82           0 :         Flow xna        = xn.getFlowInEdge(eid);
      83           0 :         Flow xAoNa      = xAoN.getFlowInEdge(eid);
      84           0 :         Flow xStarStara = xStarStar.getFlowInEdge(eid);
      85             : 
      86           0 :         top += (xStarStara - xna) * (xAoNa - xna) * e.calculateCostDerivative(xn);
      87           0 :         bot += (xStarStara - xna) * (xAoNa - xStarStara) * e.calculateCostDerivative(xn);
      88             :     }
      89             : 
      90           0 :     double a;
      91             : 
      92           0 : #pragma GCC diagnostic push
      93           0 : #pragma GCC diagnostic ignored "-Wfloat-equal"
      94           0 :     if(bot == 0.0)
      95           0 :         a = 0.0;
      96             :     else
      97           0 :         a = top / bot;
      98           0 : #pragma GCC diagnostic pop
      99             : 
     100           0 :     a = max(0.0, min(1.0 - EPSILON, a));
     101             : 
     102           0 :     xStar = Solution::interpolate(xAoN, xStarStar, a);
     103             : 
     104           0 :     xStar.materialize();
     105             : 
     106           0 :     return xStar;
     107             : }