状压dp，区间dp，矩阵快速幂

Posted 2020-11-28 king-of-dark

DP

首先先回忆一下dp，dp叫做记忆化搜索，是一种可以把暴力搜索中重复的部分重复利用，从而到达减小复杂度的目的。比如最应该熟悉的背包模型，如果你把选择的过程看成一步一步的，那么在这么多的搜索路径中一定有着很多很多的重复部分，dp就是一种把重复的部分加以利用的方法。相信大家都已经在以前的练习中已经明白了dp是什么样的思路了，接下来的两种dp会在大家已经了解经典的背包dp等模型下展开。

状态压缩dp：

首先先讲一个状压dp最最经典的模型，求哈密尔顿路径问题，也叫做旅行商问题：给你一张图，你要在每个点只能走一次的限制下把所有点都走一次。这是最经典的模型，也希望大家都能熟练的使用这个模型，和记住这个模型中代码的写法。

旅行商例题：hdu 5418

问题描述：给你n个点，m条边，找出最短的一条回路，使其经过所有的点。（n<16）

首先简化问题，我们固定从1号点出发，求然后走遍到所有点的答案，再加上从末尾点到1号点的距离即可。

首先考虑暴力搜，用dfs把每种可能性都尝试一下，可以计算一下复杂度，其实就是A（n，n），尝试把所有可行的路径都尝试一下，然后走到最后一个的时候直接返回起点。

然后你考虑如何优化这个搜索，你会发现在n的那么多种排列中，其实有非常多的部分都是重复的

比如1 2 3 4

　　1 3 2 4

　　1 2 4 3

　　1 3 4 2

　　1 4 2 3

   　1 4 3 2

　　你会发现其实有些排列末尾其实都是一样的，对答案的贡献其实都是相同的，那么前半段中答案比较大的那部分其实是没有必要再计算下去了，因为他们肯定不可能构成答案。后面都是相同的搜索，我们只需要进行一次就可以了，所以说，我们的目的就是把能够简化的搜索都要记录下来，从中选择一个最优解记录下来即可，为做到这个目的，接下来要分析的就是接下来的选择到底和什么有关。

　　现在假设你再搜索中已经选的集合设为S，没有选中的设为V，你会发现，无论你途经什么样的手段把前S个点，都走一遍，都不会对后半部分把V中的点走一遍再回到1号点。也就是说构成了dp的无后效性，无论你前面如何到达当前状态都不会对后续的选择有影响，因此前面这么多路径中，我们仅仅记录最优的一条即可。

　　接下来来看一下具体思路：

　　你需要保存的是已经走过的集合S，和当前停留的城市i，因为要决定下一步去哪个城市，所以说要保存当前所在城市，然后你考虑如何表达这个集合S呢，因为此时n是很小的，为了节省空间和时间，我们使用一个整数s来表示集合S，对于集合S中第i个城市走过的话，s中对应的第i位就是1，否则就是0.然后我们这样设立dp数组 dp[s][i],表示当前已经走过的城市集合为s，并且停留在i城市的最小代价。然后你要考虑去扩展新的城市，假设j不是S中的城市，那么走向j的代价就是dp[s][i]+dis[i][j].所以我们的状态转移方程也出来了，就是dp[s|(1<<j)][j]=max(dp[s|(1<<j)][j],dp[s][i]+dis[i][j]).

　　最后大家来看一下写法吧。

 

int M=(1<<n)-1;
        for(int i=0;i<n;i++)
            for(int j=0;j<=M;j++)
                dp[i][j]=1e9;
        dp[0][1]=0;
        for( int s=1;s<M;s++)
            for( int i=0;i<n;++i)
                for( int j=0;j<n;++j)
                    if(!(s>>j & 1)) {//如果当前城市j还没有经过,我们从i走向j
                        int next=s|(1<<j);
                        if(next==M)
                            dp[j][next]=min(dp[j][next],dp[i][s]+dis[i][j]+dis[j][0]);
                        else
                            dp[j][next]=min(dp[j][next],dp[i][s]+dis[i][j]);
                    }
        int ans=1e9;
        for(int i=1;i<n;i++)
            ans=min(ans,dp[i][M]);

 

　　现在解释一下这段代码，因为我是用的递推的写法写的，而且状态转移有一个特点，就是状态s1要是能够转移到s2，那么用整数比较的话，s1<s2是必然的。这个很显然。那么就是说如果我们从小枚举状态，是可以把所有可能转移的情况都考虑到的。

　　好的，现在咱们来总结一下状压dp：

