hdiff output

r32686/potential.f90 2017-06-01 23:30:10.680422835 +0100 r32685/potential.f90 2017-06-01 23:30:10.976426685 +0100
907:          CSMGPSAVE=CSMGP907:          CSMGPSAVE=CSMGP
908:          CSMGP=CSMGUIDEGP908:          CSMGP=CSMGUIDEGP
909:          CSMGPINDEXSAVE=CSMGPINDEX909:          CSMGPINDEXSAVE=CSMGPINDEX
910:          CSMGPINDEX=CSMGUIDEGPINDEX910:          CSMGPINDEX=CSMGUIDEGPINDEX
911:          PTGP(1:3, 1:3, 1:2*CSMGPINDEX)=PTGPGUIDE(1:3, 1:3, 1:2*CSMGPINDEX)911:          PTGP(1:3, 1:3, 1:2*CSMGPINDEX)=PTGPGUIDE(1:3, 1:3, 1:2*CSMGPINDEX)
912:       END IF912:       END IF
913:       IF (PERMDIST) THEN913:       IF (PERMDIST) THEN
914:          DO J1=1, CSMGPINDEX914:          DO J1=1, CSMGPINDEX
915:             XTEMP(1:3*NATOMS)=X(1:3*NATOMS)915:             XTEMP(1:3*NATOMS)=X(1:3*NATOMS)
916:             CALL CSMROT(XTEMP, DUMMY, 1, J1)916:             CALL CSMROT(XTEMP, DUMMY, 1, J1)
917:             CALL MINPERMDIST(XTEMP, DUMMY, NATOMS, DEBUG, BOXLX, BOXLY, BOXLZ, PERIODIC, TWOD, EREAL, DIST2, RIGID, RMAT)917:             CALL ALIGN_DECIDE(XTEMP, DUMMY, NATOMS, DEBUG, BOXLX, BOXLY, BOXLZ, PERIODIC, TWOD, EREAL, DIST2, RIGID, RMAT)
918:             CALL CSMROT(DUMMY, XTEMP,-1, J1) ! need to rotate the permuted rotated images back to the reference orientation918:             CALL CSMROT(DUMMY, XTEMP,-1, J1) ! need to rotate the permuted rotated images back to the reference orientation
919:             CSMIMAGES(1+3*NATOMS*(J1-1):3*NATOMS*J1)=XTEMP(1:3*NATOMS)919:             CSMIMAGES(1+3*NATOMS*(J1-1):3*NATOMS*J1)=XTEMP(1:3*NATOMS)
920:          END DO920:          END DO
921:       ELSE921:       ELSE
922:          DO J1=1, CSMGPINDEX922:          DO J1=1, CSMGPINDEX
923:             CSMIMAGES(1+3*NATOMS*(J1-1):3*NATOMS*J1)=X(1:3*NATOMS)923:             CSMIMAGES(1+3*NATOMS*(J1-1):3*NATOMS*J1)=X(1:3*NATOMS)
924:          END DO924:          END DO
925:       END IF925:       END IF
926:       CALL CSMMIN(X, EREAL, CSMRMS, CSMIT)926:       CALL CSMMIN(X, EREAL, CSMRMS, CSMIT)
927: 927: 
958:          END DO958:          END DO
959:          CSMAV(1:3*NATOMS)=CSMAV(1:3*NATOMS)/CSMGPINDEX959:          CSMAV(1:3*NATOMS)=CSMAV(1:3*NATOMS)/CSMGPINDEX
960: !960: !
961: !  Check the CSM for the averaged structure. It should be zero if this structure has the961: !  Check the CSM for the averaged structure. It should be zero if this structure has the
962: !  right point group. Need to reset CSMIMAGES and CSMNORM temporarily.962: !  right point group. Need to reset CSMIMAGES and CSMNORM temporarily.
963: !963: !
964:          IF (PERMDIST) THEN964:          IF (PERMDIST) THEN
965:             DO J1=1, CSMGPINDEX965:             DO J1=1, CSMGPINDEX
966:                XTEMP(1:3*NATOMS)=CSMAV(1:3*NATOMS)966:                XTEMP(1:3*NATOMS)=CSMAV(1:3*NATOMS)
967:                CALL CSMROT(XTEMP, DUMMY, 1, J1)967:                CALL CSMROT(XTEMP, DUMMY, 1, J1)
968:                CALL MINPERMDIST(XTEMP, DUMMY, NATOMS, DEBUG, BOXLX, BOXLY, BOXLZ, PERIODIC, TWOD, DUMMY2, DIST2, RIGID, RMAT)968:                CALL ALIGN_DECIDE(XTEMP, DUMMY, NATOMS, DEBUG, BOXLX, BOXLY, BOXLZ, PERIODIC, TWOD, DUMMY2, DIST2, RIGID, RMAT)
969:                CALL CSMROT(DUMMY, XTEMP,-1, J1) ! need to rotate the permuted rotated images back to the reference orientation969:                CALL CSMROT(DUMMY, XTEMP,-1, J1) ! need to rotate the permuted rotated images back to the reference orientation
970:                CSMIMAGES(1+3*NATOMS*(J1-1):3*NATOMS*J1)=XTEMP(1:3*NATOMS)970:                CSMIMAGES(1+3*NATOMS*(J1-1):3*NATOMS*J1)=XTEMP(1:3*NATOMS)
971:             END DO971:             END DO
972:          ELSE972:          ELSE
973:             DO J1=1, CSMGPINDEX973:             DO J1=1, CSMGPINDEX
974:                CSMIMAGES(1+3*NATOMS*(J1-1):3*NATOMS*J1)=CSMAV(1:3*NATOMS)974:                CSMIMAGES(1+3*NATOMS*(J1-1):3*NATOMS*J1)=CSMAV(1:3*NATOMS)
975:             END DO975:             END DO
976:          END IF976:          END IF
977:          AA(1)=0.0D0; AA(2)=0.0D0; AA(3)=6.283185307D0 ! should give an identity matrix977:          AA(1)=0.0D0; AA(2)=0.0D0; AA(3)=6.283185307D0 ! should give an identity matrix
978:          SAVECSMPMAT(1:3, 1:3)=CSMPMAT(1:3, 1:3)978:          SAVECSMPMAT(1:3, 1:3)=CSMPMAT(1:3, 1:3)
990:       END IF ! DEBUG990:       END IF ! DEBUG
991:       IF (CSMDOGUIDET) THEN ! undo guiding changes991:       IF (CSMDOGUIDET) THEN ! undo guiding changes
992:          CSMGP=CSMGPSAVE992:          CSMGP=CSMGPSAVE
993:          CSMGPINDEX=CSMGPINDEXSAVE993:          CSMGPINDEX=CSMGPINDEXSAVE
994:          PTGP(1:3, 1:3, 1:2*CSMGPINDEX)=PTGPSAVE(1:3, 1:3, 1:2*CSMGPINDEX)994:          PTGP(1:3, 1:3, 1:2*CSMGPINDEX)=PTGPSAVE(1:3, 1:3, 1:2*CSMGPINDEX)
995:          CSMNORM=CSMNORMSAVE995:          CSMNORM=CSMNORMSAVE
996:       END IF996:       END IF
997: 997: 
998:    ELSE IF (PERMOPT .OR. PERMINVOPT .OR. DISTOPT) THEN998:    ELSE IF (PERMOPT .OR. PERMINVOPT .OR. DISTOPT) THEN
999: !  EREAL is the distance in this case999: !  EREAL is the distance in this case
1000:       CALL MINPERMDIST(FIN, X, NATOMS, DEBUG, BOXLX, BOXLY, BOXLZ, PERIODIC, TWOD, EREAL, DIST2, RIGID, RMAT)1000:       CALL ALIGN_DECIDE(FIN, X, NATOMS, DEBUG, BOXLX, BOXLY, BOXLZ, PERIODIC, TWOD, EREAL, DIST2, RIGID, RMAT)
1001: 1001: 
1002:    ELSE IF (MKTRAPT) THEN1002:    ELSE IF (MKTRAPT) THEN
1003:       CALL MKTRAP(NATOMS, X, GRAD, EREAL)1003:       CALL MKTRAP(NATOMS, X, GRAD, EREAL)
1004:    ELSE IF (BLJCLUSTER) THEN1004:    ELSE IF (BLJCLUSTER) THEN
1005:       CALL RAD(X, GRAD, EREAL, GRADT)1005:       CALL RAD(X, GRAD, EREAL, GRADT)
1006:       CALL LJPSHIFTBINC(X, GRAD, EREAL, GRADT, SECT)1006:       CALL LJPSHIFTBINC(X, GRAD, EREAL, GRADT, SECT)
1007: 1007: 
1008:    ELSE IF (BLJCLUSTER_NOCUT) THEN1008:    ELSE IF (BLJCLUSTER_NOCUT) THEN
1009: ! ds656> Binary LJ without cutoff1009: ! ds656> Binary LJ without cutoff
1010:       CALL RAD(X, GRAD, EREAL, GRADT)1010:       CALL RAD(X, GRAD, EREAL, GRADT)
1373:          END IF1373:          END IF
1374:          IF (AXTELL) CALL AXT(NATOMS, X, GRAD, EREAL, GRADT, ZSTAR)1374:          IF (AXTELL) CALL AXT(NATOMS, X, GRAD, EREAL, GRADT, ZSTAR)
1375:       END IF1375:       END IF
1376: 1376: 
1377:        IF ( RIGIDINIT .AND. (.NOT. ATOMRIGIDCOORDT)) THEN1377:        IF ( RIGIDINIT .AND. (.NOT. ATOMRIGIDCOORDT)) THEN
1378: 1378: 
1379:            X(DEGFREEDOMS+1:3*NATOMS)=0.0D01379:            X(DEGFREEDOMS+1:3*NATOMS)=0.0D0
1380:            X(1:DEGFREEDOMS)=XRIGIDCOORDS(1:DEGFREEDOMS)1380:            X(1:DEGFREEDOMS)=XRIGIDCOORDS(1:DEGFREEDOMS)
1381: 1381: 
1382:            IF(GRADT) THEN1382:            IF(GRADT) THEN
1383:               GRADATOMS(:) = GRAD(1:3*NATOMS)  ! This is needed for AACONVERGENCET1383:               GRADATOMS(:) = GRAD(:)  ! This is needed for AACONVERGENCET
1384:               CALL TRANSFORMGRAD(GRAD, XRIGIDCOORDS, XRIGIDGRAD)1384:               CALL TRANSFORMGRAD(GRAD, XRIGIDCOORDS, XRIGIDGRAD)
1385:               GRAD(DEGFREEDOMS+1:3*NATOMS)=0.0D01385:               GRAD(DEGFREEDOMS+1:3*NATOMS)=0.0D0
1386:               GRAD(1:DEGFREEDOMS)=XRIGIDGRAD(1:DEGFREEDOMS)1386:               GRAD(1:DEGFREEDOMS)=XRIGIDGRAD(1:DEGFREEDOMS)
1387:            ENDIF1387:            ENDIF
1388: 1388: 
1389:            IF(SECT) THEN1389:            IF(SECT) THEN
1390:               CALL TRANSFORMHESSIAN(HESS, GRAD, XRIGIDCOORDS, XRIGIDHESS, RBAANORMALMODET)1390:               CALL TRANSFORMHESSIAN(HESS, GRAD, XRIGIDCOORDS, XRIGIDHESS, RBAANORMALMODET)
1391:               HESS(DEGFREEDOMS+1:3*NATOMS,:)=0.0D01391:               HESS(DEGFREEDOMS+1:3*NATOMS,:)=0.0D0
1392:               HESS(:, DEGFREEDOMS+1:3*NATOMS)=0.0D01392:               HESS(:, DEGFREEDOMS+1:3*NATOMS)=0.0D0
1393:               HESS(1:DEGFREEDOMS, 1:DEGFREEDOMS)=XRIGIDHESS(1:DEGFREEDOMS, 1:DEGFREEDOMS)1393:               HESS(1:DEGFREEDOMS, 1:DEGFREEDOMS)=XRIGIDHESS(1:DEGFREEDOMS, 1:DEGFREEDOMS)
1394:            ENDIF1394:            ENDIF
1395: 1395: 
1396:        END IF1396:        END IF
1397: 1397: 
 1398:        write(*,*) "Energy ", EREAL
 1399:        write(*,*) GRAD(:3*NATOMS)
 1400: 
