      FUNCTION DECIMA(ERROR)                                            MOD02917
      LOGICAL LEXP,AFTER,MINUS,ERROR                                    MOD02918
      DOUBLE PRECISION MANT,DECIMA                                      MOD02919
      INTEGER EXP                                                       MOD02920
      COMMON/FRECNM/MULTIP                                              MOD02921
      LEXP    = .FALSE.                                                 MOD02922
      MINUS   = .FALSE.                                                 MOD02923
      AFTER   = .FALSE.                                                 MOD02924
      IOFF    = 1                                                       MOD02925
      MANT    = 0.0                                                     MOD02926
      EXP     = 0                                                       MOD02927
      X       = 0.0                                                     MOD02928
      DO 100  IPOSIT=1,20                                               MOD02929
      N = IDIGIT(IPOSIT)                                                MOD02930
      GO TO (5,5,5,5,5,5,5,5,5,5,100,12,13,14,14,16,17,18),N            MOD02931
5     IF (AFTER) GO TO 6                                                MOD02932
      X = X*10.0 + N - 1.0                                              MOD02933
      GO TO 100                                                         MOD02934
6     Y = N - 1                                                         MOD02935
      X = X + Y/10.0**IOFF                                              MOD02936
      IOFF = IOFF + 1                                                   MOD02937
      GO TO 100                                                         MOD02938
12    MINUS = .TRUE.                                                    MOD02939
      GO TO 100                                                         MOD02940
13    AFTER = .TRUE.                                                    MOD02941
      IOFF = 1                                                          MOD02942
      GO TO 100                                                         MOD02943
14    MANT = X                                                          MOD02944
      IF (MINUS) MANT = -X                                              MOD02945
      MINUS = .FALSE.                                                   MOD02946
      AFTER = .FALSE.                                                   MOD02947
      LEXP  = .TRUE.                                                    MOD02948
      X     = 0.0                                                       MOD02949
      GO TO 100                                                         MOD02950
17    MULTIP=X                                                          MOD02951
      GO TO 100                                                         MOD02952
100   CONTINUE                                                          MOD02953
18    ERROR = .TRUE.                                                    MOD02954
16    IF (LEXP) EXP = X                                                 MOD02955
      IF (LEXP.AND.MINUS) EXP = -X                                      MOD02956
      IF (.NOT.LEXP) MANT = X                                           MOD02957
      IF (.NOT.LEXP.AND.MINUS) MANT = -X                                MOD02958
      DECIMA = MANT*10.0**EXP                                           MOD02959
      RETURN                                                            MOD02960
      END                                                               MOD02961
      FUNCTION ERR1(ERR001)                                             MOD02962
      INTEGER SECT,POINT,BEGIN,GETSEC,EPOINT,BLANK,AGET,AGETW           MOD02963
      LOGICAL EOL,EOS,EOF,ERROR,INDEX,GETWRD,IGET,RGET,ERR1             MOD02964
      COMMON  /FRECNT/ LINE(80),SECT,EOL,EOS,EOF,                       MOD02965
     1ERROR,BLANK,ICOMMA,POINT,BEGIN,LENGTH,                            MOD02966
     2EPOINT,NSECT,LINENM,MAXSTR,KUNIT,JUNIT,IUNIT                      MOD02967
      ERR1=ERROR                                                        MOD02968
      RETURN                                                            MOD02969
      END                                                               MOD02970
      SUBROUTINE GETNL(GET001)                                          MOD02971
      INTEGER SECT,POINT,BEGIN,GETSEC,EPOINT,BLANK,AGET,AGETW           MOD02972
      LOGICAL EOL,EOS,EOF,ERROR,INDEX,GETWRD,IGET,RGET,ERR1             MOD02973
      COMMON/UNIT/II1,II2,II3,II4,II5,ITER                              MOD02974
      COMMON/ECHOT/IECHOT                                               MOD02975
      COMMON/FRECNM/MULTIP                                              MOD02976
      COMMON  /FRECNT/ LINE(80),SECT,EOL,EOS,EOF,                       MOD02977
     1ERROR,BLANK,ICOMMA,POINT,BEGIN,LENGTH,                            MOD02978
     2EPOINT,NSECT,LINENM,MAXSTR,KUNIT,JUNIT,IUNIT                      MOD02979
      DATA    IZER/1H0/,ININE/1H9/,ICOM/1H*/,IDOL/1H$/                  MOD02980
      MULTIP=0                                                          MOD02981
      IF (.NOT.ERROR) GO TO 100                                         MOD02982
9800  WRITE (II2,9801)                                                  MOD02983
9801  FORMAT (' SYNTAX ERROR ON BELOW LINE')                            MOD02984
      WRITE(II2,101)LINE                                                MOD02985
      ERROR = .FALSE.                                                   MOD02986
      EPOINT = 1                                                        MOD02987
      STOP                                                              MOD02988
100   READ (II1,101,ERR=9800,END=910) (LINE(JJ),JJ=1,71)                MOD02989
101   FORMAT (71A1)                                                     MOD02990
      IF(IECHOT.EQ.1)WRITE(II2,102)(LINE(JJ),JJ=1,71)                   MOD02991
102   FORMAT(1X,71A1)                                                   MOD02992
      LINENM = LINENM + 1                                               MOD02993
210   CONTINUE                                                          MOD02994
      POINT = 1                                                         MOD02995
      EOL   = .FALSE.                                                   MOD02996
      EOS   = .FALSE.                                                   MOD02997
      RETURN                                                            MOD02998
910   EOF = .TRUE.                                                      MOD02999
      RETURN                                                            MOD03000
      END                                                               MOD03001
      FUNCTION GETWRD(GET001)                                           MOD03002
      INTEGER SECT,POINT,BEGIN,GETSEC,EPOINT,BLANK,AGET,AGETW           MOD03003
      LOGICAL EOL,EOS,EOF,ERROR,INDEX,GETWRD,IGET,RGET,ERR1             MOD03004
      COMMON  /FRECNT/ LINE(80),SECT,EOL,EOS,EOF,                       MOD03005
     1ERROR,BLANK,ICOMMA,POINT,BEGIN,LENGTH,                            MOD03006
     2EPOINT,NSECT,LINENM,MAXSTR,KUNIT,JUNIT,IUNIT                      MOD03007
      GETWRD = .FALSE.                                                  MOD03008
      LENGTH = 0                                                        MOD03009
      IF (EOL) RETURN                                                   MOD03010
      DO 100 BEGIN = POINT,80                                           MOD03011
      IF (LINE(BEGIN).NE.BLANK) GO TO 110                               MOD03012
100   CONTINUE                                                          MOD03013
      EOL = .TRUE.                                                      MOD03014
      POINT = 80                                                        MOD03015
      RETURN                                                            MOD03016
110   DO 170 POINT = BEGIN,80                                           MOD03017
      IF (LINE(POINT).EQ.BLANK.OR.LINE(POINT).EQ.ICOMMA)                MOD03018
     1GO TO 180                                                         MOD03019
      LENGTH = POINT - BEGIN + 1                                        MOD03020
      MAXSTR = LENGTH                                                   MOD03021
170   CONTINUE                                                          MOD03022
      GETWRD = .TRUE.                                                   MOD03023
      EOL = .TRUE.                                                      MOD03024
      RETURN                                                            MOD03025
180   IP = POINT                                                        MOD03026
      DO 200 POINT = POINT,80                                           MOD03027
      IF (LINE(POINT).EQ.ICOMMA) GO TO 210                              MOD03028
      IF (LINE(POINT).NE.BLANK) GO TO 190                               MOD03029
200   CONTINUE                                                          MOD03030
      GETWRD = .TRUE.                                                   MOD03031
      EOL =.TRUE.                                                       MOD03032
      RETURN                                                            MOD03033
190   POINT = IP                                                        MOD03034
      GETWRD = .TRUE.                                                   MOD03035
      RETURN                                                            MOD03036
210   POINT = POINT + 1                                                 MOD03037
      GETWRD = .TRUE.                                                   MOD03038
      RETURN                                                            MOD03039
      END                                                               MOD03040
      FUNCTION IDIGIT(IPOSIT)                                           MOD03041
      INTEGER AGET                                                      MOD03042
      DIMENSION IVALID(17)                                              MOD03043
      DATA IVALID/1H0,1H1,1H2,1H3,1H4,1H5,1H6,1H7,1H8,1H9,1H+,          MOD03044
     11H-,1H.,1HE,1HD,1H ,1H*/                                          MOD03045
      DO 100 IDIGIT=1,17                                                MOD03046
      IF (AGET(IPOSIT).EQ.IVALID(IDIGIT)) GO TO 110                     MOD03047
100   CONTINUE                                                          MOD03048
      IDIGIT = 18                                                       MOD03049