　　（1）首先是复杂度，状压dp把np难的问题转移成了2^n。

　　（2）其次是状态压缩的技巧，把一个状态压成一个整数，而不是使用一个数组的方式，在各个题目中都很常见。其实状态压缩的另外一个好处。是可以快速枚举每一个状态，使用for循环即可。

　　（3）对于状态压缩的题目，主要还是基于状态对后续的影响是基于前面的状态，而与如何达到状态无关，这是和基本的dp相同的。

　　先说一句题外话，哈密尔顿回路问题目前是无高效算法的，已经被证明了是一个np难的问题，只有充分条件可以证明，但是对与大数据来说是无解的。

接下来看一道例题：

　　poj-2288

　　现在你要求走一条哈密尔顿路径（不是回路），设路径是a1a2a3a4a5,那么这条路径的权值计算过程分为三部分：

　　　　　　1.所有点的点权和

　　　　　　2.所有相邻点的权值乘，比如v[a1]*v[a2]+v[a2]*v[a3]+...

　　　　　　3.如果对于相邻的三个点，a1a2a3,如果a1与a3是相连的，那么权值要加上v[a1]*v[a2]*v[a3],对于a2a3a4同理

　　点数小于等于13，你要求统计出权值最大的路径权值是多少，以及多少条路径的权值是这个（同一条路径的两种走法算一次）。

　　首先不考虑记录数量，只求值，很明显，我们选择下一点的过程中要添加的答案只与前面的状态，前面选的点是什么有关而已。

　　所以我们就记录我们所需要的东西，设dp[s][i][j],当前已经走过的点数集合的状态为s，i为上一次走过的点，j为当前停留的点。所以如果你想下一步走k号点的话，状态转移方程也是很简单的： dp[s|(1<<k)][j][k]=max(dp[s|(1<<k)][j][k],dp[s][i][j]+v[k]+v[k]*v[j]+(是否为三角形带来的权值))。

　　我们再考虑如何同时统计数量，我们考虑一个状态到另一个状态的实质，其实是把前面的搜索过程合一，而这部转移其实也是一次搜索，也就是说，如果我有x种方案到达状态s，那么状态s转移到状态s2，同样可以有x中方案走这样的搜索路径到达s2。所以说我们没次统计一下当前状态的次数，再转移的时候，如果和将要转移的状态答案相同，也就是同为最优解，直接相加即可，否则就是直接刷新之前的答案和数量。

 再来一道例题：http://codeforces.com/contest/1215/problem/E

　　给你n个带颜色的块，排成一排，n<4e5,颜色不超过20种，你可以交换相邻的两个块，让他们的颜色调换。你的目的是让所有同一颜色的块都聚到一起，问最小操作步骤。

　　首先先讲一个类似的题，将一个无序的数组用冒泡排序交换到有序，最少的交换次数是多少？

　　首先对于这个问题而言，你要知道一对相邻的数什么时候是有必要交换的，对于相邻的数a[i],a[i+1],要是a[i]>a[i+1]就一定会交换，因为最终结果对于一个点来说不存在左边的值比当前值大。所以说我们有一个猜想：就是一个点左边比他小的点，都没有必要与他交换，左边比他大的点，一定要与它交换，从而跨越到这个点的右侧。

　　所以这个点的对答案产生的贡献其实是这个点前面比他大的点的数量，我们不考虑右侧比它小的点，因为往右侧的交换，会被后面的左侧交换算到。所以你只要对每个点计算一下它前面比它大的就是答案。

　　然后回到这个题，你要做的是把所有相同的颜色都聚在一起对吧，现在考虑最终状态，是不是一定是【红红红黄黄黄绿绿绿】这种连续的，其实如果你知道最后的排列顺序，其实你就可以等效替代一下了，把红等效成0，黄是1，绿是2，转换成上面的问题就可以了。所以说对于这道题来说，你要枚举的其实是最后的状态，一共有2^n个。这个时候就应该感觉出这是一道状压的题目了。

　　现在考虑如何暴力，假如我第一个颜色选择颜色3，不贡献答案，第二个颜色选择颜色2，我现在只考虑2号点如何都移动到3号点的左侧，我只需要知道，对于所有的2号点应该交换多少次就可以了，就是对于每个2号点，左侧3号点数量之和，就是颜色2作为第二个选择的贡献。然后选择颜色1，我们需要考虑的是1号点，如何交换到2号点和3号点的左侧，对答案的贡献是，对于每个1号点，左侧2号点和3号点的和。

　　就此我们发现，对于选择下一个颜色是什么，产生的贡献只与前面选的颜色是什么有关系，而与前面的颜色顺序没有关系。所以我们就有了状态转移方程。

