[Contests]2016弱校联盟十一专场10.7

G UmBasketella

题意：

漏斗形(锥形)容器，输入其表面积(包括底部)，求其容量，高度，半径，精确到0.01

思路：

此题有两种方法，一种是推公式，比较偷懒

还有一种是三分查找确定最值，这才是正统之路

代码1（推公式）

/**************************************************************
    Problem: BNUOJ_3856
    User: soundwave
    Language: C++
    Result: Accepted
    Time: 0ms
    Memory: 692KB
****************************************************************/
//#pragma comment(linker, "/STACK:1024000000,1024000000")//手动扩栈
#include <iostream>
#include <stdio.h>
#include <math.h>

using namespace std;
const double PI = 3.1415926;
int main(){
    double s;
    double v, h, r;

    while(~scanf("%lf", &s)){
        r = sqrt(s*1.0/4/PI);
        h = sqrt((s*s)/(PI*PI*r*r) - 2*s/PI);
        v = ((s*1.0/4/PI)*PI*h) / 3;
        printf("%.2f\n%.2f\n%.2f\n", v, h, r);
    }
    return 0;
}

代码2（魔幻三分）

/**************************************************************
    Problem: BNUOJ_3856
    User: soundwave
    Language: C++
    Result: Accepted
    Time: 0ms
    Memory: 716KB
****************************************************************/
//#pragma comment(linker, "/STACK:1024000000,1024000000")
#include <iostream>
#include <stdio.h>
#include <math.h>
/*
V = (PI*r*r*h)/3;
S = PI*r(r+l);
l = sqrt(r*r+h*h);
-------------------
V = r*sqrt(s*s-2*s*PI*r)/3;
*/
using namespace std;
//acos()是反余弦函数，cosπ = -1，所以π = acos(-1)
const double PI = acos(-1.0);
const double EPS = 1e-10;
double s;
double calc(double r){
    return sqrt(s*s-2*s*PI*r*r)*r/3.0;
}
/*
double three_divide(double l, double r){
    double lmid, rmid;
    while(r-l > EPS){
        lmid = (l+r)/2.0;
        rmid = (lmid+r)/2.0;
        if(calc(lmid) > calc(rmid))
        if(x>y)
            r = rmid;
        else
            l = lmid;
    }
    return l;
}
*/
double three_divide(double l, double r){
    double lmid, rmid;
    while(r-l > EPS){
        lmid = l + (r-l)/3.0;
        rmid = r - (r-l)/3.0;
        if(calc(lmid) > calc(rmid))
            r = rmid;
        else
            l = lmid;
    }
    return (r+l)/2.0;
}
int main(){
    while(~scanf("%lf", &s)){
        double r = three_divide(0,sqrt(s*1.0/PI));
        double h = sqrt(s*s-2*s*PI*r*r)/PI/r;
        printf("%.2f\n%.2f\n%.2f\n", calc(r), h, r);
    }
    return 0;
}

反思：

不知道为什么两种三分查找的代码只能实现一个，奇怪……

D Blocks

题意：

用四种颜色去刷N个盒子，要求红色和绿色的盒子都为偶数(0也是偶数)，求方案数。

思路：

可以用矩阵乘法(矩阵快速幂优化)，或者组合数学推公式(快速幂优化)

①这里是排列组合+二次项定理推公式

来吧，推公式_(:зゝ∠)_
我们设奇数为1，偶数为0
倒着推比较简单，4^n - 红绿色(11/10/01)的情况
从n中选k(k∈[1,n])个，k为红绿总数，从k中选取不符合题意（即红绿色为11/10/01）的情况；

根据二项式定理 C(n,k)*2^(n-k) * C(2,1)*( C(k,1)+C(k,3)+C(k,5)+…… )
(1+x)^k=C(k,0)x^0+C(k,1)x^1+...+C(k,k)x^k
C(k,1)+C(k,3)+C(k,5)+…… = 2^(k-1);
C(k,0)+C(k,2)+C(k,4)+…… = 2^(k-1);

