Why does nwchem give different results from Molpro does?
Clicked A Few Times
4:44:19 PM PDT - Mon, Jul 28th 2014
I am dealing with CCSD(T) calculation for some molecules.
For some reason, I ran with both Nwchem and Molpro codes for a molecule.

While I got so different results for the same molecule with the same geometry.

Could some one tell me the reason?
By the way, I am using molpro 2012 and nwchem 6.3.

Actually, I didn't see the difference for what I did before.

Thanks a lot.

Zongtang

Gets Around
5:50:43 PM PDT - Mon, Jul 28th 2014
How badly do they diverge? Can you provide input and output files for a pair of divergent calculations? I do not have access to Molpro, but I have been in the unpleasant position of trying to harmonize an underspecified method between two different quantum chemistry packages.

Some things to consider:

-Are both packages using the exact same basis set definitions? Specify basis set definitions inline with your input file or check your basis set libraries carefully to be sure that e.g. "6-311+G*" really means what you expect in both calculations.

-Are any orbitals frozen in the calcuations? If so, are they the same orbitals? Gaussian applies orbital freezing schemes by default for post-HF calculations while NWChem does not; I don't know what Molpro does.

-Are convergence criteria sufficiently similar for both calculations? This can be a very hard question to answer a priori, as different programs specify convergence criteria in different ways and not just with different numbers. It's a question best investigated after you have eliminated more obvious culprits.

Clicked A Few Times
10:54:11 PM PDT - Mon, Jul 28th 2014
Thank you so much. Actually I ran the molpro with the same size of calculation of the same molecule to check the T1 diagnostic. The convergence criterias are similar.

I worry about that because the difference is larger than 0.1 hatree.
For one of them, the HF energy is the same and the CCSD correlation energy is very different.
so either one of calculations should be wrong.

I am thinking about the reasons such as the state of the molecule or the multi-reference character of the molecule.

Would that cause the difference?

Gets Around
2:41:20 AM PDT - Tue, Jul 29th 2014
I would guess that different programs become more sensitive to implementation choices under "difficult" conditions, meaning anything other than neutral closed shell systems near equilibrium geometries that are well-described by a single reference. A 0.1 hartree difference is indeed large. It is good to hear that you have at least one calculation where the HF energies are the same, because that clears up worries about basis sets, at least for that calculation. I would encourage you to look at the calculations where HF energies are the same more closely, to figure out where results begin to significantly diverge between programs. I would also wonder about calculations where the HF energies are not close. That will be a problem if you want to continue using both programs for a single line of research.

If you want to try to compare Molpro and NWChem results against a third program, Psi4 is pretty capable and well-documented: https://github.com/psi4/psi4public

I don't know if that well help resolve your puzzle or only deepen it, though.

Forum Vet
9:52:37 AM PDT - Tue, Jul 29th 2014
Zontang
Could you please post the input file you were using and the platform were you were running NWChem?
Thanks, Edo

Clicked A Few Times
12:08:23 AM PDT - Wed, Jul 30th 2014
I think I found the problem for one molecule.

For the nwchem job, the RHF didn't converge. (the molpro job may be not correct either)

I got the error message in the scf step.

Here, would you please have a look to help me fix this issue ?

==================
The input :

