#include <math.h>
#include <stdio.h>

  int ix=1,iy=1;

void main(){

  float  a,b;
  long int i,j,n,m;
  double urnd(void);
  double x,x0,u,w;
  double x1=0.0,x2=0.0;

  for(i=1;i<=10000;i++) urnd();

  scanf("%f%f%ld%ld",&a,&b,&m,&n);

  x=0.5;
  for(i=-m+1;i<=n;i++){
    x0=urnd();
    w=exp( (a-1.)*log(x0)+(b-1.)*log(1.-x0)
          -(a-1.)*log(x )-(b-1.)*log(1.-x ) );
    u=urnd();
    if( u<w ) x=x0;
    if( i >= 1 ){
      x1+=x/((double)n);
      x2+=x*x/((double)n);
    }
  }

  printf("# of Burn-in      = %10ld\n",m);
  printf("# of Random Draws = %10ld\n",n);
  printf("Parameters = (%7.1f,%7.1f)\n",a,b);
  printf("Mean     = %10.5lf, which should be close to %10.5f\n"
         ,x1,a/(a+b));
  printf("Variance = %10.5lf, which should be close to %10.5f\n"
         ,x2-x1*x1,a*b/((a+b)*(a+b)*(a+b+1.)));

}
/* ------------------------------------------------ */
double urnd(void)
{
   int    kx,ky;
   double rn;
/*
   Input:
     ix, iy:  Seeds
   Output:
     rn: Uniform Random Draw U(0,1)
*/
   kx=ix/53668;
   ix=40014*(ix-kx*53668)-kx*12211;

   ky=iy/52774;
   iy=40692*(iy-ky*52774)-ky*3791;

   rn=(float)(ix-iy)/2147483563.;
   rn-=(int)rn;
   if( rn<0.) rn++;

   return rn;
}