next up previous contents index
Next: Up: Nonlinear Burnett Coefficients Previous: Analysis of Data

Putting it all Together

The program is run as a series of batch jobs, each lasting about an hour. When the program starts up, it reads the file on unit 10 to obtain the run number. If it is the first run, it will generate a face-centred cubic crystal lattice, which will then be melted to provide the fluid. After that, it will read in the current state of the system on unit 11 and the contents of the shift register on unit 13. Then it performs a loop in which the system is evolved and the TTCF calculated. Every now and then, output data is written out as a check point. When the time limit has expired, all the files are updated, and the partial sum of the TTCFs is appended to unit 12.

       PROGRAM TTCF                                                             
C calculate the burnett coefficients using the new Nose-Hoover algorithm        
C Z indicates in the Nose Current-stat                                          
       INCLUDE (PARAM)                                                          
       INTEGER J,KB                                                             
       INTEGER START, ERR, IOPEN, TLIMIT, TUSED, TLEFT                          
       DOUBLE PRECISION GZ(DIMPS+2),XZ(D*NP),PZ(D*NP), CALJ                     
       DOUBLE PRECISION G2(DIMPS+2),G3(DIMPS+2),G4(DIMPS+2)                     
       DOUBLE PRECISION G1(DIMPS+2),LAMBDA                                      
       EQUIVALENCE (GZ,XZ),(GZ(D*NP+1),PZ),(GZ(DIMPS+2),LAMBDA)                 
       COMMON /GAMMAG/GZ                                                        
       COMMON /GSTORZ/G1,G2,G3,G4                                               
       INCLUDE (BUFFER)                                                         
       DOUBLE PRECISION SHIFT(NGPTS), SHIFTJ(NGPTS), SHFTJ2(NGPTS)              
       COMMON /SHIFTR/ SHIFT, SHIFTJ, SHFTJ2                                    
                                                                                
C GZ = NOSE-HOOVER PHASE SPACE POINT                                            
C G2..G5 = 2ND TO 5TH DERIVATIVE OF POSITION OF NOSE-HOOVER PS POINT            
C XZ, PZ =POSITION & MOMENTUM OF N.H. PHASE SPACE POINT                         
C KB = time step counter for NoseHoover                                         
                                                                                
       READ (10,*) START                                                        
       IF (START.EQ.0) THEN                                                     
C Initialize Nose-Hoover Phase Space                                            
         CALL FCC                                                               
         DO 100 KB=1,2000                                                       
           CALL GEAR1Z(STPSIZ,5)                                                
           CALL SHFTIN(LAMBDA,CALJ(PZ))                                         
100        CONTINUE                                                             
       ELSE IF (START.EQ.MAXRUN) THEN                                           
         STOP 999                                                               
       ELSE                                                                     
C read in previous check point                                                  
         READ(11) GZ,G1,G2,G3,G4                                                
         READ(13) SHIFT,SHIFTJ,SHFTJ2                                           
         ENDIF                                                                  
                                                                                
       KB=0                                                                     
C      WHILE CPU TIME LEFT > TIME ORIGIN TIME                                   
20     CONTINUE                                                                 
         CALL CPUTME(TLIMIT,TUSED,TLEFT)                                        
         IF (TLEFT.LE.TSTEP) GOTO 1                                             
         IF (MOD(KB,100000).EQ.99999) THEN                                      
           OPEN(UNIT=14,FILE='RKS105.TTCFSAVE.DATA',STATUS='OLD')               
           CALL FLUSH(KB)                                                       
           CLOSE (14)                                                           
           ENDIF                                                                
         KB = KB+1                                                              
         CALL GEAR1Z(STPSIZ, 5)                                                 
         CALL OUTPUT                                                            
         GOTO 20                                                                
1      CONTINUE                                                                 
       ERR=IOPEN('UNIT=12,FILE=RKS105.TTCF.OUTPUT.DATA,STATUS=MOD')             
       IF (ERR.NE.0) GOTO 11                                                    
C Error code returns to this point                                              
13     CONTINUE                                                                 
       REWIND 10                                                                
       WRITE (10,*) START+1                                                     
       REWIND 11                                                                
       WRITE (11) GZ,G1,G2,G3,G4                                                
C Flush output buffer                                                           
       WRITE (12,ERR=11) KB,NGPTS                                               
       WRITE (12,ERR=11) SLTL0                                                  
       WRITE (12,ERR=11) SLTL0J                                                 
       WRITE (12,ERR=11) LTL0J2                                                 
       REWIND 13                                                                
       WRITE (13) SHIFT,SHIFTJ,SHFTJ2                                           
       STOP                                                                     
C Error control code - this just seemed the best place to put it                
11       J=0                                                                    
C        REPEAT UNTIL a file can be opened                                      
10       CONTINUE                                                               
          CLOSE (12,ERR=12)                                                     
12        ERR=IOPEN('UNIT=12,FILE=RKS105.OUTSAVE'//CHAR(J+ICHAR('0'))//         
     &      '.DATA,STATUS=NEW')                                                 
          J=J+1                                                                 
          PRINT *,'ERR=',ERR                                                    
          IF (ERR.NE.0) GOTO 10                                                 
         GOTO 13                                                                
       END                                                                      
                                                                                
       SUBROUTINE FLUSH(KB)                                                     
       INCLUDE (PARAM)                                                          
       INTEGER J,KB                                                             
       INTEGER START, ERR, IOPEN, TLIMIT, TUSED, TLEFT                          
       DOUBLE PRECISION GZ(DIMPS+2),XZ(D*NP),PZ(D*NP), CALJ                     
       DOUBLE PRECISION G2(DIMPS+2),G3(DIMPS+2),G4(DIMPS+2)                     
       DOUBLE PRECISION G1(DIMPS+2),LAMBDA                                      
       EQUIVALENCE (GZ,XZ),(GZ(D*NP+1),PZ),(GZ(DIMPS+2),LAMBDA)                 
       COMMON /GAMMAG/GZ                                                        
       COMMON /GSTORZ/G1,G2,G3,G4                                               
       INCLUDE (BUFFER)                                                         
       DOUBLE PRECISION SHIFT(NGPTS), SHIFTJ(NGPTS), SHFTJ2(NGPTS)              
       COMMON /SHIFTR/ SHIFT, SHIFTJ, SHFTJ2                                    
                                                                                
13     CONTINUE                                                                 
       REWIND 11                                                                
       WRITE (11) GZ,G1,G2,G3,G4                                                
C Flush output buffer                                                           
       WRITE (14,ERR=11) KB,NGPTS                                               
       WRITE (14,ERR=11) SLTL0                                                  
       WRITE (14,ERR=11) SLTL0J                                                 
       WRITE (14,ERR=11) LTL0J2                                                 
       REWIND 13                                                                
       WRITE (13) SHIFT,SHIFTJ,SHFTJ2                                           
       RETURN                                                                   
C Error control code - this just seemed the best place to put it                
11       J=0                                                                    
C        REPEAT UNTIL a file can be opened                                      
10       CONTINUE                                                               
          CLOSE (12,ERR=12)                                                     
12        ERR=IOPEN('UNIT=12,FILE=RKS105.OUTSAVE'//CHAR(J+ICHAR('0'))//         
     &      '.DATA,STATUS=NEW')                                                 
          J=J+1                                                                 
          PRINT *,'ERR=',ERR                                                    
          IF (ERR.NE.0) GOTO 10                                                 
         GOTO 13                                                                
       END



Russell Standish
Thu May 18 11:43:52 EST 1995