题目:LightOJ - 1284 Lights inside 3D Grid(概率期望)
原题链接:https://vjudge.net/problem/26994/origin
题意:
一个 X ∗ Y ∗ Z X*Y*Z X∗Y∗Z的立方体,n次操作,每次找出一个子立方体,给出形式为 ( x 1 , y 1 , z 1 ) , ( x 2 , y 2 , z 1 ) (x_1,y_1,z_1),(x_2,y_2,z_1) (x1,y1,z1),(x2,y2,z1),代表立方体的对角定点,对于这个子立方体中的所有开关进行反转操作,刚开始所有开关都是关闭状态,问n次操作之后处于打开状态的开关的期望个数。
题解:
前面期望题做的大多是通过dp去做,找递推式,这个题思维惯式也打算找出递推式,但是根本找不到。
正解是考虑每个点的贡献。
将每个开关最后处于打开状态的概率加起来就是处于打开状态的期望个数。
对于每个开关,考虑对它进行一次操作,设操作之后处于打开状态的概率为 p p p,对于它的x坐标,如果 x 1 x1 x1和 x 2 x2 x2是位于x的两侧的话,那么对于x坐标来说它是位于子立方体中,它的概率可以通过容斥原理求出来,即 p x = 1 − ( X − x ) ∗ ( X − x ) + ( x − 1 ) ∗ ( x − 1 ) X ∗ X p_x=1-\frac{(X-x)*(X-x)+(x-1)*(x-1)}{X*X} px=1−X∗X(X−x)∗(X−x)+(x−1)∗(x−1),同样的,对于y和z坐标都可以通过同样的方式求出来,
最后 p = p x ∗ p y ∗ p z p=p_x*p_y*p_z p=px∗py∗pz
这是单独的一个定点被选择一次的概率,最后将每个顶点被选择奇数次的概率都加起来就是期望值。
那么对于一个顶点被选择奇数次的概率要怎么求呢?
设 f ( x ) f(x) f(x)表示某个定点被选择奇数次的概率, g ( x ) g(x) g(x)为该定点被选择偶数次的概率,可以推出以下两个式子:
① f ( x ) + g ( x ) = 1 f(x)+g(x)=1 f(x)+g(x)=1
② f ( x ) = p ∗ g ( x − 1 ) + ( 1 − p ) ∗ f ( x − 1 ) f(x)=p*g(x-1)+(1-p)*f(x-1) f(x)=p∗g(x−1)+(1−p)∗f(x−1)
其中p为该点某次被选中的概率
将①②两式联立得
f
(
x
)
=
(
1
−
2
∗
p
)
∗
f
(
x
−
1
)
+
p
f(x)=(1-2*p)*f(x-1)+p
f(x)=(1−2∗p)∗f(x−1)+p
两边同时减去
1
2
\frac{1}{2}
21并提出公因子得
f
(
x
)
−
1
2
=
(
1
−
2
∗
p
)
∗
(
f
(
x
−
1
)
−
1
2
)
f(x)-\frac{1}{2}=(1-2*p)*(f(x-1)-\frac{1}{2})
f(x)−21=(1−2∗p)∗(f(x−1)−21)
再化简得
f
(
x
)
−
1
2
f
(
x
−
1
)
−
1
2
=
1
−
2
∗
p
\frac{f(x)-\frac{1}{2}}{f(x-1)-\frac{1}{2}}=1-2*p
f(x−1)−21f(x)−21=1−2∗p
则通过递推式可以得出
f
(
n
)
−
1
2
f
(
1
)
−
1
2
=
(
1
−
2
∗
p
)
n
−
1
\frac{f(n)-\frac{1}{2}}{f(1)-\frac{1}{2}}=(1-2*p)^{n-1}
f(1)−21f(n)−21=(1−2∗p)n−1
又有
f
(
1
)
=
p
f(1)=p
f(1)=p,得
f
(
n
)
=
(
1
−
2
∗
p
)
n
+
1
2
f(n)=\frac{(1-2*p)^{n}+1}{2}
f(n)=2(1−2∗p)n+1
概率公式推出来了,最后对于每个点求和就好了。
代码:
#define push_back pb
#define make_pair mk
#define rd read()
#define mem(a,b) memset(a,b,sizeof(a))
#define bug printf("*********\n");
#define lson l,mid,rt<<1
#define rson mid+1,r,rt<<1|1
#define FIN freopen(D://code//in.txt", "r", stdin);
#define debug(x) cout<<"["<<x<<"]" <<endl;
#define IO ios::sync_with_stdio(false),cin.tie(0);
#pragma comment(linker,"/STACk:1024000000,1024000000")
//#include<bits/stdc++.h>
#include<time.h>
#include<iostream>
#include<stdlib.h>
#include<stdio.h>
#include<cmath>
#include<map>
#include<algorithm>
#include<string>
#include<string.h>
#include<set>
#include<queue>
#include<stack>
#include<functional>
using std::pair;
typedef long long ll;
typedef unsigned long long ull;
typedef pair<int,int> pii;
//const double PI=acos(-1);
const int maxn = 1e5 + 10;
const int maxm = 2e3 + 10;
const int mod = 998244353;
const int inf = 0x3f3f3f3f;
const ll INF=0x3f3f3f3f3f3f3f3f;
const double dinf=1e20;
const double eps=1e-8;
using namespace std;
ll read() {
ll X = 0, p = 1; char c = getchar();
for(; c > '9' || c < '0'; c = getchar()) if(c == '-') p = -1;
for(; c >= '0' && c <= '9'; c = getchar()) X = X * 10 + c - '0';
return X * p;
}
int sgn(double x)
{
if(fabs(x)<eps) return 0;
if(x<0) return -1;
return 1;
}
//*************************************************************
double cal(int i,int m)
{
return 1.0-1.0*((m-i)*(m-i)+(i-1)*(i-1))/(m*m);
}
double dp[maxn];
int main()
{
int t;
cin>>t;
int cas=0;
while(t--)
{
int x,y,z,n;
cin>>x>>y>>z>>n;
double res=0;
for(int i=1;i<=x;i++)
{
for(int j=1;j<=y;j++)
{
for(int k=1;k<=z;k++)
{
double p=cal(i,x)*cal(j,y)*cal(k,z);
res+=0.5-0.5*pow(1-2*p,n);
}
}
}
printf("Case %d: %.8lf\n",++cas,res);
}
}