LATTIC.f

      PROGRAM LATTIC
C
C                *** NIST*LATTICE ***
C
C     A PROGRAM TO ANALYZE LATTICE RELATIONSHIPS
C               VERSION OF SPRING 1991
C
C
C       VICKY LYNN KAREN  AND  ALAN D. MIGHELL
C
C    NATIONAL INSTITUTE OF STANDARDS AND TECHNOLOGY
C     MATERIALS SCIENCE AND ENGINEERING LABORATORY
C             GAITHERSBURG, MARYLAND 20899
C
C
C     *** PATENT PENDING ON CERTAIN ALGORITHMS ***
C     COPYRIGHT BY THE  U.S. SECRETARY OF COMMERCE
C            ON BEHALF OF THE UNITED STATES.
C
C                 ALL RIGHTS RESERVED.
C
C
C
      CHARACTER*5 AFCN,AINV,AREL,ARSS,ASYM,ATRANS
      COMMON /CELLI/ AI,BI,CI,ALPI,BETI,GAMI,VOLI
      COMMON /CONST1/ RADIAN
      COMMON /ERR1/ IERR1
      COMMON /ERR2/ IERR2
      COMMON /FCN1/ IFCN
      COMMON /HCELL1/ YA,YB,YC,YAL,YBE,YGA,YV
      COMMON /HCELL2/ ZA,ZB,ZC,ZAL,ZBE,ZGA,ZV
      COMMON /HOPT1/ IRSYM
      COMMON /MATR1/ U1,V1,W1,U2,V2,W2,U3,V3,W3
      COMMON /PRINT1/ IPRRSS
      COMMON /PROB1/ NTD,IPROB
      COMMON /PROB2/ IPCH,NUNK
      COMMON /UNIT1/ IUNITA
      COMMON /UNIT2/ IUNITB
      COMMON /UNIT4/ IUNITD
C**
      COMMON /CK02/ ICK021,ICK022,ICK023,ICK024,ICK025,ICK026,ICK027,
     $              ICK028,ICK029
      COMMON /CK05/ ICK051
C*
      DATA ARSS/'RSS  '/, ATRANS/'TRANS'/, AINV/'INV  '/
      DATA AREL/'REL  '/, ASYM/'SYM  '/
C
C
C
C     --- OPENS, CLOSES :
C         IUNITD (=10)  (ASSOCIATED FILENAME:  NIST10)
C
C     --- INITIALIZE PROGRAM CONSTANTS
      RADIAN = 180.0/3.14159265
C
C     --- INITIALIZE VARIABLES
      ICTFCN = 0
      IERRIN = 0
      IHEAD2 = 0
      IHEAD3 = 0
      IHEAD4 = 0
      IHEAD5 = 0
      IHEAD6 = 0
      IPRRSS = 0
C**
C     --- FOR CHECKING ... READ OR INITIALIZE TEST PRINT/RUN
C         VARIABLES
C     READ(IUNITA,9500)
C    $    ICK021,ICK022,ICK023,ICK024,ICK025,
C    $    ICK026,ICK027,ICK028,ICK029,
C    $    ICK051
C
      ICK021 = 0
      ICK022 = 0
      ICK023 = 0
      ICK024 = 0
      ICK025 = 0
      ICK026 = 0
      ICK027 = 0
      ICK028 = 0
      ICK029 = 0
C
      ICK051 = 0
C*
C**
C
C     --- WHEN NOT EXECUTING SPECIAL CHECK RUN, GO TO SECTION
C         OF CODE FOR NORMAL EXECUTION OF PROGRAM
      IF(ICK029.NE.1) GO TO 10
C
C        --- THIS CHECK CODE ALLOWS FOR THE REDUCTION OF A LARGE
C            NUMBER OF CELLS PER RUN (UP TO 99999).  THE REDUCED
C            CELL, TRANSFORMATION MATRIX FROM THE INITIAL CELL TO
C            THE REDUCED CELL, AND THE PROBLEM NUMBER WILL BE
C            WRITTEN TO IUNITD (ASSOCIATED FILENAME = NIST10).
C            NIST10 MUST NOT EXIST PRIOR TO PROGRAM EXECUTION AND
C            MUST BE DELETED AFTER PROGRAM EXECUTION.  THE
C            PRINT OUT IS LIMITED TO THE PROBLEM NUMBER FOLLOWED
C            BY ANY WARNING AND/OR ERROR MESSAGES.
C
C        --- INITIALIZE VARIABLES
         IFCN = 2
         IPCH = 0
         IPRRSS = 1
         NUNK = 0
