MINICERN/mathlib/gen/b/binom64.F

   1 *
   2 * $Id$
   3 *
   4 * $Log$
   5 * Revision 1.1.1.1  1996/04/01 15:01:49  mclareni
   6 * Mathlib gen
   7 *
   8 *
   9 #include "gen/pilot.h"
  10 #if defined(CERNLIB_DOUBLE)
  11       FUNCTION DBINOM(X,K)
  12 #include "gen/imp64.inc"
  13 C     PARAMETER (IMX = (2**31)-1)
  14       PARAMETER (IMX = 2147483647)
  15 #endif
  16 #if !defined(CERNLIB_DOUBLE)
  17       FUNCTION BINOM(X,K)
  18 C     PARAMETER (IMX = (2**31)-1)
  19       PARAMETER (IMX = 70368744177663)
  20 #endif
  21       CHARACTER NAME*(*)
  22       CHARACTER*80 ERRTXT
  23
  24       PARAMETER (NAME = 'KBINOM')
  25
  26       IF(K .LT. 0) THEN
  27        H=0
  28       ELSEIF(K .EQ. 0) THEN
  29        H=1
  30       ELSE
  31        H=X
  32        A=H
  33        DO 1 I = 2,K
  34        A=A-1
  35     1  H=H*(A/I)
  36       ENDIF
  37 #if defined(CERNLIB_DOUBLE)
  38       DBINOM=H
  39 #endif
  40 #if !defined(CERNLIB_DOUBLE)
  41       BINOM=H
  42 #endif
  43       RETURN
  44
  45       ENTRY KBINOM(N,K)
  46
  47     3 IF(N .GE. 0) THEN
  48        N1=N
  49        IS=1
  50       ELSE
  51        N1=K-N-1
  52        IS=(-1)**K
  53       ENDIF
  54       IF(K .LT. 0 .OR. K .GT. N1) THEN
  55        IH=0
  56       ELSEIF(K .EQ. 0 .OR. K .EQ. N1) THEN
  57        IH=IS
  58       ELSE
  59        IF(K+K .LE. N1) THEN
  60         K1=K
  61        ELSE
  62         K1=N1-K
  63        ENDIF
  64        H=N1
  65        A=H
  66        DO 2 I = 2,K1
  67        A=A-1
  68     2  H=H*(A/I)
  69        IF(H .LE. IMX) THEN
  70         IH=IS*NINT(H)
  71        ELSE
  72         IH=0
  73         WRITE(ERRTXT,101) N,K
  74         CALL MTLPRT(NAME,'B100.1',ERRTXT)
  75        ENDIF
  76       ENDIF
  77       KBINOM=IH
  78       RETURN
  79   101 FORMAT('INTEGER RANGE EXCEEDED, N =',I10,',  K = ',I10)
  80       END
  81