hiho 1166 高斯消元

一道高斯消元题目

题意：给你一个数列ai,(1<=i<=n);ai是0或1，现在你可以每次翻转一个区间，即将0变1，1变0，且每个区间被选择的概率为

n（n-1）/2,问全变成0的翻转次数期望是多少

思路：

另外，这个思路还有一些需要注意的地方，在确定这个方法是对的前提下，那么d中的1的个数必然是偶数的，这样的话直接对0列方程就要有限制条件（n+1）-i要是偶数才行，这样算出来的数是对的详见写法1；或者直接对1来列方程，只对偶数个1进行列方程，详见写法2。

写法1：
#include <bits/stdc++.h>
using namespace std;
#define rep(i,a,n) for (int i=a;i<n;i++)
#define per(i,a,n) for (int i=n-1;i>=a;i--)
typedef double LD;
const double eps = 1e-12;
int n,b[30],c;
LD g[30][30];
void gauss(int n) {
 rep(i,0,n) {
  int p=i;
  rep(j,i+1,n) if (fabs(g[j][i])>=fabs(g[p][i])) p=j;
  rep(j,i,n+1) swap(g[i][j],g[p][j]);
  if (fabs(g[i][i])<=1e-9) continue;
  rep(k,i+1,n+1) g[i][k]/=g[i][i]; g[i][i]=1;
  rep(j,i+1,n) {
   LD cof=-g[j][i];
   rep(k,i,n+1) g[j][k]+=cof*g[i][k];
  }
 }
 per(i,0,n) {
  rep(j,0,i) g[j][n]-=g[j][i]*g[i][n];
 }
}

int main() {
 scanf("%d",&n);
 memset(g,0,sizeof(g));
 g[n+1][n+1]=1;g[n+1][n+2]=0;//后面代表1的个数
 for(int i=0;i<=n;i++){
        if(i==0||i==1&&!((n+1-i)&1)){
            g[i][i]=1-2.*i*(n+1-i)/(n+1)/n;
            if(i+2<=n+1) g[i][i+2]-=1.*(n+1-i)*(n-i)/(n+1)/n;
            g[i][n+2]=1;
            continue;
        }
        else if(i==n&&!((n+1-i)&1)){
            g[i][i]=1-2.*i*(n-i+1)/(n+1)/n;
            if(i-2>=0) g[i][i-2]-=1.*(i)*(i-1)/(n+1)/n;
            g[i][n+2]=1;
            continue;
        }
        else if(!((n+1-i)&1)){
        g[i][i]=1-2.*i*(n+1-i)/(n+1)/n;
        if(i+2<=n+1) g[i][i+2]-=1.*(n+1-i)*(n-i)/(n+1)/n;
        if(i-2>=0)   g[i][i-2]-=1.*(i-1)*i/(n+1)/n;
        g[i][n+2]=1;
        }
 }
// for(int i=0;i<=n+1;i++){
//        for(int j=0;j<=n+1;j++){
//            printf("%lf   ",g[i][j]);
//        }
//        printf("%lf\n",g[i][n+2]);
// }
 gauss(n+2);
// for(int i=0;i<=n+1;i++)
//        printf("%lf\n",g[i][n+2]);
 rep(i,1,n+1) scanf("%d",&b[i]);
 rep(i,1,n+2) c+=b[i]!=b[i-1];
 printf("%.10f\n",(double)g[n+1-c][n+2]);
 return 0;
}

写法二：

#include <bits/stdc++.h>
using namespace std;
#define rep(i,a,n) for (int i=a;i<n;i++)
#define per(i,a,n) for (int i=n-1;i>=a;i--)
typedef double LD;
int n,b[30],c;
LD g[30][30];
void gauss(int n) {
	rep(i,0,n) {
		int p=i;
		rep(j,i+1,n) if (fabs(g[j][i])>=fabs(g[p][i])) p=j;
		rep(j,i,n+1) swap(g[i][j],g[p][j]);
		if (fabs(g[i][i])<=1e-9) continue;
		rep(k,i+1,n+1) g[i][k]/=g[i][i]; g[i][i]=1;
		rep(j,i+1,n) {
			LD cof=-g[j][i];
			rep(k,i,n+1) g[j][k]+=cof*g[i][k];
		}
	}
	per(i,0,n) {
		rep(j,0,i) g[j][n]-=g[j][i]*g[i][n];
	}
}

int main() {
	scanf("%d",&n);//两端为0的数列，两两异或得到的数列里1的个数为偶数
	g[0][0]=1;
	int p=(n+1)/2+1;//看的是1的个数，因为1只能是偶数个，而0的个数是没有限制的，所以这样去做，如果n+1是奇数，那么最多也只有n个1
	for (int i=2;i<=n+1;i+=2) {
		g[i/2][i/2]=1;
		g[i/2][i/2-1]-=1.*i*(i-1)/n/(n+1);
		g[i/2][i/2]-=2.*i*(n+1-i)/n/(n+1);
		g[i/2][i/2+1]-=1.*(n+1-i)*(n-i)/n/(n+1);
		g[i/2][p]=1;
	}
	gauss(p);
	rep(i,1,n+1) scanf("%d",&b[i]);
	rep(i,1,n+2) c+=b[i]!=b[i-1];
	printf("%d\n",c);
	printf("%.10f\n",(double)g[c/2][p]);
}