C
C        --- READ NUMBER OF PROBLEMS
         READ(IUNITA,9520) NTD
C
C        --- WRITE SPECIAL HEADING, NUMBER OF PROBLEMS
         CALL HEAD0
         CALL CKPT02(12)
C
C        --- CALCULATE REDUCED CELLS
         OPEN(UNIT=IUNITD,FILE='NIST10',STATUS='NEW')
         CALL REDUCE
         CLOSE(IUNITD)
         GO TO 800
   10 CONTINUE
C*
C
C     --- WRITE GENERAL HEADING FOR PROGRAM
      CALL HEAD0
C
C     --- STARTING POINT FOR EACH PROGRAM FUNCTION TO BE PROCESSED
   50 CONTINUE
      ICTFCN = ICTFCN + 1
C
C     --- READ PROGRAM LINE
C         (TYPE OF PROGRAM FUNCTION, NUMBER OF PROBLEMS)
      READ(IUNITA,6000,END=800) AFCN, NTD
C
C     --- ASSIGN A VARIABLE (IFCN) TO INDICATE THE TYPE OF
C         PROGRAM FUNCTION
      IFCN = 0
      IF(AFCN.EQ.ARSS)   IFCN = 2
      IF(AFCN.EQ.ATRANS) IFCN = 3
      IF(AFCN.EQ.AINV)   IFCN = 4
      IF(AFCN.EQ.AREL)   IFCN = 5
      IF(AFCN.EQ.ASYM)   IFCN = 6
C
C     --- CHECK FOR VALID PROGRAM FUNCTION (IF NOT VALID,
C         WRITE ERROR MESSAGE AND STOP EXECUTION OF PROGRAM)
      IF(IFCN.LT.2.OR.IFCN.GT.6) GO TO 700
C
C     --- GO TO SECTION OF CODE FOR SPECIFIED FUNCTION
      GO TO (100,200,300,400,500,600) IFCN
C
  100 CONTINUE
C
C
  200 CONTINUE
C
C
C     *** RSS FUNCTION ***
C
C     --- NEVER WRITE IUNITD, ALWAYS PRINT RSS INFORMATION
      IPCH   = 0
      IPRRSS = 0
C
C     --- CONTROL WRITE OF RSS HEADING
      IF(ICTFCN.GT.1) WRITE(IUNITB,9000)
      IF(IHEAD2.EQ.0) CALL HEAD2
C
C     --- CHECK NUMBER OF PROBLEMS INPUT ... IF NECESSARY, RESET TO 1
      IF(NTD.GE.1.AND.NTD.LE.20) GO TO 210
         IERRIN = NTD
         WRITE(IUNITB,6200)
         NTD = 1
  210 CONTINUE
C
C     --- WRITE NUMBER OF INDEPENDENT PROBLEMS TO STUDY
      CALL OUTPT1(1)
C
C     --- CALCULATE REDUCED AND DERIVATIVE CELLS
      CALL REDUCE
C
C     --- IF INPUT ERROR, SKIP EXPECTED INPUT RECORDS AND PROCEED TO
C         NEXT PROGRAM FUNCTION
      IF(IERRIN.EQ.0) GO TO 230
         DO 220 I = 1,IERRIN-1
            READ(IUNITA,6300) TEMREC
  220    CONTINUE
         IERRIN = 0
  230 CONTINUE
C
      IHEAD2 = IHEAD2 + 1
      GO TO 50
C
  300 CONTINUE
C
C
C     *** TRANS FUNCTION ***
C
C     --- CONTROL WRITE OF TRANS HEADINGS
      IF(ICTFCN.GT.1) WRITE(IUNITB,9000)
      IF(IHEAD3.EQ.0) CALL HEAD3
      CALL OUTPT1(14)
C
C     --- INITIALIZE VARIABLES
      IPCH  = 0
      IPROB = 0
C
C     --- TRANSFORM EACH INPUT CELL BY SPECIFIED MATRIX
      DO 320 I = 1,NTD
         IPROB = I
