hdiff output

r30208/dftb.f 2016-03-30 23:30:07.517389417 +0100 r30207/dftb.f 2016-03-30 23:30:07.909394574 +0100
 12: C   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the 12: C   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 13: C   GNU General Public License for more details. 13: C   GNU General Public License for more details.
 14: C 14: C
 15: C   You should have received a copy of the GNU General Public License 15: C   You should have received a copy of the GNU General Public License
 16: C   along with this program; if not, write to the Free Software 16: C   along with this program; if not, write to the Free Software
 17: C   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA 17: C   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
 18: C 18: C
 19: C 19: C
 20: C  Tiffany s TB routine for C. 20: C  Tiffany s TB routine for C.
 21: C 21: C
 22:       SUBROUTINE DFTB(N,XS,deriv1st,ENERGY,GTEST,XMULREAL,BOXLX,BOXLY,BOXLZ) 22:       SUBROUTINE DFTB(N,XS,deriv1st,ENERGY,GTEST,XMULREAL)
 23:       USE KEY 23:       USE KEY
 24:       IMPLICIT NONE 24:       IMPLICIT NONE
 25:       LOGICAL GTEST 25:       LOGICAL GTEST
 26:       INTEGER I, J, K, AI, AJ, NDIM,NMAX, I2, J2, K2, LWORK, INFO, J1, N, NDIAG 26:       INTEGER I, J, K, AI, AJ, NDIM,NMAX, I2, J2, K2, LWORK, INFO, J1, N, NDIAG
 27:       DOUBLE PRECISION DIST(N,N),XMULREAL,WORK(320*N),DPRAND,S(4*N,4*N),Sorig(4*N,4*N), 27:       DOUBLE PRECISION DIST,XMULREAL,WORK(320*N),DPRAND,S(4*N,4*N),Sorig(4*N,4*N),
 28:      &                 ENERGY,XS(3*N),HB(4*N,4*N),HA(4*N,4*N),VA(4*N,4*N),VB(4*N,4*N) 28:      &                 ENERGY,XS(3*N),HB(4*N,4*N),HA(4*N,4*N),VA(4*N,4*N),VB(4*N,4*N)
 29:  29: 
 30:       DOUBLE PRECISION R(N,N), ABSTOL,DIRCOS(N,N,3),DIAG(4*N),deriv1st(3*N),DIAGA(4*N),DIAGB(4*N) 30:       DOUBLE PRECISION R(N,N), ABSTOL,DIRCOS(N,N,3),DIAG(4*N),deriv1st(3*N),DIAGA(4*N),DIAGB(4*N)
 31:       DOUBLE PRECISION BOXLX,BOXLY,BOXLZ 
 32:       INTEGER MX(N,N), MY(N,N), MZ(N,N) 
 33:  31: 
 34:       DOUBLE PRECISION Ssssig,Sspsig,Sppsig,Spppi,EPSsac,EPSsbc,EPSpac,EPSpbc,Hsssiga,Hspsiga, 32:       DOUBLE PRECISION Ssssig,Sspsig,Sppsig,Spppi,EPSsac,EPSsbc,EPSpac,EPSpbc,Hsssiga,Hspsiga,
 35:      &        Hppsiga,Hpppia,Hsssigb,Hspsigb,Hppsigb,Hpppib,VECTOR(3) 33:      &        Hppsiga,Hpppia,Hsssigb,Hspsigb,Hppsigb,Hpppib
 36:  34: 
 37:       PARAMETER(EPSsac=-0.519289D0, EPSpac=-0.216472D0, EPSsbc=-0.434382D0, EPSpbc=-0.216472D0) 35:       PARAMETER(EPSsac=-0.519289D0, EPSpac=-0.216472D0, EPSsbc=-0.434382D0, EPSpbc=-0.216472D0)
 38:     36:    
 39:       INTEGER IWORK(20*N), IFAIL(4*N)  37:       INTEGER IWORK(20*N), IFAIL(4*N) 
 40:  38: 
 41: C     COMMON /RDIST/ R, DIRCOS, DIAGA, DIAGB 39: C     COMMON /RDIST/ R, DIRCOS, DIAGA, DIAGB
 42:  40: 
 43:       DOUBLE PRECISION UREP,REP 41:       DOUBLE PRECISION UREP,REP
 44:       LOGICAL DIAGTEST 42:       LOGICAL DIAGTEST
 45:       COMMON /FAIL/ DIAGTEST 43:       COMMON /FAIL/ DIAGTEST
 46:  44: 
 47:       LWORK=320*N 45:       LWORK=320*N
 48:       XMUL=INT(XMULREAL) 46:       XMUL=INT(XMULREAL)
 49:  47: 
 50: C   Calculate distance matrix 48: C   Calculate distance matrix
 51: C 
 52: C************ START PERIODIC BOUNDARY CONDITIONS ************************** 
 53: C  Deal with atoms leaving the box: 
 54: C 
 55:       IF (BULKT) THEN 
 56:          IF ((.NOT.FIXIMAGE).AND.(.NOT.NORESET)) THEN 
 57: C           PRINT*,'Points before:' 
 58: C           WRITE(*,'(3F15.7)') (XS(J1),J1=1,3*N) 
 59:             DO J1=1,N 
 60:                J2=3*(J1-1) 
 61:                XS(J2+1)=XS(J2+1) - BOXLX*DNINT(XS(J2+1)/BOXLX) 
 62:                XS(J2+2)=XS(J2+2) - BOXLY*DNINT(XS(J2+2)/BOXLY) 
 63:                XS(J2+3)=XS(J2+3) - BOXLZ*DNINT(XS(J2+3)/BOXLZ) 
 64:             ENDDO 
 65: C           PRINT*,'Points after:' 
 66: C           WRITE(*,'(3F15.7)') (XS(J1),J1=1,3*N) 
 67:          ENDIF 
 68: C 
 69: C  Calculation of connecting vectors; to implement the periodic 
 70: C  boundary conditions, the shortest vector between two atoms is used: 
 71: C 
 72:          DO I=1,N 
 73:             I2=3*(I-1) 
 74:             R(I,I)=0.0D0 
 75:             DO J=I+1,N 
 76:                J2=3*(J-1) 
 77:                VECTOR(1)=XS(J2+1)-XS(I2+1) 
 78:                VECTOR(2)=XS(J2+2)-XS(I2+2) 
 79:                VECTOR(3)=XS(J2+3)-XS(I2+3) 
 80:                IF (.NOT.FIXIMAGE) THEN 
 81:                   MX(J,I)=NINT(VECTOR(1)/BOXLX) 
 82:                   MY(J,I)=NINT(VECTOR(2)/BOXLY) 
 83:                   MZ(J,I)=NINT(VECTOR(3)/BOXLZ) 
 84:                ENDIF  
 85:                VECTOR(1)= VECTOR(1) - BOXLX * MX(J,I) 
 86:                VECTOR(2)= VECTOR(2) - BOXLY * MY(J,I) 
 87:                VECTOR(3)= VECTOR(3) - BOXLZ * MZ(J,I) 
 88:                DIST(I,J)=VECTOR(1)**2+VECTOR(2)**2+VECTOR(3)**2 
 89:                R(I,J)=SQRT(DIST(I,J)) 
 90:                R(J,I)=R(I,J) 
 91:                DO K=1,3 
 92:                   DIRCOS(I,J,K)=VECTOR(K)/R(I,J) 
 93:                   DIRCOS(J,I,K)=-DIRCOS(I,J,K) 
 94:                ENDDO 
 95: C              WRITE(*,'(2I4,3F15.7)') (I,J,(DIRCOS(I,J,K),K=1,3)) 
 96:             ENDDO 
 97:             DO K=1,3 
 98:                DIRCOS(I,I,K)=0.0D0 
 99:             ENDDO 