如果 k 是奇数，红绿为一奇一偶，红奇绿偶和红偶绿奇为两种情况，则
C(n,k)*2^(n-k) * C(2,1)*( C(k,1)+C(k,3)+C(k,5)+…… )
= C(n,k)*2^(n-k+1)*2^(k-1) 
= C(n,k)*2^n;
如果 k 为偶数，则
C(n,k)*2^(n-k) * ( C(k,1)+C(k,3)+C(k,5)+…… )
= C(n,k)*2^(n-k)*2^(k-1) 
= C(n,k)*2^(n-1);

于是推出总公式为 
4^n - C(n,1)*2^n - C(n,2)*2^(n-1) - C(n,3)*2^n - C(n,4)*2^(n-1) - ……
= 4^n - 2^n*( C(n,1)+C(n,3)+……) - 2^(n-1)*( C(n,2)+C(n,4)+…… )
= 4^n - 2^n*2^(n-1) - 2^(n-1)*(2^(n-1)-1);
= 2^(2n) - 2^(2n-1) - 2^(2n-2) + 2^(n-1);
= 2^(n-1) * ( 2^(n-1) + 1 );

好了_(:зゝ∠)_，终于推出来了，就是 2^(n-1) * ( 2^(n-1) + 1 ) 这个公式

②(指数型)生成函数(母函数)推公式

生成函数(母函数)是说，构造这么一个多项式函数G(x)，使得x的n次方系数为a(n)。
母函数的思想很简单，就是把离散数列和幂级数一一对应起来，把离散数列间的相互结合关系对应成为幂级数间的运算关系，最后由幂级数形式来确定离散数列的构造
注：二项展开式，即（a+b）的n次展开式

生成函数为
G(x) = (1+x^2/2!+x^4/4!+...)²*(1+x+x^2/2!+x^3/3!+...)²  
= [0.5(e^x+e^-x)]²*e^2x            
= 1/4*(e^4x+e^2x+1)  
= 1/4∑((4x)^n/n!)+1/2∑((2x)^n/n!)+1/4   n=0,1,2……

因此x^n项的系数为
a(n) = 4^(n-1) + 2^(n-1) = 2^(n-1) * ( 2^(n-1) + 1 );

如何构造生成函数真心好迷啊，估计下次见到也不会写_(:зゝ∠)_

代码1（推公式）：

/**************************************************************
    Problem: BNUOJ_3853
    User: soundwave
    Language: C++
    Result: Accepted
    Time: 0ms
    Memory: 696KB
****************************************************************/
//#pragma comment(linker, "/STACK:1024000000,1024000000")
#include <iostream>
#include <stdio.h>
/*
2^(n-1) * (2^(n-1) + 1)
*/
using namespace std;
typedef __int64 LL;
const int mod = 10007;
//二分快速幂
LL my_pow(LL x, LL c){
    LL re = 1;
    while(c>0){
        if(c&1) re = (re*x) % mod;
        x = (x*x) % mod;
        c>>=1;
    }
    return re;
}
int main(){
    int t, n;
    scanf("%d", &t);
    while(t-->0){
        scanf("%d", &n);
        LL re = my_pow(2,n-1);
        printf("%I64d\n", (re*(re+1))%mod);
    }
    return 0;
}
/*
5
1
2
3
4
5
---------
2
6
20
72
272
*/

③矩阵乘法才是王道(๑•̀ㅂ•́)و✧，推公式太麻烦了_(:зゝ∠)_

设奇为1，偶为0
设 n 块砖头的染色情况为四种abcd
a(n)：N块砖头中红色绿色都为偶数
b(n)：N块砖头中红色为偶数，绿色为奇数
c(n)：N块砖头中红色为奇数，绿色为偶数
d(n)：N块砖头中红色绿色都为奇数

可知 (n+1) 块砖头的染色情况可以由 n 块砖头的染色情况推出
a(n+1) = 2a(n) + b(n) + c(n);
b(n+1) = a(n) + 2b(n) + d(n);
c(n+1) = a(n) + 2c(n) + d(n);
d(n+1) = b(n) + c(n) + 2d(n);
a(n+1)：(黄/蓝)*a + 绿*b + 红*c + 0*d
b(n+1)：绿*a + (黄/蓝)*b + 0*c + 红*d
c(n+1)：红*a + 0*b + (黄/蓝)*c + 绿*d
d(n+1)：0*a + 红*b + 绿*c + (黄/蓝)d

