專題：動態規劃（DP）初步

內容來源：《挑戰程式設計競賽》（第2版）+《演算法競賽入門經典》（第2版）+網上資料整理彙總

一、引入

        動態規劃程式設計是對解最優化問題的一種途徑、一種方法，而不是一種特殊演算法。不像前面所述的那些搜尋或數值計算那樣，具有一個標準的數學表示式和明確清晰的解題方法。動態規劃程式設計往往是針對一種最優化問題，由於各種問題的性質不同，確定最優解的條件也互不相同，因而動態規劃的設計方法對不同的問題，有各具特色的解題方法，而不存在一種萬能的動態規劃演算法，可以解決各類最優化問題。因此除了要對基本概念和方法正確理解外，必須具體問題具體分析，以豐富的想象力去建立模型，用創造性的技巧去求解。我們也可以通過對若干有代表性的問題的動態規劃演算法進行分析、討論，逐漸學會並掌握這一設計方法。

二、DP基本模型：多階段決策過程的最優化問題

        在現實生活中，有一類活動的過程，由於它的特殊性，可將過程分成若干個互相聯絡的階段，在它的每一階段都需要作出決策，從而使整個過程達到最好的活動效果。當然，各個階段決策的選取不是任意確定的，它依賴於當前面臨的狀態，又影響以後的發展，當各個階段決策確定後，就組成一個決策序列，因而也就確定了整個過程的一條活動路線，這種把一個問題看作是一個前後關聯具有鏈狀結構的多階段過程就稱為多階段決策過程，這種問題就稱為多階段決策問題。如下圖所示：

        多階段決策過程，是指這樣的一類特殊的活動過程，問題可以按時間順序分解成若干相互聯絡的階段，在每一個階段都要做出決策，全部過程的決策是一個決策序列。

1.      基本概念

（1）階段和階段變數：

        用動態規劃求解一個問題時，需要將問題的全過程恰當地分成若干個相互聯絡的階段，以便按一定的次序去求解。描述階段的變數稱為階段變數，通常用K表示，階段的劃分一般是根據時間和空間的自然特徵來劃分，同時階段的劃分要便於把問題轉化成多階段決策過程。

（2）狀態和狀態變數：

        某一階段的出發位置稱為狀態，通常一個階段包含若干狀態。一般地，狀態可由變數來描述，用來描述狀態的變數稱為狀態變數。

（3）決策、決策變數和決策允許集合：

        在對問題的處理中作出的每種選擇性的行動就是決策。即從該階段的每一個狀態出發，通過一次選擇性的行動轉移至下一階段的相應狀態。一個實際問題可能要有多次決策和多個決策點，在每一個階段的每一個狀態中都需要有一次決策，決策也可以用變數來描述，稱這種變數為決策變數。在實際問題中，決策變數的取值往往限制在某一個範圍之內，此範圍稱為允許決策集合。

（4）策略和最優策略：

        所有階段依次排列構成問題的全過程。全過程中各階段決策變數所組成的有序總體稱為策略。在實際問題中，從決策允許集合中找出最優效果的策略成為最優策略。

（5）狀態轉移方程【分析問題的核心】

        前一階段的終點就是後一階段的起點，對前一階段的狀態作出某種決策，產生後一階段的狀態，這種關係描述了由k階段到k+1階段狀態的演變規律，稱為狀態轉移方程。

2.       兩個重要特性：最優化原理與無後效性

        一般來說，能夠採用動態規劃方法求解的問題，必須滿足最優化原理和無後效性原則：

（1）動態規劃的最優化原理：

        作為整個過程的最優策略具有：無論過去的狀態和決策如何，對前面的決策所形成的狀態而言，餘下的諸決策必須構成最優策略的性質。也可以通俗地理解為子問題的區域性最優將導致整個問題的全域性最優，即問題具有最優子結構的性質，也就是說一個問題的最優解只取決於其子問題的最優解，而非最優解對問題的求解沒有影響。