　　dp[s]代表，前面已经选了s种颜色产生的最小贡献。

　　dp[s|(1<<i)]=min(dp[s|(1<<i)],dp[s]+t)

　　t为对于i号点来说，左侧有s类的块的数量，这个是可以提前统计的。

　　这个题比较难，此处贴上我的代码，以防大家找不到题解

#include <algorithm>
#include <iostream>
#include <cstring>
#include <cstdio>
#include <cmath>
#include <vector>
#include <map>
using namespace std;
const double pi=acos(-1.0);
const int maxn=4e5+10;
const double eps=1e-6;
typedef long long ll;
const int shang=300;
ll a[30][30];
ll dp[1<<20],dp_t[20][1<<20];
int t[maxn];
int main()
{
    int n;
    cin>>n;
    for(int i=0;i<n;i++){
        scanf("%d",&t[i]);
        t[i]--;
    }
    for(int i=0;i<20;i++){
        for(int j=0;j<20;j++)
        {
            if(i==j)
                continue;
            ll count=0,ans=0;
            
            for(int k=0;k<n;k++)
                if(t[k]==i)
                    count++;
                else if(t[k]==j)
                    ans+=count;
            
            a[i][j]=ans;
//            cout<<a[i][j]<<" ";
        }
//        cout<<endl;
    }
    for(int i=0;i<20;i++)
    for(int j=0;j<(1<<20);j++)
    {
        if((j&(1<<i))!=0)
            continue;
        ll t=j,count=0,ans=0;
        while(t>0)
        {
            if(t%2==1)
            {
                ans+=a[count][i];
            }
            count++;
            t/=2;
        }
        dp_t[i][j]=ans;
    }
    dp[0]=0;
    for(int j=1;j<(1<<20);j++)
    {
        dp[j]=1e12+7;
    }
    
    for(int j=0;j<(1<<20);j++)
    {
        for(int k=0;k<20;k++)
        {
            if((j&(1<<k))==0)
            {
                dp[j|(1<<k)]=min(dp[j|(1<<k)],dp[j]+dp_t[k][j]);
//                cout<<dp[j|(1<<k)]<<endl;
            }
        }
    }
    cout<<dp[(1<<20)-1];
    
}

 

区间dp

　　 谈到区间dp一定离不开石子合并这道经典问题了。

　　对于n堆石子，每一堆都有初始值，你没次只能合并相邻的两堆，直到变为一堆，每次合并的代价是合并后的石子堆大小，问最小代价。

　　首先这个题看是能贪心，但是其实是贪心不了的，orz。

　　还是一样的套路dp先考虑暴搜，首先你要知道合成1堆肯定是从两堆和成的，而且这两堆石子一定是两个区间。而且最后一步的代价其实一定是所有的石子和，所以我们把F(i,j)表示为把i到j之间的石子和成一堆的最小代价，那么F(1,n)=F(1,k)+F(k+1,n)+sum.其中k是1和n之间的，这样我们是不是把问题分解为子问题了啊，然后我们一直递归就可以解决问题了。

　　然后就是优化过程，其实在很多时候你会，发现F(i,j)调用了很多次，而且每次其实返回结果都是一样的，所以我们在调用这个函数的时候记录一下值，如果再次调用，我们直接返回就可以了。

　　刚才我们讲的是递归的写法，但是其实比较常用的是递推的写法，先把小区间准备好，然后大区间直接用小区间的结果推出就行了。下面讲一下递推的写法。

　　

for(int i=1;i<=n;i++)
        sum[i][i]=a[i],dp[i][i]=0;
        
        for(int len=1;len<n;len++){//区间长度 
            for(int i=1;i<=n&&i+len<=n;i++){//区间起点 
                int j=i+len;//区间终点
                for(int k=i;k<=j;k++)//用k来表示分割区间 
                {
                    sum[i][j]=sum[i][k]+sum[k+1][j];
                    dp[i][j]=min(dp[i][j],dp[i][k]+dp[k+1][j]+sum[i][j]);
                 } 
            }
        }
        cout<<dp[1][n]<<endl;

————————————————
版权声明：本文为CSDN博主「紫芝」的原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接及本声明。
原文链接：https://blog.csdn.net/qq_40507857/article/details/81266843

 再来一道括号匹配，http://poj.org/problem?id=2955

有个包含（ ） [ ] 的字符串，你要找到一个最长的合法的子序列，对于一个子序列，其中的括号一定要有另一个相对应。