echo
start mno4h.1ap.cs.ccsdt.ad
permanent_dir /home/****
memory stack 1000 mb heap 200 mb global 2100 mb noverify
title "hmno4.ccsdt.ad CCSD(T)/aug-cc-pVTZ(-PP)//B3LYP/cc-pVDZ(-PP)"

charge 0

geometry
symmetry c1
 MN   -0.080460    0.007636    0.000000

 O      1.665752    0.086143    0.000000

 H      2.077332   -0.793971    0.000000

 O     -0.557994   -0.749674    1.277449

 O     -0.557994    1.488591    0.000000

 O     -0.557994   -0.749674   -1.277449
end

basis spherical
H library aug-cc-pVDZ

O library aug-cc-pVDZ



 Mn library aug-cc-pVDZ-PP  #omitted here,  for calculation, i put correct basis sets with ecp

end

set atomscf:z 1.0 -2 7.0
set atomscf:tags_z H O Mn

scf
  vectors output mno4h.1ap.cs.ccsdt.ad.movecs

  rhf

  singlet

  maxiter 200

end

task scf

ccsd
  maxiter 60

  thresh 1e-6

  freeze atomic
end
task ccsd(t)

====================
erro message:

             iter       energy          gnorm     gmax       time

            ----- ------------------- --------- --------- --------

                1     -394.7142739629  1.49D+01  4.54D+00      3.1

 Setting level-shift to  11.11 to force positive preconditioner

                2     -402.1822934857  2.64D+00  4.43D-01      3.4

                3     -402.4827505038  1.97D+00  3.83D-01      3.5

                4     -402.5946928243  8.77D-01  1.48D-01      3.7

                5     -402.6440042851  3.86D-01  6.85D-02      3.8

                6     -402.6667435139  2.21D-01  6.35D-02      4.2

                7     -402.6701957197  2.12D-02  6.03D-03      4.4

 ga_iter_lsolve: convergence stagnant ... aborting solve

                8     -402.6710875676  3.90D-02  1.07D-02      5.2

                9     -402.6712445617  1.51D-03  5.10D-04      5.7

 ga_iter_lsolve: convergence stagnant ... aborting solve



Disabled NR: increased maxiter to 210



               10     -402.6712449429  2.11D-03  8.37D-04      5.8

               11     -402.6712461508  2.14D-03  6.32D-04      5.9

               12     -402.6712468975  2.69D-03  8.66D-04      6.0

               13     -402.6712475049  2.32D-03  6.75D-04      6.1

               14     -402.6712480031  1.89D-03  3.77D-04      6.2

               15     -402.6712482699  1.31D-03  3.43D-04      6.3

               16     -402.6712484015  1.08D-03  2.57D-04      6.4

               17     -402.6712484964  7.35D-04  1.74D-04      6.5

               18     -402.6712485453  5.92D-04  1.23D-04      6.5

               19     -402.6712485693  5.14D-04  1.20D-04      6.6

               20     -402.6712485929  4.52D-04  1.04D-04      6.7

               21     -402.6712486132  4.36D-04  9.50D-05      6.7

               22     -402.6712486325  4.75D-04  1.38D-04      6.8

               23     -402.6712486526  4.77D-04  1.04D-04      6.8

               24     -402.6712486689  3.00D-04  1.07D-04      6.9

               25     -402.6712486808  2.89D-04  9.44D-05      7.0

               26     -402.6712486864  2.14D-04  5.86D-05      7.1

               27     -402.6712486913  1.68D-04  4.05D-05      7.2

               28     -402.6712486938  1.39D-04  3.12D-05      7.3

               29     -402.6712486953  1.33D-04  3.44D-05      7.4

               30     -402.6712486972  9.62D-05  1.85D-05      7.5


============

I tried to play around to get rid of symmetry or increase the iteration steps, or set atomscf, et al, the problem was still there.


Thanks a lot.

Zongtang

Forum Vet
10:37:27 AM PDT - Wed, Jul 30th 2014
basis set
Zongtang
Could you please report the basis set you are using? The following is not very useful.

Mn library aug-cc-pVDZ-PP #omitted here, for calculation, i put correct basis sets with ecp

Clicked A Few Times
11:10:20 AM PDT - Wed, Jul 30th 2014
ok, sure

the mn basis set I am using as


  1. cc-pVDZ-PP
Mn s
     2080.37   0.000159

     304.649   0.000994

     22.2505   0.051511

     7.56933   -0.393048

     1.56197   0.740590

     0.654311   0.457874

     0.107776   0.016591

     0.039409   -0.002586
Mn s
     2080.37   -0.000038

     304.649   -0.000277

     22.2505   -0.012210

     7.56933   0.101296

     1.56197   -0.238569

     0.654311   -0.258056

     0.107776   0.577934

     0.039409   0.565458
Mn s
     2080.37   0.000099

     304.649   0.000494

     22.2505   0.030238

     7.56933   -0.242022

     1.56197   0.945993

     0.654311   -0.407703

     0.107776   -1.658454

     0.039409   1.634148
Mn s
     0.039409   1.0
Mn p
     54.1365   0.003585

     11.5345   -0.068823

     2.91993   0.345476

     1.31207   0.506638

     0.565154   0.262497

     0.153142   0.021029

     0.048714   -0.002835
Mn p
     54.1365   -0.000772

     11.5345   0.016640

     2.91993   -0.096406

     1.31207   -0.166698

     0.565154   -0.014675

     0.153142   0.545351

     0.048714   0.564209
Mn p
     54.1365   -0.000978

     11.5345   0.021319

     2.91993   -0.121381

     1.31207   -0.222495

     0.565154   -0.016634

     0.153142   0.799219

     0.048714   0.296449
Mn p
     0.048714   1.0
Mn d
     34.4465   0.021926

     10.7969   0.112236

     3.83655   0.293563

     1.39838   0.427063

     0.476480   0.377381

     0.138120   0.145135
Mn d
     34.4465   -0.022540

     10.7969   -0.118668

     3.83655   -0.309567

     1.39838   -0.302288

     0.476480   0.306841

     0.138120   0.714822
Mn d
     0.138120   1.0
Mn f
     1.3394   1.0
  1. aug-cc-pVDZ-PP
Mn s
     0.0144   1.0
Mn p
     0.0155   1.0
Mn d
     0.0400   1.0
Mn f
     0.3976   1.0
end

ecp
Mn nelec 10
Mn S
2 19.77859586351460 208.36277759981999
2 8.36429231183498 40.67136312517160
Mn P
2 19.75517256031370 42.82221279051267
2 19.58138625787880 74.14325787868533
2 8.09003570055466 8.57859438555583
2 8.31459920850757 20.38995721451193
Mn D
2 24.19571159002280 -9.37867882820444
2 24.53647306721080 -14.14683917959722
2 7.95149778468298 -0.29670061844163
2 7.93161326885935 -0.44295145020000
Mn F
2 12.19160291115040 -2.56699301425514
2 14.77200557195830 -5.28909520512491

end

Forum Vet
11:25:23 AM PDT - Wed, Jul 30th 2014
Zongtang
Here is the input that will give you a converged SCF by first using the DIIS convergence out of the DFT module


echo
start mno4h.1ap.cs.ccsdt.ad
permanent_dir /home/****
memory stack 1000 mb heap 200 mb global 2100 mb noverify
title "hmno4.ccsdt.ad CCSD(T)/aug-cc-pVTZ(-PP)//B3LYP/cc-pVDZ(-PP)"

charge 0

geometry
symmetry c1
 MN   -0.080460    0.007636    0.000000
 O      1.665752    0.086143    0.000000
 H      2.077332   -0.793971    0.000000
 O     -0.557994   -0.749674    1.277449
 O     -0.557994    1.488591    0.000000
 O     -0.557994   -0.749674   -1.277449
end

basis spherical
H library aug-cc-pVDZ
O library aug-cc-pVDZ
# cc-pVDZ-PP
Mn s
      2080.37   0.000159
      304.649   0.000994
      22.2505   0.051511
      7.56933   -0.393048
      1.56197   0.740590
      0.654311   0.457874
      0.107776   0.016591
      0.039409   -0.002586
Mn s
      2080.37   -0.000038
      304.649   -0.000277
      22.2505   -0.012210
      7.56933   0.101296
      1.56197   -0.238569
      0.654311   -0.258056
      0.107776   0.577934
      0.039409   0.565458
Mn s
      2080.37   0.000099
      304.649   0.000494
      22.2505   0.030238
      7.56933   -0.242022
      1.56197   0.945993
      0.654311   -0.407703
      0.107776   -1.658454
      0.039409   1.634148
Mn s
      0.039409   1.0
Mn p
      54.1365   0.003585
      11.5345   -0.068823
      2.91993   0.345476
      1.31207   0.506638
      0.565154   0.262497
      0.153142   0.021029
      0.048714   -0.002835
Mn p
      54.1365   -0.000772
      11.5345   0.016640
      2.91993   -0.096406
      1.31207   -0.166698
      0.565154   -0.014675
      0.153142   0.545351
      0.048714   0.564209
Mn p
      54.1365   -0.000978
      11.5345   0.021319
      2.91993   -0.121381
      1.31207   -0.222495
      0.565154   -0.016634
      0.153142   0.799219
      0.048714   0.296449
Mn p
      0.048714   1.0
Mn d
      34.4465   0.021926
      10.7969   0.112236
      3.83655   0.293563
      1.39838   0.427063
      0.476480   0.377381
      0.138120   0.145135
Mn d
      34.4465   -0.022540
      10.7969   -0.118668
      3.83655   -0.309567
      1.39838   -0.302288
      0.476480   0.306841
      0.138120   0.714822
Mn d
      0.138120   1.0
Mn f
      1.3394   1.0
# aug-cc-pVDZ-PP
Mn s
      0.0144   1.0
Mn p
      0.0155   1.0
Mn d
      0.0400   1.0
Mn f
      0.3976   1.0
end

ecp
Mn nelec 10
Mn S
2    19.77859586351460   208.36277759981999
2     8.36429231183498    40.67136312517160
Mn P
2    19.75517256031370    42.82221279051267
2    19.58138625787880    74.14325787868533
2     8.09003570055466     8.57859438555583
2     8.31459920850757    20.38995721451193
Mn D
2    24.19571159002280    -9.37867882820444
2    24.53647306721080   -14.14683917959722
2     7.95149778468298    -0.29670061844163
2     7.93161326885935    -0.44295145020000
Mn F
2    12.19160291115040    -2.56699301425514
2    14.77200557195830    -5.28909520512491
end

dft
xc hfexch
convergence energy 1d-10
tolerances acccoul 12
vectors input atomic output dft.mos
maxiter 45
end
set quickguess t
task dft ignore

scf
tol2e 1d-12
vectors input dft.mos output hf.mos
end

ccsd
  maxiter 60
  thresh 1e-6
  freeze atomic
end

task ccsd(t)

Clicked A Few Times
2:08:48 PM PDT - Thu, Jul 31st 2014
Finally we found out the issue.

For our CCSD(T) calculations with ECP basis sets for metal with NWChem, the default "Freeze atomic" didn't give right number of core orbitals. This would cause the huge difference of ccsd(t) number.

Also a couple of linear dependencies were removed due to the default threshold making the different HF value. I tried tight threshold and fixed that problem.

Thanks for all your guys' responses.


Forum >> NWChem's corner >> General Topics