金山程序题2优化续 之[生成函数篇]

上篇帖子链接:金山程序题2的优化

照例先给出题目:

  给定一个数组大小m和一个数组array
  求从array中任意取得n(n<=m)个数,使得和为m,总共有多少种取法.
  例如: m = 10;  array = {1,2,3,4,5,6,7,8,9,10}  共10种取法
 1 2 3 4,  1 2 7,  1 3 6,  1 4 5,  1 9,  2 3 5,  2 8,  3 7,  4 6,  10
  方法原型 :int  getotalNum (int[] array, int m);

  这道题目我从11月11日收到起至今已经过去一个月了,处理进度分三个阶段:
  1、11.11 当天利用组合数学的算法得到初步结果,由于优化失误,当时只能正确处理正数
  2、12.06 优化了一次,时间复杂度仍然是O(sigma(1,i,n),C(n,i)),但是剪枝效果不错,比第一次效率有了百倍的提升,但是由于没有发现代码中的失误,仍旧不能正确处理负数
  3、12.13 根据腌菜同学的建议,参考了Matrix67的资料和腌菜的核心代码,进行了第二次优化。时间复杂度剧烈降低至O(m*n),而且可以正确处理负数情况,是远超上次的巨大提升。
  今天这篇日至就是12.13日优化过程中的一些情况和资料的总结。
  首先介绍了是本文的重点 生成函数,下面几段摘自Matrix67的博文:什么是生成函数? (有删改)

________________________简易分割线______________________________

  我们年级有许多漂亮的MM。一班有7个左右吧,二班大概有4个,三班最多,16个,四班最可怜,一个漂亮的MM都没有,五班据说有1个。如果用一个函数“f(班级)=漂亮MM的个数”,那么我们可以把上述信息表示成:f(1)=7,f(2)=4,f(3)=16,f(4)=0,f(5)=1,等等。
  生成函数是说,构造这么一个多项式函数g(x),使得x的n次方系数为f(n)。于是,上面的f函数的生成函数g(x)=7x+4x^2+16x^3+x^5+…。这就是传说中的生成函数了。
  生成函数最绝妙的是,某些生成函数可以化简为一个很简单的函数。也就是说,不一定每个生成函数都是用一长串多项式来表示的。比如,这个函数f(n)=1 (n当然是属于自然数的),它的生成函数就应该是g(x)=1+x+x^2+x^3+x^4+…(每一项都是一,即使n=0时也有x^0系数为1,所以有常数项)。再仔细一看,这就是一个有无穷多项的等比数列求和嘛。如果-1<x<1,那么g(x)就等于1/(1-x)了。在研究生成函数时,我们都假设级数收敛,因为生成函数的x没有实际意义,我们可以任意取值。于是,我们就说,f(n)=1的生成函数是g(x)=1/(1-x)。
  我们举一个例子说明,一些具有实际意义的组合问题也可以用像这样简单的一个函数全部表示出来。
  例1:从二班选n个MM出来有多少种选法。学过简单的排列与组合的同学都知道,答案就是C(4,n)。也就是说。从n=0开始,问题的答案分别是1,4,6,4,1,0,0,0,…(从4个MM中选出4个以上的人来方案数当然为0喽)。那么它的生成函数g(x)就应该是g(x)=1+4x+6x^2+4x^3+x^4。这不就是……二项式展开吗?于是,g(x)=(1+x)^4。
  我们再举一个例子说明一些更复杂的生成函数。例2:k=x1+x2+x3+…+xn有多少个非负整数解?这道题是学排列与组合的经典例题了。把每组解的每个数都加1,就变成n+k=x1+x2+x3+…+xn的正整数解的个数了。教材上或许会出现这么一个难听的名字叫“隔板法”:把n+k个东西排成一排,在n+k-1个空隙中插入n-1个“隔板”从而把数字分成n块,每块的东西个数为每个自变量的值。这样隔板的放置的方法和解的个数之间形成了一一对应的关系。至于隔板放置的方法我们总是知道的,就是从n+k-1个空隙中无序的找出n-1个用来放置隔板,就是C(n+k-1,n-1)。它就等于C(n+k-1,k)。而它关于n的生成函数是g(x)=C(n+0-1,0)+C(n+1-1,1)x+…+C(n+k-1,k)x^k+…. = 1/(1-x)^n 这个生成函数是如何换算来的呢?
  1/(1-x)=1+x+x^2+x^3+x^4+…是前面说过的。我们对这个式子等号两边同时求导数。于是,1/(1-x)^2=1+2x+3x^2+4x^3+5x^4+….。不断地再求导数,得到了这样一个公式:1/(1-x)^n=1+C(n,1)x^1+C(n+1,2)x^2+C(n+2,3)x^3+…+C(n+k-1,k)x^k+…。就是上面我们得到的那个公式。分析这个公式g(x)=1/(1-x)^n=(1+x+x^2+x^3+…)^n,仔细想想n个(1+x+x^2+x^3+…)相乘是什么意思。(1+x+x^2+x^3+…)^n的展开式中,k次项的系数就是我们的答案,因为它的这个系数是由原式完全展开后n个指数加起来恰好等于k的项合并起来得到的。
  所以我们总结下例2的规律,x1至xn每个数字都能取从1到k之间的值,所以每个自变量对应的生成函数 g(x)=1+x+x^2+x^3+..+x^n =1/(1-x),而整个题目的生成函数则是,G(x)=1/(1-x)^n=g(x)*g(x)*…*g(x),也就是全部自变量的生成函数的乘积。得到之后总的生成函数之后,计算出对应指数k的系数就是我们题目所要的!下面几个例题,可以理解下加深印象:
  例3:我们要从苹果、香蕉、橘子和梨中拿一些水果出来,要求苹果只能拿偶数个,香蕉的个数要是5的倍数,橘子最多拿4个,梨要么不拿,要么只能拿一个。问按这样的要求拿n个水果的方案数。
   G(x)=(1+x^2+x^4+…)*(1+x^5+x^10+…)*(1+x+x^2+x^3+x^4)*(1+x) 
    =…(划减步骤略,套等比数列通相和公式)
    =(1-x)^(-2)
     =C(1,0)+C(2,1)x+C(3,2)x^2+C(4,3)x^3…
     =1+2x+3x^2+4x^3+5x^4+….
    指数为n的系数是n+1,故n+1就是我们所求得解。

