动态 DP

例题1：最大子段和问题，区间询问版本。

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

$1 \le n, q \le 4 \times 10^5$

考虑 dp。首先考虑简单版本：令 $dp_i$ 表示以第 $i$ 个元素为结尾的最大子段和。

$$dp_{l} = a_l \\ dp_{i} = \max(dp_{i - 1} + a_i, a_i)$$

根据矩阵乘法的重新定义，我们可以把它写成矩阵，其中“乘”->“加”，“加”->“$\max$”。也就是：

$$\begin{bmatrix} a_k & a_k \\ -\inf & 0\end{bmatrix} \begin{bmatrix} dp_{k-1} \\ 0\end{bmatrix} = \begin{bmatrix} dp_k \\ 0\end{bmatrix}$$

考虑求

$$\begin{bmatrix} dp_R \\ 0\end{bmatrix}$$

有

$$\begin{bmatrix} dp_R \\ 0\end{bmatrix} = \mathbf{A_R} \begin{bmatrix} dp_{R-1} \\ 0\end{bmatrix} = ... = \mathbf{A_R}\mathbf{A_{R-1}}...\mathbf{A_{L+1}} \begin{bmatrix} a_L \\ 0\end{bmatrix}$$

也就是我们对于每次询问需要求出 $\mathbf{A_R}\mathbf{A_{R-1}}...\mathbf{A_{L+1}}$。

注意到 $\mathbf{A_i}$ 和 $a_i$ 有关，且只和 $a_i$ 有关。怎么求？考虑我们做序列连续和怎么做。前缀和/差分；树状数组；线段树...

考虑前缀积，但是由于矩阵不一定存在逆（需要高斯消元方程有解），不行。

考虑树状数组，我们复习一下树状数组的性质，发现差分需要有逆操作，但是矩阵不一定存在逆，不行。

于是我们可以用线段树。每一个节点维护一个矩阵，上传操作是矩阵乘法。不需要修改，只需要查询区间积即可。

#include<bits/stdc++.h>
using namespace std;
#define int long long
#define f(i, a, b) for(int i = (a); i <= (b); i++)
#define cl(i, n) i.clear(),i.resize(n);
#define endl '\n'
typedef long long ll;
typedef unsigned long long ull;
typedef pair<int, int> pii;
const int inf = 1e12;
void cmax(int &x, int y) {if(x < y) x = y;}
void cmin(int &x, int y) {if(x > y) x = y;}
struct mat {
    int a[5][5];
    mat() {
        memset(a, 0xcf, sizeof(a));
    }

    mat operator* (mat b) {
        mat c = mat();
        f(i, 1, 2) {
            f(k, 1, 2) {
                int r = a[i][k];
                f(j, 1, 2) {
                    c.a[i][j] = max(c.a[i][j], r + b.a[k][j]);
                }
            }
        }
        return c;
    }
    mat operator^ (int b) {
        mat c = mat(), x = *this;
        f(i, 1, 2) c.a[i][i] = 0;
        while(b){
            if(b&1)c=c*(*this);
            (*this)=(*this)*(*this);
            b>>=1;
        }
        return c;
    }
};
struct sgt {
    mat A = mat();
}t[300010];
int y[300010];
void build(int now, int l, int r) {
    if(l == r) {
        t[now].A.a[1][1] = t[now].A.a[1][2] = y[l];
        t[now].A.a[2][1] = -inf; t[now].A.a[2][2] = 0;
        return;
    }
    int mid = (l + r) >> 1;
    build(now * 2, l, mid);
    build(now * 2 + 1, mid + 1, r);
    t[now].A = t[now * 2].A * t[now * 2 + 1].A; 
}
mat query(int now, int l, int r, int x, int y) {
    if(l >= x && r <= y) {
        return t[now].A;
    }
    if(l > y || r < x) {
        mat ret = mat();
        f(i, 1, 2) ret.a[i][i] = 0;
        return ret;
    }
    int mid = (l + r) >> 1;
    return query(now * 2, l, mid, x, y) * query(now * 2 +1, mid + 1, r, x, y);
}
signed main() {
    ios::sync_with_stdio(0);
    cin.tie(NULL);
    cout.tie(NULL);
    time_t start = clock();
    //think twice,code once.
    //think once,debug forever.
    int n; cin >> n;
    f(i, 1, n) cin >> y[i];
    build(1, 1, n);
    int m; cin >> m;
    f(i,1,m){
        int l,r;cin>>l>>r;
        mat tar = mat();
        tar.a[1][1] = y[l], tar.a[2][1] = 0;
        mat rat = query(1, 1, n, l + 1, r);
        rat = rat * tar;
        cout << rat.a[1][1] << endl;
    }
    time_t finish = clock();
    //cout << "time used:" << (finish-start) * 1.0 / CLOCKS_PER_SEC <<"s"<< endl;
    return 0;
} 

现在考虑原问题怎么做。显然可以求出区间矩阵乘积，然后乘以初始矩阵即可。

例题2：带单点修改。