（2）動態規劃的無後效性原則：

        所謂無後效性原則，指的是這樣一種性質：某階段的狀態一旦確定，則此後過程的演變不再受此前各狀態及決策的影響。也就是說，“未來與過去無關”，當前的狀態是此前歷史的一個完整的總結，此前的歷史只能通過當前的狀態去影響過程未來的演變。

        即：一個問題被劃分成各個階段之後，階段K中的狀態只能由階段K+1中的狀態通過狀態轉移方程得來，與其它狀態沒有關係，特別是與未發生的狀態沒有關係。

        由此可見，對於不能劃分階段的問題，不能運用動態規劃來解；對於能劃分階段，但不符合最優化原理的，也不能用動態規劃來解；既能劃分階段，又符合最優化原理的，但不具備無後效性原則，還是不能用動態規劃來解；誤用動態規劃程式設計方法求解會導致錯誤的結果。【三個條件缺一不可】

3.       動態規劃設計方法的一般模式

        動態規劃所處理的問題是一個多階段決策問題，一般由初始狀態開始，通過對中間階段決策的選擇，達到結束狀態；或倒過來，從結束狀態開始，通過對中間階段決策的選擇，達到初始狀態。這些決策形成一個決策序列，同時確定了完成整個過程的一條活動路線，通常是求最優活動路線。

        動態規劃的設計都有著一定的模式，一般要經歷以下幾個步驟：

（1）劃分階段

        按照問題的時間或空間特徵，把問題劃分為若干個階段。在劃分階段時，注意劃分後的階段一定是有序的或者是可排序的，否則問題就無法求解。

（2）確定狀態和狀態變數

        將問題發展到各個階段時所處於的各種客觀情況用不同的狀態表示出來。當然，狀態的選擇要滿足無後效性。

（3）確定決策並寫出狀態轉移方程

        因為決策和狀態轉移有著天然的聯絡，狀態轉移就是根據上一階段的狀態和決策來匯出本階段的狀態。所以如果確定了決策，狀態轉移方程也就可以寫出。但事實上常常是反過來做，根據相鄰兩段的各個狀態之間的關係來確定決策。

（4）尋找邊界條件

        給出的狀態轉移方程是一個遞推式，需要一個遞推的終止條件或邊界條件。

四、典例

1. 數字三角形/數塔問題（DP入門題）

        有形如下圖所示的數塔，從頂部出發，在每一結點可以選擇向左走或是向右走，一起走到底層，要求找出一條路徑，使路徑上的值最大。

樣例輸入：

5

13

11 8

12 7 26

6 14 15 8

12 7 13 24 11

樣例輸出：

86（13->8->26->15->24）

【分析】這道題如果用列舉法，在數塔層數稍大的情況下（如40），則需要列舉出的路徑條數將是一個非常龐大的數目。如果用貪心法又往往得不到最優解。在用動態規劃考慮數塔問題時可以自頂向下的分析，自底向上的計算。

        從頂點出發時到底向左走還是向右走應取決於是從左走能取到最大值還是從右走能取到最大值，只要左右兩道路徑上的最大值求出來了才能作出決策。同樣的道理下一層的走向又要取決於再下一層上的最大值是否已經求出才能決策。這樣一層一層推下去，直到倒數第二層時就非常明瞭。所以實際求解時，可從底層開始，層層遞進，最後得到最大值。

        狀態轉移方程：dp[i][j]=max(dp[i+1][j],dp[i+1][j+1])+a[i][j];

