Can I use one code for geometry optimization and another code for a single point calculation on that optimized geometry? Can one publish results like that?

  • 3
    $\begingroup$ Is "code" in this context computer code? "code" is an uncountable noun in this context (and thus you can not say "a code" or "codes"). It is better to use another word and program is, I think, more fitting (e.g. "I used 7 different programs"). ("7. (programming, uncountable) Instructions for a computer, written in a programming language; the input of a translator, an interpreter or a browser, namely: source code, machine code, bytecode.") $\endgroup$ – Peter Mortensen Sep 4 '20 at 16:26
  • 1
    $\begingroup$ @PeterMortensen except, "code" can also be used as a countable noun. You also did not give any reference for it being an uncountable noun, you have just self-asserted it. I do not doubt that there are references that say it's an uncountable noun though, so please don't give us one :) $\endgroup$ – Nike Dattani Sep 4 '20 at 18:05
  • 1
    $\begingroup$ I think the discussion below is surprising. Especially the ones from solid state physics community. $\endgroup$ – Y. Zhai Sep 5 '20 at 3:33
  • 1
    $\begingroup$ @Y.Zhai I agree. I'm currently working on a big paper right now with more than 20 authors, on benchmarking GW methods. Every program gives different results for GW on the same molecule with the same geometry and even the same basis set. DFT and GW are not as simple as couple-cluster for example, but still there's no reason why you can't calculate a single-point energy at an experimental geometry or a geometry from a paper where the author's of the paper optimized it with a different code. It all depends on what the asker of this question really wants to do, but the question didn't give details $\endgroup$ – Nike Dattani Sep 5 '20 at 16:25
  • $\begingroup$ Depending on the field, it can be a very common practice. If I recall well, the early versions of Orca were encouraging to use other code (Turbomole, specifically) for geometry optimization, as Orca's main goal was to focus on properties and other software were just faster. $\endgroup$ – Greg Sep 6 '20 at 16:46


I do this all the time, and often use several codes, since each code has its own strengths and weaknesses.

I used 7 different codes in this paper where all I was doing was the ionization energy of one carbon atom (this means calculating the energy of a neutral carbon atom, and the energy of a singly ionized carbon cation and then reporting the difference between these two energies).

  • MOLPRO was used to optimize the basis sets (because not many codes can optimize basis set exponents).
  • GAUSSIAN was used to optimize the k-type functions in the basis sets (GAUSSIAN cannot automatically optimize exponents, so we had to do the optimization manually, but GAUSSIAN can treat k-type functions whereas MOLPRO cannot).
  • MOLCAS was used to calculate the integrals, because even though GAUSSIAN can use k-type and l-type functions, it cannot print the integrals in FCIDUMP format but MOLCAS can; and the integrals needed to be in FCIDUMP format because it's the format that NECI reads.
  • NECI was used to do the FCIQMC calculations (the other programs either can't do FCIQMC or can't do it with as much control as we needed).
  • CFOUR was used to calculate the DBOC and X2C corrections because MOLPRO, GAUSSIAN, MOLCAS and NECI couldn't do these at the time (and still probably cannot do DBOC).
  • MRCC was used to calculate coupled cluster calculations for the X2C and DBOC corrections beyond CCSD(T), since the other programs cannot go beyond CCSD(T).
  • Psi3 (a precursor to Psi4) was used to optimize some of the functions of the basis sets (but this part was done before we started this project, and it was the subject of previous papers, so we did not need to use Psi3 or Psi4 directly in this paper).

It is a very common thing to calculate integrals in one program, and post-SCF calculations in another program. That is the purpose of the file: To calculate integrals in one program, print them in FCIDUMP format, which is the universal format that almost all programs can read.

Why is it that MOLPRO, PySCF, (Open)MOLCAS, ARROW, DALTON, DIRAC, GAMESS, DICE, BLOCK, CheMPS2, xacc, NECI, BAGEL, HANDE, Psi4 and probably others, can all read FCIDUMP format? It's because people want to be able to calculate integrals in one program and post-SCF calculations in a different program. This is also the reason why MRCC has interfaces with CFOUR, MOLPRO, Psi4, DIRAC, ORCA, COLUMBUS, AMBER: people want to be able to calculate integrals in one program and do high-order coupled cluster in a different program. There's several other examples.

As for geometry optimization: We actually do single-point energy calculations at geometries obtained by different programs all the time when we take a geometry from the literature and do a single-point energy calculation on it. We often even take the geometry from an X-ray crystallography experiment, then do a single-point energy calculation on it. Some programs can do certain types of single-point energy calculations (maybe they have a special DFT functional implemented, for example) but cannot do geometry optimization because the optimization algorithm is simply not implemented.

But be careful: Frozen-core CCSD(T) is MOLPRO, doesn't mean the same thing as frozen-core CCSD(T) in CFOUR, for example. What a method means in one software, is not necessarily what it means in another software. In addition to making sure you are clear in your paper about which basis set, SCF, post-SCF, and other methods are used, it is a good idea to also say which program was used (for both the geometry optimization, and for the single-point energy calculation).

  • $\begingroup$ Can you explain why it is justifiable to do single point calculations on experimental geometry? It may very well not be at the minimum of the total energy in the experimental geometry. I have seen some authors try to validate XC functionals by comparing the optimized geometry to XRay crystallography data. I have done many phonon calculations using experimental lattice parameters but optimized atomic positions. But there are arguments to be made that this is a reasonable setup. $\endgroup$ – Ty Sterling Sep 4 '20 at 16:22
  • 3
    $\begingroup$ @TySterling Sounds like the first sentence of your comment can be a separate question on this site. There will be several answers, you can post a link to it here afterwards so that people that want to know the answer to your comment can see the answers over there. $\endgroup$ – Nike Dattani Sep 4 '20 at 18:41