$q$ 次修改单点 $a_i$ + 询问。

$1 \le n, q \le 2 \times 10^5$

每次在线段树中单点修改 $A_i$ 即可。

CF1814E

比赛的时候遇到的一道动态dp，场切了。

转移式子是

$$dp_i = \min (dp_{i-2}+c_{i-1}, dp_{i-3} +c_{i-1}+c_{i-2})$$

显然可以表示为矩阵 $\min, +$ 乘法。问题是要搞清楚左右乘。但是不小心写成右乘也没关系，线段树网上合并的时候改成右边乘以左边即可。

#include<bits/stdc++.h>
using namespace std;
#define int long long
//use ll instead of int.
#define f(i, a, b) for(int i = (a); i <= (b); i++)
#define cl(i, n) i.clear(),i.resize(n);
#define endl '\n'
typedef long long ll;
typedef unsigned long long ull;
typedef pair<int, int> pii;
typedef pair<ll, ll> pll;
const int inf = 1e16;
//#define cerr if(false)cerr
//#define freopen if(false)freopen
#define watch(x) cerr  << (#x) << ' '<<'i'<<'s'<<' ' << x << endl
void pofe(int number, int bitnum) {
    string s; f(i, 0, bitnum) {s += char(number & 1) + '0'; number >>= 1; } 
    reverse(s.begin(), s.end()); cerr << s << endl; 
    return;
}
template <typename TYP> void cmax(TYP &x, TYP y) {if(x < y) x = y;}
template <typename TYP> void cmin(TYP &x, TYP y) {if(x > y) x = y;}
//调不出来给我对拍！
//use std::array.
int c[200200];

struct matrix {
    //n = 3
    int a[3][3];
    matrix() {f(i,0,2)f(j,0,2)a[i][j]=inf;}
    matrix operator*(matrix op) {
        matrix res;
        f(i,0,2)f(j,0,2)f(k,0,2){
            cmin(res.a[i][j],a[i][k]+op.a[k][j]);
        }
        return res;
    }
    matrix(int cim1, int cim2) {
        a[0][0]=inf;a[0][1]=cim1;a[0][2]=cim1+cim2;
        a[1][0]=0;a[1][1]=inf;a[1][2]=inf;
        a[2][0]=inf;a[2][1]=0;a[2][2]=inf;
    }
}m[200200], it;
struct sgt {
    matrix t[4 * 200020];
    void build(int now,int l,int r) {
        // f(i,1,3)f(j,1,3)cout <<t[now].a[i][j]<<" \n"[j==3];
        // cout << "(build)\n";
        if(l==r){
            
            t[now]=m[l];// f(i,1,3)f(j,1,3)cout <<t[now].a[i][j]<<" \n"[j==3];
            // cout << "(after build)\n";
            return;
        }
        int mid=(l+r)>>1;
        build(now*2,l,mid);build(now*2+1,mid+1,r);
        t[now]=t[now*2+1]*t[now*2];
        // f(i,1,3)f(j,1,3)cout <<t[now].a[i][j]<<" \n"[j==3];
            // cout << "(after build)\n";
    }
    void modify(int now,int l,int r,int x,matrix mat) {
        if(l==r){
            t[now]=mat; return;
            
        }
        int mid=(l+r)>>1;
        if(x<=mid)modify(now*2,l,mid,x,mat);
        else modify(now*2+1,mid+1,r,x,mat);
        t[now]=t[now*2+1]*t[now*2];
    }
    // matrix query() {return t[1];}
}sgt;
signed main() {
    ios::sync_with_stdio(0);
    cin.tie(NULL);
    cout.tie(NULL);
    //freopen();
    //freopen();
    //time_t start = clock();
    //think twice,code once.
    //think once,debug forever.
    int n; cin>>n;
    c[0]=0; 
    f(i,1,n-1)cin>>c[i];
    f(i,1,n)m[i]=matrix(c[i-1],i-2<0?inf:c[i-2]);
    it.a[0][0] = 0;it.a[1][0] = inf; it.a[2][0] = inf;
    sgt.build(1,1,n);
    int q; cin>>q;
    while(q--){
        int k,x;cin>>k>>x;c[k]=x;
        m[k+1]=matrix(c[k],k-1<0?inf:c[k-1]);
        
        sgt.modify(1,1,n,k+1,m[k+1]);
        if(k<n-1){
            m[k+2]=matrix(c[k+1],k<0?inf:c[k]);
            sgt.modify(1,1,n,k+2,m[k+2]);
        }
        matrix mt = sgt.t[1]; //sgt.query();
        // f(i,1,3)f(j,1,3)cout <<mt.a[i][j]<<" \n"[j==3];
        matrix res =  mt*it;
        cout <<res.a[0][0]*2 << "\n";
    }
    //time_t finish = clock();
    //cout << "time used:" << (finish-start) * 1.0 / CLOCKS_PER_SEC <<"s"<< endl;
    return 0;
}
/*
2023/x/xx
start thinking at h:mm


start coding at h:mm
finish debugging at h:mm
*/