#include <iostream>
#include <cstdio>
#include <cstring>
#define maxn 105
using namespace std;
int n;
int a[maxn][maxn]; 
int dp[maxn][maxn];   //自底向上，記錄從點(i,j)出發到數塔底層的路徑最大和 
int main()
{
	int i,j;
	scanf("%d",&n);
	for(i=0;i<n;i++)
		for(j=0;j<=i;j++)
			scanf("%d",&a[i][j]);
	memset(dp,0,sizeof(dp));
	for(i=0;i<n;i++)  //填數塔最底層
		dp[n-1][i]=a[n-1][i];
	for(i=n-2;i>=0;i--)   //更新除數塔最底層外的各個點的路徑最大和 
		for(j=0;j<=i;j++)
			dp[i][j]=max(dp[i+1][j],dp[i+1][j+1])+a[i][j];
	printf("%d\n",dp[0][0]);
	return 0;
}

2. 序列DP

（1）最長上升子序列LIS

        輸入n及一個長度為n的數列，求出此序列的最長上升子序列長度。上升子序列指的是對於任意的i<j都滿足ai<aj的子序列。（1<=n<=1000,0<=ai<=1000000）

樣例輸入：

5

4 2 3 1 5

樣例輸出：

3（最長上升子序列為2, 3, 5）

【分析】

法一：定義dp[i]: 以ai為末尾的最長上升子序列的長度。以ai結尾的上升子序列是：

        1° 只包含ai的子序列

        2° 在滿足j<i且aj<ai的以aj結尾的上升子列末尾，追加上ai後得到的子序列

        這二者之一。這樣就能得到如下遞推關係：

        dp[i]=max{1, dp[j]+1 | j<I 且aj<ai}，複雜度為O(n²)。

#include <iostream>
#include <cstdio>
#define maxn 1005
using namespace std;
int n,a[maxn];
int dp[maxn];    //dp[i]記錄以a[i]為末尾的最長上升子序列的長度
int main()
{
	int i,j;
	int ret;
	while(scanf("%d",&n)!=EOF)
	{
	    for(i=1;i<=n;i++)
        {
            scanf("%d",&a[i]);
            dp[i]=1;      //初始化
        }
        ret=1;
        for(i=1;i<=n;i++)
        {
            for(j=1;j<i;j++)  //遍歷所有在a[i]之前的元素
            {
                if(a[j]<a[i]) //若存在aj<ai，則在以aj為結尾的上升子列末尾追加ai後得到的子序列 和 只包含ai的子序列中取長度較大者
                    dp[i]=max(dp[i],dp[j]+1);
            }
            ret=max(ret,dp[i]);  //注意隨時更新ret
        }
        printf("%d\n",ret);
	}
	return 0;
}

法二：此外，還可以定義其它的遞推關係。前面利用DP求取針對最末位的元素的最長的子序列。如果子序列的長度相同，那麼最末位的元素較小的在之後會更加有優勢，所以我們再反過來用DP針對相同長度情況下最小的末尾元素進行求解：

        dp[i]: 長度為i+1的上升子序列中末尾元素的最小值（不存在的話就是INF）。

        過程分析：最開始全部dp[i]的值都初始化為INF，然後由前到後逐個考慮數列的元素。對於每個aj，如果i=0或dp[i-1]<aj的話，就用dp[i]=min(dp[i], aj)進行更新。最終找出使得dp[i]<INF的最大的i+1就是結果。

        此DP直接實現的話，能夠與前面的方法一樣在O(n²)的時間內給出結果，但這一演算法還可以進一步優化：首先dp數列中除INF之外是單調遞增的，所以可以知道對於每個aj最多隻需要1次更新。對於這次更新究竟應在什麼位置，不必逐個遍歷，可以利用二分搜尋，這樣就可以在O(nlogn)時間內求出結果。