如果不能理解，我们可以采用dp的思想理解，然后转化为矩阵(๑•̀ㅂ•́)و✧

用dp[N][4]来表示N块砖块的染色情况，一共有四种状态：

dp[N][0]：N块砖头中红色绿色都为偶数
dp[N][1]：N块砖头中红色为偶数，绿色为奇数
dp[N][2]：N块砖头中红色为奇数，绿色为偶数
dp[N][3]：N块砖头中红色绿色都为奇数
状态转移方程如下：
dp[N+1][0] = 2 * dp[N][0] + 1 * dp[N][1] + 1 * dp[N][2] + 0 * dp[N][3]
dp[N+1][1] = 1 * dp[N][0] + 2 * dp[N][1] + 0 * dp[N][2] + 1 * dp[N][3]
dp[N+1][2] = 1 * dp[N][0] + 0 * dp[N][1] + 2 * dp[N][2] + 1 * dp[N][3]
dp[N+1][3] = 0 * dp[N][0] + 1 * dp[N][1] + 1 * dp[N][2] + 2 * dp[N][3]
可以清楚地看出来，和上一个公式的想法是一样的，(n+1) 块砖头的染色情况可以由 n 块砖头的染色情况推出
然后我们就可以转化为矩阵啦

然后求矩阵A的n次幂即可
或者，红奇绿偶和红偶绿奇合成一种情况，那么矩阵为

来，矩阵快速幂接住(๑•̀ㅂ•́)و✧

代码2（矩阵快速幂）

/**************************************************************
    Problem: BNUOJ_3853
    User: soundwave
    Language: C++
    Result: Accepted
    Time: 47ms
    Memory: 736KB
****************************************************************/
//#pragma comment(linker, "/STACK:1024000000,1024000000")
#include <iostream>
#include <stdio.h>
#include <vector>
using namespace std;

typedef vector<int> Vint;
typedef vector<Vint> VVint;
typedef __int64 LL;
const int MOD = 10007;
//矩阵乘法
VVint calc(VVint &A, VVint &B){
    VVint C(A.size(), Vint(A.size()));
    for(int i=0; i<A.size(); i++)
    for(int j=0; j<B[0].size(); j++)
    for(int k=0; k<B.size(); k++)
        C[i][j] = (C[i][j] + (A[i][k]*B[k][j])%MOD) % MOD;
    return C;
}
//二分快速幂
VVint my_pow(VVint A, LL c){
    VVint B(A.size(), Vint(A.size()));
    for(int i=0; i<A.size(); i++)
        B[i][i] = 1;
    while(c>0){
        if(c&1)
            B = calc(B,A);
        A = calc(A,A);
        c>>=1;
    }
    return B;
}
int main(){
    int t, n;
    scanf("%d", &t);
    while(t-->0){
        scanf("%d", &n);
        VVint A(4,Vint(4));
        A[0][0]=2, A[0][1]=1, A[0][2]=1, A[0][3]=0;
        A[1][0]=1, A[1][1]=2, A[1][2]=0, A[1][3]=1;
        A[2][0]=1, A[2][1]=0, A[2][2]=2, A[2][3]=1;
        A[3][0]=0, A[3][1]=1, A[3][2]=1, A[3][3]=2;
        A = my_pow(A,n);
        /*
        VVint A(3,Vint(3));
        A[0][0]=2, A[0][1]=1, A[0][2]=0;
        A[1][0]=2, A[1][1]=2, A[1][2]=2;
        A[2][0]=0, A[2][1]=1, A[2][2]=2;
        A = my_pow(A,n);
        */
        cout << A[0][0] << endl;
    }
    return 0;
}

反思：

尽管比起公式慢了，但是摆脱了推导公式的恐惧_(:зゝ∠)_