1398:    END IF1401:    END IF
1399: !1402: !
1400: !  --------------- End of possible potentials - now add fields if required ------------------------------1403: !  --------------- End of possible potentials - now add fields if required ------------------------------
1401: !1404: !
1402: 1405: 
1403:    IF (RIGIDCONTOURT) THEN1406:    IF (RIGIDCONTOURT) THEN
1404: !       sf344> sanity check, since this keyword works only for a pair of rigid particles1407: !       sf344> sanity check, since this keyword works only for a pair of rigid particles
1405:       IF (NATOMS.NE.4) THEN1408:       IF (NATOMS.NE.4) THEN
1406:          WRITE(MYUNIT,'(A)') 'potential> ERROR - RIGIDCONTOUR works only for a PAIR of rigid bodies, exiting'1409:          WRITE(MYUNIT,'(A)') 'potential> ERROR - RIGIDCONTOUR works only for a PAIR of rigid bodies, exiting'
1407:          STOP1410:          STOP
1552: !            GRAD(J2-2)=GRAD(J2-2)-XG1555: !            GRAD(J2-2)=GRAD(J2-2)-XG
1553: !            GRAD(J2-1)=GRAD(J2-1)-YG1556: !            GRAD(J2-1)=GRAD(J2-1)-YG
1554: !            GRAD(J2)=  GRAD(J2)-ZG1557: !            GRAD(J2)=  GRAD(J2)-ZG
1555: !         END DO1558: !         END DO
1556: !      END IF1559: !      END IF
1557: 1560: 
1558:       IF (CSMT .AND. (.NOT.SYMMETRIZECSM)) THEN1561:       IF (CSMT .AND. (.NOT.SYMMETRIZECSM)) THEN
1559:          DUMMY2=0.0D01562:          DUMMY2=0.0D0
1560:          RMS=0.0D01563:          RMS=0.0D0
1561:       ELSE IF (RIGIDINIT .AND. (.NOT. ATOMRIGIDCOORDT)) THEN1564:       ELSE IF (RIGIDINIT .AND. (.NOT. ATOMRIGIDCOORDT)) THEN
 1565:          write(*,*) "Rigid body RMS"