________________________简易分割线______________________________

  以上大部分内容来自于Matrix67的博文,dave为了使大家更好理解删改了一部分内容,这里是原文连接:http://www.matrix67.com/blog/archives/120

  下面就结合咱们的程序题目用生成函数的思想来解决这个问题。整理题目如下:
  无序整数数列a[0…m-1], 对于数列的任何一个数字可以标记或者不标记,使得标记的数字之和为m。
  于是:G(x)=(1+x^a[0])*(1+x^a[1])*..*(1+x^a[m-1]),这里每个多项式退化成二项式,由于我们的数列a是用户输入的,无序的,所以这个式子没有办法用数学方法划减,只能用计算机来处理。
  我们可以用一个长m+1的数数组r来存储G(x)从0到m指数的系数。把每个a[ i ]看作每个二项式的指数,依次相乘。但是问题又来了,如果指数a[ i ]是小于0的,怎么办?我们在r数组用下标表示指数,但是下标不能为负数。于是在进行计算前,我们应该把所有小于0的指数都变成非负指数。
  例如对于(1+x^-2)我可以让它乘上x^2,变成(x^2+1),可以作为正常的二项式参加运算,G(x)变成G(x)*x^2,而我们要求指数m的系数,也变成求目的指数n=m+2的系数,数组r的长度也必须增加到=n+1 = m+1+2。
  对于处理二项式的乘法,由于每个二项式中必定有一个常数项1,用1乘另外一个多项式p,多项式不变p,故可以忽略这一步 ,直接加上另外一项和p的乘积。
  例如(1+x^k)*(1+x^i)=(1+x^k)+(x^i+x^(i+k)) 如果k,i或者k+i任何一项超过了目的指数n,就可以舍弃。因为在指数都是正数的情况下,一旦超越了n不可能再变小的。在两个二项式相乘的过程中,每个r[j] (0<=j<=n)都要和a[i]相乘,故时间复杂度是O(n),对于m个二项式来言,整体复杂度就是O(m*n).
  好了原理就解释到这里,下面给出代码:

