# How do quantum software packages Work?

If one wants to calculate a moderate size Alkane (with say 10-15 carbons, assuming 100 electrons, with Restricted Hartree Fock based methods) we can simply say that the electron-electron part will be $$100^4 = 100$$ million integrals. If we remove the non-unique ones it will become $$\frac{100M}{8}= 12.5 M$$ Integrals. If we do a pre-screening it may become around 6 million. Even if we go for DFT it will be around $$\frac{100^3}{8} =125,000$$ .

I have calculated a small system like Methane with just 10 electrons (around 400 unique/pre-screened integrals) with an efficient C++ code (based on Rys Integration which is among the most efficient methods) and it took around 1 hour. When 400 cycles take 1 hour I can assume that 6 million cycles may take years! But software packages do that just in seconds!

My question is if I just want to do simple math like 2+2 for 6 million times it takes way longer than what these packages do with such complicated integrals!! Does anyone know what their trick is? Do they have a pre-calculated database or some approximation methods?

UPDATE: there are many good hints in this question and thanks to all for that. A combination of symmetry, finding negligible integrals, and using optimized methods work fine. The only unanswered part is an algorithm for finding symmetric situations.

• I think many of these packages use highly optimized numerical libraries for linear algebra like LAPACK and BLAS. They may also implement some parallel programming techniques (eg, multi-threading or MPI) to take advantage of multi-core processors. Commented Sep 23, 2014 at 1:33
• >The only unanswered part is an algorithm for finding the symmetric situations || actually, gamess puts this load onto user. Gaussian somehow can guess symmetry without it explicitly mentioning in input. It is not very hard, though. D2h and friends is easily detectable for linear molecules, and all elements of symmetry must go through geometrical center of the molecule, tranforming one atom into anouther. Since all symmetry operations are representable in matrix form, it is possible for each to points to concstuct such a matrix and try to see if it is indeed and element of symmetry. Commented Sep 28, 2014 at 10:48
• I also suggest to google for article "what every programmer should know about memory". It gives interesting insights to optimizations. Commented Sep 28, 2014 at 10:49
• @Aug You should look into the source code of free software computing packages. Psi4 is written mostly in C++ and Python. MPQC is mainly C++. There may be others available as well. There may also be some helpful information in the Q-Chem manual Commented Sep 29, 2014 at 12:17
• @Aug Doing Hartree-Fock for methane should be a few seconds if not less. The code you tested is highly inefficient for one reason or another. By the way, Rys integration is among the most efficient if implemented right. Anyway, for your small system choosing the best is not going to change much. Commented Aug 23, 2018 at 18:41

There are many tricks. Since I'm not a developer, I'm not privy to all of them, but there are some simple things to start.

• Some integrals are predictably close to zero: Even in fully first-principals calculations, there are parts of the code that prune out integrals that will clearly be negligible. This cuts out a lot of work. Consider, for example, overlap or Coulomb integrals between atoms that are very far apart.
• Gaussian Orbitals: Using Gaussian-type orbitals are incredibly fast because most quantum integrals can be derived as simple recurrence relationships. That is, you almost never need to do the integration. See, for example, Obara and Saika and the recent PRISM algorithm from Gill and Pople.

I believe there are some actual quantum developers on here and perhaps they'd be willing to share other tricks. But these are big - eliminate any integral you don't need, and don't actually take the integrals in the first place. ;-)

Yes, there are numeric tricks (using BLAS, LAPACK, etc.) and the use of multi-core processing. But these are actually minor optimizations.

UPDATE

I realized that I forgot to talk about the convergence of the self-consistent field which is also optimized. There are various convergence strategies, most notably DIIS, which help to extrapolate the final orbital coefficients and electron density.

• Yes I absolutely agree with the efficacy of rejection of negligible integrals especially in larger molecules as there will be more zero overlaps. I have used Obara-Saika for Overlap Integrals and Kinetic energy. I found that Rys iteration is way faster for e-e integrals . But I couldn't find a definite formula for selection-rejection of the integrals (I only know "Schwarz inequality" but it is not that effective). Do you know any effective selection-rejection criteria? I think the answer of this question may be such criteria.
– Aug
Commented Sep 23, 2014 at 3:09
• I'd add that many times you can exploit molecular point group symmetry. Even just having one symmetry element can reduce the storage by 2 and the work by 4. D2h, which has order eight, can reduce the work by 64! Commented Sep 23, 2014 at 19:55
• @jjgoings I think this symmetry sounds really promising but it is so difficult to find symmetric situations ( I mean by coding ) because 4 integrals can come from any different atoms but be symmetric yet. By any chance you may know a reliable algorithm to apply this symmetry?
– Aug
Commented Sep 28, 2014 at 4:37
• >> Yes, there are numeric tricks (using BLAS, LAPACK, etc.) and use of multi-core processing. But these are actually minor optimizations. || I.. wouldn't say so. use of multi-core and parallel runs may reduce job time by orders of magniture, though parallelization involves significant overhead. Also, some packages include GPGPU support. Commented Sep 28, 2014 at 8:43