#include <iostream>
#include <cstdio>
#include <algorithm>
#define maxn 1005
#define INF 99999999
using namespace std;
int n,a[maxn];
int dp[maxn];    //dp[i]:長度為i+1的上升子序列中末尾元素的最小值（不存在的話就是INF）
int main()
{
	int i,j;
	scanf("%d",&n);
	for(i=0;i<n;i++)
		scanf("%d",&a[i]);
	fill(dp,dp+n,INF);   //初始化dp陣列為INF
	for(i=0;i<n;i++)     //找到更新dp[i]的位置並用a[i]更新之
    {
        *lower_bound(dp,dp+n,a[i])=a[i];
        for(j=0;j<n;j++) //觀察dp陣列的填充過程
            printf("%d ",dp[j]);
        printf("\n");
    }
	printf("%d\n",lower_bound(dp,dp+n,INF)-dp);  //第一個INF出現的位置即為LIS長度
	return 0;
}

        拓展：使用upper_bound和lower_bound兩個函式求長度為n的有序陣列a中的k的個數

        upper_bound(a, a+n, k)-lower_bound(a, a+n, k);

        另外，求最長下降/不上升/不下降子序列思路同此題，只是判斷條件有變化。

（2）最長公共子序列LCS

        給定兩個字串s1和s2（長度均不超過1000），求出這兩個字串的最長公共子序列的長度。

【分析】定義dp[i][j]：串s1的前i個字元 和 串s2的前j個字元的最長公共子序列長度，則s1…si+1和t1…tj+1對應的公共子列可能是：

        ①si+1=tj+1時：在s1…si 和 t1…tj的公共子列末尾追加si+1（即LCS長度+1）

        ②否則可能為s1…si和t1…tj+1的公共子列長度l1 或s1…si+1和t1…tj的公共子列長度l2，二者取較大者。

        故狀態轉移方程為：

        dp[i+1][j+1]=dp[i][j]+1,                                      ①

                                max(dp[i][j+1], dp[i+1][j]),        ②

        最後dp[len1][len2]即為所求，其中len1、len2分別為串s1和s2的長度。

#include <iostream>
#include <cstdio>
#include <cstring>
using namespace std;
const int maxlen=1010;
char s1[maxlen],s2[maxlen];
int dp[maxlen][maxlen];   //dp[i][j]記錄串s1的前i個字元和串s2的前j個字元的LCS長度
int main()
{
    int i,j;
    int len1,len2;
    while(scanf("%s",s1)!=EOF)
    {
        scanf("%s",s2);
        len1=strlen(s1);
        len2=strlen(s2);
        dp[0][0]=0;         //初始化：兩串均為空時，len(LCS)=0
        for(i=1;i<=len1;i++)//s2串為空時，不論s1中有多少字元，len(LCS)=0
            dp[i][0]=0;
        for(i=1;i<=len2;i++)//s1串為空時，不論s1中有多少字元，len(LCS)=0
            dp[0][i]=0;
        for(i=0;i<len1;i++)
        {
            for(j=0;j<len2;j++)
            {
                if(s1[i]==s2[j]) //s1與s2對應位置字元相等
                    dp[i+1][j+1]=dp[i][j]+1;
                else     //其它情況：兩者取較大者
                    dp[i+1][j+1]=max(dp[i][j+1],dp[i+1][j]);
            }
        }
        /*for(i=1;i<=len1;i++)
        {
            for(j=1;j<=len2;j++)
                printf("dp[%d][%d]=%d ",i,j,dp[i][j]);
            printf("\n");
        }*/
        printf("%d\n",dp[len1][len2]);
    }
    return 0;
}

（3）最大公共子串LCS

        給定兩個字串s1和s2（長度均不超過1000），求出這兩個字串的最大公共子串的長度。

【分析】情境類似求最長公共子序列長度問題，不過需要注意的是：所求子串中的字元需要在串s1和串s2中連續出現。

        例：s1=”abcad”

                s2=”abd”

        它們的最長公共子序列長度為3（”abd”），而最大公共子串長度為2（”ab”）。

        因此，定義dp[i][j]：串s1的前i個字元 和 串s2的前j個字元的最大公共子串長度，則s1…si+1和t1…tj+1對應的公共子串可能是：

        ①si+1=tj+1時：在s1…si 和 t1…tj的公共子串末尾追加si+1（即LCS長度+1）

        ②否則dp[i][j]=0。

        分析可知狀態轉移方程：

        dp[i+1][j+1]=dp[i][j]+1,                     ①

                                0,                                   ②

