【四边形不等式优化DP】

2020-02-25 22:05:26


四边形不等式优化DP

有二元函数 $w(i,j)$,类似1D/1D动态规划中的cost函数:

$dp[i] = min(dp[j] + cost(j,i))$

  • 【原始定义】: 如对于任意 $a\le b\le c\le d$,有:$w(a,d) + w(b,c) \ge w(a,c) + w(b,d)$ 成立 就说函数 $w(x,y)$ 满足【四边形不等式】

  • 【等价定义(相交小于包含)】 对于任意 $a< a+1\le b<b+1$,有:$w(a,b+1) + w(a+1,b) \ge w(a,b) + w(a+1,b+1)$ 就说函数 $w(i,j)$ 满足 【四边形不等式】

1->2显然成立,下面证明 2->1:

step1
对于 a<c,有1式:
w(a,c+1) + w(a+1,c) >= w(a,c) + w(a+1,c+1)
对于 a+1<c,有2式:
w(a+1,c+1) + w(a+2,c) >= w(a+1,c) + w(a+2,c+1)

两式相加,有3式:
w(a,c+1) + w(a+2,c) >= w(a,c) + w(a+2,c+1)

对比1,3式,发现 (a+1) 可以扩成 (a+2),同理可以扩成 b<=c:
w(a,c+1) + w(b,c) >= w(a,c) + w(b,c+1)

step2
对于 a<c+1,有4式:
w(a,c+2) + w(a+1,c+1) >= w(a,c+1) + w(a+1,c+2)

1式+4式,有5式:
w(a,c+2) + w(a+1,c) >= w(a,c) + w(a+1,c+2)
对比1,5式: 发现 (c+1) 可以扩成 (c+2),同理可以扩成 c<=d:
w(a,d) + w(a+1,c) >= w(a,c) + w(a+1,d)

再把上式 (a+1)扩成b: 
则有对于 a<=b<=c<=d
w(a,d) + w(b,c) >= w(a,c) + w(b,d) 
成立,原始定义得证

对于1D/1D动态规划 $dp[i] = min(dp[j] + cost(j,i))$,有以下定理:

  • 【决策单调性定理】: 如果函数 $cost(j,i)$ 满足四边形不等式,则 $dp[i]$ 有决策单调性(充分条件)

证明:

dp[i]的决策点是p[i],对于j<p[i]根据最优性:
dp[p[i]] + cost(p[i],i) <= dp[j] + cost(j,i)

当i变大为i'时,j<=p[i]<=i<=i',根据四边形不等式:
cost(j,i') + cost(p[i],i) >= cost(j,i) + cost(p[i],i')
cost(p[i],i') - cost(p[i],i) <= cost(j,i') - cost(j,i)

和dp式相加得:
dp[p[i]] + cost(p[i],i') <= dp[j] + cost(j,i')

即对于i'来说,p[i]比j优,定理得证

做题时,一般是由四边形不等式的等价定义推出决策单调性。

【P3515 [POI2011]Lightning Conductor】

https://www.luogu.com.cn/problem/P3515

根据题意列出转移方程:

p = max(a[j] - a[i] + sqrt(|i-j|)), 1<=j<=N

假设j<i,设dp状态
dp[i] = min(-a[j] - sqrt(i-j)), j<i

其中
w(j,i) = -a[j] - sqrt(i-j)

四边形不等式等价定义:
w(j,i+1) + w(j+1,i) >= w(j,i) + w(j+1,i+1)
化简得:
sqrt(i-j) - sqrt(i-j-1) >= sqrt(i-j+1)- sqrt(i-j)

令 x = (i-j),
sqrt(x) - sqrt(x-1) >= sqrt(x+1)- sqrt(x)

函数 f(x) = sqrt(x) - sqrt(x-1),需要证明其单调非递增 f(x) >= f(x+1),求导:
f'(x) = (1/2)*(1/sqrt(x) - 1/sqrt(x-1))

sqrt(x) > sqrt(x-1), 1/sqrt(x) < 1/sqrt(x-1), f'(x) <= 0, 得证

推出决策单调性之后,正反各做一遍答案取max,复杂度 $O(nlogn)$。

#include<bits/stdc++.h>
#define MAXN 500005

using namespace std;
typedef long long ll;

int N,C;
int a[MAXN];

int dp[2][MAXN];
/*
dp[i] = max(a[j]-a[i]+sqrt(|i-j|))
dp0[i] = max(a[j]-a[i]+sqrt(i-j)), j<i
dp1[i] = max(a[j]-a[i]+sqrt(j-i)), j>i
*/

int x;
double f(int i, int j){
    return a[j]-a[i]+sqrt(i-j);
}

int search(int l, int r, int i, int j){
    int mid;
    while(l < r){
        mid = (l + r)/2;
        if(f(mid,i) > f(mid,j)) r = mid;
        else l = mid + 1;
    }
    return r;
}