Depending on what properties you care about.. maybe. I don't suggest this in general. I work on phonon calculations where forces must relaxed to a pretty strict tolerance. It has been my experience that the relaxed structure (using e.g. VASP) is not usually identical to what another code finds (e.g. ABINIT). Of course, the two structures are close and further relaxation with the 2nd code (in my case ABINIT) converges in only a couple steps. For something like total energy or electronic density of states, it probably doesn't really matter.

  • 1
    $\begingroup$ I think in your case, we should find out why they are different. Is there any discussion on their forums? $\endgroup$ – Y. Zhai Sep 4 '20 at 2:19
  • 1
    $\begingroup$ Checkout this link molmod.ugent.be/deltacodesdft ----------- All DFT codes are NOT created equally. For things like total energy, using the same basis set, k-point grid, level of approximation, etc in different codes should yield identical results. However, the numerical aspects (integrals in real vs reciprocal space etc.) lead to slightly different results. For strict tolerances such as for phonons, it is more likely to matter. $\endgroup$ – Ty Sterling Sep 4 '20 at 16:10
  • 1
    $\begingroup$ @Y.Zhai, TySterling, this is what I think Susi Lehtola was trying to say in this answer: mattermodeling.stackexchange.com/a/1558/5 $\endgroup$ – Nike Dattani Sep 5 '20 at 16:21

This has been done previously. In my experience the geometry is found with one code, and solvent corrections are performed with another. I have personally done this[1], and have had the results published. I performed a geometry optimization using QChem, and then solvent corrections using ORCA.


[1] T. J. Doyon, J. C. Perkins, S. A. Baker Dockrey, E. O. Romero, K. C. Skinner, P. M. Zimmerman, A. R. H. Narayan “Chemoenzymatic o-quinone methide formation” J. Am. Chem. Soc. 2019, 141, 20269–20277.


Yes, absolutely.

The reason is that different codes are required to generate exactly the same results if they refers to the same methods.

E.g., if you use Gaussian, Orca, Molpro, NWChem, Psi4, PySCF...(I cannot list the full list) to perform HF/cc-pVTZ to study H2, they should generate exactly the same results, including the energy and wave functions.

Different packages can of course output in different formats. Besides, Different packages can have different default values for some parameters (E.g., the threshold of optimization, the geminal beta of explicit correlated methods, size of integral grids, the initial orbital for SCF iteration, etc.) It is the users' responsibility to check the default value, because apparently it is impossible for the pioneer package to tell you how will the after-coming software treats the default values.

So, if you perform the exactly same computation use two different packages, here same means the same explicitly given parameters and implicit parameters, while get different results, we should consider

  • At least one of them has a bug; or
  • They should be referred as different methods, e.g., the PM7 methods in Gaussian and in MOPAC. This should be mentioned in the manual, or should be considered as a bug.

In summary, it is the methods matter, not the software. Since it is a common way to perform the geometry optimization using a cheaper method, and to do the single-point computation using a more accurate method, of course it is acceptable to use different software, considering some methods are only available in less-used packages.

  • 2
    $\begingroup$ This is true for quantum chemistry, but for solid state calculations there has been a LOT more variability. How else would a paper be published on Science that various codes actually give the same result? :P science.sciencemag.org/content/351/6280/aad3000.long $\endgroup$ – Susi Lehtola Sep 4 '20 at 12:06
  • 1
    $\begingroup$ @SusiLehtola Nice point. But what they compared are different functionals/basis sets in Fig. 4 therein. You see, in quantum chemistry, HF/cc-pVTZ will of course generate different results for B3LYP/6-311g*, that is not surprising, and is what I wrote point 2 of the unnumbered list. Maybe I missed that, were there two software packages claim that they implemented exactly the same method, while led to different results? $\endgroup$ – Y. Zhai Sep 5 '20 at 3:55
  • 1
    $\begingroup$ nope, that's not at all what they did: "In this work, we applied the original Delta protocol to 40 DFT implementations of the Perdew-Burke-Ernzerhof (PBE) functional (43)." These were all calculations with PBE. Even in quantum chemistry e.g. PBE/cc-pVTZ can mean different things depending on which algorithms are used. Just using a different linear dependence threshold will affect absolute energies by quite a lot, even if relative energies are reproduced to similar accuracy. Then there's more choices of using RI, looser thresholds, etc. It's not a given that you get the same answer. $\endgroup$ – Susi Lehtola Sep 5 '20 at 12:29
  • $\begingroup$ You often have to spend a considerable amount of time to make calculations match exactly! $\endgroup$ – Susi Lehtola Sep 5 '20 at 12:29

Sure. But you have to be aware that the referees can question about that.

If there is no physical reason to mix codes, why mixing?


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.