First things first, I would check your code: there is no reason 400 unique integrals should take an hour, even with the most naive implementation. You are computing at a rate of an integral every ten seconds. While I agree with all the comments pointing to increasing parallelization, I think efficient scalar code is more important. Why scale up inefficient code? Considering GPUs is a whole 'nother beast, and from the problems, my coworkers have expressed trying to program them I'd recommend avoiding them for the time being.

The best thing you can do is to exploit molecular symmetry. Heads up though: this won't fix everything, because many interesting molecules lack symmetry. But for methane, as you tried, this should help a ton. Of course, this isn't as simple as just looking up some algorithm, because depending on your program you need to rotate the nuclear orientations to be symmetric, or develop routines to detect point group symmetry, etc. That being said, there are --- as far as I know --- two primary ways of exploiting point group symmetry for integrals. The first is to create symmetry-adapted basis functions, so your basis set is symmetric, which means your wave function/integrals/Fock matrix/whatever will be symmetric as well. Ernie Davidson did this in the seventies. A slightly more modern reference by Helgaker and Taylor is given here. Daniel Crawford's electronic structure tutorial explores this as well (it's really well written --- a great introduction).

The other way you can do it, and I believe this is more popular, is to create symmetry unique integrals and then build a partial Fock matrix only to take care of symmetrization later. If you see something like "generating petite list" in the output of your favorite electronic structure program, this is probably what it's doing. This method allows for integral screening more efficiently, I believe. Dupuis and King did this in the late seventies, based on earlier work by Dacre and Elder. These are the standard references you'll find. There may be some more recent work that I'm unaware of. The problem of point group symmetry had been effectively solved by 1980.

Of course, many interesting molecules don't have symmetry to take advantage of, so symmetry won't fix all your problems. The Schwarz equality is well addressed in Molecular Electronic Structure Theory (actually, most electronic structure theory is well addressed in that lovely pink book). If you haven't checked that out, please do. I think it will answer many of your questions about how modern codes work.

I don't want to repeat the other answers (Geoff gave you the papers I'd have suggested as well), so I'll finish with two lesser-considered references. First, there is an awesome (but tiny!) book called Parallel Computing in Quantum Chemistry that lays out the basics for efficient code including integral evaluation. It's not too hard to read, and I found the work on local MP2 to be pretty interesting. I'm often surprised more people haven't heard of it. The other place I would check out is Ed Valeev's LibInt package. He's a pretty talented programmer, and it's used in many modern packages like Psi and MPQC and ORCA, not to mention free for anyone to download and use. Some of my coworkers have played around with it for their own "homemade" electronic structure code, and speak highly of it. The user's manual is actually quite helpful in understanding integral evaluation. He has also written several exercises and tutorials, which you may find useful.

• Thank you for the good suggestions. I wrote some "Homemade codes" to find the symmetry ( it is not very efficient as it can not find signed symmetries/ plane symmetries and also is rather slow) but I cannot say how faster the whole calculation is now ( maybe 10 times faster )! I got surprised and now I know that the main magic is the symmetry and negligible integral rejection. I'm gonna dig into the suggested sources and will post the results here.
– Aug
Commented Sep 29, 2014 at 19:52

To further reduce storage requirements, integrals may be calculated on the fly, employing different caching strategies. For some methods storing all the integrals is plainly impossible.

There are two things here:

1. Quantum chemistry problem
2. Optimization problem

Quantum chemistry or the Schrodinger equation is simplified using DFT and some other techniques:

a) ECP or effective core potentials and pseudopotentials: Instead of solving for all electrons just the outer shell electrons are solved and core electrons are considered using pre-calculated integrals
b) ri approximation: Resolution of identity, it treats two electron terms as a product of single-particle orbital function thus greatly reducing computational cost

And as mentioned above there are some other good techniques are used to deal with quantum chemistry part.

Optimization problem: In optimization usually second derivative (Hessian matrix) is used to find minima or maxima but in computational chemistry code only the first derivative or jacobian is used, thus it is highly sensitive to initial guess and computational setup.

Calculating Hessian is the most costly operation of computational chemistry code, either approximate hessian is used or no hessian is used in this type of optimization problem. If you want to make sure that you got the right geometry you can calculate Hessian with your optimized geometry and for geometry optimization, all the frequencies should be positive and for the transition state, only one frequency should be negative but you will find that to meet this requirement is highly demanding.

Moreover to reduce computational cost we can apply periodic boundary conditions. It enables us to solve the problem for minimum nonrepetitive unit cell ( Irreducible Kpoints). Finding symmetry is also a crucial part of this type of problem.

Cut-off energy is used to reduce the number of basis sets considered. And also optimization techniques also have a great effect on computational cost, you may use a fast but aggressive algorithm but that may not be able to find the minima or maxima.

You can read VASP or TURBOMOLE manual (available freely online), and there you will find some good techniques for reducing the computational cost. And above all, as it is mentioned above efficient implementation of parallelization can help to reduce computational costs dramatically. One another potential method is to use GPU, but it is still in its primary stage.