#include <iostream>
#include <cstdio>
#include <cstring>
using namespace std;
const int maxlen=1010;
char s1[maxlen],s2[maxlen];
int dp[maxlen][maxlen];    //dp[i][j]為串s1的前i個字元和串s2的前j個字元的最大公共子串長度
int main()
{
    int i,j;
    int len1,len2,ret;     //ret記錄結果
    while(scanf("%s",s1)!=EOF)
    {
        scanf("%s",s2);
        memset(dp,0,sizeof(dp));   //初始化：開始LCS長度均為0
        len1=strlen(s1);
        len2=strlen(s2);
        ret=0;
        for(i=0;i<len1;i++)
        {
            for(j=0;j<len2;j++)
            {
                if(s1[i]==s2[j])
                    dp[i+1][j+1]=dp[i][j]+1;
                else
                    dp[i+1][j+1]=0;
                ret=max(ret,dp[i+1][j+1]);    //隨時更新最大值
            }
        }
        printf("%d\n",ret);
    }
    return 0;
}

1. 攔截導彈(Noip1999)

　　某國為了防禦敵國的導彈襲擊，發展出一種導彈攔截系統。但是這種導彈攔截系統有一個缺陷：雖然它的第一發炮彈能夠到達任意的高度，但是以後每一發炮彈都不能高於前一發的高度。某天，雷達捕捉到敵國的導彈來襲。由於該系統還在試用階段，所以只有一套系統，因此有可能不能攔截所有的導彈。

　　輸入導彈數n及n顆導彈依次飛來的高度（雷達給出的高度資料是不大於30000的正整數，導彈數不超過1000），計算這套系統最多能攔截多少導彈，如果要攔截所有導彈最少要配備多少套這種導彈攔截系統。

樣例輸入：

8

389 207155 300 299 170 158 65

樣例輸出：

6（最多能攔截的導彈數）

2（要攔截所有導彈最少要配備的系統數）

【分析】DP+貪心法

        第一問即經典的最長不下降子序列問題，可以用一般的DP演算法，也可以用高效演算法（求LIS問題的法二），但這個題的資料規模不需要。

　　第二問用貪心法即可。每顆導彈來襲時，使用能攔截這顆導彈的防禦系統中上一次攔截導彈高度最低的那一套來攔截。若不存在符合這一條件的系統，則使用一套新系統。

#include <iostream>
#include <cstdio>
#include <cstring>
using namespace std;
const int maxn=1010;
int n;
int dp[maxn];   //dp[i]記錄前i發導彈中最多攔截的顆數
struct shell        //炮彈
{
    int height;     //高度
    int is_catched; //是否被攔截
} s[maxn];
int main()
{
    int i,j;
    int cur;        //cur記錄當前導彈高度
    int maxans,minret;  //最多攔截的導彈數.攔截導彈最少配備的系統數
    while(scanf("%d",&n)!=EOF)
    {
        for(i=1;i<=n;i++)
        {
            scanf("%d",&s[i].height);
            s[i].is_catched=0;
            dp[i]=1;
        }
        maxans=1,minret=0;
        for(i=1;i<=n;i++)
        {
            for(j=1;j<i;j++)
            {
                if(s[j].height>=s[i].height)
                    dp[i]=max(dp[i],dp[j]+1);
            }
            maxans=max(maxans,dp[i]);
            if(s[i].is_catched==1)
                continue;
            cur=s[i].height;
            minret++;
            for(j=i+1;j<=n;j++)
            {
                if(s[j].height<=cur)
                {
                    cur=s[j].height;
                    s[j].is_catched=1;
                }
            }
        }
        printf("%d\n",maxans);
        printf("%d\n",minret);
    }
    return 0;
}

