feat: max-plus convolution

content/convolutions/MaxPlusConvolution.hpp

@@ -0,0 +1,104 @@
+#pragma once
+
+/**
+ * Author: Teetat T.
+ * Date: 2024-09-01
+ * Description: Max Plus Convolution. Find $C[k] = max_{i+j=k}\{A[i] + B[j]\}$.
+ * Time: $O(N)$.
+ */
+
+// SMAWCK algorithm for finding row-wise maxima.
+// f(i,j,k) checks if M[i][j] <= M[i][k].
+// f(i,j,k) checks if M[i][k] is at least as good as M[i][j].
+// higher is better.
+template<class F>
+vector<int> smawck(const F &f,const vector<int> &rows,const vector<int> &cols){
+    int n=(int)rows.size(),m=(int)cols.size();
+    if(max(n,m)<=2){
+        vector<int> ans(n,-1);
+        for(int i=0;i<n;i++){
+            for(int j:cols){
+                if(ans[i]==-1||f(rows[i],ans[i],j)){
+                    ans[i]=j;
+                }
+            }
+        }
+        return ans;
+    }
+    if(n<m){
+        // reduce
+        vector<int> st;
+        for(int j:cols){
+            while(true){
+                if(st.empty()){
+                    st.emplace_back(j);
+                    break;
+                }else if(f(rows[(int)st.size()-1],st.back(),j)){
+                    st.pop_back();
+                }else if(st.size()<n){
+                    st.emplace_back(j);
+                    break;
+                }else{
+                    break;
+                }
+            }
+        }
+        return smawck(f,rows,st);
+    }   
+    vector<int> ans(n,-1);
+    vector<int> new_rows;
+    for(int i=1;i<n;i+=2){
+        new_rows.emplace_back(rows[i]);
+    }
+    auto res=smawck(f,new_rows,cols);
+    for(int i=0;i<new_rows.size();i++){
+        ans[2*i+1]=res[i];
+    }
+    for(int i=0,l=0,r=0;i<n;i+=2){
+        if(i+1==n)r=m;
+        while(r<m&&cols[r]<=ans[i+1])r++;
+        ans[i]=cols[l++];
+        for(;l<r;l++){
+            if(f(rows[i],ans[i],cols[l])){
+                ans[i]=cols[l];
+            }
+        }
+        l--;
+    }
+    return ans;
+}
+
+template<class F>
+vector<int> smawck(const F &f,int n,int m){
+    vector<int> rows(n),cols(m);
+    iota(rows.begin(),rows.end(),0);
+    iota(cols.begin(),cols.end(),0);
+    return smawck(f,rows,cols);
+}
+
+// Max Plus Convolution.
+// b must be convex, i.e. b[i]-b[i-1]>=b[i+1]-b[i].
+template<class T>
+vector<T> max_plus_convolution_arbitary_convex(vector<T> a,const vector<T> &b){
+    if(a.empty()||b.empty())return {};
+    if((int)b.size()==1){
+        for(auto &x:a)x+=b[0];
+        return a;
+    }
+    int n=(int)a.size(),m=(int)b.size();
+    auto f=[&](int i,int j){
+        return a[j]+b[i-j];
+    };
+    auto cmp=[&](int i,int j,int k){
+        if(i<k)return false;
+        if(i-j>=m)return true;
+        return f(i,j)<=f(i,k);
+    };
+    auto best=smawck(cmp,n+m-1,n);
+    vector<T> ans(n+m-1);
+    for(int i=0;i<n+m-1;i++){
+        ans[i]=f(i,best[i]);
+    }
+    return ans;
+}
+

content/convolutions/chapter.tex

@@ -4,4 +4,5 @@ \chapter{Convolutions}
 \kactlimport{GCDConvolution.hpp}
 \kactlimport{LCMConvolution.hpp}
 \kactlimport{ORConvolution.hpp}
-\kactlimport{XORConvolution.hpp}
+\kactlimport{XORConvolution.hpp}
+\kactlimport{MaxPlusConvolution.hpp}