/************************************************************
*名称:    problem02.c                                        *
*描述:    整形数组中取n个元素和等于元素个数m(m&gt;=n)		     *
*假设:    最大数-sum(负数) 与 取法总数 均不超过INT_MAX           *
*思路:    利用生成多项式的原理来求目的指数的系数                      *
*环境:    Code::Blocks & Windows 7 & x86                     *
*备注:    davelv于09-12-14                                   *
*************************************************************/
#include &lt;stdio.h&gt;
#include &lt;stdlib.h&gt;
#include &lt;time.h&gt;

//求有几种取法的函数
int  getotalNum (int array[], int m);
//多项式相乘
int multiplyPolynomial(const int array[], int result[], const int m, const int n);

int main(void)
{
    int *buf, m, i;

    printf("Please input NO. of integer(s)");
    scanf("%d",&m);
    printf("Please input integer(s):");
    buf = (int*) malloc( m * sizeof(int) );
    if (buf == NULL)
        return -1;
    for (i = 0; i &lt; m; ++i)
    {
        scanf("%d", buf + i);
    }
    printf("Totle %d way(s)/n", getotalNum(buf, m));
    return 0;
}

//返回取法的总数,如果是-1则表示函数失败
//array[]([0,m-1]),待处理数列
//m 数列中数字个数
//原理为利用生成函数得到生成多项式,求目的指数的系数即为取法总数
int  getotalNum (int array[], int m)
{
    int i, n, totle;//累加变量,目的指数,目的系数
    int *r;//计算缓存指针
    time_t t = clock();
    //把负数处理成正数
    for (i=n=0; i &lt; m; ++i)
    {
        if (array[i] &lt; 0)
        {
            array[i] = -array[i];
            n += array[i];//目的指数随之增加
        }
    }
    n += m;//加上原指数m,得到目的指数
    r = (int *)calloc( n+1, sizeof(int) );//分配并清零缓存
    if (r == NULL)
    {
        totle = -1;
    }
    else
    {
        totle = multiplyPolynomial(array, r, m, n);
        free(r);
    }
    printf("%ldms/n",clock()-t);
    return totle;
}
//返回 指数为n的多项式的系数
//m是多项式个数
//n最高指数/待求系数所在的变量指数
//array[i]([0.m-1])待乘多项式的指数向量
//result[i]([0,n]):result存放多项式乘积结果,i为当前项指数,整个数组提前清零。
int multiplyPolynomial(const int array[], int result[], const int m, const int n)
{
    int i , j , t;//循环变量,循环变量,中间变量

    result[0] = 1;//常数相的系数为1
    for (i = 0; i &lt; m; ++i)  // 对 array 中每个指数/多项式进行操作
    {
        for (j = n; j &gt;= 0; --j) //当前多项式/指数与上次结果求笛卡尔积
        {
            t = j + array[i];
            if ( t &lt;= n)
                result[t] += result[j];  //只有指数和小于最高指数时才会保存
        }
    }
    return result[n];//返回结果--指数为n的系数
}

 

  新的函数代码量降低了近一半(除去注释),更简洁,更优雅,更高效。
  使用0至n-1共n个数据作为一个测试项,测试数据测试如下(由于数据量问题,仅测试运行时间,对于数据溢出没有作处理):
 

N		old(ms)			new(ms)  
100		560			0
150		14490			0
200		402009			0
500		/			1
1000		/			6
2000		/			25

  可见时间复杂度从O(sigma(1,i,n),C(n,i))提升到O(m*n)的效果是极其显著的,。
  结论:
  1、时间复杂度的下降是提高运行效率的最有效手段
  2、数学方法非常重要,尤其对于这种数值处理的问题。
  3、PKU和TSU的同学果然都是无比强大,再次感谢Martix67和腌菜同学,并致以崇高的敬意!
 

同时欢迎各位同学前来批评指正,只有讨论的越深入,我们才越能了解事物的本质。

This entry was posted in 程序设计 and tagged . Bookmark the permalink.