设dp[i][j],表示从i到j能构成的最长合法子序列的长度。

首先还是考虑如何暴搜，想一想最后的答案是怎么形成的，会发现其实还是只有两种形成方式

　　1.首尾相同，计算一次答案，然后和内部的答案加在一起

　　2.分成两部分，两部分的答案加在一起

　　这两种方式对应着两种可能，一种是这个合法的序列是一个括号括起来的，另一种是多个分割的括号括起来的。

　　对应于（（）（（）））   和 （）（（））（（））（（））

　　所以我们就有了状态转移方程，dp[i][j]=max(dp[i][k]+dp[k+1][j],dp[i][j]) (i<=k<=j),当a[i]和a[j]互补时，外加dp[i][j]=max(dp[i][j],dp[i+1][j-1]+1).

　　这样就把所有的情况遍历到了。

矩阵快速幂

　　相信大家都已经会了快速幂了，如果没会也没有关系。

　　首先还是讲经典例题求斐波那契数列，斐波那契数列是前两项是1，以为每一项都是前两项和的一个数列，比如1，1，2，3，5，8...也就是说i>2时，F[i]=F[i-1]+F[i-2].

　　现在考虑如何求这么一个矩阵，如果朴素的求时o（n）的，但是一般这种题，n都是1e18左右的范围，肯定是不行的，现在我们假设两个矩阵M，T：

　　F[2] F[1] 　*　1   1　　　　F[3] F[2]   

　　0    0　　　　 1   0  　　　0　　0

　　可以发现这两个特殊的矩阵乘在一起还是很神奇的，你会发现矩阵M乘矩阵T一次，就把矩阵M中的所有值都推进了一次，所以说我只要求F[n]，我只需要把矩阵M乘（n-1）次矩阵T即可。但是矩阵乘法时满足结合律的，你可以先把T^(n-1)先算出来，然后在与M相乘。于是乎我们就用快速幂的套路去加速。

　　我们先回顾一下快速幂的操作，求a的b次方，我们假设b为11，二进制为1011，是不是a^11=a^8+a^2+a^1. 我们进行二进制拆分后，a一次一次进行倍增，推出a^1,a^2,a^4,a^8的值，就可以把a^11在logn的时间内计算出来了。上代码：

ll qpow(ll a,ll b,ll m){
    ll ans=1;
    ll k=a;
    while(b){
        if(b&1)ans=ans*k%m;
        k=k*k%m;
        b>>=1;
    }
    return ans;
}

其实矩阵套进去也是很好写的，仅仅改一下参数，返回值，重构一下乘法而已，复杂度为m^3*logn：

struct mat
{
    ll m[MAXN][MAXN];//矩阵结构体
}unit;//unit为单位矩阵，即主对角线全部为1，这样任何矩阵与单位矩阵相乘都为它本身

mat msub(mat a,mat b)//矩阵相乘函数
{
    mat ret;
    ll x;
    for(int i=0;i<n;i++)
    {
        for(int j=0;j<n;j++)
        {
            x=0;
            for(int k=0;k<n;k++)
            {
                x+=((a.m[i][k]*b.m[k][j])%mod);//取余
            }
            ret.m[i][j]=x%mod;//取余
        }
    }
    return ret;
}


void init_unit()//初始化单位矩阵
{
    for(int i=0;i<MAXN;i++)
    {
        unit.m[i][i]=1;
    }
}

mat qpow(mat a,ll x)//快速幂
{
    mat ans=unit;
    while(x)
    {
        if(x&1) ans=msub(ans,a);
        a=msub(a,a);
        x>>=1;
    }
    return ans;
}

 

以上就是矩阵快速幂的基础了，接下来讲一讲如何变形。

F[i]=a*F[i-1]+b*F[i-2].如何构造矩阵？

F[i]=a*F[i-1]+b*F[i-2]+c*F[i-3].如何构造矩阵？

F[i]=a*F[i-1]+b*F[i-2]+i^2.如何构造矩阵？

推荐例题：hdu5950

黑科技：更快的矩阵快速幂：常系数齐次线性递推（有需要的同学自行百度，是个板子，会用就好），用这个方法，可以省掉矩阵快速幂多余的部分，因为我们只想求F[n]，但是我们在过程中也求出了F[n-1],F[n-2]...,所以我们有一定的浪费，用这个方法，可以做到m^2*logn的复杂度，做到求一个数。

　　杜教BM：一个板子，如果把一个线性递推式的前2*n项，这里的n指开始出现规律的第一个下标，导入这个类中，你可以以很快的速度求出第x项。和上面的方法速度差不多。