struct Node{
    int l,r,j;
    Node(int l=0, int r=0, int j=0):l(l), r(r), j(j){}
};

Node q[MAXN];
void work(int x){
    int head = 0, tail = -1;
    q[++tail] = Node(1,N,1);
    dp[x][1] = 0;
    for(int i=2;i<=N;i++){
        while(head<=tail && q[head].r<i) ++head;
        int j = q[head].j;
        dp[x][i] = ceil(f(i,j));
        int pos = N + 1;
        while(head<=tail){
            int l = q[tail].l, r = q[tail].r, j = q[tail].j;
            if(f(l,i) > f(l,j)){
                pos = l;
                --tail;
                continue;
            }
            if(f(r,i) > f(r,j)) pos = search(l,r,i,j);
            q[tail].r = pos - 1;
            if(pos <= N) q[++tail] = Node(pos,N,i);
            break;
        }
    }
}

int main(){
    ios::sync_with_stdio(0);
    cin>>N;
    for(int i=1;i<=N;i++){
        cin>>a[i];
    }
    work(0);
    reverse(a+1, a+1+N);
    work(1);

    int ans = 0;
    for(int i=1;i<=N;i++){
        ans = max(dp[0][i], dp[1][N-i+1]);
        ans = max(ans, 0);
        cout<<ans<<endl;
    }
    return 0;
}

【P1912 [NOI2009]诗人小G】

https://www.luogu.com.cn/problem/P1912

和上题同样的推导方式,整理化简后,需要证明:

函数 f(x) = |x|^P - |x+c|^P,P>0, c>=0 单调非递增

需要分成3个定义域考虑:

case1. x<=-c

f(x) = (-x)^P - (-x-c)^P,求导得:
f'(x) = P*(-1)^P*[(-x)^(P-1) - (-x-c)^(P-1)]

当P为偶数时,(-x)^(P-1)为奇函数,有[(-x)^(P-1) - (-x-c)^(P-1)]<=0
当P为奇数时,(-x)^(P-1)为偶函数,有(-1)^P*[(-x)^(P-1) - (-x-c)^(P-1)]<=0

case2 和 case3 同理可证

复杂度 $O(nlogn)$

#include <bits/stdc++.h>
#define MAXN 100005
using namespace std;

typedef long double ll;
const ll inf = 1e18;

int T,N,L,P;
ll s[MAXN], dp[MAXN];
string t;
/*
dp[i] = min{dp[j] + (|sum[i]-sum[j]+(i-j-1)-L|)^P, 0<=j<i}
*/

inline ll qpow(ll x){
    ll res = 1;
    int n = P;
    while(n){
        if(n&1) res = res * x;
        x = x * x;
        n /= 2;
    }
    return res;
}

ll f(int i, int j){
    return dp[j] + qpow(abs(s[i]-s[j]+(i-j-1)-L));
}

int search(int l, int r, int i, int j){
    int mid;
    while(l < r){
        mid = (l + r)/2;
        if(f(mid,i) < f(mid,j)) r = mid;
        else l = mid + 1;
    }
    return r;
}

struct Node{
    int l,r,j;
    Node(int l=0, int r=0, int j=0):l(l), r(r), j(j){}
};

Node q[MAXN];

int main(){
    ios::sync_with_stdio(0);

    cin>>T;

    while(T--){
        cin>>N>>L>>P;
        for(int i=1;i<=N;i++){
            cin>>t;
            s[i] = s[i-1] + t.length();
        }

        int head = 0, tail = -1;
        q[++tail] = Node(1,N,0);

        bool ok = 1;
        for(int i=1;i<=N;i++){
            while(head<=tail && q[head].r<i) ++head;

            int j = q[head].j;
            dp[i] = f(i,j);

            int pos = N + 1;
            while(head<=tail){
                int l = q[tail].l, r = q[tail].r, j = q[tail].j;

                if(f(l,i) < f(l,j)){
                    pos = l;
                    --tail;
                    continue;
                }

                if(f(r,i) < f(r,j)) pos = search(l,r,i,j);
                q[tail].r = pos - 1;

                if(pos <= N) q[++tail] = Node(pos,N,i);
                break;
            }
        }

        if(dp[N] > inf || dp[N] < 0) cout<<"Too hard to arrange"<<endl;
        else cout<<(long long)dp[N]<<endl;
        cout<<"--------------------"<<endl;
    }
    return 0;
}

【2D/1D 决策单调性定理】:

对于DP:

dp[i][j] = min{dp[i][k] + dp[k+1][j] + w(i,j), i<=k<j}
或
dp[i][j] = min{dp[k][j-1] + w(k,i), k<i}
决策点记为 = p[i][j]

  1. 函数w满足四边形不等式
  2. 对于任意 $a\le b\le c\le d$,有$w(a,d) \ge w(b,c)$

