max_plus_convolution.hpp¶
计算两个凹序列的 max-plus 卷积；用于具有四边形不等式/凹性的 DP 合并。
Complexity: Time: O(n + m) for the supported convex/concave cases. Space: O(n + m) output and monotone-optimum workspace.
AC 记录：min_plus_convolution_concave_arbitrary, min_plus_convolution_convex_arbitrary, min_plus_convolution_convex_convex。
Implementation¶
当前头文件，省略 include guard；依赖见 #include。
/// @complexity Time: O(n + m) for the supported convex/concave cases.
/// Space: O(n + m) output and monotone-optimum workspace.
#include "noya/smawk.hpp"
#include <algorithm>
#include <cassert>
#include <cstdint>
#include <functional>
#include <limits>
#include <vector>
namespace noya {
/// @brief Max-plus convolution of two concave sequences.
template <class T>
std::vector<T> two_concave_maxplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
if (a.empty())
return b;
if (b.empty())
return a;
const int n = int(a.size());
const int m = int(b.size());
int p = 0, q = 0;
std::vector<T> c(n + m - 1);
c[0] = a[0] + b[0];
for (int i = 1; i < n + m - 1; i++) {
if (p + 1 == n) {
q++;
} else if (q + 1 == m) {
p++;
} else {
if (a[p + 1] - a[p] > b[q + 1] - b[q]) {
p++;
} else {
q++;
}
}
c[i] = a[p] + b[q];
}
return c;
}
/// @brief Min-plus convolution of two convex sequences.
/// The first differences of a convex sequence are nondecreasing. Merging the
/// two difference sequences therefore describes, in order, every step of the
/// lower boundary of their min-plus convolution.
template <class T>
std::vector<T> convex_convex_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
if (a.empty())
return b;
if (b.empty())
return a;
const int n = int(a.size());
const int m = int(b.size());
std::vector<T> na(n);
std::vector<T> nb(m);
for (int i = 0; i < n; i++) {
na[i] = -a[i];
}
for (int i = 0; i < m; i++) {
nb[i] = -b[i];
}
auto nc = two_concave_maxplus_convolution(na, nb);
std::vector<T> c(n + m - 1);
for (int i = 0; i < n + m - 1; i++)
c[i] = -nc[i];
return c;
}
/// @brief Backward-compatible name for convex-convex min-plus convolution.
template <class T>
std::vector<T> two_concave_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
return convex_convex_minplus_convolution(a, b);
}
/// @brief Max-plus convolution where b is concave.
/// Concavity makes the implicit matrix a[i] + b[row-i] totally monotone, so
/// SMAWK finds all row maxima after only O(n + m) value comparisons.
template <class T>
std::vector<T> concave_maxplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
if (a.empty())
return b;
if (b.empty())
return a;
const int n = int(a.size());
const int m = int(b.size());
const auto get = [&](const int &i, const int &j) -> T {
return a[j] + b[i - j];
};
const auto sel = [&](const int &i, const int &j, const int &k) -> bool {
if (i < k)
return false;
if (i - j >= m)
return true;
return get(i, j) <= get(i, k);
};
const auto amx = smawk(n + m - 1, n, sel);
std::vector<T> c(n + m - 1);
for (int i = 0; i < n + m - 1; i++)
c[i] = get(i, amx[i]);
return c;
}
/// @brief Min-plus convolution of an arbitrary sequence and a convex sequence.
/// Convexity makes the implicit matrix a[i] + b[row-i] totally monotone.
/// Negating both inputs turns row minima into row maxima, which SMAWK finds in
/// linear time without materializing the matrix.
template <class T>
std::vector<T> arbitrary_convex_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
if (a.empty())
return b;
if (b.empty())
return a;
const int n = int(a.size());
const int m = int(b.size());
std::vector<T> na(n);
std::vector<T> nb(m);
for (int i = 0; i < n; i++) {
na[i] = -a[i];
}
for (int i = 0; i < m; i++) {
nb[i] = -b[i];
}
auto nc = concave_maxplus_convolution(na, nb);
std::vector<T> c(n + m - 1);
for (int i = 0; i < n + m - 1; i++)
c[i] = -nc[i];
return c;
}
/// @brief Min-plus convolution of a convex sequence and an arbitrary sequence.
template <class T>
std::vector<T> convex_arbitrary_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
return arbitrary_convex_minplus_convolution(b, a);
}
namespace max_plus_convolution_internal {
template <class Value>
std::vector<int> monotone_row_minima(int rs, int cs, Value val) {
std::vector<int> res(rs);
int str = 1;
while (str < rs) {
str <<= 1;
}
for (; str > 0; str >>= 1) {
for (int row = str - 1; row < rs; row += 2 * str) {
int lhs = row >= str ? res[row - str] : 0;
int las = row + str < rs ? res[row + str] : cs - 1;
res[row] = lhs;
for (int col = lhs + 1; col <= las; col++) {
if (val(row, col) < val(row, res[row])) {
res[row] = col;
}
}
}
}
return res;
}
template <class T>
std::vector<T>
arbitrary_concave_minplus_convolution(const std::vector<T> &arb,
const std::vector<T> &con) {
if (arb.empty()) {
return con;
}
if (con.empty()) {
return arb;
}
for (int i = 0; i + 2 < int(con.size()); i++) {
assert(con[i + 1] - con[i] >= con[i + 2] - con[i + 1]);
}
const int wid = int(arb.size());
const int sz = int(con.size());
const int hei = wid + sz - 1;
std::vector<int> lc(hei), mc(hei, wid - 1);
for (int row = sz; row < hei; row++) {
lc[row] = row - sz + 1;
}
for (int row = 0; row <= hei - sz; row++) {
mc[row] = row;
}
std::vector<int> lr(wid), mr(wid);
for (int col = 0; col < wid; col++) {
lr[col] = col;
mr[col] = sz - 1 + col;
}
std::vector<T> res(hei, std::numeric_limits<T>::max());
std::function<void(int, int, int, int)> sol = [&](int rl, int lro, int cl,
int cr) {
if (mc[rl] >= cr && cl >= lc[lro]) {
auto val = [&](int ri, int rc) {
int col = cr - rc;
int row = rl + ri;
return arb[col] + con[row - col];
};
auto mnm = monotone_row_minima(lro - rl + 1, cr - cl + 1, val);
for (int row = rl; row <= lro; row++) {
res[row] = std::min(res[row], val(row - rl, mnm[row - rl]));
}
return;
}
if (std::int64_t(lro - rl) * (cr - cl) < 1024) {
for (int row = rl; row <= lro; row++) {
int fro = std::max(lc[row], cl);
int to = std::min(mc[row], cr);
for (int col = fro; col <= to; col++) {
res[row] = std::min(res[row], arb[col] + con[row - col]);
}
}
return;
}
if (lro - rl > cr - cl) {
int mid = (rl + lro) / 2;
int nl = std::min(mc[mid], cr);
if (cl <= nl) {
sol(rl, mid, cl, nl);
}
int fst = std::max(lc[mid], cl);
if (fst <= cr) {
sol(mid + 1, lro, fst, cr);
}
} else {
int mid = (cl + cr) / 2;
int nl = std::min(mr[mid], lro);
if (rl <= nl) {
sol(rl, nl, cl, mid);
}
int fst = std::max(lr[mid], rl);
if (fst <= lro) {
sol(fst, lro, mid + 1, cr);
}
}
};
sol(0, hei - 1, 0, wid - 1);
return res;
}
} // namespace max_plus_convolution_internal
/// @brief Min-plus convolution of a concave sequence and an arbitrary
/// sequence. Valid pairs form a diagonal staircase rather than one rectangular
/// Monge matrix. Recursively splitting that staircase produces fully valid
/// rectangles; after reversing their columns, concavity makes row minima
/// monotone and they are found together. Small boundary rectangles are scanned
/// directly.
template <class T>
std::vector<T> concave_arbitrary_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
return max_plus_convolution_internal::arbitrary_concave_minplus_convolution(
b, a);
}
/// @brief Backward-compatible name for arbitrary-convex min-plus convolution.
template <class T>
std::vector<T> concave_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
return arbitrary_convex_minplus_convolution(a, b);
}
} // namespace noya
#ifndef NOYA_MAXPLUS_CONVOLUTION_HPP
#define NOYA_MAXPLUS_CONVOLUTION_HPP 1
/// @complexity Time: O(n + m) for the supported convex/concave cases.
/// Space: O(n + m) output and monotone-optimum workspace.
#include "noya/smawk.hpp"
#include <algorithm>
#include <cassert>
#include <cstdint>
#include <functional>
#include <limits>
#include <vector>
namespace noya {
/// @brief Max-plus convolution of two concave sequences.
template <class T>
std::vector<T> two_concave_maxplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
if (a.empty())
return b;
if (b.empty())
return a;
const int n = int(a.size());
const int m = int(b.size());
int p = 0, q = 0;
std::vector<T> c(n + m - 1);
c[0] = a[0] + b[0];
for (int i = 1; i < n + m - 1; i++) {
if (p + 1 == n) {
q++;
} else if (q + 1 == m) {
p++;
} else {
if (a[p + 1] - a[p] > b[q + 1] - b[q]) {
p++;
} else {
q++;
}
}
c[i] = a[p] + b[q];
}
return c;
}
/// @brief Min-plus convolution of two convex sequences.
/// The first differences of a convex sequence are nondecreasing. Merging the
/// two difference sequences therefore describes, in order, every step of the
/// lower boundary of their min-plus convolution.
template <class T>
std::vector<T> convex_convex_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
if (a.empty())
return b;
if (b.empty())
return a;
const int n = int(a.size());
const int m = int(b.size());
std::vector<T> na(n);
std::vector<T> nb(m);
for (int i = 0; i < n; i++) {
na[i] = -a[i];
}
for (int i = 0; i < m; i++) {
nb[i] = -b[i];
}
auto nc = two_concave_maxplus_convolution(na, nb);
std::vector<T> c(n + m - 1);
for (int i = 0; i < n + m - 1; i++)
c[i] = -nc[i];
return c;
}
/// @brief Backward-compatible name for convex-convex min-plus convolution.
template <class T>
std::vector<T> two_concave_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
return convex_convex_minplus_convolution(a, b);
}
/// @brief Max-plus convolution where b is concave.
/// Concavity makes the implicit matrix a[i] + b[row-i] totally monotone, so
/// SMAWK finds all row maxima after only O(n + m) value comparisons.
template <class T>
std::vector<T> concave_maxplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
if (a.empty())
return b;
if (b.empty())
return a;
const int n = int(a.size());
const int m = int(b.size());
const auto get = [&](const int &i, const int &j) -> T {
return a[j] + b[i - j];
};
const auto sel = [&](const int &i, const int &j, const int &k) -> bool {
if (i < k)
return false;
if (i - j >= m)
return true;
return get(i, j) <= get(i, k);
};
const auto amx = smawk(n + m - 1, n, sel);
std::vector<T> c(n + m - 1);
for (int i = 0; i < n + m - 1; i++)
c[i] = get(i, amx[i]);
return c;
}
/// @brief Min-plus convolution of an arbitrary sequence and a convex sequence.
/// Convexity makes the implicit matrix a[i] + b[row-i] totally monotone.
/// Negating both inputs turns row minima into row maxima, which SMAWK finds in
/// linear time without materializing the matrix.
template <class T>
std::vector<T> arbitrary_convex_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
if (a.empty())
return b;
if (b.empty())
return a;
const int n = int(a.size());
const int m = int(b.size());
std::vector<T> na(n);
std::vector<T> nb(m);
for (int i = 0; i < n; i++) {
na[i] = -a[i];
}
for (int i = 0; i < m; i++) {
nb[i] = -b[i];
}
auto nc = concave_maxplus_convolution(na, nb);
std::vector<T> c(n + m - 1);
for (int i = 0; i < n + m - 1; i++)
c[i] = -nc[i];
return c;
}
/// @brief Min-plus convolution of a convex sequence and an arbitrary sequence.
template <class T>
std::vector<T> convex_arbitrary_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
return arbitrary_convex_minplus_convolution(b, a);
}
namespace max_plus_convolution_internal {
template <class Value>
std::vector<int> monotone_row_minima(int rs, int cs, Value val) {
std::vector<int> res(rs);
int str = 1;
while (str < rs) {
str <<= 1;
}
for (; str > 0; str >>= 1) {
for (int row = str - 1; row < rs; row += 2 * str) {
int lhs = row >= str ? res[row - str] : 0;
int las = row + str < rs ? res[row + str] : cs - 1;
res[row] = lhs;
for (int col = lhs + 1; col <= las; col++) {
if (val(row, col) < val(row, res[row])) {
res[row] = col;
}
}
}
}
return res;
}
template <class T>
std::vector<T>
arbitrary_concave_minplus_convolution(const std::vector<T> &arb,
const std::vector<T> &con) {
if (arb.empty()) {
return con;
}
if (con.empty()) {
return arb;
}
for (int i = 0; i + 2 < int(con.size()); i++) {
assert(con[i + 1] - con[i] >= con[i + 2] - con[i + 1]);
}
const int wid = int(arb.size());
const int sz = int(con.size());
const int hei = wid + sz - 1;
std::vector<int> lc(hei), mc(hei, wid - 1);
for (int row = sz; row < hei; row++) {
lc[row] = row - sz + 1;
}
for (int row = 0; row <= hei - sz; row++) {
mc[row] = row;
}
std::vector<int> lr(wid), mr(wid);
for (int col = 0; col < wid; col++) {
lr[col] = col;
mr[col] = sz - 1 + col;
}
std::vector<T> res(hei, std::numeric_limits<T>::max());
std::function<void(int, int, int, int)> sol = [&](int rl, int lro, int cl,
int cr) {
if (mc[rl] >= cr && cl >= lc[lro]) {
auto val = [&](int ri, int rc) {
int col = cr - rc;
int row = rl + ri;
return arb[col] + con[row - col];
};
auto mnm = monotone_row_minima(lro - rl + 1, cr - cl + 1, val);
for (int row = rl; row <= lro; row++) {
res[row] = std::min(res[row], val(row - rl, mnm[row - rl]));
}
return;
}
if (std::int64_t(lro - rl) * (cr - cl) < 1024) {
for (int row = rl; row <= lro; row++) {
int fro = std::max(lc[row], cl);
int to = std::min(mc[row], cr);
for (int col = fro; col <= to; col++) {
res[row] = std::min(res[row], arb[col] + con[row - col]);
}
}
return;
}
if (lro - rl > cr - cl) {
int mid = (rl + lro) / 2;
int nl = std::min(mc[mid], cr);
if (cl <= nl) {
sol(rl, mid, cl, nl);
}
int fst = std::max(lc[mid], cl);
if (fst <= cr) {
sol(mid + 1, lro, fst, cr);
}
} else {
int mid = (cl + cr) / 2;
int nl = std::min(mr[mid], lro);
if (rl <= nl) {
sol(rl, nl, cl, mid);
}
int fst = std::max(lr[mid], rl);
if (fst <= lro) {
sol(fst, lro, mid + 1, cr);
}
}
};
sol(0, hei - 1, 0, wid - 1);
return res;
}
} // namespace max_plus_convolution_internal
/// @brief Min-plus convolution of a concave sequence and an arbitrary
/// sequence. Valid pairs form a diagonal staircase rather than one rectangular
/// Monge matrix. Recursively splitting that staircase produces fully valid
/// rectangles; after reversing their columns, concavity makes row minima
/// monotone and they are found together. Small boundary rectangles are scanned
/// directly.
template <class T>
std::vector<T> concave_arbitrary_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
return max_plus_convolution_internal::arbitrary_concave_minplus_convolution(
b, a);
}
/// @brief Backward-compatible name for arbitrary-convex min-plus convolution.
template <class T>
std::vector<T> concave_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
return arbitrary_convex_minplus_convolution(a, b);
}
} // namespace noya
#endif // NOYA_MAXPLUS_CONVOLUTION_HPP
#include <algorithm>
#include <cassert>
#include <cstdint>
#include <functional>
#include <limits>
#include <numeric>
#include <vector>
/// @complexity Time: O(n + m) for the supported convex/concave cases.
/// Space: O(n + m) output and monotone-optimum workspace.
/// @complexity Time: O(rows + columns) matrix probes.
/// Space: O(rows + columns).
namespace noya {
/// @brief SMAWK algorithm: compute row minima of a totally monotone matrix.
/// A stack reduction leaves at most one candidate column per row, recursion
/// solves the odd rows, and monotone argmins bound the scan that interpolates
/// each even row. Every row and column is discarded or scanned only O(1) times.
/// @return Vector where ans[i] is the column index of the minimum in row i.
template <class Select>
std::vector<int> smawk(const int nr, const int nc, const Select &sel) {
const std::function<std::vector<int>(const std::vector<int> &,
const std::vector<int> &)>
sol = [&](const std::vector<int> &row,
const std::vector<int> &col) -> std::vector<int> {
const int n = int(row.size());
if (n == 0)
return {};
std::vector<int> c2;
for (const int i : col) {
while (!c2.empty() && sel(row[c2.size() - 1], c2.back(), i))
c2.pop_back();
if (c2.size() < n)
c2.push_back(i);
}
std::vector<int> r2;
for (int i = 1; i < n; i += 2)
r2.push_back(row[i]);
const std::vector<int> a2 = sol(r2, c2);
std::vector<int> ans(n);
for (int i = 0; i != a2.size(); i += 1)
ans[i * 2 + 1] = a2[i];
int j = 0;
for (int i = 0; i < n; i += 2) {
ans[i] = c2[j];
const int end = i + 1 == n ? c2.back() : ans[i + 1];
while (c2[j] != end) {
j += 1;
if (sel(row[i], ans[i], c2[j]))
ans[i] = c2[j];
}
}
return ans;
};
std::vector<int> row(nr);
std::iota(row.begin(), row.end(), 0);
std::vector<int> col(nc);
std::iota(col.begin(), col.end(), 0);
return sol(row, col);
}
} // namespace noya
namespace noya {
/// @brief Max-plus convolution of two concave sequences.
template <class T>
std::vector<T> two_concave_maxplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
if (a.empty())
return b;
if (b.empty())
return a;
const int n = int(a.size());
const int m = int(b.size());
int p = 0, q = 0;
std::vector<T> c(n + m - 1);
c[0] = a[0] + b[0];
for (int i = 1; i < n + m - 1; i++) {
if (p + 1 == n) {
q++;
} else if (q + 1 == m) {
p++;
} else {
if (a[p + 1] - a[p] > b[q + 1] - b[q]) {
p++;
} else {
q++;
}
}
c[i] = a[p] + b[q];
}
return c;
}
/// @brief Min-plus convolution of two convex sequences.
/// The first differences of a convex sequence are nondecreasing. Merging the
/// two difference sequences therefore describes, in order, every step of the
/// lower boundary of their min-plus convolution.
template <class T>
std::vector<T> convex_convex_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
if (a.empty())
return b;
if (b.empty())
return a;
const int n = int(a.size());
const int m = int(b.size());
std::vector<T> na(n);
std::vector<T> nb(m);
for (int i = 0; i < n; i++) {
na[i] = -a[i];
}
for (int i = 0; i < m; i++) {
nb[i] = -b[i];
}
auto nc = two_concave_maxplus_convolution(na, nb);
std::vector<T> c(n + m - 1);
for (int i = 0; i < n + m - 1; i++)
c[i] = -nc[i];
return c;
}
/// @brief Backward-compatible name for convex-convex min-plus convolution.
template <class T>
std::vector<T> two_concave_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
return convex_convex_minplus_convolution(a, b);
}
/// @brief Max-plus convolution where b is concave.
/// Concavity makes the implicit matrix a[i] + b[row-i] totally monotone, so
/// SMAWK finds all row maxima after only O(n + m) value comparisons.
template <class T>
std::vector<T> concave_maxplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
if (a.empty())
return b;
if (b.empty())
return a;
const int n = int(a.size());
const int m = int(b.size());
const auto get = [&](const int &i, const int &j) -> T {
return a[j] + b[i - j];
};
const auto sel = [&](const int &i, const int &j, const int &k) -> bool {
if (i < k)
return false;
if (i - j >= m)
return true;
return get(i, j) <= get(i, k);
};
const auto amx = smawk(n + m - 1, n, sel);
std::vector<T> c(n + m - 1);
for (int i = 0; i < n + m - 1; i++)
c[i] = get(i, amx[i]);
return c;
}
/// @brief Min-plus convolution of an arbitrary sequence and a convex sequence.
/// Convexity makes the implicit matrix a[i] + b[row-i] totally monotone.
/// Negating both inputs turns row minima into row maxima, which SMAWK finds in
/// linear time without materializing the matrix.
template <class T>
std::vector<T> arbitrary_convex_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
if (a.empty())
return b;
if (b.empty())
return a;
const int n = int(a.size());
const int m = int(b.size());
std::vector<T> na(n);
std::vector<T> nb(m);
for (int i = 0; i < n; i++) {
na[i] = -a[i];
}
for (int i = 0; i < m; i++) {
nb[i] = -b[i];
}
auto nc = concave_maxplus_convolution(na, nb);
std::vector<T> c(n + m - 1);
for (int i = 0; i < n + m - 1; i++)
c[i] = -nc[i];
return c;
}
/// @brief Min-plus convolution of a convex sequence and an arbitrary sequence.
template <class T>
std::vector<T> convex_arbitrary_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
return arbitrary_convex_minplus_convolution(b, a);
}
namespace max_plus_convolution_internal {
template <class Value>
std::vector<int> monotone_row_minima(int rs, int cs, Value val) {
std::vector<int> res(rs);
int str = 1;
while (str < rs) {
str <<= 1;
}
for (; str > 0; str >>= 1) {
for (int row = str - 1; row < rs; row += 2 * str) {
int lhs = row >= str ? res[row - str] : 0;
int las = row + str < rs ? res[row + str] : cs - 1;
res[row] = lhs;
for (int col = lhs + 1; col <= las; col++) {
if (val(row, col) < val(row, res[row])) {
res[row] = col;
}
}
}
}
return res;
}
template <class T>
std::vector<T>
arbitrary_concave_minplus_convolution(const std::vector<T> &arb,
const std::vector<T> &con) {
if (arb.empty()) {
return con;
}
if (con.empty()) {
return arb;
}
for (int i = 0; i + 2 < int(con.size()); i++) {
assert(con[i + 1] - con[i] >= con[i + 2] - con[i + 1]);
}
const int wid = int(arb.size());
const int sz = int(con.size());
const int hei = wid + sz - 1;
std::vector<int> lc(hei), mc(hei, wid - 1);
for (int row = sz; row < hei; row++) {
lc[row] = row - sz + 1;
}
for (int row = 0; row <= hei - sz; row++) {
mc[row] = row;
}
std::vector<int> lr(wid), mr(wid);
for (int col = 0; col < wid; col++) {
lr[col] = col;
mr[col] = sz - 1 + col;
}
std::vector<T> res(hei, std::numeric_limits<T>::max());
std::function<void(int, int, int, int)> sol = [&](int rl, int lro, int cl,
int cr) {
if (mc[rl] >= cr && cl >= lc[lro]) {
auto val = [&](int ri, int rc) {
int col = cr - rc;
int row = rl + ri;
return arb[col] + con[row - col];
};
auto mnm = monotone_row_minima(lro - rl + 1, cr - cl + 1, val);
for (int row = rl; row <= lro; row++) {
res[row] = std::min(res[row], val(row - rl, mnm[row - rl]));
}
return;
}
if (std::int64_t(lro - rl) * (cr - cl) < 1024) {
for (int row = rl; row <= lro; row++) {
int fro = std::max(lc[row], cl);
int to = std::min(mc[row], cr);
for (int col = fro; col <= to; col++) {
res[row] = std::min(res[row], arb[col] + con[row - col]);
}
}
return;
}
if (lro - rl > cr - cl) {
int mid = (rl + lro) / 2;
int nl = std::min(mc[mid], cr);
if (cl <= nl) {
sol(rl, mid, cl, nl);
}
int fst = std::max(lc[mid], cl);
if (fst <= cr) {
sol(mid + 1, lro, fst, cr);
}
} else {
int mid = (cl + cr) / 2;
int nl = std::min(mr[mid], lro);
if (rl <= nl) {
sol(rl, nl, cl, mid);
}
int fst = std::max(lr[mid], rl);
if (fst <= lro) {
sol(fst, lro, mid + 1, cr);
}
}
};
sol(0, hei - 1, 0, wid - 1);
return res;
}
} // namespace max_plus_convolution_internal
/// @brief Min-plus convolution of a concave sequence and an arbitrary
/// sequence. Valid pairs form a diagonal staircase rather than one rectangular
/// Monge matrix. Recursively splitting that staircase produces fully valid
/// rectangles; after reversing their columns, concavity makes row minima
/// monotone and they are found together. Small boundary rectangles are scanned
/// directly.
template <class T>
std::vector<T> concave_arbitrary_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
return max_plus_convolution_internal::arbitrary_concave_minplus_convolution(
b, a);
}
/// @brief Backward-compatible name for arbitrary-convex min-plus convolution.
template <class T>
std::vector<T> concave_minplus_convolution(const std::vector<T> &a,
const std::vector<T> &b) {
return arbitrary_convex_minplus_convolution(a, b);
}
} // namespace noya