C
C        --- WRITE SUB-HEADING
         CALL OUTPT1(15)
C
C        --- READ TRANSFORMATION MATRIX AND CELL
         READ(IUNITA,6400) U1, V1, W1, U2, V2, W2, U3, V3, W3
         READ(IUNITA,6500) AI, BI, CI, ALPI, BETI, GAMI
C
C        --- CALCULATE VOLUME OF THE INPUT CELL
C            (AND CHECK FOR VALID INPUT DATA)
         CALL VOLUME(AI,BI,CI,ALPI,BETI,GAMI,VOLI)
C
C        --- CHECK THAT INPUT CELL IS VALID
C            (INDICATED BY ERROR FLAG SET IN *VOLUME*) ...
C            IF NOT, WRITE ERROR MESSAGE AND GO TO NEXT PROBLEM
         IF(IERR1.EQ.0) GO TO 310
            WRITE(IUNITB,6600)
            GO TO 320
  310    CONTINUE
C
C        --- TRANSFORM CELL AND WRITE
         CALL TRANS(1)
  320 CONTINUE
      IHEAD3 = IHEAD3 + 1
      GO TO 50
C
  400 CONTINUE
C
C
C     *** INV FUNCTION ***
C
C     --- CONTROL WRITE OF INV HEADINGS
      IF(ICTFCN.GT.1) WRITE(IUNITB,9000)
      IF(IHEAD4.EQ.0) CALL HEAD4
      CALL OUTPT1(16)
C
C      --- INITIALIZE VARIABLES
      IPROB = 0
C
C     --- CALCULATE INVERSE FOR EACH 3X3 MATRIX
      DO 420 I = 1,NTD
         IPROB = I
         READ(IUNITA,6400) U1, V1, W1, U2, V2, W2, U3, V3, W3
         CALL OUTPT1(17)
         CALL OUTPT1(11)
C
C        --- CALCULATE DETERMINANT AND CHECK FOR VALID INPUT MATRIX
         CALL DETERM
         IF(IERR2.EQ.0) GO TO 410
C
C           --- INVALID MATRIX (DETERMINANT IS ZERO, >= 100, OR
C               <= -100), WRITE ERROR MESSAGE AND GO TO NEXT PROBLEM
            WRITE(IUNITB,6700)
            GO TO 420
  410    CONTINUE
         CALL INVERS(1)
  420 CONTINUE
      IHEAD4 = IHEAD4 + 1
      GO TO 50
C
  500 CONTINUE
C
C
C     *** REL FUNCTION ***
C
C     --- CONTROL WRITE OF REL HEADINGS
      IF(ICTFCN.GT.1) WRITE(IUNITB,9000)
      IF(IHEAD5.EQ.0) CALL HEAD5
      CALL OUTPT2(1)
C
C     --- READ AND INTERPRET REL LINE
C         (GENERATES OR READS POSSIBLE MATRIX ELEMENTS)
      CALL RD051
C
C     --- WRITE INPUT TOLERANCES AND POSSIBLE MATRIX ELEMENTS
      CALL OUTPT2(2)
C
C      --- INITIALIZE VARIABLE
      IPROB = 0
C
C     --- GENERATE H MATRICES RELATING TWO LATTICES
      DO 520 I = 1,NTD
         IPROB = I
C
C        --- WRITE SUB-HEADING FOR LATTICE GROUP
         CALL OUTPT2(3)
C
C        --- READ INPUT CELLS
         READ(IUNITA,6500) YA,YB,YC,YAL,YBE,YGA
         READ(IUNITA,6500) ZA,ZB,ZC,ZAL,ZBE,ZGA
C
C        --- CALCULATE VOLUME OF THE INPUT CELL
C            (AND CHECK FOR VALID INPUT DATA)
         IERRY = 0
         CALL VOLUME(YA,YB,YC,YAL,YBE,YGA,YV)
         IERRY = IERR1
         CALL VOLUME(ZA,ZB,ZC,ZAL,ZBE,ZGA,ZV)
