convolution/min_plus_convolution.hpp

• View this file on GitHub
• Last update: 2025-09-02 16:52:07+09:00
• Include: #include "convolution/min_plus_convolution.hpp"

Depends on

• opt/monotone_minima.hpp

Verified with

• convolution/test/min_plus_convolution_convex_arbitrary.test.cpp
• convolution/test/min_plus_convolution_convex_convex.test.cpp

Code

#pragma once
#include "../opt/monotone_minima.hpp"
#include <cassert>

template <typename T>
std::vector<T> min_plus_convolution_convex_convex(const std::vector<T> &a, const std::vector<T> &b) {
    if (a.empty() || b.empty()) {
        return std::vector<T>(0);
    }
    int n = (int)a.size();
    int m = (int)b.size();
    std::vector<T> c(n + m - 1);
    c[0] = a[0] + b[0];
    int ita = 0, itb = 0;
    for (int i = 0; i < n + m - 2; ++i) {
        if (itb == m - 1 || (ita != n - 1 && a[ita + 1] - a[ita] < b[itb + 1] - b[itb])) {
            c[i + 1] = c[i] + (a[ita + 1] - a[ita]);
            ++ita;
        } else {
            c[i + 1] = c[i] + (b[itb + 1] - b[itb]);
            ++itb;
        }
    }
    return c;
}

template <typename T>
std::vector<T> min_plus_convolution_convex_arbitrary(const std::vector<T> &a, const std::vector<T> &b) {
    if (a.empty() || b.empty()) {
        return std::vector<T>(0);
    }
    int n = (int)a.size();
    int m = (int)b.size();
    auto select = [&](int i, int j, int k) -> bool {
        if (i < k) {
            return false;
        }
        if (i - j >= n) {
            return true;
        }
        return a[i - j] + b[j] > a[i - k] + b[k];
    };
    std::vector<int> argmin = monotone_minima(n + m - 1, m, select);
    std::vector<T> c(n + m - 1);
    for (int i = 0; i < n + m - 1; ++i) {
        c[i] = a[i - argmin[i]] + b[argmin[i]];
    }
    return c;
}

template <typename T>
bool is_convex(const std::vector<T> &a) {
    int n = (int)a.size();
    for (int i = 0; i < n - 2; ++i) {
        if (a[i + 1] - a[i] > a[i + 2] - a[i + 1]) {
            return false;
        }
    }
    return true;
}

// is_convex(a) || is_convex(b)
template <typename T>
std::vector<T> min_plus_convolution(const std::vector<T> &a, const std::vector<T> &b) {
    bool ica = is_convex(a);
    bool icb = is_convex(b);
    if (ica && icb) {
        return min_plus_convolution_convex_convex(a, b);
    } else if (ica) {
        return min_plus_convolution_convex_arbitrary(a, b);
    } else if (icb) {
        return min_plus_convolution_convex_arbitrary(b, a);
    } else {
        assert(false);
    }
}

#line 2 "opt/monotone_minima.hpp"
#include <vector>

// f(i, j, k) = true <=> (i, j) -> (i, k)
// j < k.
template <typename F>
std::vector<int> monotone_minima(int h, int w, F f) {
    std::vector<int> argmin(h);
    auto rec = [&](auto rec, int xl, int xr, int yl, int yr) -> void {
        if (xl == xr) {
            return;
        }
        int xm = (xl + xr) / 2;
        argmin[xm] = yl;
        for (int i = yl + 1; i < yr; ++i) {
            if (f(xm, argmin[xm], i)) {
                argmin[xm] = i;
            }
        }
        rec(rec, xl, xm, yl, argmin[xm] + 1);
        rec(rec, xm + 1, xr, argmin[xm], yr);
    };
    rec(rec, 0, h, 0, w);
    return argmin;
}
#line 3 "convolution/min_plus_convolution.hpp"
#include <cassert>

template <typename T>
std::vector<T> min_plus_convolution_convex_convex(const std::vector<T> &a, const std::vector<T> &b) {
    if (a.empty() || b.empty()) {
        return std::vector<T>(0);
    }
    int n = (int)a.size();
    int m = (int)b.size();
    std::vector<T> c(n + m - 1);
    c[0] = a[0] + b[0];
    int ita = 0, itb = 0;
    for (int i = 0; i < n + m - 2; ++i) {
        if (itb == m - 1 || (ita != n - 1 && a[ita + 1] - a[ita] < b[itb + 1] - b[itb])) {
            c[i + 1] = c[i] + (a[ita + 1] - a[ita]);
            ++ita;
        } else {
            c[i + 1] = c[i] + (b[itb + 1] - b[itb]);
            ++itb;
        }
    }
    return c;
}

template <typename T>
std::vector<T> min_plus_convolution_convex_arbitrary(const std::vector<T> &a, const std::vector<T> &b) {
    if (a.empty() || b.empty()) {
        return std::vector<T>(0);
    }
    int n = (int)a.size();
    int m = (int)b.size();
    auto select = [&](int i, int j, int k) -> bool {
        if (i < k) {
            return false;
        }
        if (i - j >= n) {
            return true;
        }
        return a[i - j] + b[j] > a[i - k] + b[k];
    };
    std::vector<int> argmin = monotone_minima(n + m - 1, m, select);
    std::vector<T> c(n + m - 1);
    for (int i = 0; i < n + m - 1; ++i) {
        c[i] = a[i - argmin[i]] + b[argmin[i]];
    }
    return c;
}

template <typename T>
bool is_convex(const std::vector<T> &a) {
    int n = (int)a.size();
    for (int i = 0; i < n - 2; ++i) {
        if (a[i + 1] - a[i] > a[i + 2] - a[i + 1]) {
            return false;
        }
    }
    return true;
}

// is_convex(a) || is_convex(b)
template <typename T>
std::vector<T> min_plus_convolution(const std::vector<T> &a, const std::vector<T> &b) {
    bool ica = is_convex(a);
    bool icb = is_convex(b);
    if (ica && icb) {
        return min_plus_convolution_convex_convex(a, b);
    } else if (ica) {
        return min_plus_convolution_convex_arbitrary(a, b);
    } else if (icb) {
        return min_plus_convolution_convex_arbitrary(b, a);
    } else {
        assert(false);
    }
}

Back to top page