110   RETURN                                                            MOD03050
      END                                                               MOD03051
      SUBROUTINE JACOB2(NTOT,XYZ)                                       MOD03553
      IMPLICIT REAL*8 (A-H,O-Z)                                         MOD03554
      COMMON/UNIT/INN,IOUT,IP                                           MOD03555
      COMMON/CGELE/N                                                    MOD03556
      DIMENSION XYZ(NTOT,3),NOD(21),XM(21)                              MOD03557
      DO 300 I = 1, NTOT                                                MOD03558
      READ(IP,301) N,(XYZ(I,J),J=1,3)                                   MOD03559
300   CONTINUE                                                          MOD03560
301   FORMAT(I10,/,3E12.5)                                              MOD03561
      READ(IP,302) IDUMMY                                               MOD03562
302   FORMAT(I10)                                                       MOD03563
      READ(IP,401) NEL                                                  MOD03564
401   FORMAT(I10)                                                       MOD03565
      DO 2000 I = 1, NEL                                                MOD03566
      IPL = 0                                                           MOD03567
      READ(IP,501) KIND,NGROUP,DEN,AREA                                 MOD03568
501   FORMAT(I10,I11,2E12.5)                                            MOD03569
      READ (IP,502) N,(NOD(J),J=1,21)                                   MOD03570
502   FORMAT(I10,21I5)                                                  MOD03571
      GO TO (1100,1100,1200,1200,1300,1200,1500,1200,                   MOD03572
     *1400,1300,1200,1200,1200,1500,1500,1200,2000,2000),KIND           MOD03573
1100  CALL DIM1(XYZ,DEN,AREA,XM1,NOD,NTOT)                              MOD03574
      GO TO 2000                                                        MOD03575
1200  CALL DIM2(XYZ,DEN,AREA,XM,NOD,NNOD,KIND,NTOT,N)                   MOD03576
      GO TO 2000                                                        MOD03577
1300  CALL DIM3(XYZ,DEN,XM,NOD,NNOD,NTOT,N)                             MOD03578
      GO TO 2000                                                        MOD03579
1400  CALL CURVE(XYZ,DEN,AREA,XM1,NOD,NTOT)                             MOD03580
      GO TO 2000                                                        MOD03581
1500  CONTINUE                                                          MOD03582
2000  CONTINUE                                                          MOD03583
      RETURN                                                            MOD03584
      END                                                               MOD03585
      SUBROUTINE SURINT (LIMIT,NUM,SEL,COOR,INDEX,GLOB,IER,GAUSS)       MOD06868
      INTEGER SEL                                                       MOD06869
      DIMENSION COOR(LIMIT,3),INDEX(NUM),GLOB(3),GAUSS(2)               MOD06870
      DIMENSION A(2,2),DER(2),LC(2)                                     MOD06871
      EQUIVALENCE (A(1,1),R),(A(2,2),S),(A(1,2),RB),(A(2,1),SB),        MOD06872
     1(DER(1),D1),(DER(2),D2)                                           MOD06873
      DIMENSION PHI(8),PHID(2,8),XIA(3,2),E(3)                          MOD06874
      LOGICAL SKIP                                                      MOD06875
      GAUSS(1)=0.D0                                                     MOD06876
      GAUSS(2)=0.D0                                                     MOD06877
      SKIP=.FALSE.                                                      MOD06878
      IER=0                                                             MOD06879
211   CONTINUE                                                          MOD06880
      R=GAUSS(1)                                                        MOD06881
      S=GAUSS(2)                                                        MOD06882
      IF (NUM.NE.3) GO TO 311                                           MOD06883
      PHI(1)=R                                                          MOD06884
      PHI(2)=S                                                          MOD06885
      PHI(3)=1.D0-R-S                                                   MOD06886
      PHID(1,1)=1.D0                                                    MOD06887
      PHID(2,2)=1.D0                                                    MOD06888
      PHID(1,3)=-1.D0                                                   MOD06889
      PHID(2,3)=-1.D0                                                   MOD06890
      PHID(2,1)=0.D0                                                    MOD06891
      PHID(1,2)=0.D0                                                    MOD06892
      GO TO 322                                                         MOD06893
311   CONTINUE                                                          MOD06894
      RB=1.D0-R                                                         MOD06895
      SB=1.D0-S                                                         MOD06896
      RD=-R-R                                                           MOD06897
      SD=-S-S                                                           MOD06898
      R=1.D0+R                                                          MOD06899
      S=1.D0+S                                                          MOD06900
      RR=R*RB                                                           MOD06901
      SS=S*SB                                                           MOD06902
      DO 11 I=1,4                                                       MOD06903
      GO TO (1,2,3,4),I                                                 MOD06904
1     D=R*S                                                             MOD06905
      D1=S                                                              MOD06906
      D2=R                                                              MOD06907
      GO TO 22                                                          MOD06908
2     D=RB*S                                                            MOD06909
      D1=-S                                                             MOD06910
      D2=RB                                                             MOD06911
      GO TO 22                                                          MOD06912
3     D=RB*SB                                                           MOD06913
      D1=-SB                                                            MOD06914
      D2=-RB                                                            MOD06915
      GO TO 22                                                          MOD06916
4     D=R*SB                                                            MOD06917
      D1=SB                                                             MOD06918
      D2=-R                                                             MOD06919
22    PHI(I)=.25D0*D                                                    MOD06920
      PHID(1,I)=.25D0*D1                                                MOD06921
      PHID(2,I)=.25D0*D2                                                MOD06922
11    CONTINUE                                                          MOD06923
      JB=SEL                                                            MOD06924
      K=4                                                               MOD06925
      DO 44 IND=1,4                                                     MOD06926
      IF (JB.EQ.0) GO TO 100                                            MOD06927
      L=JB/2                                                            MOD06928
      IF (JB.EQ.L+L) GO TO 44                                           MOD06929
      K=K+1                                                             MOD06930
      GO TO (101,102,103,104),IND                                       MOD06931
101   D=RR*S                                                            MOD06932
      D1=RD*S                                                           MOD06933
      D2=RR                                                             MOD06934
      GO TO 55                                                          MOD06935
102   D=RB*SS                                                           MOD06936
      D1=-SS                                                            MOD06937
      D2=RB*SD                                                          MOD06938
      GO TO 55                                                          MOD06939
103   D=RR*SB                                                           MOD06940
      D1=RD*SB                                                          MOD06941
      D2=-RR                                                            MOD06942
      GO TO 55                                                          MOD06943
104   D=R*SS                                                            MOD06944
      D1=SS                                                             MOD06945
      D2=R*SD                                                           MOD06946
55    LC(1)=IND                                                         MOD06947
      LC(2)=IND+1                                                       MOD06948
      IF (LC(2).GT.4) LC(2)=1                                           MOD06949
      PHI(K)=.5D0*D                                                     MOD06950
      D=.25D0*D                                                         MOD06951
      DO 66 IA=1,2                                                      MOD06952
      PHID(IA,K)=DER(IA)*.5D0                                           MOD06953
66    DER(IA)=DER(IA)*.25D0                                             MOD06954
      DO 77 I=1,2                                                       MOD06955
      J=LC(I)                                                           MOD06956
      PHI(J)=PHI(J)-D                                                   MOD06957
      DO 77 IA=1,2                                                      MOD06958
77    PHID(IA,J)=PHID(IA,J)-DER(IA)                                     MOD06959
44    JB=L                                                              MOD06960
100   CONTINUE                                                          MOD06961
322   DO 333 I=1,3                                                      MOD06962
      D=GLOB(I)                                                         MOD06963
      DO 344 K=1,NUM                                                    MOD06964
      JB=INDEX(K)                                                       MOD06965
344   D=D-PHI(K)*COOR(JB,I)                                             MOD06966
333   E(I)=D                                                            MOD06967
      IF (SKIP) GO TO 222                                               MOD06968
      DO 111 IA=1,2                                                     MOD06969
      DO 111 I=1,3                                                      MOD06970
      D=0.D0                                                            MOD06971
      DO 122 K=1,NUM                                                    MOD06972
      JB=INDEX(K)                                                       MOD06973
122   D=D+PHID(IA,K)*COOR(JB,I)                                         MOD06974
111   XIA(I,IA)=D                                                       MOD06975
      DO 133 IA=1,2                                                     MOD06976
      DO 133 JB=1,IA                                                    MOD06977
      D=0.D0                                                            MOD06978
      DO 144 K=1,3                                                      MOD06979
144   D=D+XIA(K,IA)*XIA(K,JB)                                           MOD06980
133   A(IA,JB)=D                                                        MOD06981
      D=A(1,1)*A(2,2)-A(2,1)**2                                         MOD06982
      IF (D.GT.0.) GO TO 155                                            MOD06983
      IER=1                                                             MOD06984
      RETURN                                                            MOD06985