20 Responses to 金山程序题2优化续 之[生成函数篇]

  1. zydj_2006 says:

    感谢分享,你的文章已经推荐到学生大本营首页。希望你的经验能够帮助更多的同学。

  2. 不错,我也要多学习啊,算法的力量在次证明

    算法就是王道,我喜欢。。。呵呵~!

    谢谢你的改进。

  3. NMJKTAO says:

    呵呵,高手,真正的高手哦!

  4. davelv says:

    <div class="quote"><span class="q"><b>CSDN学生大本营辅导员 张媛(职业指导爱好者)<img src=image/group/admin.gif border=0></b>: 感谢分享,你的文章已经推荐到学生大本营首页。希望你的经验能够帮助更多的同学。</span></div>谢谢辅导员的推荐,我会继续努力的

  5. davelv says:

    <div class="quote"><span class="q"><b>中国地质大学(武汉) 杨辉鱼(C/C++学生)</b>: 不错,我也要多学习啊,算法的力量在次证明

    算法就是王道,我喜欢。。。呵呵~!

    谢谢你的改进。</span></div>当初用数学验证了这个方法的正确性之后,我还是真的大为震惊呢.

    数学真是个好东西啊…

    以后得看看数值分析和组合数学方面的内容

  6. 我想起来了,这个题实际上是个动态规划题,哎哟!

    现在才想想到,01背包问题,就是dp,复杂度就为o(m^2)。

    哎哟 ,实在是…

    我可以把代码写写。

  7. davelv says:

    <div class="quote"><span class="q"><b>中国地质大学(武汉) 杨辉鱼(C/C++学生)</b>: 我想起来了,这个题实际上是个动态规划题,哎哟!

    现在才想想到,01背包问题,就是dp,复杂度就为o(m^2)。

    哎哟 ,实在是…

    我可以把代码写写。</span></div>DP问题么?我怎么感觉某些地方不适合DP.

    代码写出来偶要看..THX!

  8. 今天可能不行了

  9. 等等在说吧。

    我一定给你写个。

    主要好久没写算法题了,连dp都反应一半天。

    不过我很确定,就是01背包问题。

  10. davelv says:

    <div class="quote"><span class="q"><b>中国地质大学(武汉) 杨辉鱼(C/C++学生)</b>: 等等在说吧。

    我一定给你写个。

    主要好久没写算法题了,连dp都反应一半天。

    不过我很确定,就是01背包问题。</span></div>Matrix67也说跟DP有关系.

    时间长没有关系,把核心算法写出来先看看.

  11. /*

    数据规范:

    输入:

    第一行为数组大小

    第二行为元素

    输出:

    方法数

    input sample

    10

    1 2 3 4 5 6 7 8 9 10

    output sample

    10

    */

    #include <iostream>

    using namespace std;

    //数组做多为100,000个数,数组元素为非负数.

    const int N = 100000;

    const int M = 100000;

    int dp[N+1];

    int a[M],m;

    void DP()

    {

            memset(dp,0,sizeof(dp));

            dp[a[0]] = 1;

            for(int i =1; i < m;i++)

            {

                    for(int j = N ; j >= 0; j–)

                    {

                            if(dp[j])

                            {

                                    dp[j+a<i>] += dp[j];

                            }

                    }

                    dp[a<i>] ++;

         

  12. 上面程序那个a是a<i>,不知怎么会没显示。

    我只写了个处理非负数的,处理负数应该差不多,还要改改。

    不过就是dp了,你找几组数据测试一下。

  13. davelv says:

    <div class="quote"><span class="q"><b>中国地质大学(武汉) 杨辉鱼(C/C++学生)</b>: 上面程序那个a是a<i>,不知怎么会没显示。

    我只写了个处理非负数的,处理负数应该差不多,还要改改。

    不过就是dp了,你找几组数据测试一下。</span></div>跟dz!论坛一样,i下标作为斜体处理了,我看下.谢谢你咯

  14. davelv says:

    我用了几个测试数据程序都正确.

    DP()函数有这行代码:for(j = N ; j >= 0; j–) 这里的N可以用m代替,以获得更好的效率

    DP跟生成函数代码出奇的相似,莫非他们之间还有某种神秘的关系不成?我好好研究下

  15. franklong18 says:

    这是面试题?

  16. davelv says:

    <div class="quote"><span class="q"><b>中国科学技术大学 龙柏(C/C++学生)</b>: 这是面试题?</span></div>是面试后,技术人员给的,让你在几个小时之内完成交给他..

  17. 对的,,继续加油

  18. howlet2 says:

    虽然没看懂,但能将这个算法复杂度降到O(m*n) 一定是好算法

回复 zydj_2006 取消回复

您的电子邮箱地址不会被公开。 必填项已用*标注