hdiff output

r31544/Cv.HSA.f90 2016-11-29 12:30:11.395716626 +0000 r31543/Cv.HSA.f90 2016-11-29 12:30:11.883723136 +0000
  1: !  1: svn: E195012: Unable to find repository location for 'svn+ssh://svn.ch.private.cam.ac.uk/groups/wales/trunk/UTILS/GMIN/Cv.HSA.f90' in revision 31543
  2: !  Calculate Cv curve from weights and visit statistics 
  3: !  for replicas at different temperatures REPT(J1). 
  4: !  Minimise chi^2 statistic to extract best probabilities. 
  5: ! 
  6: PROGRAM CALCCV 
  7: IMPLICIT NONE 
  8: DOUBLE PRECISION TMIN, TMAX, TINT, Z0, Z1, Z2, TEMPERATURE, DUMMY, DELTA, ONEMEXP, PSUM, DP, PSUM2 
  9: DOUBLE PRECISION, ALLOCATABLE :: ENERGY(:), LNFRQS(:), MEANE(:), ZSAVE(:) 
 10: INTEGER, ALLOCATABLE :: HORDER(:) 
 11: INTEGER J1, NTEMP, J2, KAPPA, NMIN 
 12:  
 13: OPEN(UNIT=10,FILE='Cv.data',STATUS='OLD') 
 14: READ(10,*) TMIN, TMAX, NTEMP, KAPPA, NMIN 
 15: CLOSE(10) 
 16: ALLOCATE(MEANE(NTEMP),ZSAVE(NTEMP)) 
 17:  
 18: PRINT '(3(A,I4))','minima=',NMIN 
 19: PRINT '(A,2G15.5,A,I8)','minimum and maximum T for Cv calculation=',TMIN,TMAX,' number of Cv data points=',NTEMP 
 20: PRINT '(A,I8)','number of degrees of freedom=',KAPPA 
 21:  
 22: ALLOCATE(ENERGY(NMIN),HORDER(NMIN),LNFRQS(NMIN)) 
 23:  
 24: OPEN(UNIT=10,FILE='min.data',STATUS='OLD') 
 25: DO J1=1,NMIN 
 26:    READ(10,*) ENERGY(J1),LNFRQS(J1),HORDER(J1) 
 27: ENDDO 
 28: CLOSE(10) 
 29:  
 30: TINT=(TMAX-TMIN)/(NTEMP-1) 
 31: ! 
 32: !  Calculate Z0, Z1 and Z2 over the required T range. Omit factors of (kT/h)^kappa, 
 33: !  which cancel. 
 34: ! 
 35: OPEN(UNIT=1,FILE='Cv.out',STATUS='UNKNOWN') 
 36: DO J1=1,NTEMP 
 37:    Z0=0.0D0 
 38:    Z1=0.0D0 
 39:    Z2=0.0D0 
 40:    TEMPERATURE=TMIN+(J1-1)*TINT 
 41:    DO J2=1,NMIN 
 42:       DUMMY=EXP(-(ENERGY(J2)-ENERGY(1))/TEMPERATURE+(LNFRQS(1)-LNFRQS(J2))/2.0D0)/HORDER(J2) 
 43:       Z0=Z0+DUMMY 
 44:       Z1=Z1+DUMMY*(ENERGY(J2)-ENERGY(1)) 
 45:       Z2=Z2+DUMMY*(ENERGY(J2)-ENERGY(1))**2 
 46:    ENDDO 
 47:    MEANE(J1)=Z1/Z0 
 48:    ZSAVE(J1)=Z0 
 49:    WRITE(1,'(7G20.10)') TEMPERATURE, Z0, Z1, Z2, & 
 50:   &                     KAPPA - (Z1/(Z0*TEMPERATURE))**2 + Z2/(Z0*TEMPERATURE**2), MEANE(J1), ZSAVE(J1) 
 51: ENDDO 
 52: CLOSE(1) 
 53:  
 54: OPEN(UNIT=1,FILE='prob.out',STATUS='UNKNOWN') 
 55: DO J1=1,NTEMP 
 56:    TEMPERATURE=TMIN+(J1-1)*TINT 
 57:    PSUM=0.0D0 
 58:    PSUM2=0.0D0 
 59:    DO J2=1,NMIN 
 60:       DUMMY=EXP(-(ENERGY(J2)-ENERGY(1))/TEMPERATURE+(LNFRQS(1)-LNFRQS(J2))/2.0D0)/HORDER(J2) 
 61:       DP=DUMMY*(ENERGY(J2)-ENERGY(1)-MEANE(J1))/(ZSAVE(J1)*TEMPERATURE**2) 
 62:       PSUM=PSUM+DP*(ENERGY(J2)-ENERGY(1)-MEANE(J1)) 
 63: !     PRINT '(A,2I6,4G20.10)','J1,J2,T,ZSAVE,DUMMY,DP=',J1,J2,TEMPERATURE,ZSAVE(J1),DUMMY,DP 
 64:       PSUM2=PSUM2+(ZSAVE(J1)/DUMMY)*(TEMPERATURE*DP)**2 
 65:    ENDDO 
 66:    WRITE(1,'(5G20.10)') TEMPERATURE, KAPPA+PSUM, KAPPA+PSUM2 
 67: ENDDO 
 68: CLOSE(1) 
 69:  
 70: END PROGRAM CALCCV 
 71:  


legend
Lines Added 
Lines changed
 Lines Removed

hdiff - version: 2.1.0