C
C        --- CHECK THAT INPUT CELL IS VALID
C            (INDICATED BY ERROR FLAG SET IN *VOLUME*) ...
C            IF NOT, WRITE ERROR MESSAGE AND GO TO NEXT PROBLEM
         IF(IERRY.EQ.0.AND.IERR1.EQ.0) GO TO 510
            WRITE(IUNITB,6600)
            GO TO 520
  510    CONTINUE
C
C        --- WRITE INPUT CELLS
         CALL OUTPT2(5)
C
C        --- WRITE COLUMN HEADING FOR OUTPUT
         CALL OUTPT2(6)
C
C        --- GENERATE H MATRICES
         CALL HMATRX
C
C        --- WRITE OUTPUT INDICATING END OF PROBLEM
         CALL OUTPT2(8)
  520 CONTINUE
C
      IHEAD5 = IHEAD5 + 1
      GO TO 50
C
  600 CONTINUE
C
C
C     *** SYM FUNCTION ***
C
      IF(ICTFCN.GT.1.OR.IRSYM.NE.1) WRITE(IUNITB,9000)
C
C     --- READ AND INTERPRET SYM LINE
C         (GENERATES POSSIBLE MATRIX ELEMENTS)
      CALL RD061
C
C     --- CHECK FOR R,R1 OPTION
      IF(IRSYM.EQ.1) GO TO 630
C
C        --- I OPTION HAS BEEN SPECIFIED
C
C        --- CONTROL WRITE OF SYM HEADINGS
         IF(IHEAD6.EQ.0) CALL HEAD6
C
C        --- WRITE HEADING FOR SYM DET THROUGH CT ANALYSIS
         CALL OUTPT2(9)
C
C        --- WRITE INPUT TOLERANCES AND POSSIBLE MATRIX ELEMENTS
         CALL OUTPT2(2)
C
C         --- INITIALIZE VARIABLE
         IPROB = 0
C
C        --- GENERATE H MATRICES RELATING LATTICE TO ITSELF
         DO 620 I = 1,NTD
            IPROB = I
C
C           --- WRITE SUB-HEADING FOR SYMMETRY
            CALL OUTPT2(10)
C
C           --- READ INPUT CELL
            READ(IUNITA,6500) YA,YB,YC,YAL,YBE,YGA
C
C           --- CALCULATE VOLUME OF THE INPUT CELL
C               (AND CHECK FOR VALID INPUT DATA)
            CALL VOLUME(YA,YB,YC,YAL,YBE,YGA,YV)
C
C           --- CHECK THAT INPUT CELL IS VALID
C               (INDICATED BY ERROR FLAG SET IN *VOLUME*) ...
C               IF NOT, WRITE ERROR MESSAGE AND GO TO NEXT PROBLEM
            IF(IERR1.EQ.0) GO TO 610
               WRITE(IUNITB,6600)
               GO TO 620
  610       CONTINUE
C
C           --- WRITE INPUT CELL
            CALL OUTPT2(11)
C
C           --- WRITE COLUMN HEADING FOR OUTPUT
            CALL OUTPT2(6)
C
C           --- GENERATE H MATRICES
            ZA  = YA
            ZB  = YB
            ZC  = YC
            ZAL = YAL
            ZBE = YBE
            ZGA = YGA
            CALL HMATRX
C
C           --- WRITE OUTPUT INDICATING END OF PROBLEM
            CALL OUTPT2(8)
  620    CONTINUE
C
         IHEAD6 = IHEAD6 + 1
         GO TO 50
C
  630 CONTINUE
C
C     --- R,R1 OPTION HAS BEEN SPECIFIED
C
C     --- ALWAYS WRITE IUNITD, CONTROL WRITE OF RSS HEADING
      IPCH = 1
      IF(IPRRSS.EQ.0) CALL HEAD2
C
C     --- CHECK NUMBER OF PROBLEMS INPUT ... IF NECESSARY, RESET TO 1
      IF(NTD.GE.1.AND.NTD.LE.20) GO TO 640
         IERRIN = NTD
         WRITE(IUNITB,6200)
         NTD = 1
  640 CONTINUE
C
C     --- WRITE NUMBER OF INDEPENDENT PROBLEMS TO STUDY
      CALL OUTPT1(1)
C
C     --- WRITE HEADING FOR SYM DET THROUGH CT ANALYSIS
      CALL OUTPT2(9)
C
C     --- INITIALIZE VARIABLE
      NUNK = 0