100:          ENDDO 
101: C************ END PERIODIC BOUNDARY CONDITIONS **************************** 
102:       ELSE 
103:  49: 
104:       DO I=1,N 50:       DO I=1,N
105:        I2=3*(I-1) 51:        I2=3*(I-1)
106:        R(I,I)=0.0D0 52:        R(I,I)=0.0D0
107:           DO J=I+1,N 53:           DO J=I+1,N
108:             J2=3*(J-1) 54:             J2=3*(J-1)
109:             DIST(I,J)=(XS(I2+1)-XS(J2+1))**2+(XS(I2+2)-XS(J2+2))**2+(XS(I2+3)-XS(J2+3))**2 55:             DIST=(XS(I2+1)-XS(J2+1))**2+(XS(I2+2)-XS(J2+2))**2+(XS(I2+3)-XS(J2+3))**2
110:              R(I,J)=SQRT(DIST(I,J)) 56:              R(I,J)=SQRT(DIST)
111:              R(J,I)=R(I,J) 57:              R(J,I)=R(I,J)
112:  58: 
113: C   Now the direction cosines 59: C   Now the direction cosines
114: C   Which incidentally are the projections of a unit vector onto  60: C   Which incidentally are the projections of a unit vector onto 
115: C   the axes and NOT the cosines. 61: C   the axes and NOT the cosines.
116: C   Note we have not calculated DIRCOS(I,I) 62: C   Note we have not calculated DIRCOS(I,I)
117:  63: 
118:             DO K=1,3 64:             DO K=1,3
119:                DIRCOS(I,J,K)=(XS(J2+K)-XS(I2+K))/R(I,J) 65:                DIRCOS(I,J,K)=(XS(J2+K)-XS(I2+K))/R(I,J)
120:                DIRCOS(J,I,K)=-DIRCOS(I,J,K) 66:                DIRCOS(J,I,K)=-DIRCOS(I,J,K)
121:                DIRCOS(I,I,K)=0.0D0 67:                DIRCOS(I,I,K)=0.0D0
122:             END DO 68:             END DO
123: C            PRINT*,'D1 IS',DIRCOS(I,J,1) 69: C            PRINT*,'D1 IS',DIRCOS(I,J,1)
124: C            PRINT*,'D2 IS',DIRCOS(I,J,2) 70: C            PRINT*,'D2 IS',DIRCOS(I,J,2)
125: C            PRINT*,'D3 IS',DIRCOS(I,J,3) 71: C            PRINT*,'D3 IS',DIRCOS(I,J,3)
126:  72: 
127:          END DO 73:          END DO
128:       END DO  74:       END DO 
129:       ENDIF 
130:  75: 
131:  76: 
132: C   Calculate Repulsive Free Energy 77: C   Calculate Repulsive Free Energy
133:  78: 
134:       UREP=0.0D0 79:       UREP=0.0D0
135:  80: 
136:       DO I=1,N 81:       DO I=1,N
137:          DO J=I+1,N 82:          DO J=I+1,N
138:  83: 
139: CC          PRINT*,'js and atnos are',I,J,IATNUM(I),IATNUM(J) 84: CC          PRINT*,'js and atnos are',I,J,IATNUM(I),IATNUM(J)
1845:          ENDIF1790:          ENDIF
1846: 1791: 
1847: 1792: 
1848: C       PRINT*,'in drepcc, rep is',rep1793: C       PRINT*,'in drepcc, rep is',rep
1849:       RETURN1794:       RETURN
1850:       END1795:       END
1851: C1796: C
1852: C**************************************************************************1797: C**************************************************************************
1853: C1798: C
1854: C  Subroutine diff forms the cartesian second derivative matrix using1799: C  Subroutine diff forms the cartesian second derivative matrix using
1855: C  numerical differentiation. Single-sided!1800: C  numerical differentiation.
1856: C1801: C
1857: C**************************************************************************1802: C**************************************************************************
1858: C1803: C
1859:       SUBROUTINE SECDFTB(N,X,VNEW,XMULREAL,BOXLX,BOXLY,BOXLZ)1804:       SUBROUTINE SECDFTB(N,X,VNEW,XMULREAL)
1860:       USE MODHESS1805:       USE MODHESS
1861:       IMPLICIT NONE1806:       IMPLICIT NONE
1862:       INTEGER N, J1, J21807:       INTEGER N, J1, J2
1863:       DOUBLE PRECISION X(3*N), 1808:       DOUBLE PRECISION X(3*N), 
1864:      1                 DIF, FG1(3*N),V1,VNEW(3*N),XMULREAL,BOXLX,BOXLY,BOXLZ1809:      1                 DIF, FG1(3*N),V1,VNEW(3*N),XMULREAL
1865: 1810: 
1866:       DIF=1.0D-61811:       DIF=1.0D-6
1867: 1812: 
1868:       DO J1=1,3*N1813:       DO J1=1,3*N
1869:          X(J1)=X(J1)+DIF1814:          X(J1)=X(J1)+DIF
1870:          CALL DFTB(N,X,FG1,V1,.TRUE.,XMULREAL,BOXLX,BOXLY,BOXLZ)1815:          CALL DFTB(N,X,FG1,V1,.TRUE.,XMULREAL)
1871:          X(J1)=X(J1)-DIF1816:          X(J1)=X(J1)-DIF
1872:          DO J2=J1,3*N1817:          DO J2=J1,3*N
1873:             HESS(J2,J1)=(FG1(J2)-VNEW(J2))/DIF1818:             HESS(J2,J1)=(FG1(J2)-VNEW(J2))/DIF
1874:             HESS(J1,J2)=HESS(J2,J1)1819:             HESS(J1,J2)=HESS(J2,J1)
1875:          ENDDO1820:          ENDDO
1876:       ENDDO1821:       ENDDO
1877: 1822: 
1878:       RETURN1823:       RETURN
1879:       END1824:       END
1880: 1825: 


