最大子段和

参考资料：《计算机算法设计与分析》第3版：Page59-64，王晓东 著

一、问题描述：

给定由n个整数（可能为负整数）组成的序列a₁, a₂, ..., a_n，求该序列的形如(s(i,j)=a_i + a_i+1 + ... a_j, 1<=i<=j<=n)的子段和的最大值。当所有整数均为负整数时定义其最大子段和为0。最优值为：max{0, max{s(i,j), 1<=i<=j<=n}}。

二、简单算法：

用数组a[]，存储给定的n个整数a₁, a₂, ..., a_n。

int max_sum(int n,int *a,int &besti,int &bestj){int sum=0;besti=bestj=-1;for(int i=1;i<=n;++i){for(int j=i;j<=n;++j){int this_sum=0;for(int k=i;k<=j;++k){this_sum+=a[k];}if(this_sum>sum){sum=this_sum;besti=i;bestj=j;}}}return sum;}

其中算法时间复杂度为O(n^3)。事实上，s(i,j)=s(i,j-1)+a_j，故可以将最内层循环省去，避免重复计算，改进算法：

int max_sum_optimize1(int n,int *a,int &besti,int &bestj){int sum=0;besti=bestj=-1;for(int i=1;i<=n;++i){int this_sum=0;for(int j=i;j<=n;++j){this_sum+=a[j];if(this_sum>sum){sum=this_sum;besti=i;bestj=j;}}}return sum;}

改进后算法时间复杂度为O(n^2)：算法设计技巧上的改进，充分利用已经得到的结果，避免重复计算，节省时间。

三、分治算法：

对最大子段和问题的本身结构，还可以从算法设计的策略上对上述O(n^2)计算时间的算法加以更深刻的改进。从问题的解的结果上可以看出，它适合用分治法求解。将序列a[1:n]，分为长度相当的两段a[1:n/2]和a[n/2+1:n]，分别求出这两个子段的最大子段和，则序列a[1:n]的最大子段和有三种情形：

(1). a[1:n]的最大子段和与a[1:n/2]的最大子段和相同；

(2). a[1:n]的最大子段和与a[n/2+1:n]的最大子段和相同；

(3). a[1:n]的最大子段和为s(i,j), 1<=i<=n/2, n/2+1<=j<=n；

对(1)、(2)可递归求解；对(3)必有a[n/2]和a[n/2+1]在最优子序列中，故在a[1:n/2]中可计算s1=max{0, max{s(i,n/2), 1<=i<=n/2}}，并在在a[n/2+1:n]中可计算s2=max{0, max{s(n/2+1,j), n/2+1<=j<=n}}，则s1+s2即为(3)时的最优值。故最大子段和的分治算法如下：

int max_sum_dac(int *a,int left,int right){int sum=0;if(left==right){sum=a[left]>0?a[left]:0;}else{int center=(left+right)/2;int left_sum=max_sum_dac(a,left,center);int right_sum=max_sum_dac(a,center+1,right);int s1=0,s2=0,temp=0;for(int i=center;i>=left;--i){temp+=a[i];if(temp>s1){s1=temp;}}temp=0;for(int i=center+1;i<=right;++i){temp+=a[i];if(temp>s2){s2=temp;}}sum=s1+s2;if(sum<left_sum){sum=left_sum;}if(sum<right_sum){sum=right_sum;}}return sum;}int max_sum_optimize2(int n,int *a){return max_sum_dac(a,1,n);}

此算法时间复杂度为O(nlogn)。

四、动态规划算法：

若记b[j]=max{s(i,j), 1<=i<=j}，1<=j<=n，则所求最大子段和为：max{b[j], 1<=j<=n}。由b[j]的定义，当b[j-1]>0时b[j]=b[j-1]+a[j]，否则b[j]=a[j]。由此，可得：b[j]=max{b[j-1]+a[j], a[j]}, 1<=j<=n。

int max_sum_optimize3(int n,int *a){int sum=0,b=0;for(int i=1;i<=n;++i){if(b>0){b+=a[i];}else{b=a[i];}if(b>sum){sum=b;}}return sum;}

上述算法时间复杂度为O(n)，空间复杂度为O(n)（存放数组a[]）。

五、推广1-最大子矩阵和：

