Gas Station Problem - cheapest and least amount of stations

Question

I am working on a problem that consists of the following: You are driving a car with a certain fuel usage m (in our example we will take 8l/100km) and you are driving a straight of length x (example: 1000km). The car starts out with a certain amount of fuel f (example: 22l). Your car has a fuel tank of size g (example: 55l) and there are gas stations (which have a price per liter) plotted around the straight (e.g 100km (1.45$/l), 400km(1.40$/l) and 900km(1.22$/l). The aim of the algorithm I'm having a hard time creating is: With the least amount of stops (so not the cheapest route, but the one with the least stops at gas stations) find the cheapest way and tell the user how much liters he has to tank at which gas station and how much it will cost.

Currently using recursion and for loops (presumably O(n^2)) I've created an algorithm that can solve some problems up to a certain complexity, it starts to struggle when there are about 100 gas stations.

How my algorithm works:

• Go to the tank stations that are available from the start (the 22l in the example)
• Go to each tank station and list the tank stations (or the end) in range by having a full tank (since the car can fuel up at the tank station you can have a full tank) I do save this in a list so it's not calculated twice.
• Then for loop every one of those tank stations that can be reached and recursion occurs, pretty much I then save the smallest amount of stops required and rinse and repeat and voila I get the smallest answer which is (in our example) stopping at 100 fueling 10.00 liters, paying 14.5$ and then stopping at 400 fueling 48 liters and paying 67.20$

The problems I still have:

• How (Is it even possible)can I reduce the complexity to O(N log N) or linear, so that all (even if it's 100+ gas stations) can be checked. At the moment the recursive methods go down to 10+ recursions in recursions which makes anything above 100 gas stations pretty much unsolvable for this algorithm.

• At the moment my algorithm only fuels up as much as it requires to reach the next gas station or the end: What's the best way to handle the problem if the first gas station is cheaper than the second and you can fuel up "n liters more" so you buy less liters at the second gas station. Since in the ideal case you have 0 liters left at the end of the trip.

Extra notes:

• Arriving at a gas station with 0 liters of fuel counts as having arrived.
• EDIT: All paths of the same price and least amount of stations must be found.

My current code (snippet) I think that the methods names are self explanitory, add comment if they are not.,

    void findRoutes(List<GasStation> reachableStations, List<GasStation> previousStations) {
        int currentSteps = previousStations.size();
        if (currentSteps > leastSteps) {
            return;
        }
        // Reached the end (reachableStations will be null if it can reach the end)
        if (reachableStations == null) {
            // less steps
            if (currentSteps < leastSteps) {
                routes.clear();
                routes.add(previousStations);
                leastSteps = previousStations.size();
                return;
            } else {
                // same amount of steps
                routes.add(previousStations);
                return;
            }
        }
        // would be too many steps
        if (currentSteps + 1 > leastSteps) {
            return;
        }
        // Check those further away so we get a smaller step amount quicker
        Collections.reverse(reachableStations);
        for (GasStation reachableStation : reachableStations) {
            List<GasStation> newPrevious = new LinkedList<>(previousStations);
            newPrevious.add(reachableStation);
            findRoutes(reachableStation.getReachableGasStations(), newPrevious);
        }
    }

Answer 1

tl;dr: By following the paper mentioned in the comments, a C# implementation of a solver (below) handles the case of 500 randomly distributed stations in about 14 ms on an aging laptop, so in particular, this handles the 100 station case easily, and is orders of magnitude faster than using a MIP solver, as suggested in the comments.

Typically, the gas station problem (which we should really start calling the charging station problem, but that's a different story) assumes that the starting fuel amount is 0: The more general case may be reduced to the 0 case by adding a new starting station with free fuel at a distance to your initial starting point that causes the car to reach your initial starting point with a tank containing your given amount of fuel. This can be done without ruining the general complexity of the solution below.

Noting this, the problem boils down to that described in To Fill or not to Fill: The Gas Station Problem, as noted by @PySeeker in the comments. In particular, $O(N \log N)$ seems optimistic. In the paper, the relevant theorem handles your case in $O(\Delta N^2 \log N)$, where $\Delta$ is minimum number of stops required (which you can easily precompute in linear time if necessary). Another paper, A fast algorithm for the gas station problem, describes how to solve the problem in $O(\Delta N^2 + N^2 \log N)$, but let's just focus on the former paper here.

Its Theorem 2.2 solves the problem for a fixed $\Delta$, where you're really only interested in the lowest possible one. Since their recurrence is set up to solve the problem for increasing $\Delta$, this amounts to simply halt the algorithm once $A(s, \Delta, 0)$, in the notation of the paper, becomes finite.

Note also that compared to the general problem which handles general graphs whose edge weights form a metric (a requirement noted in the second of the above-mentioned papers but for some reason not in the first one), your situation is simpler with vertices $0, \dots, N - 1$ and distances $d_{uv} = d[v] - d[u]$.

One thing to be a bit careful when implementing the algorithm is that while the description in the paper is fine, the pseudo-code is rather buggy/lacking (cf. e.g. this question). Below we implement the various fixes needed to get the algorithm up and running, as well as a some amount of indexing to help improve performance.

You mention in your edit that besides the value of the optimal solution, you would also like to get the actual paths taken. The algorithm below only outputs the value, i.e. $A(0, \Delta, 0)$, but by keeping track of the argmin in a separate table whenever the table of values is updated, you'll immediately get the desired path as well. This is completely analogous to how you would implement e.g. Dijkstra's algorithm.

You don't specify a language in the question, so I took the liberty of writing this in C#; the code is very C'y, so it should be straightforward to port it to Java if necessary (s/List/ArrayList/g). The notation follows the paper, so let me simply refer to that for comments/documentation (but let me also apologize that without the paper at hand, the implementation is likely impossible to read).

The solution doesn't go all the way: As mentioned above, a different algorithm with better complexity exists, and it would be natural to try that one as well; it's not particularly complicated. Moreover, the implementation at hand has a natural performance optimization that's not included: It's not necessary to grow the table for all $q$; for example, the source vertex $u = 0$ will depend only on the previous row through R[0], so by precomputing the minimal value of $\Delta$, we can avoid some redundant calculation.

private static double Solve(double[] c, double[] d, double U)
{
    int n = d.Length;
    int t = n - 1;
    var R = new int[n][];
    var indep = new double[n][];
    var GV = new List<List<double>>();
    var GVar = new List<Dictionary<int, int>>();
    for (int u = 0; u < n; u++)
    {
        var l = new List<int>();
        for (int v = u + 1; v < n; v++)
        {
            if (d[v] - d[u] <= U)
                l.Add(v);
            else break;
        }

        R[u] = l.ToArray();
        indep[u] = new double[l.Count];
    }

    for (int u = 0; u < n; u++)
    {
        var l = new List<double> { 0 };
        var gvar = new Dictionary<int, int>();
        int i = 1;
        for (int w = 0; w < u; w++)
        {
            if (c[w] < c[u] && d[u] - d[w] <= U)
            {
                l.Add(U - (d[u] - d[w]));
                gvar[w] = i++;
            }
        }

        GV.Add(l);
        GVar.Add(gvar);
    }

    int q = 0;
    var Cq1 = new double[n][];
    var Cq2 = new double[n][];
    for (int i = 0; i < n; i++)
    {
        Cq1[i] = new double[GV[i].Count];
        Cq2[i] = new double[GV[i].Count];
        for (int j = 0; j < GV[i].Count; j++)
        {
            Cq1[i][j] = double.PositiveInfinity;
            Cq2[i][j] = double.PositiveInfinity;
        }
    }

    var toggle = true;
    while (true)
    {
        var Cq = toggle ? Cq1 : Cq2;
        var prev = !toggle ? Cq1 : Cq2;
        toggle = !toggle;
        for (int i = 0; i < n; i++)
        {
            for (int j = 0; j < GV[i].Count; j++)
                Cq[i][j] = double.PositiveInfinity;
        }
        for (int u = 0; u < n; u++)
        {
            Grow(u, q, t, c, d, U, R[u], GV[u], q == 0 ? null : prev, Cq, indep[u], GVar);
            if (u == 0 && !double.IsPositiveInfinity(Cq[u][0]))
                return Cq[u][0];
        }
        q++;
    }
}

private static void Grow(int u, int q, int t, double[] c, double[] d, double U,
        int[] r, List<double> gv, double[][] prev, double[][] ret, double[] indep,
        List<Dictionary<int, int>> GVar)
{
    double cost = c[u];
    if (q == 0 || u == t)
    {
        for (int i = 0; i < gv.Count; i++)
        {
            var g = gv[i];
            if (q == 0 && g <= d[t] - d[u] && d[t] - d[u] <= U)
                ret[u][i] = (d[t] - d[u] - g) * cost;
        }
        return;
    }

    for (var i = 0; i < r.Length; i++)
    {
        var v = r[i];
        indep[i] = c[v] <= cost ? prev[v][0] + (d[v] - d[u]) * cost : prev[v][GVar[v][u]] + U * cost;
    }

    Array.Sort(indep, r);
    var j = 0;
    var w = r[j];
    for (int gi = 0; gi < gv.Count; gi++)
    {
        var g = gv[gi];
        while (g > d[w] - d[u] && c[w] <= cost)
        {
            j++;
            if (j == r.Length) return;
            w = r[j];
        }

        ret[u][gi] = indep[j] - g * cost;
    }
}

Example usage, and benchmark on a 500 station case:

static void Main(string[] args)
{
    var rng = new Random();
    var sw = new Stopwatch();
    for (int k = 0; k < 100; k++)
    {
        int n = 500;
        var prices = Enumerable.Range(1, n).Select(_ => rng.NextDouble()).ToArray();
        var distances = Enumerable.Range(1, n).Select(_ => rng.NextDouble() * n).OrderBy(x => x).ToArray();
        var capacity = 15;
        sw.Start();
        var result = Solve(prices, distances, capacity);
        sw.Stop();
        var time = sw.Elapsed;
        Console.WriteLine($"{time} {result}");
        sw.Reset();
    }
}