C
C Program INV4 to subtract two sets of grid values of a function
C describing the model:  GRID1-GRID2=GRID3.
C
C Date: 1996, September 30
C Coded by Ludek Klimes
C
C ......................................................................
C
C                                                    
C Description of the data files:
C
C Main input data file read from the interactive device (*):
C (1) 'GRID1','GRID2','GRID3',/
C     'GRID1','GRID2','GRID3'... Names of the input and output
C             files described below.
C     /...    Obligatory slash to enable future compatible extensions.
C     Default: 'GRID1'='inv3.dat', 'GRID2'='inv3.out',
C             'GRID3'='dif.out'.
C
C Input files 'GRID1','GRID2':
C (1) TEXTG,IGROUP
C     Identification of the group.
C     TEXTG...A string.  If its first character is 'I' or 'S', the
C             computed function describes a surface.  Otherwise, it
C             describes a material parameter.
C     IGROUP..Index of a surface, or of a complex block.
C (2) TEXTF
C     This input is not performed for a surface, i.e. if the first
C     character of textg (see above) is 'I' or 'S'.
C     TEXTF...String identifying a material parameter.
C (3) K1,K2,K3
C     K1,K2,K3... Indices of coordinates.
C (4) N1,N2,N3
C     N1,N2,N3... Numbers of grid lines.
C (5) X1(1),...,X1(N1)
C     The grid coordinates corresponding to the first independent
C     variable.
C (6) X2(1),...,X2(N2)
C     The grid coordinates corresponding to the second independent
C     variable.
C (7) X3(1),...,X3(N3)
C     The grid coordinates corresponding to the third independent
C     variable.
C (8) (((W(I1,I2,I3),I1=1,N1),I2=1,N2),I3=1,N3)
C     The values of function W at grid points. Function value
C     W(I1,I2,I3) corresponds to point (X1(I1),X2(I2),X3(I3)).
C
C Output file 'GRID3':
C (1) 'DIFFERENCES ............................' (a string)
C (2) (((D(I1,I2,I3),I1=1,N1),I2=1,N2),I3=1,N3)
C     Differences of function values.
C
C-----------------------------------------------------------------------
C
C     Filenames:
      CHARACTER*80 FILE1,FILE2,FILE3
C
C     Logical unit numbers:
      INTEGER LU1,LU2,LU3
      PARAMETER (LU1=11)
      PARAMETER (LU2=12)
      PARAMETER (LU3=13)
C
C     Auxiliary storage locations:
      CHARACTER*3 TEXT
      INTEGER MX
      PARAMETER (MX=100)
      INTEGER ICB,K1,K2,K3,N1,N2,N3,I1,I2,I3
      REAL X(MX,3),V1(MX),V2(MX)
C
C.......................................................................
C
C     Opening data files and reading the input data:
      WRITE(*,'(A)') ' Enter names of input and output files: '
      FILE1='inv3.dat'
      FILE2='inv3.out'
      FILE3='dif.out'
      READ(*,*) FILE1,FILE2,FILE3
      OPEN(LU1,FILE=FILE1,STATUS='OLD')
      OPEN(LU2,FILE=FILE2,STATUS='OLD')
      OPEN(LU3,FILE=FILE3)
      WRITE(*,'(A)') '+                                       '
C
      READ(LU1,*) TEXT,ICB
      READ(LU2,*) TEXT,I2
      IF(ICB.NE.I2) THEN
C       INV4-01
        PAUSE 'Error INV4-01: Different index of a function'
        STOP
      END IF
      IF(TEXT(1:1).NE.'I'.AND.TEXT(1:1).NE.'S') THEN
        READ(LU1,*) TEXT
        READ(LU2,*) TEXT
      END IF
      READ(LU1,*) K1,K2,K3
      READ(LU2,*) I1,I2,I3
      IF(K1.NE.I1.OR.K2.NE.I2.OR.K3.NE.I3) THEN
C       INV4-02
        PAUSE 'Error INV4-02: Different indices of coordinates'
        STOP
      END IF
      READ(LU1,*) N1,N2,N3
      READ(LU2,*) I1,I2,I3
      IF(N1.NE.I1.OR.N2.NE.I2.OR.N3.NE.I3) THEN
C       INV4-03
        PAUSE 'Error INV4-03: Different numbers of grid lines'
        STOP
      END IF
      IF(MAX0(N1,N2,N3).GT.MX) THEN
C       INV4-04
        PAUSE 'Error INV4-04: Too many grid lines'
        STOP
      END IF
      READ(LU1,*) (X(I1,K1),I1=1,N1)
      READ(LU2,*) (X(I1,K1),I1=1,N1)
      READ(LU1,*) (X(I2,K2),I2=1,N2)
      READ(LU2,*) (X(I2,K2),I2=1,N2)
      READ(LU1,*) (X(I3,K3),I3=1,N3)
      READ(LU2,*) (X(I3,K3),I3=1,N3)
      WRITE(LU3,'(A)') '''DIFFERENCES ............................'''
      DO 23 I3=1,N3
        DO 22 I2=1,N2
          READ(LU1,*) (V1(I1),I1=1,N1)
          READ(LU2,*) (V2(I1),I1=1,N1)
          DO 21 I1=1,N1
            V1(I1)=V1(I1)-V2(I1)
   21     CONTINUE
          WRITE(LU3,'(8F10.6)') (V1(I1),I1=1,N1)
   22   CONTINUE
        IF(N1.NE.1.AND.N2.NE.1) WRITE(LU3,*)
   23 CONTINUE
C
      STOP
      END
C
C=======================================================================
C