则决策点满足:$p[i][j-1] \le p[i][j] \le p[i+1][j]$。证明略(见 https://oi-wiki.org/dp/opt/quadrangle/ )。

转移时直接使用上述结论,复杂度 $O(n^3) \to O(n^2)$

【AcWing 282. 石子合并】

https://www.acwing.com/problem/content/description/284/

经典区间DP:

$dp[l][r] = min(dp[l][k] + dp[k+1][r] + sum(l,r), l\le k < r)$

此处 $w(l,r) = sum(l,r)$ 是区间和且 $a[i]>0$ ,满足四边形不等式且取到等号:

  1. 对于任意 $a \le b\le c\le d$, 有 $w(a,d) + w(b,c) = w(a,c) + w(b,d)$
  2. 对于任意 $a \le b\le c\le d$,有 $w(a,d)\ge w(b,c)$

可以使用四边形不等式优化,$O(n^2)$。 注意若求max,则无法优化

#include<bits/stdc++.h>
#define MAXN 305
#define sum(i,j) (s[j] - s[i-1])
using namespace std;

int N;
int a[MAXN], s[MAXN];
int dp[MAXN][MAXN], p[MAXN][MAXN];
int main(){

    cin>>N;
    for(int i=1;i<=N;i++){
        cin>>a[i];
        s[i] = s[i-1] + a[i];
        p[i][i] = i;
    }

    for(int len=2;len<=N;len++){
        for(int i=1;i+len-1<=N;i++){
            int j = i+len-1;
            dp[i][j] = 1e9;
            for(int k=p[i][j-1];k<=p[i+1][j];k++){
                if(dp[i][k] + dp[k+1][j] + sum(i,j) < dp[i][j]){
                    dp[i][j] = dp[i][k] + dp[k+1][j] + sum(i,j);
                    p[i][j] = k;
                }
            }
        }
    }

    cout<<dp[1][N]<<endl;

    return 0;
}

【AcWing 336. 邮局】

https://www.acwing.com/problem/content/338/

数轴上 $n$ 个村庄建 $p$ 个小学,使得上学总距离最小。

n 个村庄坐标 x[1...n],前缀和 s[1...n]
dp[i][j] = min{dp[k][j-1] + w(k+1,i), k<i}

m = (j+i)/2, len = (i-j+1)
w(j,i) = s[i] + s[j-1] - 2*s[m] + x[m], len is odd
     = s[i] + s[j-1] - 2*s[m], len is even

对于函数 $w(i,j)$, 不妨设 len is odd (even 同理)

w(j,i+1) + w(j+1,i) >= w(j,i) + w(j+1,i+1)

m1 = (j+i+1)/2, m2 = (j+i)/2 = m1, m3 = (j+i+2)/2 = m1+1
s[i+1]+s[j-1]-2*s[m1] + s[i]+s[j]-2*s[m1] >= s[i]+s[j-1]-2*s[m1]+a[m1] + s[i+1]+s[j]-2*s[m3]+a[m3]

-2*s[m1] >= +a[m1]-2*s[m3]+a[m3]

-2*s[m1] >= +a[m1]-2*(s[m1]+a[m3])+a[m3]

a[m1+1] >= a[m1], 得证

可以使用四边形不等式优化,决策点满足:$p[i][j-1] \le p[i][j] \le p[i+1][j]$,注意递推顺序。

#include<bits/stdc++.h>
#define MAXN 305
#define sum(i,j) (s[j] - s[i-1])
using namespace std;

int N,P;
int a[MAXN], s[MAXN], w[MAXN][MAXN];
int dp[MAXN][MAXN], p[MAXN][MAXN];
int main(){

    cin>>N>>P;
    for(int i=1;i<=N;i++){
        cin>>a[i];
        s[i] = s[i-1] + a[i];
        p[i][i] = i;
    }

    for(int l=1;l<=N;l++){
        for(int r=l;r<=N;r++){
            int m = (l+r)/2;
            w[l][r] = (m-l+1)*a[m] - (s[m] - s[l-1]) + (s[r] - s[m]) - (r-m)*a[m];
        }
    }

    for(int i=N;i>=1;i--){
        dp[i][1] = w[1][i]; 
    }

    for(int j=2;j<=P;j++){
        p[N+1][j] = N;
        for(int i=N;i>=1;i--){
            dp[i][j] = 1e9;
            for(int k=p[i][j-1];k<=p[i+1][j];k++){
                if(dp[k][j-1] + w[k+1][i] < dp[i][j]){
                    dp[i][j] = dp[k][j-1] + w[k+1][i];
                    p[i][j] = k;
                }
            }
        }
    }

    cout<<dp[N][P]<<endl;

    return 0;
}