Energy of water

72 views
Skip to first unread message

sumit agrawal

unread,
Oct 6, 2021, 3:32:57 PMOct 6
to cp...@googlegroups.com
Hii,

Why is the energy of water coming around 17 a.u in cp2k. In gaussian it is around 76 a.u. I use the same functional (b3lyp) and same basis set (aug-cc-pvdz) in both the software. Am I doing something wrong in cp2k?

My input and output file for cp2k are attached herewith.


Please help me with this.

Thanks,
Sumit
water.inp
water.out

Krack Matthias (PSI)

unread,
Oct 6, 2021, 5:13:10 PMOct 6
to cp...@googlegroups.com

Hi Sumit

 

Your CP2K calculation employs pseudopotentials (the O 1s electrons are not included) whereas your Gaussian calculation is most likely an all-electron calculation. Just check the electron count in the outputs. CP2K and Gaussian should list 8 and 10 electrons, respectively.

 

Matthias

--
You received this message because you are subscribed to the Google Groups "cp2k" group.
To unsubscribe from this group and stop receiving emails from it, send an email to cp2k+uns...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/cp2k/CALSDoYY1E7w%2Bq9XVZfMPetT%2B%2B0Hhi8a3e6jRSFS9dBMKnZL0tg%40mail.gmail.com.

sumit agrawal

unread,
Oct 7, 2021, 6:17:18 AMOct 7
to cp...@googlegroups.com
Dear Matthias,

Yes, the number of electrons is different. But in cp2k I am using an all-electron method (GAPW). I thought it would take care of core electrons. But this is not the case. Can you suggest to me how I can introduce core elcetons in my cp2k calculation?

Thanks,
Sumit

Marcella Iannuzzi

unread,
Oct 7, 2021, 6:27:44 AMOct 7
to cp2k

Dear Sumit, 

The GAPW method can be used with pseudo potentials and with all electrons descriptions of the electronic structure. 
As Matthias writes, your input indicates the use of pseudo potentials in the KIND sections. 
These should be replaced with  ALL
Regards
Marcella

sumit agrawal

unread,
Oct 8, 2021, 3:11:53 AMOct 8
to cp...@googlegroups.com
Hii Marcella,

Thanks for the suggestion. Now I have it. But if I compare the energy in both the calculations( cp2k vs gaussian) the difference in energy is almost 0.38 a.u. In gaussian energy is roughly 76.4446 a.u and in cp2k it is 76.0632 a.u. Can this difference be reduced ? (Since i have used the same functional and basis set in both). I also check the frequency in both the calculations. There is roughly 20 cm-1 difference in the frequency. 

Since I just want to match my gaussian calculations with cp2k to know more about cp2k. Is it possible that gaussian calculations and cp2k calculation both are different  (if yes then what are those) and you can not get the same data in both the calculations?  Or my input having some problem. My input and outfiles are attached.

Please let me know.

Thanks,
Sumit

water.inp
water.out

Krack Matthias (PSI)

unread,
Oct 8, 2021, 5:08:32 AMOct 8
to cp...@googlegroups.com

Hi Sumit

 

Add ‘PERIODIC none’ to the &CELL section and add the section

  &POISSON

   PERIODIC NONE

   PSOLVER MT

  &END POISSON

Moreover, drop the sections &INTERACTION_POTENTIAL and &VDW_POTENTIAL and set EPSFIT 1.0E-8 and EPS_DEFAULT 1.0E-12 in the &QS section.

That should result in an energy closer to the Gaussian.

 

HTH

 

Matthias

 

hut...@chem.uzh.ch

unread,
Oct 8, 2021, 8:51:44 AMOct 8
to cp...@googlegroups.com
Hi

two additional possible sources of differences between Gaussian
and CP2K:

1) you need to use the 5d option in Gaussian in order to use spherical
d functions (not cartesian) as in CP2K

2) Gaussian has a slightly different definition of B3LYP, see
tests/QS/regtest-hybrid-2/H2O-hybrid-b3lyp-g03.inp
on how to reproduce that definition with CP2K.

regards

Juerg Hutter
--------------------------------------------------------------
Juerg Hutter Phone : ++41 44 635 4491
Institut für Chemie E-mail: hut...@chem.uzh.ch
Universität Zürich
Winterthurerstrasse 190
CH-8057 Zürich, Switzerland
---------------------------------------------------------------

