anchor.diff:  !                | A41 A42 A43 A44 A45 B46 B47 |
   !                | A51 A52 A53 A54 A55 B56 B57 |
   !
  -! This procedure performs Gauss-Jordan elimination with repeated calls
  -! to MGJECOL to zero out non-diagonal elements, and normalize to 1.0 
  -! resulting in the following:
  -!
  -!                | 1.0 0.0 0.0 0.0 0.0 X16 X17 |
  -!                | 0.0 1.0 0.0 0.0 0.0 X26 X27 |
  -!                | 0.0 0.0 1.0 0.0 0.0 X36 X37 |
  -!                | 0.0 0.0 0.0 1.0 0.0 X46 X47 |
  -!                | 0.0 0.0 0.0 0.0 1.0 X56 X57 |
  +! This procedure performs Gauss elimination with repeated calls
  +! to MGECOL leave an upper triangle matrix.  Junk values are left in
  +! place below the diagonal to save operations.  This results in the
  +! following:
  +!
  +!                | A11 A12 A13 A14 A15 B16 B17 |
  +!                | JNK A22 A23 A24 A25 B26 B27 |
  +!                | JNK JNK A33 A34 A35 B36 B37 |
  +!                | JNK JNK JNK A44 A45 B46 B47 |
  +!                | JNK JNK JNK JNK A55 B56 B57 |
  +!
  +! The final step is to perform back substitution callin MBSUB resulting
  +! in the following:
  +!
  +!                | A11 A12 A13 A14 A15 X16 X17 |
  +!                | JNK A22 A23 A24 A25 X26 X27 |
  +!                | JNK JNK A33 A34 A35 X36 X37 |
  +!                | JNK JNK JNK A44 A45 X46 X47 |
  +!                | JNK JNK JNK JNK A55 X56 X57 |
   !
   ! X values are where the solution matrix is stored.
   !
   !------------------------------------------------------------------------------+
  -      SUBROUTINE MGAUSJ(M, ROWS, COLS)
  -      IMPLICIT NONE
  -      DOUBLE PRECISION M(ROWS,COLS),TMP
  -      INTEGER ROWS, COLS
  -      INTEGER I,J
  -!
  -      IF (ROWS.GT.COLS) THEN
  -        WRITE(0,910), ROWS, COLS
  -        STOP
  -      END IF
  -!
  -      DO 110 I = 1, ROWS
  -      CALL MPPVT(M,ROWS,COLS,I)
  -      CALL MGJECOL(M,ROWS,COLS,I)
  -      TMP = M(I,I)
  -      DO 110 J = I, COLS
  -      M(I,J) = M(I,J) / TMP
  -  110 CONTINUE
  -!
  -  910 FORMAT('*** FATAL ERRROR: MGAUSJ: MATRIX ROWS LARGER THAN'       &
  +  SUBROUTINE MGAUSB(M)
  +    IMPLICIT NONE
  +    REAL(KIND=8),INTENT(INOUT) :: M(:,:)
  +    INTEGER(KIND=4)            :: SH(2)
  +    INTEGER I
  +!
  +    SH = SHAPE(M)
  +    IF (SH(1).GT.SH(2)) THEN
  +      WRITE(0,910), SH(1), SH(2)
  +      STOP
  +    END IF
  +!
  +    CALL MGAUS(M)
  +    CALL MBSUB(M)
  +!
  +  910 FORMAT('*** FATAL ERRROR: MGAUSB: MATRIX ROWS LARGER THAN'   &
                ' COLUMNS',/,18X,I8, '   .GT.   ',I8)
     999 RETURN
  -      END SUBROUTINE
  -!------------------------------------------------------------------------------+
  -! MGAUSB: Matrix GAUSs elimination with Back substitution
  +  END SUBROUTINE
  +!------------------------------------------------------------------------------+
  +! MGAUS: Matrix GAUSs elimination
   !------------------------------------------------------------------------------+
   !    M : Matrix
   ! ROWS : Number of rows
  
positive.code:         ELSE
            M(ROW,I) = M(ROW,I) - RATIO * M(COL,I)
          ENDIF
    110   CONTINUE
        ELSE
          M(ROW,COL) = 0D0
        END IF
    999 RETURN
        END SUBROUTINE MGEROW
  !------------------------------------------------------------------------------+
  ! MBSUB: Matrix Back SUBstitution
  !------------------------------------------------------------------------------+
  
negative.code: !
  !------------------------------------------------------------------------------+
  !    M : Matrix
  !  ROW : Matrix row on which to perform elimination
  !  COL : Associated column to be zeroed for elimination.
  ! SCOL : Column to begin elimination
  !------------------------------------------------------------------------------+
    SUBROUTINE MGEROW(M,ROW,COL,SCOL)
      IMPLICIT NONE
      REAL(KIND=8),INTENT(INOUT) :: M(:,:)
      INTEGER(KIND=4),INTENT(IN) :: ROW,COL,SCOL
      DOUBLE PRECISION TOL,RATIO
      PARAMETER (TOL=1.0D-12)
      INTEGER I,J,SH(2)
  !
        SH = SHAPE(M)
        IF (ABS(M(ROW,COL)).GE.TOL) THEN
          RATIO = M(ROW,COL) / M(COL,COL)
          DO 110 I = SCOL, SH(2)
          IF (I.EQ.COL) THEN
  !
  !  By definition this location must be zero.  Instead of doing the
  