home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 9 / 09.iso / e / e003 / 12.ddi / TAP6M.FOR < prev    next >
Encoding:
Text File  |  1988-01-04  |  37.8 KB  |  475 lines
  1.         PROGRAM TAP6
  2. C       ROCKETDYNE THERMAL ANALYZER PROGRAM   VERSION TAP-6  AUG.1969           
  3. C       THIS IS THE MAIN CONTROL PROGRAM FOR TAP-6 WHICH WILL SOLVE             
  4. C       PARABOLIC DIFFERENTIAL EQUATIONS USING A THREE TIME LEVEL METHOD        
  5. C       DEVELOPED AT ROCKETDYNE. THIS PROGRAM CAN SOLVE PROBLEMS                
  6. C       INVOLVING 999 TEMPERATURE NODES AND 2999 ADMITTANCES INCLUDING          
  7. C       THE EFFECTS OF ENERGY TRANSPORT THROUGH THE MATERIAL.                   
  8. C                ROCKETDYNE ADVANCED SYSTEMS  DEPT. 589                         
  9. C                THIS IS A REDUCED MEMORY REQUIREMENT DECK                      
  10. C                MODIFIED FOR FLOW NOV.1972                                     
  11. C                 MODIFIED 3/16/82                                              
  12.       INTEGER NY1(3000),NY2(1000),NZA(1000),INXP(1000),INWP(1000)               
  13.      X,NXP(1000),KEY(1000),IZ(1000),IPR(500),IPH(100)                           
  14.      X,NTHEAD(100),NTTAIL(100),NFLAG(1000)                                      
  15.       DIMENSION  ZZ(16000),DYTCQ(11030),Y(3000),T(1000)                         
  16.      X,C(1000)        ,Q(1000)        ,GAIN(1000) ,A(1000),B(1000)              
  17.      X,D(1000),SUMY(1000),SUMYT(1000), TOLD(1000),L(1000)                       
  18.      X,PHT1(100),PHH1(100),PHT2(100),PHH2(100),QCH1(100)                        
  19.      X,QCH2(100),CON(20),PRTIM(10),DELTIM(10),TITLE(17,10),TABL(1000)           
  20.      X,TPROG(2),AA(7)                                                           
  21.        COMMON/I2BYT/NY1,NY2,NZA,NXP,INXP,INWP,KEY,IPR,IPH                       
  22.      X,NTHEAD,NTTAIL,NFLAG,IZ                                                   
  23.       COMMON /LAB1/ PLOT(100,10),IPLOT(10),TMAX(1000)                           
  24.       COMMON /LAB2/DATE(3)                                                      
  25.       COMMON /LAB3/ICASE,XLIM(2),YLIM(2),CT(5),XLAB(5),YLAB(5)                  
  26.       COMMON /LAB4/ SPROG(2)                                                    
  27.       COMMON /LAB5/ TITLE                                                       
  28.         COMMON ZZ                                                               
  29.       EQUIVALENCE (ZZ(001),TIME),         (ZZ(002),Y(1),DYTCQ(1))               
  30.      X,          (ZZ(3002),T(1)),                 (ZZ(4002),C(1))               
  31.      X,          (ZZ(5002),Q(1)),                 (ZZ(6002),D(1) )              
  32.      X,          (ZZ(7002),A(1) ),              (ZZ(8002), B(1) )               
  33.      X,          (ZZ(9002),L(1) ),                (ZZ(10002),GAIN(1) )          
  34.      X,     (ZZ(10001),LYMIN),     (ZZ(10501),LYMAX)                            
  35.      X,          (ZZ(11002),CON(1)),              (ZZ(11010),CON9,SIGMA)        
  36.      X,          (ZZ(11011),CON10,TZERO), (ZZ(11032),PRTIM(1))                  
  37.      X,          (ZZ(11042),DELTIM(1)),           (ZZ(11053),IREAD)             
  38.      X,          (ZZ(11056),NYMIN),               (ZZ(11057),MAX2T)             
  39.      X,          (ZZ(11058),MIN2Y),               (ZZ(11059),MAX2Y)             
  40.       EQUIVALENCE(ZZ(11060),MIN1T),               (ZZ(11061),MAX1T)             
  41.      X,          (ZZ(11062),N2)                                                 
  42.       EQUIVALENCE(ZZ(11069), LFLAG),              (ZZ(11070), IERFL)            
  43.      X,          (ZZ(11071),PRTMX),               (ZZ(11072),DTHETA)            
  44.      X,          (ZZ(11073),RCMIN),               (ZZ(11074), JSAVE)            
  45.      X,          (ZZ(11075),DELTAT),              (ZZ(11076),DTMAX )            
  46.      X,          (ZZ(11077), IMIN),               (ZZ(11078), IMAX )            
  47.      X,          (ZZ(11079), IXMIN),              (ZZ(11080), IXMAX)            
  48.      X,          (ZZ(11081),IPHMAX),              (ZZ(11082),IPRMAX)            
  49.      X,          (ZZ(11083),IPFLAG),              (ZZ(11084),ISFLAG)            
  50.      X,          (ZZ(11085),MAXNOT),              (ZZ(11086),MAXNOY)            
  51.      X,          (ZZ(11087),NARITH),              (ZZ(11088),NCRIT )            
  52.       EQUIVALENCE (ZZ(11089),NTIME),              (ZZ(11090), NT )              
  53.      X,          (ZZ(11091),NTLAST),              (ZZ(11092),NWW)               
  54.      X,          (ZZ(11093),NXX),                 (ZZ(11094),NZZ)               
  55.      X,          (ZZ(11095), W ),                 (ZZ(11096), X )               
  56.      X,          (ZZ(11097), Z ) ,                (ZZ(11098), NTAB)             
  57.      X,          (ZZ(11099),OTIME2),              (ZZ(11100),ODTIME)            
  58.      X,          (ZZ(11101),DTIME),               (ZZ(11102),PHT1(1))           
  59.      X,          (ZZ(11202),PHH1(1)),             (ZZ(11302),PHT2(1))           
  60.      X,          (ZZ(11402),PHH2(1)),             (ZZ(11502),QCH1(1))           
  61.      X,          (ZZ(11602),QCH2(1))                                            
  62.       EQUIVALENCE(ZZ(12002),SUMY(1)),             (ZZ(13002),TOLD(1))           
  63.      X,          (ZZ(14002),SUMYT(1)),            (ZZ(15002),TABL(1))           
  64.       DATA  TPROG/4H TAP,4H-6F /                                                
  65. C                                                                               
  66. C                                                                               
  67.       DATA BLANK/4H    /     
  68. C                                                                               
  69.       REWIND 8                                                                  
  70.       NODE=1000                                                                 
  71.         MEGA=1000000                                                            
  72.       CALL CDATEV(DATE)                                                         
  73.       SPROG(1)= TPROG(1)                                                        
  74.       SPROG(2)= TPROG(2)                                                        
  75. C                                                                               
  76.       CALL HEADER                                                               
  77. C     CLEAR CORE LOCATIONS AND EQUATE PROGRAM VARIABLES                         
  78.       DO 102 I=1,9                                                              
  79.   102 XLIM(I)=0.0                                                               
  80.   105 DO 110 I=1,16000                                                          
  81.   110 ZZ(I) = 0.0                                                               
  82.       DO 111 I=1,170                                                            
  83.   111 TITLE(I,1)=BLANK                                                          
  84.       DO 112 I=1,11400                                                          
  85.   112 NY1(I)=0                                                                  
  86.       NSTART = -98                                                              
  87.       IERFL  = 1                                                                
  88.       IMIN   = 1                                                                
  89.       MAXNOY = 0                                                                
  90.       LYMAX = 0                                                                 
  91.       IPHK   = 0                                                                
  92.       IF1=-3                                                                    
  93.       IF2=-2                                                                    
  94.       IF3= 1                                                                    
  95.         IF4=-4                                                                  
  96.       JSTART=1                                                                  
  97.       ICASE=0                                                                   
  98.       ISTART = 0                                                                
  99.       NTLAST = -1                                                               
  100.       NCZERO = 1                                                                
  101.       SIGMA=0.3304E-14                                                          
  102.       TZERO  = 459.69                                                           
  103.       CON(1) = 1.0E35                                                           
  104.       CON(3) = 1.0                                                              
  105.       CON(4)= .001                                                              
  106.       CON17 = 8.                                                                
  107.   115 IERFL = 1                                                                 
  108.       ICASE= ICASE+1                                                            
  109.       CFLAG = 0                                                                 
  110.       ICZERO = 0                                                                
  111.       DELTAT = 0.0                                                              
  112.       ODTIME = 0.0                                                              
  113.       OTIME2 = 0.                                                               
  114.       DTMAX = 0.                                                                
  115.       DTHETA= 0.0                                                               
  116.       DTIME= 0.0                                                                
  117.       ATIME = 0.                                                                
  118.       IF (NSTART+97) 145,135,125                                                
  119.   125 WRITE (6,130)                                                             
  120.       IMIN=1                                                                    
  121.   130 FORMAT (1H1,75H THIS RUN USES THE SAME INITIAL DATA AS THE PREVIOU        
  122.      1S RUN WITH SOME NEW DATA)                                                 
  123.       READ(8)(ZZ(I),I=1,16000)                                                  
  124.       REWIND 8                                                                  
  125.       GO TO 200                                                                 
  126.   135 WRITE (6,140)                                                             
  127.   140 FORMAT (1H1,66H THIS RUN IS A CONTINUATION OF THE PREVIOUS RUN WIT        
  128.      1H SOME NEW DATA)                                                          
  129.       GO TO 200                                                                 
  130.   145 WRITE(6,150) DATE,SPROG                                                   
  131.       IMIN=1                                                                    
  132.       ISTART= 0                                                                 
  133.       IMAX = 0                                                                  
  134.       DO 146 I=1,1010                                                           
  135.   146 PLOT(I,1)= 0.0                                                            
  136.   150 FORMAT(1H1,' ALL NEW INPUT DATA FOR THIS CASE',51X,3A4,1X,2A4)            
  137. C                                                                               
  138. C     INPUT                                                                     
  139. C                                                                               
  140.   200 CALL COUNTV                                                               
  141.       IF(ICASE.GT.1)WRITE(6,210)(TITLE(J,1),J=1,15),ICASE,SPROG,DATE            
  142.   210 FORMAT(1H /,1X,15A4,'  CASE NO.',I4,2X,2A4,1X,2A4)                        
  143.       CALL TAPIN(ISTART,NSTART)                                                 
  144.       IF(CON(15).GT.0.)IMIN=1                                                   
  145.                                                                                 
  146.       IF(CON(5).GT.0.0.AND.CON(5).LT.1.0)CON17=CON(5)*8.0                       
  147.       IRAD = CON(16)                                                            
  148.       ICON12 = CON(12)                                                          
  149.       CALL TIMEV(TIME2)                                                         
  150.       WRITE(6,806)DATE,TIME2                                                    
  151.       IF (IERFL-2) 220,105,999                                                  
  152.   220 IF(JSTART.EQ.1)GO TO 230                                                  
  153.       IF(NSTART.NE.-96) GO TO 240                                               
  154.   230 WRITE(8)(ZZ(I),I=1,16000)                                                 
  155.       REWIND 8                                                                  
  156.   240 TIME= PRTIM(1)                                                            
  157.       NTIME= 1                                                                  
  158.       IPFLAG = 1                                                                
  159.       ISFLAG = 0                                                                
  160.       NARITH = 1                                                                
  161.       IXMIN = 499                                                               
  162.       IXMAX = 0                                                                 
  163.       TZERO = ABS(TZERO)                                                        
  164.       WRITE(6,242) IMAX,IMIN                                                    
  165.   242 FORMAT(1H ,19X,'THE INITIAL *030 DO LOOP RANGE IS IMAX=',I4,              
  166.      1'IMIN=',I4)                                                               
  167.       CALL ARITH                                                                
  168.       N2=1                                                                      
  169.       MINZC = 999                                                               
  170.       MAXZC = 1                                                                 
  171.       MINZY= 2999                                                               
  172.       MAXZY = 1                                                                 
  173.       MAX1T= 1                                                                  
  174.       MIN1T= 999                                                                
  175.       MIN1Y = 2999                                                              
  176.       MAX1Y = 1                                                                 
  177.       MIN3T=999                                                                 
  178.       MAX3T=1                                                                   
  179.       N2= 1+ICON12                                                              
  180.       IF(CON(5).GT.0.) N2= 4                                                    
  181.       ICFLAG = 0                                                                
  182.       IF(N2.EQ.4) DELTAT = 1.                                                   
  183. C     PREPARE FOR ZEROTH ITERATION                                              
  184.       IF(LYMAX.LE.LYMIN) GO TO 244                                              
  185.       DO 243 JJ=1,LYMAX                                                         
  186.       J=NY2(JJ)                                                                 
  187.         J1=MOD(NY1(J),MEGA)/NODE                                                
  188.       J2= MOD(NY1(J),NODE)                                                      
  189.       IF(NFLAG(J1).NE.IF1.AND.NFLAG(J2).NE.IF1) GO TO 243                       
  190.       MIN3T= MIN0(J1,J2,MIN3T)                                                  
  191.       MAX3T= MAX0(J1,J2,MAX3T)                                                  
  192.   243 CONTINUE                                                                  
  193.         J=NY2(1)                                                                
  194.         J1=MOD(NY1(J),MEGA)/NODE                                                
  195.         NYMIN=LYMAX                                                             
  196.         IF(NFLAG(J1).EQ.-1)NYMIN=0                                              
  197.   244 CONTINUE                                                                  
  198.       DO 250 JJ=1,MAXNOY                                                        
  199.       IF(NY1(JJ).LE.0) GO TO 250                                                
  200.       JT= NY1(JJ)/MEGA                                                          
  201.       IF(JT.GT.0) GO TO 250                                                     
  202.        J1=MOD(NY1(JJ),MEGA)/NODE                                                
  203.       J2= MOD(NY1(JJ),NODE)                                                     
  204.         IF(NFLAG(J1).GE.1.OR.NFLAG(J2).GE.1)MAX1Y=MAX0(MAX1Y,JJ)                
  205.         IF(NFLAG(J1).GE.1.OR.NFLAG(J2).GE.1)MIN1Y=MIN0(MIN1Y,JJ)                
  206.       IF(NFLAG(J1).NE.IF4.AND.NFLAG(J2).NE.IF4) GO TO 250                       
  207.       MAXZY= MAX0(JJ,MAXZY)                                                     
  208.       MINZY= MIN0(JJ,MINZY)                                                     
  209.       ICZERO = 1                                                                
  210.   250 CONTINUE                                                                  
  211.       CALL SETUP(JSTART)                                                        
  212.       IF(IERFL.GT.1) GO TO 830                                                  
  213.       DO 260 J=1,999                                                            
  214.       IF(NFLAG(J).EQ.0) GO TO 260                                               
  215.       IF(N2.NE.4)TMAX(J)= AMAX1(TMAX(J),T(J))                                   
  216.       IF(NFLAG(J).EQ.IF4)MAXZC=MAX0(MAXZC,J)                                    
  217.       IF(NFLAG(J).EQ.IF4)MINZC=MIN0(MINZC,J)                                    
  218.         IF(NFLAG(J).GE.1)MIN1T=MIN0(J,MIN1T)                                    
  219.         IF(NFLAG(J).GE.1)MAX1T=MAX0(J,MAX1T)                                    
  220.   260 TOLD(J)= T(J)                                                             
  221.       WRITE(6,262) IMAX,IMIN                                                    
  222.   262 FORMAT(1H ,19X,'THE *030 DO LOOP RANGE IS IMAX=',I4,',IMIN=',I4)          
  223.       WRITE(6,265)                                                              
  224.   265 FORMAT(20X,'THE DO LOOP RANGES FOR THIS PROBLEM ARE AS FOLLOWS'/          
  225.      120X,      '  NFLAG       MIN Y      MAXY       MINT       MAXT ' )        
  226.       WRITE(6,270)IF3,MIN1Y,MAX1Y,MIN1T,MAX1T                                   
  227.       IF(ICZERO.EQ.1)WRITE(6,270)IF4,MINZY,MAXZY,MINZC,MAXZC                    
  228.       IF(MIN3T.LT.MAX3T)WRITE(6,270)IF1,NY2(1),NY2(LYMAX),MIN3T,MAX3T           
  229.   270 FORMAT(20X,5(3X,I5,3X))                                                   
  230.       CALL TIMEV(TIME2)                                                         
  231.       WRITE(6,806) DATE,TIME2                                                   
  232.       WRITE(6,275)                                                              
  233.   275 FORMAT(1H1,'  INITIAL CONDITIONS FOR THIS RUN')                           
  234.       GO TO 310                                                                 
  235.   300 CONTINUE                                                                  
  236.       CALL ARITH                                                                
  237.       ATIME= TIME + CON(13)                                                     
  238. C                                                                               
  239. C                RELAXATION OF CONNECTED ZERO CAPACITANCE NODES                 
  240. C                                                                               
  241. C                CALCULATION OF SUMY AND SUMYT                                  
  242. C                                                                               
  243.   310 IF(IPFLAG.GT.0)CALL TAPOUT                                                
  244.       DO 311 J=1,MAXNOT                                                         
  245.       SUMY(J) = 0.                                                              
  246.   311 SUMYT(J)= Q(J)                                                            
  247.   312 DTMAX2= 0.                                                                
  248.       DTMAX = 0.                                                                
  249.       IF(ICZERO.LE.0) GO TO 340                                                 
  250.       DO 325 J= MINZY,MAXZY                                                     
  251.        JT=NY1(J)/MEGA                                                           
  252.       IF(JT.GT.0.OR.Y(J).LE.0.0) GO TO 325                                      
  253.       J1= MOD(NY1(J),MEGA)/NODE                                                 
  254.       J2= MOD(NY1(J),NODE)                                                      
  255.       IF(NFLAG(J1).NE.IF4.AND.NFLAG(J2).NE.IF4)GO TO 325                        
  256.       YJ= Y(J)                                                                  
  257.       IF(CON(16).LE.0.) GO TO 315                                               
  258.       IF(IRAD.EQ.2.AND.NFLAG(J1).NE.IF4)GO TO 315                               
  259.       IF(IRAD.EQ.2.AND.NFLAG(J2).NE.IF4)GO TO 315                               
  260.       TJ1= T(J1)+CON10                                                          
  261.       TJ2= T(J2)+CON10                                                          
  262.       TCO= CON9*(TJ1**2+TJ2**2)*(TJ1+TJ2)                                       
  263.       YJ= TCO*YJ                                                                
  264.   315 IF(NFLAG(J1).NE.IF4)GO TO 320                                             
  265.       SUMY(J1)=SUMY(J1)+YJ                                                      
  266.         SUMYT(J1)=SUMYT(J1)+YJ*T(J2)                                            
  267.   320 IF(NFLAG(J2).NE.IF4)GO TO 325                                             
  268.       SUMY(J2)=SUMY(J2)+ YJ                                                     
  269.       SUMYT(J2)=SUMYT(J2)+YJ*T(J1)                                              
  270.   325 CONTINUE                                                                  
  271. C                                                                               
  272. C                CALCULATION OF CONNECTED ZERO CAPACITANCE NODE TEMP.           
  273. C                                                                               
  274.       DO 335 J= MINZC,MAXZC                                                     
  275.       IF(NFLAG(J).NE.IF4) GO TO 335                                             
  276.       IF(SUMY(J).LE.0.) GO TO 335                                               
  277.       TNEW= SUMYT(J)/SUMY(J)+(T(J)-TOLD(J))/CON17                               
  278.       DTMAX2= AMAX1(DTMAX2,ABS(TNEW-T(J)))                                      
  279.       TOLD(J)=T(J)                                                              
  280.       T(J)=TNEW                                                                 
  281.       SUMY(J)= 0.                                                               
  282.       SUMYT(J)= Q(J)                                                            
  283.   335 CONTINUE                                                                  
  284.       IF(DTMAX2.GT.CON(4))GO TO 312                                             
  285.   340 CONTINUE                                                                  
  286. C                                                                               
  287. C                CALCULATION OF SUMY AND SUMYT FOR INTEGRATION OF T             
  288. C                                                                               
  289.         DO 350 J=1,MAXNOY                                                       
  290.       JT= NY1(J)/MEGA                                                           
  291.       IF(JT.GT.0.OR.Y(J).LE.0.0) GO TO 350                                      
  292.       J1= MOD(NY1(J),MEGA)/NODE                                                 
  293.       J2= MOD(NY1(J),NODE)                                                      
  294.       YJ= Y(J)                                                                  
  295.       IF(IRAD   .NE.1 ) GO TO 345                                               
  296.       IF(NFLAG(J1).NE.IF4.AND.NFLAG(J2).NE.IF4) GO TO 345                       
  297.       TJ1= T(J1) +CON10                                                         
  298.       TJ2= T(J2)+CON10                                                          
  299.       TCO= CON9*(TJ1**2+TJ2**2)*(TJ1+TJ2)                                       
  300.       YJ= TCO*YJ                                                                
  301.   345 CONTINUE                                                                  
  302.       IF(NFLAG(J1).EQ.0.OR.NFLAG(J1).EQ.-1) GO TO 346                           
  303.       IF(NFLAG(J1).EQ.IF4)GO TO 346                                             
  304.       SUMY(J1)=SUMY(J1)+ YJ                                                     
  305.       SUMYT(J1)= SUMYT(J1)+ YJ*T(J2)                                            
  306.   346 IF(NFLAG(J2).EQ.0.OR.NFLAG(J2).EQ.-1) GO TO 347                           
  307.       IF(NFLAG(J2).EQ.IF4)GO TO 347                                             
  308.       SUMY(J2)=SUMY(J2)+ YJ                                                     
  309.       SUMYT(J2)=SUMYT(J2)+YJ*T(J1)                                              
  310.   347 CONTINUE                                                                  
  311.   350 CONTINUE                                                                  
  312.       IF(LYMAX.LE.LYMIN) GO TO 361                                              
  313. C                                                                               
  314. C                ADDITION OF FIRST DERIVATIVE TERMS TO SUMYT                    
  315. C                                                                               
  316.       DO 360 JJ=1,LYMAX                                                         
  317.       J=NY2(JJ)                                                                 
  318.       IF(J.LE.0.OR.Y(J).EQ.0.) GO TO 360                                        
  319.       JP= J                                                                     
  320.       IF(JJ.GT.1) JP= NY2(JJ-1)                                                 
  321.       J1= MOD(NY1(J),MEGA)/NODE                                                 
  322.       JTP= MOD(NY1(JP),NODE)                                                    
  323.       J2= MOD(NY1(JJ),NODE)                                                     
  324.       YJ= .5*Y(J)                                                               
  325.        IF(NFLAG(J).LT.1)GO TO 360                                               
  326.       IF(JTP.NE.J1.OR.NFLAG(J2).EQ.-1) GO TO 355                                
  327.       JTP= MOD(NY1(JP),MEGA)/NODE                                               
  328.       TERM=-YJ*(T(JTP)-T(J2))                                                   
  329.       SUMYT(J1)= SUMYT(J1)+ TERM                                                
  330.       GO TO 360                                                                 
  331.   355 SUMYT(J1)= SUMYT(J1)+YJ*T(J2)                                             
  332.       SUMY(J1)= SUMY(J1)+YJ                                                     
  333.   360 CONTINUE                                                                  
  334.   361 CONTINUE                                                                  
  335.       CALL TIMER                                                                
  336.       IF(IPHMAX.GT.0) CALL LATENT(1)                                            
  337.       DTIME= DTHETA                                                             
  338. C  FLOW NODE CALCULATIONS                                                       
  339.       DO 420 JJ=1,LYMAX                                                         
  340.         IF(NYMIN.EQ.0)GO TO 380                                                 
  341.       K=LYMAX+1-JJ                                                              
  342.       J=NY2(K)                                                                  
  343.       J1= MOD(NY1(J),MEGA)/NODE                                                 
  344.       J2=MOD(NY1(J),NODE)                                                       
  345.         GO TO 385                                                               
  346. 380   J=NY2(JJ)                                                                 
  347.         J2=MOD(NY1(J),MEGA)/NODE                                                
  348.         J1=MOD(NY1(J),NODE)                                                     
  349. 385     IF(NFLAG(J1).NE.IF1)GO TO 420                                           
  350.        DT=1.                                                                    
  351.       TOLD(J1)=T(J1)                                                            
  352.         IF(C(J1).GT.0.0)DT=Y(J)*DTIME/C(J1)                                     
  353.         CJ=AMAX1(0.0,C(J1)*(1.-DT))                                             
  354.         IF(N2.EQ.4)CJ=0.0                                                       
  355.   400 DENOM=SUMY(J1)+CJ/DTIME +Y(J)                                             
  356.       IF(DENOM.LE.0.0)GO TO 420                                                 
  357.         IF(DT.GE.1.)TDUM=TOLD(J2)/DT+(DT-1.)*T(J2)/DT                           
  358.         IF(DT.LT.1.)TDUM=DT*TOLD(J2)+(1.-DT)*TOLD(J1)                           
  359.       T(J1)=(SUMYT(J1)+CJ*TDUM/DTIME +Y(J)*TDUM)/DENOM                          
  360.   420 CONTINUE                                                                  
  361. C                                                                               
  362. 450     GO TO (500,1000,1000,605),N2                                            
  363. C                                                                               
  364. C     RELAX ALL FINITE CAPACITY NODES                                           
  365. C                NEWTON INTEGRATION METHOD (FORWARD DIFFERENCE)                 
  366. C                                                                               
  367.   500 CONTINUE                                                                  
  368.       DO 540 J=MIN1T,MAX1T                                                      
  369.         IF(NFLAG(J).LT.1)GO TO 540                                              
  370.   510 IF(C(J).GT.0.) GO TO 530                                                  
  371.   520 IF(SUMY(J).GT.0.0)T(J)=SUMYT(J)/SUMY(J)                                   
  372.       GO TO 540                                                                 
  373.   530 TNEW=DTIME* SUMYT(J)/C(J)+(1.-DTIME*SUMY(J)/C(J))*T(J)                    
  374.       DT= ABS(TNEW-T(J))                                                        
  375.       TOLD(J)= T(J)                                                             
  376.       T(J)= TNEW                                                                
  377.       TMAX(J)=AMAX1(TMAX(J),T(J))                                               
  378.       IF(DT.LE.DTMAX)GO TO 540                                                  
  379.       DTMAX = DT                                                                
  380.       JSAVE = J                                                                 
  381.   540 CONTINUE                                                                  
  382.       ODTIME= DTIME                                                             
  383.   545 CONTINUE                                                                  
  384.   560 GO TO (700,830),IERFL                                                     
  385.   605 CONTINUE                                                                  
  386. C                                                                               
  387. C     STEADY-STATE                                                              
  388. C                                                                               
  389.   610 DO 650 J= MIN1T,MAX1T                                                     
  390.   620  IF(NFLAG(J).GE.1.AND.SUMY(J).GT.0.0)GO TO 630                            
  391.       GO TO 650                                                                 
  392.   630 DT= SUMYT(J)/SUMY(J)-T(J)+(T(J)-TOLD(J))/CON17                            
  393.       IF(CON(5).LE.1.) GO TO 631                                                
  394.       DTO = T(J)-TOLD(J)                                                        
  395.       IF((DTO*DT).GT.0.0)DT=SUMYT(J)/SUMY(J) -TOLD(J)                           
  396.   631 CONTINUE                                                                  
  397.       TOLD(J)= T(J)                                                             
  398.       T(J)     = T(J)+DT                                                        
  399.       TMAX(J)= DT                                                               
  400.       DT= ABS(DT)                                                               
  401.       IF (DT .LE. DTMAX) GO TO 650                                              
  402.   640 DTMAX  = DT                                                               
  403.       JSAVE  = J                                                                
  404.   650 CONTINUE                                                                  
  405. C     TEST FOR PRINT-OUT SIGNAL                                                 
  406.   700 IF(IERFL.GT.1) GO TO 830                                                  
  407.       IF(IPHMAX.GT.0) CALL LATENT(2)                                            
  408.       IF(ISFLAG.GT.0) GO TO 770                                                 
  409.       DTMAX= DTMAX/DTHETA                                                       
  410.       IF(DTMAX.LT.CON(4)) GO TO 740                                             
  411.       IF(IPFLAG.GT.0) GO TO 300                                                 
  412.       IF(TIME.LT.ATIME) GO TO 310                                               
  413.   710 GO TO 300                                                                 
  414.   740 ICFLAG = ICFLAG+1                                                         
  415.       IF(ICFLAG.LT.2) GO TO 300                                                 
  416.   750 WRITE (6,760) DTMAX, CON(4)                                               
  417.   760 FORMAT (30H0 GREATEST TEMPERATURE CHANGE ,E12.5, 36H IS LESS THAN      
  418.      1CONVERGENCE TEST VALUE,E12.5)                                             
  419.       ICFLAG = 0                                                                
  420.   770 ISFLAG=1                                                                  
  421.       CALL ARITH                                                                
  422.   780 CALL TAPOUT                                                               
  423. C     TEST FOR PLOTTING                                                         
  424.   800 ISFP1=0                                                                   
  425.   804 CONTINUE                                                                  
  426.         IF(CON(6).GT.0.0)CALLSPILL                                              
  427.       DO 805 J=1,10                                                             
  428.   805 ISFP1=ISFP1+IPLOT(J)                                                      
  429.       IF(ISFP1.GT.0)CALL PLOTPT                                                 
  430. C     TEST FOR FINAL DUMP SIGNAL                                                
  431.   820 NS     = NSTART+100                                                       
  432.       CALL TIMEV(TIME2)                                                         
  433.       WRITE(6,806) DATE,TIME2                                                   
  434.   806 FORMAT(1H0,31X,3A4,'   ELAPSED TIME= ',F12.3,' SECONDS ')                
  435.       GO TO (860,105,115,115),NS                                                
  436.   830 WRITE (6,840)                                                             
  437.   840 FORMAT (52H0 THIS RUN HAS BEEN DISCONTINUED, FINAL VALUES ARE -)          
  438.       CALL TAPOUT                                                               
  439.       GO TO 800                                                                 
  440.   860 CONTINUE                                                                  
  441.   999 CALL EXIT
  442.       STOP                                                                      
  443. C                                                                               
  444.                                                                                 
  445. C      THREE TIME LEVEL INTEGRATION METHOD                                      
  446. C                                                                               
  447.  1000 CONTINUE                                                                  
  448.       IF(ODTIME.LE.0.) ODTIME= DTIME                                            
  449.       TR= 1.+ ODTIME/DTIME                                                      
  450.  1015 DTT= ODTIME + DTIME                                                       
  451.       DO 1040 J= MIN1T,MAX1T                                                    
  452.         IF(NFLAG(J).LT.1) GO TO 1040                                            
  453.       IF(SUMY(J).LE.0.) GO TO 1040                                              
  454.  1017 CONTINUE                                                                  
  455.       TNEW=TR*C(J)*T(J)+DTT*SUMYT(J)-DTIME*SUMY(J)*TOLD(J)                      
  456.       DENOM = TR*C(J)+ODTIME*SUMY(J)                                            
  457.         IF(N2.GT.2)GOTO 1020                                                    
  458.         IF(C(J).GT.0.0)GOTO 1025                                                
  459.         T(J)=SUMYT(J)/SUMY(J)                                                   
  460. 1020    TNEW=TNEW+2.*L(J)*(DTT*T(J)/DTIME-TOLD(J))/ODTIME                       
  461.                                                                                 
  462.         RC=AMIN1(RCMIN,SQRT(L(J)/SUMY(J))/2.)                                   
  463. 1025   TOLD(J)=T(J)                                                             
  464.         T(J)=TNEW/DENOM                                                         
  465.       DT= ABS(T(J)-TOLD(J))                                                     
  466.       TMAX(J)=AMAX1(TMAX(J),T(J))                                               
  467.       IF(DT.LE.DTMAX)GO TO 1040                                                 
  468.       DTMAX= DT                                                                 
  469.       JSAVE = J                                                                 
  470.  1040 CONTINUE                                                                  
  471.       OTIME2=ODTIME                                                             
  472.       ODTIME = DTIME                                                            
  473.  1050 GO TO 700                                                                 
  474.       END                                                                       
  475.