1562:          DUMMY2=SUM(GRAD(1:DEGFREEDOMS)**2)1566:          DUMMY2=SUM(GRAD(1:DEGFREEDOMS)**2)
1563:          RMS=MAX(SQRT(DUMMY2/DEGFREEDOMS), 1.0D-100)1567:          RMS=MAX(SQRT(DUMMY2/DEGFREEDOMS), 1.0D-100)
 1568:          write(*,*) "RMS:", RMS
1564:          IF ((RMS < 5.0D0 * BQMAX) .AND. AACONVERGENCET .AND. (.NOT. COMPRESSRIGIDT) .AND. (.NOT. COMPON)) THEN1569:          IF ((RMS < 5.0D0 * BQMAX) .AND. AACONVERGENCET .AND. (.NOT. COMPRESSRIGIDT) .AND. (.NOT. COMPON)) THEN
 1570:             write(*,*) "Using AACONVERGENCE"
1565:             CALL AACONVERGENCE (GRADATOMS, XRIGIDCOORDS, XRIGIDGRAD, RMS)1571:             CALL AACONVERGENCE (GRADATOMS, XRIGIDCOORDS, XRIGIDGRAD, RMS)
 1572:             write(*,*) "RMS", RMS
1566:          END IF1573:          END IF
1567:       ELSE IF (MLP3T.OR.MLPB3T.OR.MLQT) THEN1574:       ELSE IF (MLP3T.OR.MLPB3T.OR.MLQT) THEN
1568:          DUMMY2=SUM(GRAD(1:NATOMS)**2)1575:          DUMMY2=SUM(GRAD(1:NATOMS)**2)
1569:          RMS=MAX(DSQRT(DUMMY2/(NATOMS)), 1.0D-100)1576:          RMS=MAX(DSQRT(DUMMY2/(NATOMS)), 1.0D-100)
1570:       ELSE IF (.NOT.THOMSONT) THEN1577:       ELSE IF (.NOT.THOMSONT) THEN
1571:          DUMMY2=SUM(GRAD(1:3*NATOMS)**2)1578:          DUMMY2=SUM(GRAD(1:3*NATOMS)**2)
1572:          RMS=MAX(DSQRT(DUMMY2/(3*NATOMS)), 1.0D-100)1579:          RMS=MAX(DSQRT(DUMMY2/(3*NATOMS)), 1.0D-100)
 1580:          write(*,*) "Using regular RMS", RMS
