grimme-lab / xtb

Semiempirical Extended Tight-Binding Program Package
https://xtb-docs.readthedocs.io/
GNU Lesser General Public License v3.0
556 stars 139 forks source link

Thermostatistical corrections for T=0 #66

Closed awvwgk closed 4 years ago

awvwgk commented 4 years ago

Describe the bug In case you try to calculate thermostatistical corrections at T=0K interesting things happen.

To Reproduce Tested with version 6.2.1 (36dbe8b). xtb will complain for good reasons

printf "\$thermo\ntemp=0.0\n\$end\n" > input
xtb -I input <geometry> --define
... skip ...
# WARNING! Your input has following faults:
#  - A temperature of 0.0 K sounds strange to me

Ignoring the warning and still running the program with this input will produce meaningless and wrong values.

$ printf "\$thermo\ntemp=0.0\n\$end\n" > input
$ xtb -I input <geometry> --ohess
... skip ...
       T/K    H(0)-H(T)+PV         H(T)/Eh          T*S/Eh         G(T)/Eh
   ------------------------------------------------------------------------
   ------------------------------------------------------------------------

         :::::::::::::::::::::::::::::::::::::::::::::::::::::
         ::                  THERMODYNAMIC                  ::
         :::::::::::::::::::::::::::::::::::::::::::::::::::::
         :: total free energy         -42.153937359591 Eh   ::
         ::.................................................::
         :: total energy              -42.153937359591 Eh   ::
         :: zero point energy           0.182091755134 Eh   ::
         :: G(RRHO) w/o ZPVE           -0.182091755134 Eh   ::
         :: G(RRHO) contrib.            0.000000000000 Eh   ::
         :::::::::::::::::::::::::::::::::::::::::::::::::::::
... skip ...

Notice that no thermodynamic function have been evaluated, since there is no printout regarding the partition function or temperature, so the whole output is obviously useless and wrong here.

Using a very small number instead will give meaningful results (the one the user expected for T=0K):

$ printf "\$thermo\ntemp=1.0e-10\n\$end\n" > input
$ xtb -I input <geometry> --ohess
... skip ...
       T/K    H(0)-H(T)+PV         H(T)/Eh          T*S/Eh         G(T)/Eh
   ------------------------------------------------------------------------
      0.00    0.126676E-14    0.182087E+00   -0.289393E-13    0.182087E+00
   ------------------------------------------------------------------------

         :::::::::::::::::::::::::::::::::::::::::::::::::::::
         ::                  THERMODYNAMIC                  ::
         :::::::::::::::::::::::::::::::::::::::::::::::::::::
         :: total free energy         -41.971849822766 Eh   ::
         ::.................................................::
         :: total energy              -42.153937303642 Eh   ::
         :: zero point energy           0.182087480876 Eh   ::
         :: G(RRHO) w/o ZPVE            0.000000000000 Eh   ::
         :: G(RRHO) contrib.            0.182087480876 Eh   ::
         :::::::::::::::::::::::::::::::::::::::::::::::::::::
... skip ...

Run with a coffeine molecule, but this is not the point here.

Expected behaviour Maybe kill the calculation after a full hessian run, maybe skip the thermo part, maybe a special case for T=0K, maybe fix the user input by adding some noise.

Additional context Thanks to @fabothch for reporting and discussion.

fabothch commented 4 years ago

@awvwgk Thx for adding the issue!

awvwgk commented 4 years ago

Maybe I will just make it an input error if the set_thermo reader does not find at least one correct temperature. Would certainly solve this issue by deferring it.