　　N位同學站成一排，音樂老師要請其中的(N-K)位同學出列，使得剩下的K位同學排成合唱隊形。

　　合唱隊形是指這樣的一種隊形：設K位同學從左到右依次編號為1, 2, …, K，他們的身高分別為T₁, T₂, …, T_K，則他們的身高滿足T₁< T₂ < … < T_i , T_i > T_i+1> … > T_K (1≤i≤K)。

　　你的任務是，已知所有N位同學的身高，計算最少需要幾位同學出列，可以使得剩下的同學排成合唱隊形。

輸入：

　　輸入的第一行是一個整數N（2 ≤ N ≤ 100），表示同學的總數。第二行有N個整數，用空格分隔，第i個整數Ti（130 ≤ T_i ≤ 230）是第i位同學的身高（釐米）。

輸出：

　　一行，這一行只包含一個整數，就是最少需要幾位同學出列。

樣例輸入：

8

186 186 150 200 160 130 197 220

樣例輸出：

4

【分析】最長上升子序列+最長下降子序列綜合

        我們按照由左而右和由右而左的順序，將n個同學的身高排成數列。如何分別在這兩個數列中尋求遞增的、未必連續的最長子序列，就成為問題的關鍵。設：

　　a為身高序列，其中a[i]為同學i的身高；

　　b為由左而右身高遞增的人數序列，其中 b[i]為同學1‥同學i間（包括同學i）身高滿足遞增順序的最多人數。顯然b[i]=max{b[j]|同學j的身高<同學i的身高}+1；

　　c為由右而左身高遞增的人數序列，其中c[i]為同學n‥同學i間（包括同學i）身高滿足遞增順序的最多人數。顯然c[i]=max{c[j]|同學j的身高<同學i的身高}+1；

　　由上述狀態轉移方程可知，計算合唱隊形的問題具備了最優子結構性質(要使b[i]和c[i]最大，子問題的解b[j]和c[k]必須最大(1≤j≤i-1，i+1≤k≤n))和重迭子問題的性質(為求得b[i]和c[i]，必須一一查閱子問題的解b[1]‥b[i-1]和c[i+1]‥c[n])，因此可採用動態程式設計的方法求解。

　　顯然，合唱隊的人數為max{b[i]+c[i]}-1(公式中同學i被重複計算，因此減1)，n減去合唱隊人數即為解。

#include <iostream>
#include <cstdio>
#include <cstring>
using namespace std;
const int maxn=105;
int N,a[maxn];
//idp[i]:以ai為末尾的最長上升子序列長度 ddp[i]:以ai為開頭的最長下降子序列長度
int idp[maxn],ddp[maxn];
int main()
{
    int i,j;
    int ret;     //ret記錄ai左側LIS和右側LDS的長度和
    while(scanf("%d",&N)!=EOF)
    {
        for(i=1;i<=N;i++)
        {
            scanf("%d",&a[i]);
            idp[i]=ddp[i]=1;
        }
        ret=0;
        for(i=1;i<=N;i++)    //順推，求maxlis
        {
            for(j=1;j<i;j++)
            {
                if(a[j]<a[i])
                    idp[i]=max(idp[i],idp[j]+1);
            }
        }
        for(i=N;i>=1;i--)    //逆推，求maxlds
        {
            for(j=i+1;j<=N;j++)
            {
                if(a[j]<a[i])
                    ddp[i]=max(ddp[i],ddp[j]+1);
            }
        }
        for(i=1;i<=N;i++)    //注意每個a[i]計算了兩次
        {
            printf("%d %d\n",idp[i],ddp[i]);
            if(idp[i]+ddp[i]>ret)
                ret=idp[i]+ddp[i];
        }
        printf("%d\n",N-ret+1);  //這裡N-ret+1即為結果，注意加1
    }
    return 0;
}