r30208/potential.f 2016-03-30 23:30:07.717392046 +0100 r30207/potential.f 2016-03-30 23:30:08.109397207 +0100
3024:             WRITE(15,*) NATOMS3024:             WRITE(15,*) NATOMS
3025:             WRITE(15,*) "AFTER STEP",ITG03,"Energy=",ENERGY3025:             WRITE(15,*) "AFTER STEP",ITG03,"Energy=",ENERGY
3026:             OPEN(7,file="gaussian.name")3026:             OPEN(7,file="gaussian.name")
3027:             DO J1=1,NATOMS3027:             DO J1=1,NATOMS
3028:                READ(7,*) STRING3028:                READ(7,*) STRING
3029:                WRITE(15,"(A2,3F20.10)") ADJUSTL(TRIM(STRING)),COORDS(3*J1-2:3*J1)3029:                WRITE(15,"(A2,3F20.10)") ADJUSTL(TRIM(STRING)),COORDS(3*J1-2:3*J1)
3030:             ENDDO3030:             ENDDO
3031:             CLOSE(7)3031:             CLOSE(7)
3032:             CLOSE(15)3032:             CLOSE(15)
3033:             ELSE IF (DFTBT) THEN3033:             ELSE IF (DFTBT) THEN
3034:             CALL DFTB(NATOMS,COORDS,VNEW,ENERGY,GTEST,PARAM1,PARAM2,PARAM3,PARAM4)3034:             CALL DFTB(NATOMS,COORDS,VNEW,ENERGY,GTEST,PARAM1)
3035:             IF (SSTEST) CALL SECDFTB(NATOMS,COORDS,VNEW,PARAM1,PARAM2,PARAM3,PARAM4)3035:             IF (SSTEST) CALL SECDFTB(NATOMS,COORDS,VNEW,PARAM1)
3036:             IF (PTEST) THEN3036:             IF (PTEST) THEN
3037:                WRITE(*,10) ' potential> Energy for last cycle=',ENERGY,' eV'3037:                WRITE(*,10) ' potential> Energy for last cycle=',ENERGY,' eV'
3038:                WRITE(ESTRING,10) 'Energy for last cycle=',ENERGY,' eV'3038:                WRITE(ESTRING,10) 'Energy for last cycle=',ENERGY,' eV'
3039:             ENDIF3039:             ENDIF
3040:          ELSE IF (AMBER12T) THEN3040:          ELSE IF (AMBER12T) THEN
3041:             VNEW = 0.0D03041:             VNEW = 0.0D0
3042:             IF (RIGIDINIT .AND. (.NOT. ATOMRIGIDCOORDT)) THEN3042:             IF (RIGIDINIT .AND. (.NOT. ATOMRIGIDCOORDT)) THEN
3043:                XRIGIDCOORDS(1:DEGFREEDOMS) = COORDS(1:DEGFREEDOMS)3043:                XRIGIDCOORDS(1:DEGFREEDOMS) = COORDS(1:DEGFREEDOMS)
3044:                CALL TRANSFORMRIGIDTOC(1, NRIGIDBODY, COORDS, XRIGIDCOORDS)3044:                CALL TRANSFORMRIGIDTOC(1, NRIGIDBODY, COORDS, XRIGIDCOORDS)
3045:             ENDIF3045:             ENDIF


legend
Lines Added 
Lines changed
 Lines Removed

hdiff - version: 2.1.0