share/odepack/fortran/jgroup.f

   1 *DECK JGROUP
   2       SUBROUTINE JGROUP (N,IA,JA,MAXG,NGRP,IGP,JGP,INCL,JDONE,IER)
   3       INTEGER N, IA, JA, MAXG, NGRP, IGP, JGP, INCL, JDONE, IER
   4       DIMENSION IA(*), JA(*), IGP(*), JGP(*), INCL(*), JDONE(*)
   5 C-----------------------------------------------------------------------
   6 C This subroutine constructs groupings of the column indices of
   7 C the Jacobian matrix, used in the numerical evaluation of the
   8 C Jacobian by finite differences.
   9 C
  10 C Input:
  11 C N      = the order of the matrix.
  12 C IA,JA  = sparse structure descriptors of the matrix by rows.
  13 C MAXG   = length of available storage in the IGP array.
  14 C
  15 C Output:
  16 C NGRP   = number of groups.
  17 C JGP    = array of length N containing the column indices by groups.
  18 C IGP    = pointer array of length NGRP + 1 to the locations in JGP
  19 C          of the beginning of each group.
  20 C IER    = error indicator.  IER = 0 if no error occurred, or 1 if
  21 C          MAXG was insufficient.
  22 C
  23 C INCL and JDONE are working arrays of length N.
  24 C-----------------------------------------------------------------------
  25       INTEGER I, J, K, KMIN, KMAX, NCOL, NG
  26 C
  27       IER = 0
  28       DO 10 J = 1,N
  29  10     JDONE(J) = 0
  30       NCOL = 1
  31       DO 60 NG = 1,MAXG
  32         IGP(NG) = NCOL
  33         DO 20 I = 1,N
  34  20       INCL(I) = 0
  35         DO 50 J = 1,N
  36 C Reject column J if it is already in a group.--------------------------
  37           IF (JDONE(J) .EQ. 1) GO TO 50
  38           KMIN = IA(J)
  39           KMAX = IA(J+1) - 1
  40           DO 30 K = KMIN,KMAX
  41 C Reject column J if it overlaps any column already in this group.------
  42             I = JA(K)
  43             IF (INCL(I) .EQ. 1) GO TO 50
  44  30         CONTINUE
  45 C Accept column J into group NG.----------------------------------------
  46           JGP(NCOL) = J
  47           NCOL = NCOL + 1
  48           JDONE(J) = 1
  49           DO 40 K = KMIN,KMAX
  50             I = JA(K)
  51  40         INCL(I) = 1
  52  50       CONTINUE
  53 C Stop if this group is empty (grouping is complete).-------------------
  54         IF (NCOL .EQ. IGP(NG)) GO TO 70
  55  60     CONTINUE
  56 C Error return if not all columns were chosen (MAXG too small).---------
  57       IF (NCOL .LE. N) GO TO 80
  58       NG = MAXG
  59  70   NGRP = NG - 1
  60       RETURN
  61  80   IER = 1
  62       RETURN
  63 C----------------------- End of Subroutine JGROUP ----------------------
  64       END