155   CONTINUE                                                          MOD06986
      DO 255 IA=1,2                                                     MOD06987
      RR=0.D0                                                           MOD06988
      DO 266 I=1,3                                                      MOD06989
266   RR=RR+E(I)*XIA(I,IA)                                              MOD06990
255   DER(IA)=RR                                                        MOD06991
      D=1.D0/D                                                          MOD06992
      RR=A(1,1)*D                                                       MOD06993
      A(1,1)=A(2,2)*D                                                   MOD06994
      A(2,2)=RR                                                         MOD06995
      RR=-A(2,1)*D                                                      MOD06996
      A(2,1)=RR                                                         MOD06997
      A(1,2)=RR                                                         MOD06998
      RR=0.D0                                                           MOD06999
      DO 277 IA=1,2                                                     MOD07000
      D=0.D0                                                            MOD07001
      DO 288 JB=1,2                                                     MOD07002
288   D=D+A(IA,JB)*DER(JB)                                              MOD07003
      RR=RR+D**2                                                        MOD07004
277   GAUSS(IA)=GAUSS(IA)+D                                             MOD07005
      SKIP=RR.LE.1.E-6                                                  MOD07006
      GO TO 211                                                         MOD07007
222   DO 166 I=1,3                                                      MOD07008
166   GLOB(I)=GLOB(I)-E(I)                                              MOD07009
      RETURN                                                            MOD07010
      END                                                               MOD07011
      SUBROUTINE PTS1 (P1,P2,R,P3)                                      MOD05219
      COMMON/UNIT/INN,IO                                                MOD05220
      DIMENSION P1(3),P2(3),P3(3)                                       MOD05221
      DX=P2(1)-P1(1)                                                    MOD05222
      DY=P2(2)-P1(2)                                                    MOD05223
      DZ=P2(3)-P1(3)                                                    MOD05224
      DEN=SQRT(DX*DX+DY*DY+DZ*DZ)                                       MOD05225
      IF(DEN.EQ.0.0)WRITE(IO,10)                                        MOD05226