C
C     --- CALCULATE REDUCED AND DERIVATIVE CELLS
      OPEN(UNIT=IUNITD,FILE='NIST10',STATUS='NEW')
      CALL REDUCE
      CLOSE(IUNITD)
C
C     --- IF NECESSARY, RESET NUMBER OF UNKNOWNS AND
C         WRITE ERROR/WARNING MESSAGE
      IRSET1 = 0
      IF(NUNK.LT.1.OR.NUNK.GT.900) IRSET1 = 1
      IF(NUNK.GT.900) NUNK = 900
      IF(IRSET1.EQ.1) WRITE(IUNITB,6800)
C
C     --- CONTROL WRITE OF SYM HEADING
      WRITE(IUNITB,9000)
      IF(IHEAD6.EQ.0) CALL HEAD6
C
C     --- WRITE HEADING FOR SYM DET THROUGH CT ANALYSIS
      CALL OUTPT2(9)
C
C     --- WRITE INPUT TOLERANCES AND POSSIBLE MATRIX ELEMENTS
      CALL OUTPT2(2)
C
C     --- OPEN FILE TO READ GENERATED INPUT CELLS (IUNITD)
      OPEN(UNIT=IUNITD,FILE='NIST10',STATUS='OLD')
C
C**   FOR CDC COMPUTER ONLY
C     REWIND IUNITD
C
C      --- INITIALIZE VARIABLE
      IPROB = 0
C
C     --- GENERATE H MATRICES RELATING LATTICE TO ITSELF
      DO 650 I = 1,NUNK
         IPROB = I
C
C        --- WRITE SUB-HEADING FOR SYMMETRY
         CALL OUTPT2(10)
C
C        --- READ INPUT CELL
         READ(IUNITD,6900) YA,YB,YC,YAL,YBE,YGA,YV
C
C        --- WRITE INPUT CELL
         CALL OUTPT2(11)
C
C        --- WRITE COLUMN HEADING FOR OUTPUT
         CALL OUTPT2(6)
C
C        --- GENERATE H MATRICES
         ZA  = YA
         ZB  = YB
         ZC  = YC
         ZAL = YAL
         ZBE = YBE
         ZGA = YGA
         CALL HMATRX
C
C        --- WRITE OUTPUT INDICATING END OF PROBLEM
         CALL OUTPT2(8)
  650 CONTINUE
C
C     --- IF INPUT ERROR, SKIP EXPECTED INPUT RECORDS AND PROCEED TO
C         NEXT PROGRAM FUNCTION
      IF(IERRIN.EQ.0) GO TO 670
         DO 660 I = 1,IERRIN-1
            READ(IUNITA,6300) TEMREC
  660    CONTINUE
         IERRIN = 0
  670 CONTINUE
C
      CLOSE(UNIT=IUNITD,STATUS='DELETE')
C
      IPRRSS = 0
      IPCH   = 0
      NUNK   = 0
C
      IHEAD6 = IHEAD6 + 1
      GO TO 50
C
C
  700 CONTINUE
C
C     --- ERROR IN SPECIFIED PROGRAM FUNCTION
      WRITE(IUNITB,8500)
  800 CONTINUE
      STOP
 6000 FORMAT(A5,3X,I2)
 6200 FORMAT(///1X,'*LATTIC* ERROR ... Invalid number of problems input.
     $'/1X,19X,'Number of independent problems set to 1.')
 6300 FORMAT(A80)
 6400 FORMAT(9F8.2)
 6500 FORMAT(10X,6F10.5)
 6600 FORMAT(/1X,'*LATTIC* ERROR ... Input cell has illegal cell paramet
     $er(s) and/or cell volume.'/)
 6700 FORMAT(/1X,'*LATTIC* ERROR ... Invalid matrix, check determinant.'
     $/)
 6800 FORMAT(1H1////1X,'*LATTIC* ERROR ... Check number of lattices gene
     $rated (input to unit 10).'/)
 6900 FORMAT(6F10.5,F10.2)
 8500 FORMAT(1H1////1X,'*LATTIC* ERROR ... Invalid program function was
     $specified.')
 9000 FORMAT(1H1)
 9500 FORMAT(10X,9I1,1X,I1,9X)
 9520 FORMAT(5X,I5)
       END