1573:       ELSE1581:       ELSE
1574:          DUMMY2=SUM(GRAD(1:2*NATOMS)**2)1582:          DUMMY2=SUM(GRAD(1:2*NATOMS)**2)
1575:          RMS=MAX(DSQRT(DUMMY2/(2*NATOMS)), 1.0D-100)1583:          RMS=MAX(DSQRT(DUMMY2/(2*NATOMS)), 1.0D-100)
1576:       END IF1584:       END IF
1577:       IF(DEBUG.AND.(RMS.NE.RMS)) THEN1585:       IF(DEBUG.AND.(RMS.NE.RMS)) THEN
1578:          WRITE(MYUNIT,'(A)' ) 'potential> WARNING - RMS force is NaN - if using AMBER igb=1, can be due to negative Born radii'1586:          WRITE(MYUNIT,'(A)' ) 'potential> WARNING - RMS force is NaN - if using AMBER igb=1, can be due to negative Born radii'
1579:       END IF1587:       END IF
1580: 1588: 
1581: !   IF ((CSMRMS < GUIDECUT).AND.(CSMDOGUIDET)) THEN1589: !   IF ((CSMRMS < GUIDECUT).AND.(CSMDOGUIDET)) THEN
1582:       IF (CSMDOGUIDET) THEN1590:       IF (CSMDOGUIDET) THEN


legend
Lines Added 
Lines changed
 Lines Removed

hdiff - version: 2.1.0