10    FORMAT('   +++ F +++ A POINT FIND OUT ON THE CENTER OF            MOD05227
     1THE SPHERE .',/,10X,'CHANGE THE CENTER OF THE SPHERE')            MOD05228
      COSA=DX/DEN                                                       MOD05229
      COSB=DY/DEN                                                       MOD05230
      COSC=DZ/DEN                                                       MOD05231
      P3(1)=R*COSA+P1(1)                                                MOD05232
      P3(2)=R*COSB+P1(2)                                                MOD05233
      P3(3)=R*COSC+P1(3)                                                MOD05234
      RETURN                                                            MOD05235
      END                                                               MOD05236
      SUBROUTINE PTS2 (P1,P2,P3,A,B,C,D)                                MOD05237
      DIMENSION P1(3),P2(3),P3(3)                                       MOD05238
      A=P2(1)-P1(1)                                                     MOD05239
      B=P2(2)-P1(2)                                                     MOD05240
      C=P2(3)-P1(3)                                                     MOD05241
      D=A*P3(1)+B*P3(2)+C*P3(3)                                         MOD05242
      RETURN                                                            MOD05243
      END                                                               MOD05244
      SUBROUTINE PTS3 (P1,P2,A,B,C,D,P3)                                MOD05245
      DIMENSION P1(3),P2(3),P3(3)                                       MOD05246
      REAL NX,NY,NZ                                                     MOD05247
      NX=P2(1)-P1(1)                                                    MOD05248
      NY=P2(2)-P1(2)                                                    MOD05249
      NZ=P2(3)-P1(3)                                                    MOD05250
      IF (NX.EQ.0.)  GO TO 10                                           MOD05251
      P3(1)= (D+B*NY/NX*P1(1)-B*P1(2)+C*NZ/NX*P1(1)-C*P1(3))            MOD05252
     $/ (A+B*NY/NX+C*NZ/NX)                                             MOD05253
      P3(2)= NY*(P3(1)-P1(1))/NX+P1(2)                                  MOD05254
      P3(3)= NZ*(P3(1)-P1(1))/NX+P1(3)                                  MOD05255
      GO TO 30                                                          MOD05256
10    CONTINUE                                                          MOD05257
      IF (NZ.EQ.0.)  GO TO 20                                           MOD05258
      P3(3)= (D-A*P1(1)+B*NY/NZ*P1(3)-B*P1(2))/(B*NY/NZ+C)              MOD05259
      P3(1)= P1(1)                                                      MOD05260
      P3(2)= NY/NZ*(P3(3)-P1(3))+P1(2)                                  MOD05261
      GO TO 30                                                          MOD05262
20    CONTINUE                                                          MOD05263
      P3(1)= P1(1)                                                      MOD05264
      P3(3)= P1(3)                                                      MOD05265
      P3(2)= (D-A*P1(1)-C*P1(3))/B                                      MOD05266
30    CONTINUE                                                          MOD05267
      RETURN                                                            MOD05268
      END                                                               MOD05269
      SUBROUTINE VOLM (N3D,INCR,INCS,INCT,NUMDR,NUMDS,NUMDT             MOD07199
     1,KS,LIMIT,NSTOR,X,Y,Z)                                            MOD07200
      LOGICAL IGET,RGET                                                 MOD07201
      DIMENSION N3D(20),NIB(6),NSTOR(LIMIT),X(LIMIT),Y(LIMIT),Z(LIMIT)  MOD07202
      DIMENSION P1(3),P2(3),P3(3),P4(3),P5(3),P6(3),P7(3),P8(3)         MOD07203
      DATA NIB/1,4,2,1,5,1/                                             MOD07204
      IF (KS.LT.0 .OR. KS.GT.6)  GO TO 600                              MOD07205
      DO 500 N=1,KS                                                     MOD07206
5     CALL GETNL(GET001)                                                MOD07207
      IF(IGET(NS))GO TO 5                                               MOD07208
      IF(IGET(NTYP))GO TO 5                                             MOD07209
      IF (NS.LT.1 .OR. NS.GT.6)   GO TO 500                             MOD07210
      IF (NTYP.LT.1 .OR. NTYP.GT.3)   GO TO 500                         MOD07211
      M=NIB(NS)                                                         MOD07212
      NODE=N3D(M)                                                       MOD07213
      GO TO (10,10,20,20,30,30), NS                                     MOD07214
10    NUMD1=NUMDT                                                       MOD07215
      NUMD2=NUMDS                                                       MOD07216
      INC1=INCS                                                         MOD07217
      INC2=INCT                                                         MOD07218
      GO TO 70                                                          MOD07219
20    NUMD1=NUMDT                                                       MOD07220
      NUMD2=NUMDR                                                       MOD07221
      INC1=INCR                                                         MOD07222
      INC2=INCT                                                         MOD07223
      GO TO 70                                                          MOD07224
30    NUMD1=NUMDS                                                       MOD07225
      NUMD2=NUMDR                                                       MOD07226
      INC1=INCR                                                         MOD07227
      INC2=INCS                                                         MOD07228
70    CONTINUE                                                          MOD07229
      GO TO (100,200,300), NTYP                                         MOD07230
100   CONTINUE                                                          MOD07231
      DO 102 K=1,3                                                      MOD07232
      IF(RGET(P1(K)))GO TO 5                                            MOD07233
102   CONTINUE                                                          MOD07234
      IF(RGET(RAD))GO TO 5                                              MOD07235
      DO 150 I=1,NUMD1                                                  MOD07236
      DO 140 J=1,NUMD2                                                  MOD07237
      IF(NODE.GT.LIMIT)CALL PRTERR(7)                                   MOD07238
      P2(1)=X(NODE)                                                     MOD07239
      P2(2)=Y(NODE)                                                     MOD07240
      P2(3)=Z(NODE)                                                     MOD07241
      CALL PTS1 (P1,P2,RAD,P3)                                          MOD07242
      NSTOR(NODE)=NODE                                                  MOD07243
      X(NODE)=P3(1)                                                     MOD07244
      Y(NODE)=P3(2)                                                     MOD07245
      Z(NODE)=P3(3)                                                     MOD07246
      NODE=NODE+INC1                                                    MOD07247
140   CONTINUE                                                          MOD07248
      NODE=NODE-NUMD2*INC1+INC2                                         MOD07249
150   CONTINUE                                                          MOD07250
      GO TO 500                                                         MOD07251
200   CONTINUE                                                          MOD07252
      DO 202 K=1,3                                                      MOD07253
      IF(RGET(P1(K)))GO TO 5                                            MOD07254
202   CONTINUE                                                          MOD07255
      DO 203 K=1,3                                                      MOD07256
      IF(RGET(P2(K)))GO TO 5                                            MOD07257
203   CONTINUE                                                          MOD07258
      IF(RGET(RAD))GO TO 5                                              MOD07259
      DO 250 I=1,NUMD1                                                  MOD07260
      DO 240 J=1,NUMD2                                                  MOD07261
      IF(NODE.GT.LIMIT)CALL PRTERR(7)                                   MOD07262
      P3(1)=X(NODE)                                                     MOD07263
      P3(2)=Y(NODE)                                                     MOD07264
      P3(3)=Z(NODE)                                                     MOD07265
      CALL PTS2 (P1,P2,P3,A,B,C,D)                                      MOD07266
      CALL PTS3 (P1,P2,A,B,C,D,P4)                                      MOD07267
      CALL PTS1 (P4,P3,RAD,P5)                                          MOD07268
      NSTOR(NODE)=NODE                                                  MOD07269
      X(NODE)=P5(1)                                                     MOD07270
      Y(NODE)=P5(2)                                                     MOD07271
      Z(NODE)=P5(3)                                                     MOD07272
      NODE=NODE+INC1                                                    MOD07273
240   CONTINUE                                                          MOD07274
      NODE=NODE-NUMD2*INC1+INC2                                         MOD07275
250   CONTINUE                                                          MOD07276
      GO TO 500                                                         MOD07277
300   CONTINUE                                                          MOD07278
      DO 302 K=1,3                                                      MOD07279
      IF(RGET(P1(K)))GO TO 5                                            MOD07280
302   CONTINUE                                                          MOD07281
      DO 303 K=1,3                                                      MOD07282
      IF(RGET(P2(K)))GO TO 5                                            MOD07283
303   CONTINUE                                                          MOD07284
      DO 304 K=1,3                                                      MOD07285
      IF(RGET(P3(K)))GO TO 5                                            MOD07286
304   CONTINUE                                                          MOD07287
      DO 305 K=1,3                                                      MOD07288
      IF(RGET(P4(K)))GO TO 5                                            MOD07289
305   CONTINUE                                                          MOD07290
      DO 350 I=1,NUMD1                                                  MOD07291
      DO 340 J=1,NUMD2                                                  MOD07292
      IF(NODE.GT.LIMIT)CALL PRTERR(7)                                   MOD07293
      P5(1)=X(NODE)                                                     MOD07294
      P5(2)=Y(NODE)                                                     MOD07295
      P5(3)=Z(NODE)                                                     MOD07296
      CALL PTS2 (P1,P2,P5,A,B,C,D)                                      MOD07297
      CALL PTS3 (P1,P2,A,B,C,D,P6)                                      MOD07298
      CALL PTS3 (P3,P4,A,B,C,D,P7)                                      MOD07299
      RAD=SQRT((P7(1)-P6(1))**2+(P7(2)-P6(2))**2+(P7(3)-P6(3))**2)      MOD07300
      CALL PTS1 (P6,P5,RAD,P8)                                          MOD07301
      NSTOR(NODE)=NODE                                                  MOD07302
      X(NODE)=P8(1)                                                     MOD07303
      Y(NODE)=P8(2)                                                     MOD07304
      Z(NODE)=P8(3)                                                     MOD07305
      NODE=NODE+INC1                                                    MOD07306
340   CONTINUE                                                          MOD07307
      NODE=NODE-NUMD2*INC1+INC2                                         MOD07308
350   CONTINUE                                                          MOD07309
500   CONTINUE                                                          MOD07310
600   CONTINUE                                                          MOD07311
      RETURN                                                            MOD07312
      END                                                               MOD07313
      SUBROUTINE SRFC (NOD,INCS,INCT,NUMDS,NUMDT,KS,LIMIT,NSTOR,X,Y,Z)  MOD06709
      LOGICAL IGET,RGET                                                 MOD06710
      DIMENSION NOD(8),NIB(4),NSTOR(LIMIT),X(LIMIT),Y(LIMIT),Z(LIMIT)   MOD06711
      DIMENSION P1(3),P2(3),P3(3),P4(3),P5(3),P6(3),P7(3),P8(3)         MOD06712
     1,ID(2,2),CORD(20,3),NN(4)                                         MOD06713
      DATA NIB/1,2,4,1/                                                 MOD06714
      IF (KS.GT.3 .OR. KS.LT.-4)   GO TO 500                            MOD06715
      IF (KS.LT.0)   GO TO 300                                          MOD06716
      GO TO (22,100,200), KS                                            MOD06717
22    CONTINUE                                                          MOD06718
10    CALL GETNL(GET001)                                                MOD06719
      DO 25 K=1,3                                                       MOD06720
      IF(RGET(P1(K))) GO TO 10                                          MOD06721
25    CONTINUE                                                          MOD06722
      IF(RGET(RAD)) GO TO 10                                            MOD06723
      NODE=NOD(1)                                                       MOD06724
      IF(NODE.GT.LIMIT)CALL PRTERR(7)                                   MOD06725
      DO 50 I=1,NUMDT                                                   MOD06726
      DO 40 J=1,NUMDS                                                   MOD06727
      P2(1)=X(NODE)                                                     MOD06728
      P2(2)=Y(NODE)                                                     MOD06729
      P2(3)=Z(NODE)                                                     MOD06730
      CALL PTS1 (P1,P2,RAD,P3)                                          MOD06731
      NSTOR(NODE)=NODE                                                  MOD06732
      X(NODE)=P3(1)                                                     MOD06733
      Y(NODE)=P3(2)                                                     MOD06734
      Z(NODE)=P3(3)                                                     MOD06735
      NODE=NODE+INCS                                                    MOD06736
      IF(NODE.GT.LIMIT)CALL PRTERR(7)                                   MOD06737
40    CONTINUE                                                          MOD06738
      NODE=NODE-NUMDS*INCS+INCT                                         MOD06739
      IF(NODE.GT.LIMIT)CALL PRTERR(7)                                   MOD06740
50    CONTINUE                                                          MOD06741
      GO TO 500                                                         MOD06742
100   CONTINUE                                                          MOD06743
110   CALL GETNL(GET001)                                                MOD06744
      DO 125 K=1,3                                                      MOD06745
      IF(RGET(P1(K)))GO TO 110                                          MOD06746
125   CONTINUE                                                          MOD06747
      DO 126 K=1,3                                                      MOD06748
      IF(RGET(P2(K)))GO TO 110                                          MOD06749
126   CONTINUE                                                          MOD06750
      IF(RGET(RAD))GO TO 100                                            MOD06751
      NODE=NOD(1)                                                       MOD06752
      DO 150 I=1,NUMDT                                                  MOD06753
      DO 140 J=1,NUMDS                                                  MOD06754
      IF(NODE.GT.LIMIT)CALL PRTERR(7)                                   MOD06755
      P3(1)=X(NODE)                                                     MOD06756
      P3(2)=Y(NODE)                                                     MOD06757
      P3(3)=Z(NODE)                                                     MOD06758
      CALL PTS2 (P1,P2,P3,A,B,C,D)                                      MOD06759
      CALL PTS3 (P1,P2,A,B,C,D,P4)                                      MOD06760
      CALL PTS1 (P4,P3,RAD,P5)                                          MOD06761
      NSTOR(NODE)=NODE                                                  MOD06762
      X(NODE)=P5(1)                                                     MOD06763
      Y(NODE)=P5(2)                                                     MOD06764
      Z(NODE)=P5(3)                                                     MOD06765
      NODE=NODE+INCS                                                    MOD06766
140   CONTINUE                                                          MOD06767
      NODE=NODE-NUMDS*INCS+INCT                                         MOD06768
150   CONTINUE                                                          MOD06769
      GO TO 500                                                         MOD06770
200   CONTINUE                                                          MOD06771
201   CALL GETNL(GET001)                                                MOD06772
      DO 202 K=1,3                                                      MOD06773
      IF(RGET(P1(K)))GO TO 201                                          MOD06774
202   CONTINUE                                                          MOD06775
      DO 203 K=1,3                                                      MOD06776
      IF(RGET(P2(K)))GO TO 201                                          MOD06777
203   CONTINUE                                                          MOD06778
      DO 204 K=1,3                                                      MOD06779
      IF(RGET(P3(K)))GO TO 201                                          MOD06780
204   CONTINUE                                                          MOD06781
      DO 205 K=1,3                                                      MOD06782
      IF(RGET(P4(K)))GO TO 201                                          MOD06783
205   CONTINUE                                                          MOD06784
      NODE=NOD(1)                                                       MOD06785
      DO 250 I=1,NUMDT                                                  MOD06786
      DO 240 J=1,NUMDS                                                  MOD06787
      IF(NODE.GT.LIMIT)CALL PRTERR(7)                                   MOD06788
      P5(1)=X(NODE)                                                     MOD06789
      P5(2)=Y(NODE)                                                     MOD06790
      P5(3)=Z(NODE)                                                     MOD06791
      CALL PTS2 (P1,P2,P5,A,B,C,D)                                      MOD06792
      CALL PTS3 (P1,P2,A,B,C,D,P6)                                      MOD06793
      CALL PTS3 (P3,P4,A,B,C,D,P7)                                      MOD06794
      RAD=SQRT((P7(1)-P6(1))**2+(P7(2)-P6(2))**2+(P7(3)-P6(3))**2)      MOD06795
      CALL PTS1 (P6,P5,RAD,P8)                                          MOD06796
      NSTOR(NODE)=NODE                                                  MOD06797
      X(NODE)=P8(1)                                                     MOD06798
      Y(NODE)=P8(2)                                                     MOD06799
      Z(NODE)=P8(3)                                                     MOD06800
      NODE=NODE+INCS                                                    MOD06801
240   CONTINUE                                                          MOD06802
      NODE=NODE-NUMDS*INCS+INCT                                         MOD06803
250   CONTINUE                                                          MOD06804
      GO TO 500                                                         MOD06805
300   CONTINUE                                                          MOD06806
      KS=IABS(KS)                                                       MOD06807
      DO 400 N=1,KS                                                     MOD06808
301   CALL GETNL(GET01)                                                 MOD06809
      IF(IGET(NS))GO TO 301                                             MOD06810
      DO 320 K=1,3                                                      MOD06811
      IF(RGET(P1(K)))GO TO 301                                          MOD06812
320   CONTINUE                                                          MOD06813
      DO 321 K=1,3                                                      MOD06814
      IF(RGET(P2(K)))GO TO 301                                          MOD06815
321   CONTINUE                                                          MOD06816
      IF(RGET(RAD))GO TO 301                                            MOD06817
      IF (NS.LT.1 .OR. NS.GT.4)   GO TO 400                             MOD06818
      M=NIB(NS)                                                         MOD06819
      NODE=NOD(M)                                                       MOD06820
      INC=INCS                                                          MOD06821
      IF (NS.EQ.2 .OR. NS.EQ.4)   INC=INCT                              MOD06822
      NUMD=NUMDS                                                        MOD06823
      IF (NS.EQ.2 .OR. NS.EQ.4)   NUMD=NUMDT                            MOD06824
      DO 380 M=1,NUMD                                                   MOD06825
      IF(NODE.GT.LIMIT)CALL PRTERR(7)                                   MOD06826
      P3(1)=X(NODE)                                                     MOD06827
      P3(2)=Y(NODE)                                                     MOD06828
      P3(3)=Z(NODE)                                                     MOD06829
      CALL PTS2 (P1,P2,P3,A,B,C,D)                                      MOD06830
      CALL PTS3 (P1,P2,A,B,C,D,P4)                                      MOD06831
      CALL PTS1 (P4,P3,RAD,P5)                                          MOD06832
      NSTOR(NODE)=NODE                                                  MOD06833
      X(NODE)=P5(1)                                                     MOD06834
      Y(NODE)=P5(2)                                                     MOD06835
      Z(NODE)=P5(3)                                                     MOD06836
      NODE=NODE+INC                                                     MOD06837
380   CONTINUE                                                          MOD06838
400   CONTINUE                                                          MOD06839
500   RETURN                                                            MOD06840
      END                                                               MOD06841
      SUBROUTINE MESHGE(XN,INCS,INCT,INCR,PERS,PERT,PERR,ID,KS1         MOD04275
     1,LIMIT,NSTOR,X,Y,Z)                                               MOD04276
      DIMENSION XN(20,4),CORD(20,3),H(20),CZ(3),NOD(8),N3D(20)          MOD04277
      DIMENSION XC(3),XI(3),XX(3)                                       MOD04278
      DIMENSION NSTOR(LIMIT),X(LIMIT),Y(LIMIT),Z(LIMIT)                 MOD04279
210   KO=1                                                              MOD04280
      IF(ID.NE.2) GO TO 600                                             MOD04281
410   CONTINUE                                                          MOD04282
      DO 407 I=1,8                                                      MOD04283
      NOD(I)=XN(I,1)                                                    MOD04284
407   CONTINUE                                                          MOD04285
      IF(PERS.LE.0.0) PERS=100.                                         MOD04286
      IF(PERT.LE.0.0) PERT=100.                                         MOD04287
      IF(INCS.EQ.0) INCS=1                                              MOD04288
      IF(INCT.EQ.0) INCT=NOD(2)-NOD(1)+1                                MOD04289
      DO 480 I=1,8                                                      MOD04290
      J=NOD(I)                                                          MOD04291
450   DO 460 K=1,3                                                      MOD04292
460   CORD(I,K)=0.0                                                     MOD04293
      IF(J.EQ.0) GO TO 480                                              MOD04294
      DO 470 K=1,3                                                      MOD04295
      KK=K+1                                                            MOD04296
470   CORD(I,K)=XN(I,KK)                                                MOD04297
480   CONTINUE                                                          MOD04298
      NI=4                                                              MOD04299
      DO 490 I=5,8                                                      MOD04300
      IF(NOD(I).GT.0) NI=8                                              MOD04301
      H(I)=0.0                                                          MOD04302
490   CONTINUE                                                          MOD04303
      DO 500 I=1,4                                                      MOD04304
      IF(NOD(I).EQ.0) CALL PRTERR(2)                                    MOD04305
      IF(NOD(I).EQ.0) RETURN                                            MOD04306
500   CONTINUE                                                          MOD04307
      NODE=NOD(1)                                                       MOD04308
      NUMDS=(NOD(2)-NOD(1))/INCS+1                                      MOD04309
      NUMDT=(NOD(3)-NOD(2))/INCT+1                                      MOD04310
      NX=NOD(1)-INCS+(NUMDS*NUMDT*INCS)+(INCT-NUMDS*INCS)*(NUMDT-1)     MOD04311
      IF(NX.EQ.NOD(3)) GO TO 515                                        MOD04312
      CALL PRTERR(3)                                                    MOD04313
      RETURN                                                            MOD04314
515   CONTINUE                                                          MOD04315
      T=-1.                                                             MOD04316
      PERS=PERS/100.0                                                   MOD04317
      PERT=PERT/100.0                                                   MOD04318
      LLL=NUMDS-2                                                       MOD04319
      SUMS=1.0                                                          MOD04320
      DO 520 I=1,LLL                                                    MOD04321
520   SUMS=SUMS+PERS**I                                                 MOD04322
      LLL=NUMDT-2                                                       MOD04323
      SUMT=1.0                                                          MOD04324
      DO 530 I=1,LLL                                                    MOD04325
530   SUMT=SUMT+PERT**I                                                 MOD04326
      FACTT=1.0/PERT                                                    MOD04327
      IF(PERS.EQ.1.0)SUMS=NUMDS-1                                       MOD04328
      IF(PERT.EQ.1.0)SUMT=NUMDT-1                                       MOD04329
      DS=2.0/SUMS                                                       MOD04330
      DT=2.0/SUMT                                                       MOD04331
      DO 590 I=1,NUMDT                                                  MOD04332
      S=-1.                                                             MOD04333
      FACTS=1.0/PERS                                                    MOD04334
      DO 580 J=1,NUMDS                                                  MOD04335
      IF(NI.NE.8) GO TO 540                                             MOD04336
      IF(NOD(5).GT.0)                                                   MOD04337
     $H(5)=(1.-S**2)*(1.-T)*0.5                                         MOD04338
      IF(NOD(6).GT.0)                                                   MOD04339
     $H(6)=(1.-T**2)*(1.+S)*0.5                                         MOD04340
      IF(NOD(7).GT.0)                                                   MOD04341
     $H(7)=(1.-S**2)*(1.+T)*0.5                                         MOD04342
      IF(NOD(8).GT.0)                                                   MOD04343
     $H(8)=(1.-T**2)*(1.-S)*0.5                                         MOD04344
540   H(1)=(1.-S)*(1.-T)*   0.25-(H(5)+H(8)) *0.5                       MOD04345
      H(2)=(1.+S)*(1.-T)*   0.25-(H(5)+H(6)) *0.5                       MOD04346
      H(3)=(1.+S)*(1.+T)*   0.25-(H(6)+H(7)) *0.5                       MOD04347
      H(4)=(1.-S)*(1.+T)*   0.25-(H(7)+H(8)) *0.5                       MOD04348
      DO 550 JJ=1,3                                                     MOD04349
550   CZ(JJ)=0.0                                                        MOD04350
      DO 560 II=1,NI                                                    MOD04351
      DO 560 JJ=1,3                                                     MOD04352
560   CZ(JJ)=CZ(JJ)+H(II)*CORD(II,JJ)                                   MOD04353
      IF(NODE.GT.LIMIT)CALL PRTERR(7)                                   MOD04354
      NSTOR(NODE)=NODE                                                  MOD04355
      X(NODE)=CZ(1)                                                     MOD04356
      Y(NODE)=CZ(2)                                                     MOD04357
      Z(NODE)=CZ(3)                                                     MOD04358
      NODE=NODE+INCS                                                    MOD04359
      FACTS=FACTS*PERS                                                  MOD04360
580   S=S+DS*FACTS                                                      MOD04361
      NODE=NODE-(NUMDS)*INCS+INCT                                       MOD04362
      FACTT=FACTT*PERT                                                  MOD04363
590   T=T+DT*FACTT                                                      MOD04364
      IF(KS1.NE.0)CALL SRFC(NOD,INCS,INCT,NUMDS,NUMDT,KS1               MOD04365
     1,LIMIT,NSTOR,X,Y,Z)                                               MOD04366
      RETURN                                                            MOD04367
600   IF(ID.NE.3)GO TO 1000                                             MOD04368
      DO 601 M=1,20                                                     MOD04369
      N3D(M)=XN(M,1)                                                    MOD04370
601   CONTINUE                                                          MOD04371
      IF(INCR.LE.0) INCR=1                                              MOD04372
      NUMDR=(N3D(4)-N3D(1))/INCR+1                                      MOD04373
      IF(INCS.LE.0)                                                     MOD04374
     $INCS=N3D(4)-N3D(1)+1                                              MOD04375
      NUMDS=(N3D(3)-N3D(4))/INCS+1                                      MOD04376
      IF(INCT.LE.0)                                                     MOD04377
     $INCT=N3D(3)-N3D(1)+1                                              MOD04378
      NUMDT=(N3D(7)-N3D(3))/INCT+1                                      MOD04379
      IF(PERR.LE.0.0)PERR=100.0                                         MOD04380
      IF(PERS.LE.0.0)PERS=100.0                                         MOD04381
      IF(PERT.LE.0.0)PERT=100.0                                         MOD04382
      DO 680 I=1,20                                                     MOD04383
      J=N3D(I)                                                          MOD04384
      IF(J.GE.0)GO TO 650                                               MOD04385
      CALL PRTERR(4)                                                    MOD04386
      RETURN                                                            MOD04387
650   DO 660 K=1,3                                                      MOD04388
660   CORD(I,K)=0.0                                                     MOD04389
      IF(J.EQ.0) GO TO 680                                              MOD04390
      DO 670 K=1,3                                                      MOD04391
      KK=K+1                                                            MOD04392
670   CORD(I,K)=XN(I,KK)                                                MOD04393
680   CONTINUE                                                          MOD04394
      NI=8                                                              MOD04395
      DO 690 I=9,20                                                     MOD04396
      IF(N3D(I).GT.0) NI=20                                             MOD04397
690   H(I)=0.0                                                          MOD04398
      DO 700 I=1,8                                                      MOD04399
      IF(N3D(I).EQ.0) CALL PRTERR(2)                                    MOD04400
      IF(N3D(I).EQ.0) RETURN                                            MOD04401
700   CONTINUE                                                          MOD04402
      NNN=NUMDR*INCR*NUMDS+(INCS-NUMDR*INCR)*(NUMDS-1)                  MOD04403
      NNN=N3D(1)+NNN*NUMDT+(INCT-NNN)*(NUMDT-1)-INCR                    MOD04404
      IF(NNN.NE.N3D(7)) CALL PRTERR(5)                                  MOD04405
      IF(NNN.NE.N3D(7))RETURN                                           MOD04406
      PERR=PERR/100.0                                                   MOD04407
      PERS=PERS/100.0                                                   MOD04408
      PERT=PERT/100.0                                                   MOD04409
      T=-1.                                                             MOD04410
      LLL=NUMDR-2                                                       MOD04411
      SUMR=1.0                                                          MOD04412
      DO 730 I=1,LLL                                                    MOD04413
730   SUMR=SUMR+PERR**I                                                 MOD04414
      LLL=NUMDS-2                                                       MOD04415
      SUMS=1.0                                                          MOD04416
      DO 740 I=1,LLL                                                    MOD04417
740   SUMS=SUMS+PERS**I                                                 MOD04418
      LLL=NUMDT-2                                                       MOD04419
      SUMT=1.0                                                          MOD04420
      DO 750 I=1,LLL                                                    MOD04421
750   SUMT=SUMT+PERT**I                                                 MOD04422
      IF(PERR.EQ.1.0) SUMR=NUMDR-1                                      MOD04423
      IF(PERS.EQ.1.0) SUMS=NUMDS-1                                      MOD04424
      IF(PERT.EQ.1.0) SUMT=NUMDT-1                                      MOD04425
      FACTT=1.0/PERT                                                    MOD04426
      NODE=N3D(1)                                                       MOD04427
      DR=2.0/SUMR                                                       MOD04428
      DS=2.0/SUMS                                                       MOD04429
      DT=2.0/SUMT                                                       MOD04430
      DO 820 MMM=1,NUMDT                                                MOD04431
      NODEI=NODE                                                        MOD04432
      S=-1.                                                             MOD04433
      FACTS=1.0/PERS                                                    MOD04434
      DO 810 I=1,NUMDS                                                  MOD04435
      R=-1.                                                             MOD04436
      FACTR=1.0/PERR                                                    MOD04437
      DO 800 J=1,NUMDR                                                  MOD04438
      RP=1.+R                                                           MOD04439
      SP=1.+S                                                           MOD04440
      TP=1.+T                                                           MOD04441
      RM=1.-R                                                           MOD04442
      SM=1.-S                                                           MOD04443
      TM=1.-T                                                           MOD04444
      IF(NI.NE.20) GO TO 760                                            MOD04445
      RR=1.-R*R                                                         MOD04446
      SS=1.-S*S                                                         MOD04447
      TT=1.-T*T                                                         MOD04448
      IF(N3D(9).GT.0)                                                   MOD04449
     $H(9)=RM*SS*TM*0.25                                                MOD04450
      IF(N3D(10).GT.0)                                                  MOD04451
     $H(10)=RR*SP*TM*0.25                                               MOD04452
      IF(N3D(11).GT.0)                                                  MOD04453
     $H(11)=RP*SS*TM*0.25                                               MOD04454
      IF(N3D(12).GT.0)                                                  MOD04455
     $H(12)=RR*SM*TM *0.25                                              MOD04456
      IF(N3D(13).GT.0)                                                  MOD04457
     $H(13)=RM*SS*TP*0.25                                               MOD04458
      IF(N3D(14).GT.0)                                                  MOD04459
     $H(14)=RR*SP*TP*0.25                                               MOD04460
      IF(N3D(15).GT.0)                                                  MOD04461
     $H(15)=RP*SS*TP*0.25                                               MOD04462
      IF(N3D(16).GT.0)                                                  MOD04463
     $H(16)=RR*SM*TP*0.25                                               MOD04464
      IF(N3D(17).GT.0)                                                  MOD04465
     $H(17)=RM*SM*TT*0.25                                               MOD04466
      IF(N3D(18).GT.0)                                                  MOD04467
     $H(18)=RM*SP*TT *0.25                                              MOD04468
      IF(N3D(19).GT.0)                                                  MOD04469
     $H(19)=RP*SP*TT*0.25                                               MOD04470
      IF(N3D(20).GT.0)                                                  MOD04471
     $H(20)=RP*SM*TT*0.25                                               MOD04472
760   TM=0.125*TM                                                       MOD04473
      TP=0.125*TP                                                       MOD04474
      H(1)=RM*SM*TM     -0.5*(H( 9)+H(17)+H(12))                        MOD04475
      H(2)=RM*SP*TM     -0.5*(H( 9)+H(18)+H(10))                        MOD04476
      H(3)=RP*SP*TM     -0.5*(H(10)+H(19)+H(11))                        MOD04477
      H(4)=RP*SM*TM     -0.5*(H(11)+H(20)+H(12))                        MOD04478
      H(5)=RM*SM*TP     -0.5*(H(13)+H(17)+H(16))                        MOD04479
      H(6)=RM*SP*TP     -0.5*(H(13)+H(18)+H(14))                        MOD04480
      H(7)=RP*SP*TP     -0.5*(H(14)+H(19)+H(15))                        MOD04481
      H(8)=RP*SM*TP     -0.5*(H(15)+H(20)+H(16))                        MOD04482
      DO 770 JJ=1,3                                                     MOD04483
770   CZ(JJ)=0.0                                                        MOD04484
      DO 780 II=1,NI                                                    MOD04485
      DO 780 JJ=1,3                                                     MOD04486
780   CZ(JJ)=CZ(JJ)+H(II)*CORD(II,JJ)                                   MOD04487
      IF(NODE.GT.LIMIT)CALL PRTERR(7)                                   MOD04488
      NSTOR(NODE)=NODE                                                  MOD04489
      X(NODE)=CZ(1)                                                     MOD04490
      Y(NODE)=CZ(2)                                                     MOD04491
      Z(NODE)=CZ(3)                                                     MOD04492
      NODE=NODE+INCR                                                    MOD04493
      FACTR=FACTR*PERR                                                  MOD04494
800   R=R+DR*FACTR                                                      MOD04495
      NODE=NODE-NUMDR*INCR+INCS                                         MOD04496
      FACTS=FACTS*PERS                                                  MOD04497
810   S=S+DS*FACTS                                                      MOD04498
      NODE=NODEI+INCT                                                   MOD04499
      FACTT=FACTT*PERT                                                  MOD04500
820   T=T+DT*FACTT                                                      MOD04501
      IF(KS1.NE.0)CALL VOLM(N3D,INCR,INCS,INCT,NUMDR,NUMDS,NUMDT        MOD04502
     1,KS1,LIMIT,NSTOR,X,Y,Z)                                           MOD04503
      RETURN                                                            MOD04504
1000  CALL PRTERR(6)                                                    MOD04505
      RETURN                                                            MOD04506
      END                                                               MOD04507
	FUNCTION ARSIN(A)
	IMPLICIT REAL*4(A-H,O-Z)
	ARSIN=ASIN(A)
	RETURN
	END
	FUNCTION ARCOS(A)
	IMPLICIT REAL*4(A-H,O-Z)
	ARCOS=ACOS(A)
	RETURN
	END
	FUNCTION DFLOAT(A)
	IMPLICIT REAL*8(A-H,O-Z)
	DFLOAT=DBLE(A)
	RETURN
	END
	FUNCTION DARSIN(A)
	IMPLICIT REAL*8(A-H,O-Z)
	DARSIN=DASIN(A)
	RETURN
	END
	FUNCTION DARCOS(A)
	IMPLICIT REAL*8(A-H,O-Z)
	DARCOS=DACOS(A)
	RETURN
	END
      FUNCTION AGET(IIPOS)                                              MOD02895
      INTEGER SECT,POINT,BEGIN,GETSEC,EPOINT,BLANK,AGET,AGETW           MOD02896
      LOGICAL EOL,EOS,EOF,ERROR,INDEX,GETWRD                            MOD02897
      COMMON  /FRECNT/ LINE(80),SECT,EOL,EOS,EOF,                       MOD02898
     1ERROR,BLANK,ICOMMA,POINT,BEGIN,LENGTH,                            MOD02899
     2EPOINT,NSECT,LINENM,MAXSTR,KUNIT,JUNIT,IUNIT                      MOD02900
      AGET = BLANK                                                      MOD02901
      IPOSIT = IIPOS                                                    MOD02902
500   IF (IPOSIT.GT.MAXSTR.OR.IPOSIT.LE.0) RETURN                       MOD02912
      IF ((BEGIN+IPOSIT-1).LE.80) AGET = LINE (BEGIN+IPOSIT-1)          MOD02913
      AGETW = AGET                                                      MOD02914
      RETURN                                                            MOD02915
      END                                                               MOD02916
      FUNCTION IGET(II)                                                 MOD03052
      INTEGER SECT,POINT,BEGIN,GETSEC,EPOINT,BLANK,AGET,AGETW           MOD03053
      LOGICAL EOL,EOS,EOF,ERROR,INDEX,GETWRD,IGET,RGET,ERR1             MOD03054
      COMMON  /FRECNT/ LINE(80),SECT,EOL,EOS,EOF,                       MOD03055
     1ERROR,BLANK,ICOMMA,POINT,BEGIN,LENGTH,                            MOD03056
     2EPOINT,NSECT,LINENM,MAXSTR,KUNIT,JUNIT,IUNIT                      MOD03057
      COMMON/FRECNM/MULTIP                                              MOD03058
      LOGICAL MODE                                                      MOD03059
      DOUBLE PRECISION RNUMBR,DECIMA                                    MOD03060
      DATA    IPER/1H./,LETE/1HE/,LETD/1HD/                             MOD03061
      II=0                                                              MOD03062
      MODE = .FALSE.                                                    MOD03063
      IGET = .FALSE.                                                    MOD03064
100   IF(MULTIP.GT.0)GO TO 200                                          MOD03070
      IF (.NOT.GETWRD(GET001))RETURN
      IF (LENGTH.EQ.0) RETURN
      RNUMBR = DECIMA(ERROR)                                            MOD03072
      IF(MULTIP.GT.0)GO TO 200                                          MOD03073
      IF (ERROR.AND.EPOINT.EQ.1) EPOINT = BEGIN                         MOD03074
      IF (MODE)XX = RNUMBR                                              MOD03075
      IF (.NOT.MODE) II = RNUMBR                                        MOD03076
      IGET=ERR1(ERR001)                                                 MOD03077
      RGET=ERR1(ERR001)                                                 MOD03078
      RETURN                                                            MOD03079
200   MULTIP=MULTIP-1                                                   MOD03080
      RETURN                                                            MOD03081
      END                                                               MOD03082
      FUNCTION AGETW(AGE001)                                            MOD02895
      INTEGER SECT,POINT,BEGIN,GETSEC,EPOINT,BLANK,AGET,AGETW           MOD02896
      LOGICAL EOL,EOS,EOF,ERROR,INDEX,GETWRD                            MOD02897
      COMMON  /FRECNT/ LINE(80),SECT,EOL,EOS,EOF,                       MOD02898
     1ERROR,BLANK,ICOMMA,POINT,BEGIN,LENGTH,                            MOD02899
     2EPOINT,NSECT,LINENM,MAXSTR,KUNIT,JUNIT,IUNIT                      MOD02900
      AGETW = BLANK                                                     MOD02905
      IPOSIT = 1                                                        MOD02906
      IF (.NOT.GETWRD(GET001))RETURN                                   
       IF (LENGTH.EQ.0) RETURN
      DO 450 ILOOP=BEGIN,80                                             MOD02908
      IF (LINE(ILOOP).EQ.ICOMMA) GO TO 460                              MOD02909
450   CONTINUE                                                          MOD02910
460   MAXSTR = ILOOP - BEGIN                                            MOD02911
500   IF (IPOSIT.GT.MAXSTR.OR.IPOSIT.LE.0) RETURN                       MOD02912
      IF ((BEGIN+IPOSIT-1).LE.80) AGET = LINE (BEGIN+IPOSIT-1)          MOD02913
      AGETW = AGET                                                      MOD02914
      RETURN                                                            MOD02915
      END                                                               MOD02916
      FUNCTION RGET(XX)                                                 MOD03052
      INTEGER SECT,POINT,BEGIN,GETSEC,EPOINT,BLANK,AGET,AGETW           MOD03053
      LOGICAL EOL,EOS,EOF,ERROR,INDEX,GETWRD,IGET,RGET,ERR1             MOD03054
      COMMON  /FRECNT/ LINE(80),SECT,EOL,EOS,EOF,                       MOD03055
     1ERROR,BLANK,ICOMMA,POINT,BEGIN,LENGTH,                            MOD03056
     2EPOINT,NSECT,LINENM,MAXSTR,KUNIT,JUNIT,IUNIT                      MOD03057
      COMMON/FRECNM/MULTIP                                              MOD03058
      LOGICAL MODE                                                      MOD03059
      DOUBLE PRECISION RNUMBR,DECIMA                                    MOD03060
      DATA    IPER/1H./,LETE/1HE/,LETD/1HD/                             MOD03061
      MODE = .TRUE.                                                     MOD03067
      RGET = .FALSE.                                                    MOD03068
      XX=0.0                                                            MOD03069
100   IF(MULTIP.GT.0)GO TO 200                                          MOD03070
      IF (.NOT.GETWRD(GET001))RETURN
      IF (LENGTH.EQ.0) RETURN
      RNUMBR = DECIMA(ERROR)                                            MOD03072
      IF(MULTIP.GT.0)GO TO 200                                          MOD03073
      IF (ERROR.AND.EPOINT.EQ.1) EPOINT = BEGIN                         MOD03074
      IF (MODE)XX = RNUMBR                                              MOD03075
      IF (.NOT.MODE) II = RNUMBR                                        MOD03076
      IGET=ERR1(ERR001)                                                 MOD03077
      RGET=ERR1(ERR001)                                                 MOD03078
      RETURN                                                            MOD03079
200   MULTIP=MULTIP-1                                                   MOD03080
      RETURN                                                            MOD03081
      END                                                               MOD03082
      SUBROUTINE CENT(XYZ,AM,IECHO,IB,NTOT)                             MOD00451
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC               MOD00452
      IMPLICIT REAL*8 (A-H,O-Z)                                         MOD00453
      COMMON /UNIT/ INN,IOUT,IP                                         MOD00454
      COMMON/CGELE/N                                                    MOD00455
      DIMENSION XYZ(NTOT,3),AM(1),                                      MOD00456
     *NBON(3),JBON(3),XD(3),TOTMOM(3),                                  MOD00457
     *NOD(21),XM(21)                                                    MOD00458
      DO 150 I = 1, 3                                                   MOD00459
      NBON(I) = 1                                                       MOD00460
150   CONTINUE                                                          MOD00461
      DO 200 I = 1, NTOT                                                MOD00462
      READ(IB,201) N,(JBON(J),J=1,3)                                    MOD00463
      DO 200 K = 1, 3                                                   MOD00464
      IF (JBON(K) .EQ. 0) NBON(K)=0                                     MOD00465
200   CONTINUE                                                          MOD00466
201   FORMAT(4I5)                                                       MOD00467
      NN = 0                                                            MOD00468
      DO 250 I = 1, 3                                                   MOD00469
      IF (NBON(I) .EQ. 0) NN = NN + 1                                   MOD00470
250   CONTINUE                                                          MOD00471
      IF (NN .GT. 0) GO TO 255                                          MOD00472
      WRITE(IOUT,251)                                                   MOD00473
251   FORMAT(1H1,////,10X,'*** ZERO DEGREE OF FREEDOM SYSTEM ***')      MOD00474
      STOP                                                              MOD00475
255   NNT = NN * NTOT                                                   MOD00476
      DO 260 L = 1, NNT                                                 MOD00477
      AM(L) = 0.0                                                       MOD00478
260   CONTINUE                                                          MOD00479
      DO 300 I = 1, NTOT                                                MOD00480
      READ(IP,301) N,(XYZ(I,J),J=1,3)                                   MOD00481
300   CONTINUE                                                          MOD00482
301   FORMAT(I10,/,3E12.5)                                              MOD00483
      READ(IP,302) IDUMMY                                               MOD00484
302   FORMAT(I10)                                                       MOD00485
      READ(IP,401) NEL                                                  MOD00486
401   FORMAT(I10)                                                       MOD00487
      DO 2000 I = 1, NEL                                                MOD00488
      IPL = 0                                                           MOD00489
      READ(IP,501) KIND,NGROUP,DEN,AREA                                 MOD00490
501   FORMAT(I10,I11,2E12.5)                                            MOD00491
      READ (IP,502) N,(NOD(J),J=1,21)                                   MOD00492
502   FORMAT(I10,21I5)                                                  MOD00493
      IF (AREA .LT. 0.0) KIND = 15                                      MOD00494
      GO TO (1100,1100,1200,1200,1300,1200,1500,1200,                   MOD00495
     *1400,1300,1200,1200,1200,1500,1500,1200,2000,2000),KIND           MOD00496
1100  CALL DIM1(XYZ,DEN,AREA,XM1,NOD,NTOT)                              MOD00497
      GO TO 1700                                                        MOD00498
1200  CALL DIM2(XYZ,DEN,AREA,XM,NOD,NNOD,KIND,NTOT,N)                   MOD00499
      GO TO 1800                                                        MOD00500
1300  CALL DIM3(XYZ,DEN,XM,NOD,NNOD,NTOT,N)                             MOD00501
      GO TO 1800                                                        MOD00502
1400  CALL CURVE(XYZ,DEN,AREA,XM1,NOD,NTOT)                             MOD00503
      GO TO 1700                                                        MOD00504
1500  GO TO 2000                                                        MOD00505
1700  DO 1750 K1 = 1, 2                                                 MOD00506
      INOD = NOD(K1) * NN                                               MOD00507
      KK = 0                                                            MOD00508
      DO 1750 K2 = 1, 3                                                 MOD00509
      IF (NBON(K2) .EQ. 1) GO TO 1750                                   MOD00510
      AM(INOD+KK-NN+1) = XM1 + AM(INOD+KK-NN+1)                         MOD00511
      KK = KK + 1                                                       MOD00512
1750  CONTINUE                                                          MOD00513
      IF (IECHO .GT. 0)   WRITE (IOUT,8888) I, XM1                      MOD00514
8888  FORMAT(10X,'***EL,MASS***',I5,10E10.3,/,28X,10E10.3)              MOD00515
      GO TO 2000                                                        MOD00516
1800  DO 1850 K1 = 1, NNOD                                              MOD00517
      INOD = NOD(K1) * NN                                               MOD00518
      KK = 0                                                            MOD00519
      DO 1850 K2 = 1, 3                                                 MOD00520
      IF (NBON(K2) .EQ. 1) GO TO 1850                                   MOD00521
      AM(INOD+KK+1-NN) = XM(K1) + AM(INOD+KK+1-NN)                      MOD00522
      KK = KK + 1                                                       MOD00523
1850  CONTINUE                                                          MOD00524
      IF (IECHO .GT. 0)   WRITE(IOUT,8888) I,(XM(J),J=1,NNOD)           MOD00525
2000  CONTINUE                                                          MOD00526
      IF (IECHO .GT. 0)  WRITE(IOUT,9998)                               MOD00527
      DO 2001 I = 1, NTOT                                               MOD00528
      IF (IECHO .GT. 0)  WRITE(IOUT,9999) I, (AM(I*NN-NN+J),J=1,NN)     MOD00529
9998  FORMAT(///,5X,'NODE  LUMP-MASS  X          Y           Z',//)     MOD00530
9999  FORMAT(I8,5X,3E12.5)                                              MOD00531
2001  CONTINUE                                                          MOD00532
      TOTDEN = 0.0                                                      MOD00533
      IDTOT  = NN * NTOT                                                MOD00534
      DO 2010 I = 1, 3                                                  MOD00535
      TOTMOM(I) = 0.0                                                   MOD00536
2010  CONTINUE                                                          MOD00537
      DO 2100 I = 1, IDTOT, NN                                          MOD00538
      TOTDEN = TOTDEN + AM(I)                                           MOD00539
2100  CONTINUE                                                          MOD00540
      KK = 0                                                            MOD00541
      DO 2700 I = 1, 3                                                  MOD00542
      IF (NBON(I) .EQ. 1) GO TO 2800                                    MOD00543
      KK = KK + 1                                                       MOD00544
      DO 2500 J = 1, NTOT                                               MOD00545
      J1 = J * NN                                                       MOD00546
      TOTMOM(I) = TOTMOM(I) + AM(J1+KK-NN)*XYZ(J,I)                     MOD00547
2500  CONTINUE                                                          MOD00548
2800  CONTINUE                                                          MOD00549
2700  CONTINUE                                                          MOD00550
                                                                        MOD00551
      DO 2600 I = 1, 3                                                  MOD00552
      XD(I) = TOTMOM(I) / TOTDEN                                        MOD00553
2600  CONTINUE                                                          MOD00554
      GRAV = 32.0                                                       MOD00555
      TOTWGT = TOTDEN * GRAV                                            MOD00556
      WRITE (6,3001) (XD(I),I=1,3),TOTDEN,TOTWGT                        MOD00557
3001  FORMAT(1H1,/////,                                                 MOD00558
     *10X,'***  CENTRAL GRAVITY ***',//,                                MOD00559
     *10X,'X - DIRECTION =',E12.5,/,                                    MOD00560
     *10X,'Y - DIRECTION =',E12.5,/,                                    MOD00561
     *10X,'Z - DIRECTION =',E12.5,/,                                    MOD00562
     *10X,'TOTAL DENSITY =',E12.5,/,                                    MOD00563
     *10X,'TOTAL WEIGHT  =',E12.5//)                                    MOD00564
      RETURN                                                            MOD00565
      END                                                               MOD00566
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC                 MOD00567
      FUNCTION ICRSEL (POINT,NUM,COMP)                                  MOD03449
CC..        ADAPTED FOR MULTIPLE POINT CONSTRAINTS BY                   MOD03450
CC..    DR. ELIEZER MENDELSSOHN                                         MOD03451
CC..    SAP USERS' GROUP                                                MOD03452
CC..    DEPT. OF CIVIL ENG.                                             MOD03453
CC..    UNIVERSITY OF SOUTHERN CALIFORNIA                               MOD03454
CC..    UNIVERSITY PARK                                                 MOD03455
CC..    LOS ANGELES, CA. 90007                                          MOD03456
CC..    PHONE (213)743-5508                                             MOD03457
      INTEGER POINT(8),COMP(NUM)                                        MOD03458
      DO 11 I=1,4                                                       MOD03459
11    COMP(I)=POINT(I)                                                  MOD03460
      K=0                                                               MOD03461
      L=4                                                               MOD03462
      INC=1                                                             MOD03463
      DO 22 IND=5,8                                                     MOD03464
      IF (L.GE.NUM) GO TO 999                                           MOD03465
      J=POINT(IND)                                                      MOD03466
      IF (J.LT.1) GO TO 22                                              MOD03467
      L=L+1                                                             MOD03468
      COMP(L)=J                                                         MOD03469
      K=K+INC                                                           MOD03470
22    INC=INC+INC                                                       MOD03471
999   ICRSEL=K                                                          MOD03472
      RETURN                                                            MOD03473
      END                                                               MOD03474
