LightOJ 1284 Lights inside 3D Grid （贡献求期望）

$给你一个x*y*z的长方体，然后让你先选A点再选B点，然后把AB为两个顶点之间的翻转$
$初始都是0，问你这样操作k次，最后是1的个数的期望$
$这题不能正着算期望，因为没法找到转移$
$所以我们想到求贡献，每个位置，对最后答案的贡献$
$就是求这个位置被翻了奇数次的概率$
$假设这个位置被翻的概率为p，被翻奇数次的概率为\sum_{i=1}^kC(k,i)p^i(1-p)^{k-i},i为奇数$
$怎么求这个概率呢，二项式定理偶数部分\sum_{i=0}^kC(k,i)p^i(1-p)^{k-i},i为偶数$
$p_奇+p_偶=(1-p+p)^{2k}$
$p_偶-p_奇=(1-p-p)^{2k}$
$p_奇=0.5-0.5*(1-2p)^{2k}$
$接下来就是这个p怎么求的问题了$
$开头我想的是直接求，三个方位，分别选在点(x,y,z)的两边，并且按照顺序乘以8，但是发现这样有重复计算$
$因为如果两个点都选在(x,y,z)上，那么乘以8没有意义$
$所以要翻过来求，p=(1-p_{x轴两点同向})(1-p_{y轴两点同向})(1-p_{z轴两点同向})$
$1-p_{x轴两点同向}直接就是算出两点在不同两侧，并且避免了重复计算$

---

代码

#include <map>
#include <set>
#include <stack>
#include <queue>
#include <cmath>
#include <string>
#include <vector>
#include <cstdio>
#include <cctype>
#include <cstring>
#include <sstream>
#include <cstdlib>
#include <iostream>
#include <algorithm>
#pragma comment(linker,"/STACK:102400000,102400000")

using namespace std;
#define   MAX           10005
#define   MAXN          1000005
#define   maxnode       2005
#define   sigma_size    2
#define   lson          l,m,rt<<1
#define   rson          m+1,r,rt<<1|1
#define   lrt           rt<<1
#define   rrt           rt<<1|1
#define   middle        int m=(r+l)>>1
#define   LL            long long
#define   ull           unsigned long long
#define   mem(x,v)      memset(x,v,sizeof(x))
#define   lowbit(x)     (x&-x)
#define   pii           pair<int,int>
#define   bits(a)       __builtin_popcount(a)
#define   mk            make_pair
#define   limit         10000

//const int    prime = 999983;
const int    INF   = 0x3f3f3f3f;
const LL     INFF  = 0x3f3f;
//const double pi    = acos(-1.0);
const double inf   = 1e18;
//const double eps   = 1e-9;
const LL     mod   = 1e9+7;
const ull    mx    = 133333331;

/*****************************************************/
inline void RI(int &x) {
      char c;
      while((c=getchar())<'0' || c>'9');
      x=c-'0';
      while((c=getchar())>='0' && c<='9') x=(x<<3)+(x<<1)+c-'0';
 }
/*****************************************************/

double cal(double a,double b){
    return 1-(a-1)*(a-1)/(b*b)-(b-a)*(b-a)/(b*b);
}

double qpow(double a,int n){
    double ans=1;
    while(n){
        if(n&1) ans=ans*a;
        a=a*a;
        n>>=1;
    }
    return ans;
}
int main(){
    //freopen("in.txt","r",stdin);
    int t,kase=0;
    cin>>t;
    while(t--){
        kase++;
        double x,y,z;
        int k;
        cin>>x>>y>>z>>k;
        double ans=0;
        for(int i=1;i<=x;i++){
            for(int j=1;j<=y;j++){
                for(int h=1;h<=z;h++){
                    double p=cal(i,x)*cal(j,y)*cal(h,z); 
                    p=1-2*p;
                    ans+=0.5-0.5*pow(p,k);
                }
            }
        }
        printf("Case %d: %.6f\n",kase,ans);
    }
    return 0;
}