-----cp...@googlegroups.com wrote: -----
To: cp...@googlegroups.com
From: "sumit agrawal"
Sent by: cp...@googlegroups.com
Date: 10/08/2021 09:12AM
Subject: Re: [CP2K:15999] Energy of water
To view this discussion on the web visit https://groups.google.com/d/msgid/cp2k/CALSDoYZ7j1uuB8hTkb%3DibP7AtGWomxGVZd9BgQ1CCL-tzAgs-A%40mail.gmail.com.


[attachment "water.inp" removed by Jürg Hutter/at/UZH]
[attachment "water.out" removed by Jürg Hutter/at/UZH]

Mozhdeh Mohammadpour

unread,
Oct 8, 2021, 9:01:29 AMOct 8
to cp2k
Thanks Mattias, I tested these suggestions, nothing changed for my system.

Krack Matthias (PSI)

unread,
Oct 8, 2021, 9:20:16 AMOct 8
to cp...@googlegroups.com

Hi Sumit

 

That’s weird, because I get with my changes already a much lower energy

-76.4075988 Ha

and after including Jürg’s suggestions for the &XC section I get a total energy of

-76.4446450 Ha

 

Matthias

 

sumit agrawal

unread,
Oct 8, 2021, 10:02:56 AMOct 8
to cp...@googlegroups.com
Hii Matthias,

I tested your suggestions. I got the water energy -76.407537347201654 a.u. Now I am closer to gaussian calculation. Thank you.

Dear Prof. Juerg Hutter, 

I changed my XC section but after that my cp2k energy is even less than gaussian energy. It was roughly -83 a.u. My input and out files are attached.

Please help me with this.

Thanks,
Sumit
water.out
water.inp

Krack Matthias (PSI)

unread,
Oct 8, 2021, 10:36:47 AMOct 8
to cp...@googlegroups.com

Hi

 

Use the &XC section as suggested by Jürg:

 

#These are the coefficients used for B3LYP using VWN3, this is the default Gaussian definition, but not recommended

    &XC

      &XC_FUNCTIONAL

       &LYP

         SCALE_C 0.81

       &END

       &BECKE88

         SCALE_X 0.72

       &END

       &VWN

         FUNCTIONAL_TYPE VWN3

         SCALE_C 0.19

       &END

       &XALPHA

         SCALE_X 0.08

       &END

      &END XC_FUNCTIONAL

      &HF

        &SCREENING

          EPS_SCHWARZ 1.0E-12

        &END

        &MEMORY

          MAX_MEMORY 100

        &END

        FRACTION 0.20

      &END

    &END XC

 

Use VWN3 instead of VWN5 and define the HF fraction of 0.2.

 

Matthias

 

sumit agrawal

unread,
Oct 8, 2021, 3:49:17 PMOct 8
to cp...@googlegroups.com
Dear Matthias,

Now it's fine. The energy and frequency in both the calculations are very close. Thanks for your all suggestions.

Few more questions: When do we have to use PBC and when not? I mean in single molecule calculation, can we use PBC ? And when I go to a higher cluster do I need to define PBC or can I do calculation without PBC? 

Second is, What should be the cell size? For single molecule calculations and some Lattice structure optimization? 

Third, How to play with the EPS_DEFAULT, CUTT_OFF, etc values in my cp2k calculations.

And last one is, When i do vibrational analysis, should i need to define & MOTION section ? 


Thanks,
Sumit


sumit agrawal

unread,
Oct 9, 2021, 4:01:59 AMOct 9
to cp...@googlegroups.com
Hii,

One more thing, when I checked the angle and OH distance in cp2k and gaussian, there was a difference in both the parameters. In cp2k the angle and OH distance are 102 and 0.78 A but in gaussian it was 104 and 0.96 A. Why is it? Energy and frequencies are almost the same in both the calculations. 

Thanks,
Sumit

Krack Matthias (PSI)

unread,
Oct 9, 2021, 6:10:51 AMOct 9
to cp...@googlegroups.com

Hi Sumit

 

Molecules and clusters in gas phase are usually treated as isolated systems without PBC independent of their size. To that end, you have to choose a sufficiently large simulation cell using “PERIODIC none” (add at least 8 A to the diameter of the molecule/cluster in each Cartesian direction). EPS_DEFAULT 1.0E-12 is a reasonable choice for most cases. The density cut-off depends on the employed (pseudo)potentials. The electron count printed at the end of the SCF run is a good indicator, if you need to further decrease EPS_DEFAULT or increase the cut-off. No &MOTION but a &VIBRATIONAL_ANALYSIS section is needed for RUN_TYPE VIBRATIONAL_ANALYSIS.

Have a look at the tutorials or search this forum for further information.

sumit agrawal

unread,
Oct 9, 2021, 6:48:19 AMOct 9
to cp...@googlegroups.com
Hii Matthias,