给定一个m行n列的整数矩阵A，试求矩阵A的一个子矩阵，使其各元素之和为最大。最大子矩阵和是最大子段和问题向二维的推广。用二维数组a[1:m][1:n]表示矩阵A。子数组a[i1:i2][j1:j2]表示左上角和右下角行列左边分别为(i1,j1)、(i2,j2)的子矩阵，其各元素之和记为：

s(i1,i2,j1,j2)=(a[i1][j1]+a[i1+1][j1]+...+a[i2][j1])+(a[i1][j1+1]+a[i1+1][j1+1]+...+a[i2][j1+1])+...+(a[i1][j2]+a[i1+1][j2]+...+a[i2][j2])

最大子矩阵和问题的最优值为：max{s(i1,i2,j1,j2), 1<=i1<=i2<=m, 1<=j1<=j2<=n}。直接枚举法时间复杂度为O(m^2 * n^2)。但有：

max{s(i1,i2,j1,j2), 1<=i1<=i2<=m, 1<=j1<=j2<=n}=max{max{s(i1,i2,j1,j2), 1<=j1<=j2<=n}, 1<=i1<=i2<=m}}=max{t(i1,i2), 1<=i1<=i2<=m}

其中：t(i1,i2)=max{s(i1,i2,j1,j2), 1<=j1<=j2<=n}。设b[j]=a[i1][j]+a[i1+1][j]+...+a[i2][j]，则t(i1,i2)=max{b[j1]+b[j1+1]+...+b[j2], 1<=j1<=j2<=n}，此问题即为一维最大子段和问题。故借助一维最大子段和问题的动态规划算法max_sum_optimize3，可设计出最大子矩阵和的动态规划算法max_sum2如下：

int max_sum2(int m,int n,int **a){int sum=0;int *b=new int[n+1];for(int i=1;i<=m;++i){for(int k=1;k<=n;++k){b[k]=0;}for(int j=i;j<=m;++j){for(int k=1;k<=n;++k){b[k]+=a[j][k];}int temp=max_sum_optimize3(n,b);if(temp>sum){sum=temp;}}}return sum;}

此算法的时间复杂度为O(m^2 * n)。

六、推广2-最大m子段和问题：

给定由n个整数（可能为负整数）组成的序列a₁, a₂, ..., a_n，以及一个正整数m且1<=m<=n，要求确定序列a₁, a₂, ..., a_n的m个不相交子段，使这m个子段的总和达到最大。最大m子段和是最大子段和问题在子段个数上的推广。

设b(i,j)表示数组a的前j项中i个子段和的最大值，且第i个子段含a[j]，1<=i<=m, i<=j<=n，则所求的最优值为：max{b(i,j), m<=j<=n}。计算b(i,j)的递归式为：b(i,j)=max{b(i,j-1)+a[j], max{b(i-1,t)+a[j], i-1<=t<j}}, 1<=i<=m, i<=j<=n。其中，b(i,j-1)+a[j]表示第i个子段含a[j-1]，max{b(i-1,t)+a[j], i-1<=t<j}表示第i个子段仅含a[j]。初始值，b(0,j)=0, 1<=j<=n; b(i,0)=0, 1<=i<=m。最大m子段和问题的动态规划算法如下：

int max_sum_m(int m,int n,int *a){if(n<m||m<1){return 0;}int **b=new int*[m+1];for(int i=0;i<=m;++i){b[i]=new int[n+1];}for(int i=0;i<=m;++i){b[i][0]=0;}for(int j=1;j<=n;++j){b[0][j]=0;}for(int i=1;i<=m;++i){b[i][i]=b[i-1][i-1]+a[i];for(int j=i+1;j<=n;++j){b[i][j]=b[i][j-1]+a[j];for(int k=i-1;k<j;++k){if(b[i][j]<b[i-1][k]+a[j]){b[i][j]=b[i-1][k]+a[j];}}}}int sum=0;for(int j=m;j<=n;++j){if(sum<b[m][j]){sum=b[m][j];}}for(int i=0;i<=m;++i){delete[] b[i];}delete[] b;return sum;}

