Bug fixes for TB Stress Tensor (DFTB and xTB), two new regtests, (#392)

some regtests have been updated
This commit is contained in:
Juerg Hutter 2019-05-29 11:17:46 +02:00 committed by GitHub
parent 56af05a0b8
commit 89a19eabde
No known key found for this signature in database
GPG key ID: 4AEE18F83AFDEB23
12 changed files with 218 additions and 49 deletions

View file

@ -225,10 +225,12 @@ CONTAINS
fij(i) = -gmij*SUM(TRANSPOSE(pblock)*dsint(:, :, i))
END IF
END DO
CALL virial_pair_force(virial%pv_virial, -1._dp, fij, rij)
fi = 1.0_dp
IF (iatom == jatom) fi = 0.5_dp
CALL virial_pair_force(virial%pv_virial, fi, fij, rij)
IF (atprop%stress) THEN
CALL virial_pair_force(atprop%atstress(:, :, irow), -0.5_dp, fij, rij)
CALL virial_pair_force(atprop%atstress(:, :, icol), -0.5_dp, fij, rij)
CALL virial_pair_force(atprop%atstress(:, :, irow), fi*0.5_dp, fij, rij)
CALL virial_pair_force(atprop%atstress(:, :, icol), fi*0.5_dp, fij, rij)
END IF
END IF
END DO
@ -281,10 +283,12 @@ CONTAINS
fij(i) = fi
END DO
IF (use_virial) THEN
CALL virial_pair_force(virial%pv_virial, 1._dp, fij, rij)
fi = 1.0_dp
IF (iatom == jatom) fi = 0.5_dp
CALL virial_pair_force(virial%pv_virial, fi, fij, rij)
IF (atprop%stress) THEN
CALL virial_pair_force(atprop%atstress(:, :, iatom), 0.5_dp, fij, rij)
CALL virial_pair_force(atprop%atstress(:, :, jatom), 0.5_dp, fij, rij)
CALL virial_pair_force(atprop%atstress(:, :, iatom), fi*0.5_dp, fij, rij)
CALL virial_pair_force(atprop%atstress(:, :, jatom), fi*0.5_dp, fij, rij)
END IF
END IF

View file

@ -231,10 +231,12 @@ CONTAINS
END IF
END DO
IF (use_virial) THEN
CALL virial_pair_force(virial%pv_virial, 1._dp, fij, rij)
fi = 1.0_dp
IF (iatom == jatom) fi = 0.5_dp
CALL virial_pair_force(virial%pv_virial, fi, fij, rij)
IF (atprop%stress) THEN
CALL virial_pair_force(atprop%atstress(:, :, iatom), 0.5_dp, fij, rij)
CALL virial_pair_force(atprop%atstress(:, :, jatom), 0.5_dp, fij, rij)
CALL virial_pair_force(atprop%atstress(:, :, iatom), fi*0.5_dp, fij, rij)
CALL virial_pair_force(atprop%atstress(:, :, jatom), fi*0.5_dp, fij, rij)
END IF
END IF
END IF
@ -433,10 +435,12 @@ CONTAINS
fij(i) = -2.0_dp*gmij*SUM(TRANSPOSE(pblock)*dsint(:, :, i))
END IF
END DO
CALL virial_pair_force(virial%pv_virial, -1._dp, fij, rij)
fi = 1.0_dp
IF (iatom == jatom) fi = 0.5_dp
CALL virial_pair_force(virial%pv_virial, fi, fij, rij)
IF (atprop%stress) THEN
CALL virial_pair_force(atprop%atstress(:, :, iatom), -0.5_dp, fij, rij)
CALL virial_pair_force(atprop%atstress(:, :, jatom), -0.5_dp, fij, rij)
CALL virial_pair_force(atprop%atstress(:, :, iatom), 0.5_dp*fi, fij, rij)
CALL virial_pair_force(atprop%atstress(:, :, jatom), 0.5_dp*fi, fij, rij)
END IF
END IF
END DO
@ -492,10 +496,12 @@ CONTAINS
fij(i) = fi
END DO
IF (use_virial) THEN
CALL virial_pair_force(virial%pv_virial, 1._dp, fij, rij)
fi = 1.0_dp
IF (iatom == jatom) fi = 0.5_dp
CALL virial_pair_force(virial%pv_virial, fi, fij, rij)
IF (atprop%stress) THEN
CALL virial_pair_force(atprop%atstress(:, :, iatom), 0.5_dp, fij, rij)
CALL virial_pair_force(atprop%atstress(:, :, jatom), 0.5_dp, fij, rij)
CALL virial_pair_force(atprop%atstress(:, :, iatom), fi*0.5_dp, fij, rij)
CALL virial_pair_force(atprop%atstress(:, :, jatom), fi*0.5_dp, fij, rij)
END IF
END IF
END IF
@ -737,7 +743,7 @@ CONTAINS
dsblock = dsblock-dsblockm
dsblock = dsblock/(2.0_dp*ddr)
clist%acint(1:n1, 1:n2, i) = dsblock(1:n1, 1:n2)
clist%acint(1:n1, 1:n2, i) = -dsblock(1:n1, 1:n2)
ENDDO
DEALLOCATE (dsblock, dsblockm)
END IF

View file

@ -368,6 +368,9 @@ CONTAINS
END IF
ENDDO
IF (use_virial) THEN
!deb iatom=jatom f0=0.5*f0
IF (iatom == jatom) f0 = 0.5_dp*f0
!deb
CALL virial_pair_force(virial%pv_virial, -f0, force_ab, rij)
CALL virial_pair_force(virial%pv_virial, -f0, force_w, rij)
IF (atprop%stress) THEN
@ -413,10 +416,12 @@ CONTAINS
force(jkind)%repulsive(:, atom_b) = &
force(jkind)%repulsive(:, atom_b)+force_rr(:)
IF (use_virial) THEN
CALL virial_pair_force(virial%pv_virial, -1._dp, force_rr, rij)
f0 = -1.0_dp
IF (iatom == jatom) f0 = -0.5_dp
CALL virial_pair_force(virial%pv_virial, f0, force_rr, rij)
IF (atprop%stress) THEN
CALL virial_pair_force(atprop%atstress(:, :, iatom), -0.5_dp, force_rr, rij)
CALL virial_pair_force(atprop%atstress(:, :, jatom), -0.5_dp, force_rr, rij)
CALL virial_pair_force(atprop%atstress(:, :, iatom), f0*0.5_dp, force_rr, rij)
CALL virial_pair_force(atprop%atstress(:, :, jatom), f0*0.5_dp, force_rr, rij)
END IF
END IF
END IF

View file

@ -486,10 +486,10 @@ CONTAINS
DO j = 1, nj
la = laoa(i)+1
lb = laob(j)+1
gcij(i, j) = 0.5_dp*(gchrg(irow, la, 1)+gchrg(icol, lb, 1))
gcij(i, j) = 0.5_dp*(gchrg(iatom, la, 1)+gchrg(jatom, lb, 1))
END DO
END DO
gmij = 0.5_dp*(gmcharge(irow, 1)+gmcharge(icol, 1))
gmij = 0.5_dp*(gmcharge(iatom, 1)+gmcharge(jatom, 1))
icol = MAX(iatom, jatom)
irow = MIN(iatom, jatom)
NULLIFY (pblock)
@ -506,7 +506,7 @@ CONTAINS
f1 = -2.0_dp*SUM(TRANSPOSE(pblock)*dsint(:, :, i)*gcij)
f2 = -2.0_dp*gmij*SUM(TRANSPOSE(pblock)*dsint(:, :, i))
END IF
fij(i) = fij(i)+f1+f2
fij(i) = f1+f2
END DO
DEALLOCATE (gcij)
fi = 1.0_dp
@ -576,7 +576,7 @@ CONTAINS
IF (calculate_forces) THEN
atom_i = atom_of_kind(iatom)
atom_j = atom_of_kind(jatom)
IF (irow == jatom) THEN
IF (irow /= iatom) THEN
gmij = -gmij
gcij = -gcij
END IF

View file

@ -4,7 +4,7 @@
base1.inp 1 1.0E-12 -4.08611632954554
base2.inp 1 1.0E-12 -4.08611630192434
base3.inp 1 1.0E-12 -4.08519715577319
str1.inp 1 1.0E-12 -32.67545045484233
str1.inp 1 1.0E-12 -4.08612900850098
str2.inp 1 1.0E-12 -32.67544573650541
str3.inp 1 1.0E-12 -4.08521550626912
#EOF

View file

@ -1,4 +1,4 @@
@SET NREP 2
@SET NREP 1
&FORCE_EVAL
STRESS_TENSOR ANALYTICAL
&DFT

View file

@ -47,11 +47,11 @@ h2o_disp2.inp 1 4e-09 -
h2o_disp3.inp 1 3e-09 -130.76542980403056
# DFTB3
h2o-6.inp 1 4e-09 -131.25337238175840
h2o-7.inp 1 1e-08 -131.26536731088510
h2o-7.inp 1 1e-08 -131.26536755467674
# Properties stress
h2o-atprop1.inp 1 4e-09 -130.77763673429416
h2o-atprop2.inp 1 8e-09 -131.68099551186123
h2o-atprop3.inp 1 8e-09 -131.33465952358338
h2o-atprop2.inp 1 8e-09 -131.68099550861825
h2o-atprop3.inp 1 8e-09 -131.33465894123017
#
c_kp1.inp 1 1.0E-14 -13.76936383448730
c_kp2.inp 1 7e-14 -13.76786043591555

View file

@ -10,11 +10,13 @@ ch2o_smear.inp 1 1.0E-12 -7.19650182
tmol.inp 1 1.0E-12 -41.90853885813245
h2.inp 1 1.0E-12 -1.03458111733093
h2o-md.inp 1 1.0E-12 -185.15731412756392
h2o_str.inp 1 1.0E-12 -5.76530887622206
h2o-atprop.inp 1 4.0E-08 -185.17649920422642
h2o_str.inp 1 1.0E-12 -46.12964016870387
h2o-atprop.inp 1 1.0E-10 -185.17668389962344
h2o-atprop0.inp 1 1.0E-12 -187.45499307089042
si_geo.inp 1 1.0E-12 -14.55678348854957
si_kp.inp 1 1.0E-12 -14.75431646498735
h2o_dimer.inp 1 1.0E-12 -11.54506384130837
AdeThyvdW.inp 1 1.0E-12 -58.89501875211175
ice.inp 1 1.0E-10 -370.501223771625689
ice2.inp 1 1.0E-10 -46.30633640965694
#EOF

View file

@ -1,3 +1,4 @@
@SET NREP 2
&FORCE_EVAL
STRESS_TENSOR ANALYTICAL
&DFT
@ -5,38 +6,34 @@
METHOD xTB
&xTB
DO_EWALD T
COULOMB_INTERACTION T
TB3_INTERACTION T
&PARAMETER
DISPERSION_PARAMETER_FILE dftd3.dat
&END PARAMETER
&END
EPS_DEFAULT 1.E-14
&END QS
&SCF
SCF_GUESS MOPAC
&MIXING
METHOD DIRECT_P_MIXING
ALPHA 0.10
ALPHA 0.25
&END
MAX_SCF 2000
MAX_SCF 200
EPS_SCF 1.e-10
&END SCF
&KPOINTS
SCHEME NONE
#SCHEME GAMMA
#SCHEME MONKHORST-PACK 1 1 1
SCHEME NONE
&END KPOINTS
&POISSON
&EWALD
ALPHA 2.0
EWALD_TYPE SPME
GMAX 50
O_SPLINE 6
ALPHA 3.0
&END EWALD
&END POISSON
&END DFT
&SUBSYS
&TOPOLOGY
MULTIPLE_UNIT_CELL ${NREP} ${NREP} ${NREP}
&END
&CELL
ABC 7.0 7.0 7.0
MULTIPLE_UNIT_CELL ${NREP} ${NREP} ${NREP}
&END CELL
&COORD
O -4.583 5.333 1.560
@ -48,14 +45,13 @@
&GLOBAL
PROJECT h2o
RUN_TYPE DEBUG
# RUN_TYPE ENERGY_FORCE
PRINT_LEVEL LOW
&END GLOBAL
&DEBUG
DEBUG_FORCES F
DEBUG_STRESS_TENSOR T
DX 0.001
EPS_NO_ERROR_CHECK 0.000001
STOP_ON_MISMATCH F
DX 0.0002
EPS_NO_ERROR_CHECK 0.0001
STOP_ON_MISMATCH T
&END DEBUG

View file

@ -62,5 +62,5 @@
DEBUG_STRESS_TENSOR T
DX 0.002
EPS_NO_ERROR_CHECK 0.000001
STOP_ON_MISMATCH F
STOP_ON_MISMATCH T
&END DEBUG

View file

@ -0,0 +1,85 @@
@SET NREP 2
&GLOBAL
PRINT_LEVEL low
RUN_TYPE DEBUG
&END GLOBAL
&DEBUG
DEBUG_FORCES no
DEBUG_STRESS_TENSOR yes
EPS_NO_ERROR_CHECK 1.0E-4
DX 2.E-4
STOP_ON_MISMATCH yes
&END DEBUG
&FORCE_EVAL
METHOD Quickstep
STRESS_TENSOR analytical
&PRINT
&STRESS_TENSOR ON
&END
&END
&DFT
&KPOINTS
SCHEME NONE
&END
&QS
EPS_DEFAULT 1.0E-12
EXTRAPOLATION ASPC
EXTRAPOLATION_ORDER 3
METHOD xTB
&xTB
DO_EWALD T
&END xTB
&END QS
&SCF
EPS_SCF 1.0E-6
MAX_SCF 10
SCF_GUESS MOPAC
&OT on
MINIMIZER DIIS
PRECONDITIONER FULL_SINGLE_INVERSE
&END OT
&OUTER_SCF on
EPS_SCF 1.0E-6
MAX_SCF 10
&END OUTER_SCF
&END SCF
&END DFT
&SUBSYS
&TOPOLOGY
MULTIPLE_UNIT_CELL ${NREP} ${NREP} ${NREP}
&END
&CELL
ABC 6.358 6.358 6.358
MULTIPLE_UNIT_CELL ${NREP} ${NREP} ${NREP}
&END CELL
&COORD
O 3.4613786876 5.6912049480 0.4543727289 H2O
H 2.8183735096 5.2012524560 -0.0760381983 H2O
H 2.9114406133 6.4013882133 0.9299789968 H2O
O 0.2726839173 5.8131769039 3.6524426114 H2O
H -0.3305134720 5.2542646859 3.0350424519 H2O
H 0.8608072875 6.3736221304 3.0591507235 H2O
O 0.3733209371 2.7092162634 0.4720541405 H2O
H -0.1034698773 2.0668534417 -0.1013800786 H2O
H -0.2883205804 3.2313214887 0.9861387459 H2O
O 5.1679371793 4.2946104941 2.0104883056 H2O
H 4.5349860728 3.7358082901 2.5469828636 H2O
H 4.6439911293 4.8921754874 1.4197090565 H2O
O 1.9804166604 1.0500358045 2.0298628573 H2O
H 2.5408437547 1.7234516156 2.4956106850 H2O
H 1.3203099947 1.5937364881 1.5079432584 H2O
O 1.8865307172 4.1962271939 5.1674455966 H2O
H 1.1978491058 4.7123566273 4.6845538753 H2O
H 1.4552978348 3.6225372508 5.7958754719 H2O
O 3.5550539674 2.8142815783 3.4327698835 H2O
H 2.9552982643 3.3909531502 4.0032656392 H2O
H 3.9895494155 2.1557772081 4.0448159552 H2O
O 5.1227218630 1.1112395839 5.1795427775 H2O
H 5.7513589539 0.5921262315 4.6175072767 H2O
H 4.6567143202 0.4540581626 5.7714857593 H2O
&END COORD
&END SUBSYS
&END FORCE_EVAL

View file

@ -0,0 +1,71 @@
&GLOBAL
PRINT_LEVEL low
RUN_TYPE DEBUG
&END GLOBAL
&DEBUG
DEBUG_FORCES no
DEBUG_STRESS_TENSOR yes
EPS_NO_ERROR_CHECK 1.0E-4
DX 1.E-3
STOP_ON_MISMATCH yes
&END DEBUG
&FORCE_EVAL
METHOD Quickstep
STRESS_TENSOR analytical
&PRINT
&STRESS_TENSOR ON
&END
&END
&DFT
&KPOINTS
SCHEME GAMMA
&END
&QS
EPS_DEFAULT 1.0E-14
EXTRAPOLATION ASPC
EXTRAPOLATION_ORDER 3
METHOD xTB
&xTB
DO_EWALD T
&END xTB
&END QS
&SCF
EPS_SCF 1.0E-9
MAX_SCF 100
&END SCF
&END DFT
&SUBSYS
&CELL
ABC 6.358 6.358 6.358
&END CELL
&COORD
O 3.4613786876 5.6912049480 0.4543727289 H2O
H 2.8183735096 5.2012524560 -0.0760381983 H2O
H 2.9114406133 6.4013882133 0.9299789968 H2O
O 0.2726839173 5.8131769039 3.6524426114 H2O
H -0.3305134720 5.2542646859 3.0350424519 H2O
H 0.8608072875 6.3736221304 3.0591507235 H2O
O 0.3733209371 2.7092162634 0.4720541405 H2O
H -0.1034698773 2.0668534417 -0.1013800786 H2O
H -0.2883205804 3.2313214887 0.9861387459 H2O
O 5.1679371793 4.2946104941 2.0104883056 H2O
H 4.5349860728 3.7358082901 2.5469828636 H2O
H 4.6439911293 4.8921754874 1.4197090565 H2O
O 1.9804166604 1.0500358045 2.0298628573 H2O
H 2.5408437547 1.7234516156 2.4956106850 H2O
H 1.3203099947 1.5937364881 1.5079432584 H2O
O 1.8865307172 4.1962271939 5.1674455966 H2O
H 1.1978491058 4.7123566273 4.6845538753 H2O
H 1.4552978348 3.6225372508 5.7958754719 H2O
O 3.5550539674 2.8142815783 3.4327698835 H2O
H 2.9552982643 3.3909531502 4.0032656392 H2O
H 3.9895494155 2.1557772081 4.0448159552 H2O
O 5.1227218630 1.1112395839 5.1795427775 H2O
H 5.7513589539 0.5921262315 4.6175072767 H2O
H 4.6567143202 0.4540581626 5.7714857593 H2O
&END COORD
&END SUBSYS
&END FORCE_EVAL