MINICERN/packlib/kernlib/kernnum/f011fort/fact.inc

   1 *
   2 * $Id$
   3 *
   4 * $Log$
   5 * Revision 1.1.1.1  1996/02/15 17:49:03  mclareni
   6 * Kernlib
   7 *
   8 *
   9 *
  10 * fact.inc
  11 *
  12           IF(IDIM .GE. N  .AND.  N .GT. 0)  GOTO 110
  13              CALL TMPRNT(HNAME,N,IDIM,0)
  14              RETURN
  15  110      IFAIL  =  NORMAL
  16           JFAIL  =  JRANGE
  17           NXCH   =  0
  18           DET    =  ONE
  19           DO 144    J  =  1, N
  20  120         K  =  J
  21              P  =  PIVOTF(A(J,J))
  22              IF(J .EQ. N)  GOTO 122
  23              JP1  =  J+1
  24              DO 121    I  =  JP1, N
  25                 Q  =  PIVOTF(A(I,J))
  26                 IF(Q .LE. P)  GOTO 121
  27                    K  =  I
  28                    P  =  Q
  29  121            CONTINUE
  30              IF(K .NE. J)  GOTO 123
  31  122         IF(P .GT. 0.)  GOTO 130
  32                 DET    =  ZERO
  33                 IFAIL  =  IMPOSS
  34                 JFAIL  =  JRANGE
  35                 RETURN
  36  123         DO 124    L  =  1, N
  37                 TF      =  A(J,L)
  38                 A(J,L)  =  A(K,L)
  39                 A(K,L)  =  TF
  40  124            CONTINUE
  41              NXCH      =  NXCH + 1
  42              IR(NXCH)  =  IPAIRF(J,K)
  43  130         DET     =  DET * A(J,J)
  44              A(J,J)  =  ONE / A(J,J)
  45              T  =  SIZEF(DET)
  46              IF(T .LT. G1)  THEN
  47                 DET    =  ZERO
  48                 IF(JFAIL .EQ. JRANGE)  JFAIL  =  JUNDER
  49              ELSEIF(T .GT. G2)  THEN
  50                 DET    =  ONE
  51                 IF(JFAIL .EQ. JRANGE)  JFAIL  =  JOVER
  52              ENDIF
  53              IF(J .EQ. N)  GOTO 144
  54              JM1  =  J-1
  55              JP1  =  J+1
  56              DO 143   K  =  JP1, N
  57                 S11  =  -A(J,K)
  58                 S12  =  -A(K,J+1)
  59                 IF(J .EQ. 1)  GOTO 142
  60                 DO 141  I  =  1, JM1
  61                    S11  =  DOTF(A(I,K),A(J,I),S11)
  62                    S12  =  DOTF(A(I,J+1),A(K,I),S12)
  63  141               CONTINUE
  64  142            A(J,K)    =  -S11 * A(J,J)
  65                 A(K,J+1)  =  -DOTF(A(J,J+1),A(K,J),S12)
  66  143            CONTINUE
  67  144         CONTINUE
  68  150      IF(MOD(NXCH,2) .NE. 0)  DET  =  -DET
  69           IF(JFAIL .NE. JRANGE)   DET  =  ZERO
  70           IR(N)  =  NXCH