Thanks for all the suggestions. Very helpful.

Now I am very close to energy, frequency and OH distance. All the parameters are well matched with the gaussian calculation. But the angle is not matching. In cp2k it is coming 102 and in gaussian it is coming around 104. Can you comment on this ? Because after that I will do the same calculations for water dimer, trimer or higher clusters.



Krack Matthias (PSI)

unread,
Oct 9, 2021, 6:50:08 AMOct 9
to cp...@googlegroups.com

Hi Sumit

 

Then there is still something wrong with your input. I get for an isolated water molecule with B3LYP (Gaussian setup) and aug-cc-pVDZ basis set: r(OH) = 0.965 A, a(HOH) = 104.75°, Energy = -76.444645 a.u.

sumit agrawal

unread,
Oct 9, 2021, 8:16:11 AMOct 9
to cp...@googlegroups.com
Dear Matthias,

My input and out files are below. Can you point out the mistakes in the input file? My outputs are r(OH) = 0.9642 A, a(HOH) = 102.1°, Energy = -76.444583 a.u.

Thanks,
Sumit



h2o-pos-1.xyz
water.inp
water.out

Krack Matthias (PSI)

unread,
Oct 9, 2021, 11:46:47 AMOct 9
to cp...@googlegroups.com

Try to decrease the four convergence threshold values in the &MOTION section by at least an order of magnitude.

 

Matthias

 

sumit agrawal

unread,
Oct 16, 2021, 8:05:33 AMOct 16
to cp...@googlegroups.com
Dear matthias,

Thanks for all the suggestions. Now it's fine. I also tried with some other organic molecules and this is working. 

Now I did the water dimer in cp2k as well as in gaussian. I got the one negative frequency in cp2k and frequencies are not matching at all with the gaussian calculation. I think it is not considering H-Bond interaction. I used the same inputs as I used for water monomer except the coordinates. The comparison is below and my input file is attached. Can you please give some suggestions on " how to optimed clusters ( dimer, trimer etc..) in cp2k ".

H2o-dimer(cp2k)

H2O-dimer(gaussian)

Energy = -152.89098920

Energy = -152.8967947

Frequencies

Frequencies

-44

128

15

156

95

157

98

184

136

360

224

634

1619

1617

1622

1637

3787

3672

3792

3789

3898

3874

3900

3895


dimer.inp

Victor Volkov

unread,
Oct 16, 2021, 4:03:16 PMOct 16
to cp...@googlegroups.com
Dear Sumit:
the imaginary frequency ("negative" is a jargon) we must expect for a transition state.
If I remember, there are two types of transition states for a selected "reactive degree of freedom": with 1 or 2 imaginary frequencies.

If you observe imaginary frequencies, then, formally, you structure is not in a global minimum.

Using Gaussian I observe such often when I optimize large metallo-organic complexes.
Sometimes, in some polynuclear complexes, low frequency imaginary frequencies
should be present, naturally. 

If the case of insufficient optimization,
I increase the grid and the accuracy for numerical evaluation of two-electron integrals and their derivatives.
The default accuracy in g16  is 10-12. When insufficient increase accuracy to 10^-14 or higher.
The default FineGrid in g16 is set to 75 radial shells and 302 angular points per shell,
resulting in about 7000 points per atom.
UltraFine requests a is set to 99 and 590, respectively.
SuperFineGrid is pruned to 175,974 for first-row atoms and 250,974 for atoms in the second and later rows.

Analogous challenges are present in cp2k.

Often, results with imaginary frequencies are sufficient for a comparative tentative discussion
of the nature normal modes, the frequency of which are about 10-20 times larger than the negative one.

This means that if you observe an imaginary frequency 40cm-1
(what is a very poor performance for cp2k)
then the nature of the computed normal modes above 800 cm-1 is alright.

However, such results are for a tentative comparison:
overall - the results you got for water are trash, not publishable.

I may be wrong,
but this is what I think, for the moment.
Victor









sumit agrawal

unread,
Oct 18, 2021, 10:59:29 PMOct 18
to cp...@googlegroups.com
Dear Victor,

I understand that my structure is not a global minimum. I am quite comfortable in gaussian calculations. I am trying to learn cp2k and comparing the gaussian data to cp2k calculations. For the water monomer I got almost the same energy, frequencies and structure data. But when I tried a water dimer with the same input that I used for the water monomer, it gives very different results in cp2k. I just want to know what I have done wrong in my input file (which I had attached in my previous email).  Can you please suggest me something on this?

Thanks,
Sumit

Reply all
Reply to author
Forward
0 new messages