此时算法时间复杂度为O(mn^2)，空间复杂度为O(mn)。由b(i,j)的定义式可知，在计算b(i,j)时只用到b(i-1,t), i-1<=t<j和b(i,j-1)：故算法中只需保存数组b的当前行，而不必存储整个数组；max{b(i-1,t)+a[j], i-1<=t<j的值可以在计算第b(i-1,j)时预先计算并保存起来，计算b(i,j)时不必重新计算，节省了计算时间和空间；另外，最后要求max{b(m,j), m<=j<=n}，故在求b(i,j), 1<=i<=m时，只需求b(i,t), i<=t<=n-m+i，而b(i,t), n-m+i+1<=t<=n不必计算。改进后算法如下：

int max_sum_m_optimize(int m,int n,int *a){if(n<m||m<1){return 0;}int *b=new int[n+1];int *c=new int[n+1];for(int j=0;j<=n-m;++j){c[j]=0;}b[0]=0;for(int i=1;i<=m;++i){b[i]=b[i-1]+a[i];int max=b[i];for(int j=i+1;j<=n-m+i;++j){b[j]=b[j-1]>c[j-1]?b[j-1]+a[j]:c[j-1]+a[j];c[j-1]=max;if(max<b[j]){max=b[j];}}c[n-m+i]=max;}int sum=0;for(int j=m;j<=n;++j){if(sum<b[j]){sum=b[j];}}return sum;}

此时，算法复杂度为O(m(n-m))，空间复杂度为O(n)，当m或(n-m)为常数时，算法复杂度为O(n)，空间复杂度为O(n)。

七、测试代码：

#include<cstdio>#include<cstdlib>#include<ctime>#define random(x) (rand()%(x)-(x)/2)int max_sum(int n,int *a,int &besti,int &bestj);int max_sum_optimize1(int n,int *a,int &besti,int &bestj);int max_sum_dac(int *a,int left,int right);int max_sum_optimize2(int n,int *a);int max_sum_optimize3(int n,int *a);int max_sum2(int m,int n,int **a);int max_sum_m(int m,int n,int *a);int max_sum_m_optimize(int m,int n,int *a);int main(int argc,char *argv[]){const int SIZE=20;srand(time(0));int a1[SIZE+1];for(int i=1;i<=SIZE;++i){a1[i]=random(101);printf("%d,",a1[i]);}printf("\n");int besti,bestj,sum;sum=max_sum(SIZE,a1,besti,bestj);printf("%d,%d,%d\n",sum,besti,bestj);sum=max_sum_optimize1(SIZE,a1,besti,bestj);printf("%d,%d,%d\n",sum,besti,bestj);sum=max_sum_optimize2(SIZE,a1);printf("%d\n",sum);sum=max_sum_optimize3(SIZE,a1);printf("%d\n\n",sum);sum=max_sum_m(1,SIZE,a1);printf("%d\n",sum);sum=max_sum_m(SIZE/3,SIZE,a1);printf("%d\n",sum);sum=max_sum_m(SIZE/2,SIZE,a1);printf("%d\n",sum);sum=max_sum_m(SIZE,SIZE,a1);printf("%d\n\n",sum);sum=max_sum_m_optimize(1,SIZE,a1);printf("%d\n",sum);sum=max_sum_m_optimize(SIZE/3,SIZE,a1);printf("%d\n",sum);sum=max_sum_m_optimize(SIZE/2,SIZE,a1);printf("%d\n",sum);sum=max_sum_m_optimize(SIZE,SIZE,a1);printf("%d\n\n",sum);int **a2;a2=new int*[SIZE+1];for(int i=1;i<=SIZE;++i){a2[i]=new int[SIZE+1];}for(int i=1;i<=SIZE;++i){for(int j=1;j<=SIZE;++j){a2[i][j]=random(101);}}sum=max_sum_m_optimize(1,SIZE,a2[1]);printf("%d\n",sum);sum=max_sum_optimize3(SIZE,a2[1]);printf("%d\n\n",sum);sum=max_sum2(1,SIZE,a2);printf("%d\n",sum);sum=max_sum2(SIZE/3,SIZE,a2);printf("%d\n",sum);sum=max_sum2(SIZE/2,SIZE,a2);printf("%d\n",sum);sum=max_sum2(SIZE,SIZE,a2);printf("%d\n",sum);for(int i=1;i<=SIZE;++i){delete[] a2[i];}delete[] a2;return 0;}