binomial.C File Reference

Detailed Description

tutorial illustrating the use of TMath::Binomial can be run with:

root > .x binomial.C
root > .x binomial.C+ with ACLIC

0.0268068313599
4.31177592278

Processing /mnt/vdb/lsf/workspace/root-makedoc-v608/rootspi/rdoc/src/v6-08-00-patches/tutorials/math/binomial.C...

TMath::Binomial simple test
Build the Tartaglia triangle
============================
n= 0                             1
n= 1                           1   1
n= 2                         1   2   1
n= 3                       1   3   3   1
n= 4                     1   4   6   4   1
n= 5                   1   5  10  10   5   1
n= 6                 1   6  15  20  15   6   1
n= 7               1   7  21  35  35  21   7   1
n= 8             1   8  28  56  70  56  28   8   1
n= 9           1   9  36  84 126 126  84  36   9   1
n=10         1  10  45 120 210 252 210 120  45  10   1
n=11       1  11  55 165 330 462 462 330 165  55  11   1
n=12     1  12  66 220 495 792 924 792 495 220  66  12   1

TMath::Binomial fancy test
Verify Newton formula for (x+y)^n
x,y in [-2,2] and n from 0 to 9  
=================================
Average Error = 1.064988e-13

#include <TMath.h>
#include <TRandom.h>

void binomialSimple() {
  //
  // Simple test for the binomial distribution
  //
  printf("\nTMath::Binomial simple test\n");
  printf("Build the Tartaglia triangle\n");
  printf("============================\n");
  const Int_t max=13;
  Int_t j;
  for(Int_t i=0;i<max;i++) {
    printf("n=%2d",i);
    for(j=0;j<(max-i);j++) printf("  ");
    for(j=0;j<i+1;j++) printf("%4d",TMath::Nint(TMath::Binomial(i,j)));
    printf("\n");
  }
}

void binomialFancy() {
  Double_t x;
  Double_t y;
  Double_t res1;
  Double_t res2;
  Double_t err;
  Double_t serr=0;
  const Int_t nmax=10000;
  printf("\nTMath::Binomial fancy test\n");
  printf("Verify Newton formula for (x+y)^n\n");
  printf("x,y in [-2,2] and n from 0 to 9  \n");
  printf("=================================\n");
  TRandom r;
  for(Int_t i=0; i<nmax; i++) {
    do {
        x=2*(1-2*r.Rndm());
        y=2*(1-2*r.Rndm());
    } while (TMath::Abs(x+y)<0.75); //Avoid large cancellations
    for(Int_t j=0; j<10; j++) {
       res1=TMath::Power(x+y,j);
       res2=0;
       for(Int_t k=0; k<=j; k++)
          res2+=TMath::Power(x,k)*TMath::Power(y,j-k)*TMath::Binomial(j,k);
       if((err=TMath::Abs(res1-res2)/TMath::Abs(res1))>1e-10)
          printf("res1=%e res2=%e x=%e y=%e err=%e j=%d\n",res1,res2,x,y,err,j);
       serr +=err;
     }
  }
  printf("Average Error = %e\n",serr/nmax);
}

void binomial () {
   binomialSimple();
   binomialFancy();
}

Author

Federico Carminati

Definition in file binomial.C.