A function that is already there

The variational principle gives a guarantee that is easy to state and easy to over-read. Enlarging a basis set can only lower the computed energy, because the old basis is contained in the new one and its answer is still available. So a bigger basis is never worse.

That is true, and it is a statement about the energy alone. This essay is about what else the same act does.

What a basis has to be, before it is any good #

A basis is a set of functions that spans a space, and spanning is a property with two halves. The set has to reach everywhere the answer might be, which is what more functions buys. And its members have to be distinguishable from one another, which is what having many similar functions in a small space destroys.

The second condition is not a matter of taste or numerical hygiene. A set of functions that is linearly dependent is not a basis at all — there is no unique expansion of anything in it — and a set that is nearly dependent is one where the expansion exists and means very little.

Crowding is inevitable rather than a failure of the search, and the reason is a count: the functions a basis has to cover span a factor of several in every measure of size, and there is only so much room between the tightest and the most diffuse. The more functions are added, the closer together two of them must be.

Hybrids are a basis makes the corresponding point about a set that is exactly independent and is a rotation of another: what a basis is for is to give a unique set of coefficients, and near-dependence is the state in which it barely does.

The quantity nobody reports #

A basis of Gaussians is a set of functions, and the matrix that says how much they resemble one another is the overlap matrix. Normalise every function and its diagonal is all ones; the off-diagonal entries are the overlaps between pairs.

Its eigenvalues are what matter, and the object itself is the one hybrids that were never orthogonal is about in a different setting: two functions that are supposed to be a basis and are not independent enough to be one. They sum to the number of functions, so their average is one, and the smallest of them measures how nearly some combination of the set vanishes: a value of one means orthogonal, and a value approaching zero means a combination that is very nearly nothing at all, whose coefficients therefore mean very little.

Every quantum-chemical calculation forms that matrix and inverts a square root of it. Almost no published calculation reports its smallest eigenvalue.

It is worth saying why. A calculation that runs to completion and returns a lower energy has, by the only test most users apply, succeeded. The eigenvalue is a diagnostic of a failure mode that does not announce itself, and the discipline of looking at it belongs to the people who write the programs rather than to the people who use them.

What growing a basis does to it #

The exponents in that figure are searched for rather than tabulated — the same descent the atom does not bring its own orbital uses to let a molecule choose how large its functions should be — a Gaussian is the wrong shape explains why quoting somebody else’s would put the interesting part behind a table — so every basis in the sequence is as good as a basis of that size can be.

Even so, the eigenvalues run 1.000, 0.445, 0.220, 0.114, 0.0615, 0.0342. A basis gets better and less independent by the same act, and the reason is not a defect of the search. There is only so much room: functions spanning one atom’s worth of space, all of the same shape, get more alike the more of them there are.

Extrapolating that halving is instructive. Thirty functions would put the smallest eigenvalue near 10⁻⁸ and forty-five near 10⁻¹², which is where double-precision arithmetic gives out. That is not a hypothetical — it is why large basis sets in real calculations are handled with a threshold that discards the smallest eigenvectors, and why the threshold has to exist at all.

What the functions in that sequence look like is an optimised set of Gaussians and the exponential they are trying to be, with each new function placed where the fit is worst. That is the right thing for an optimiser to do and it is exactly what drives two exponents together — the worst-fitted region shrinks, so the new function lands near an old one.

The mechanism is visible in the picture. The first function covers the bulk of the density; the second is added tight, to attend to the cusp; the third and fourth fill in between. By the sixth there is no region of the radial coordinate that is not already covered by two functions, and a new one is necessarily a partial duplicate of its neighbours.

That is a general property of fitting a smooth function with a set of similar shapes, and it has nothing to do with Gaussians in particular. It is why every basis of any kind eventually runs into this, and why the eigenvalue’s decay is smooth rather than sudden.

The tight functions are for the cusp at the nucleus and the diffuse ones for the tail; everything between is interpolation. Neither end is a Gaussian’s shape, which is why a basis of them needs several functions to do what one exponential does, and why adding more of them has to run out of room eventually.

What a duplicate does #

The other way to lose independence is not to grow the basis but to place a function badly, and the extreme case makes the mechanism visible.

A well-placed fifth function is worth 5.3 × 10⁻⁴ hartree — about 1.4 kilojoules a mole, which is a real quantity on a chemical scale. The near-duplicate at 1.01 is worth 3.4 × 10⁻⁹, which is nothing at all.

So the two failures are separable and they fail at different rates. The energy stops improving long before the matrix stops being invertible, which means there is a wide region in which a function is contributing nothing and costing conditioning — and no warning is issued anywhere in that region, because the energy is still going down.

The coefficients do not blow up, and that is the trap #

There is a version of this warning that is often given and is wrong in a way worth correcting, because getting it wrong is what makes the real problem invisible.

The usual account says that near-dependence makes the coefficients blow up: two nearly identical functions acquire large coefficients of opposite sign that nearly cancel, and the wavefunction becomes a difference of large numbers. That is what happens if the equations are solved by a method that has no defence against it.

It is not what happens here. The largest coefficient in the five-function basis is 0.186 at a ratio of three and 0.186 at a ratio of 1.001 — unchanged to three figures across five orders of magnitude in the smallest eigenvalue. The reason is that the calculation orthogonalises the basis symmetrically, which is the procedure that produces the minimum-norm solution among all the equivalent ones, and a minimum-norm solution does not have large cancelling coefficients by construction.

So the standard symptom is absent, and the standard cure — watch for large coefficients — detects nothing. What the ill-conditioning costs is instead paid in the arithmetic that formed the orthogonalisation, silently, and the first sign of it is the refusal at the end.

Where the refusal comes from #

The calculation is refused at a ratio of 1.00003, and it is worth saying exactly what refuses and why, because it is not an arbitrary threshold.

Solving a variational problem in a non-orthogonal basis requires the inverse square root of the overlap matrix, which is formed from the eigenvalues by taking each to the power of −½. An eigenvalue of 4 × 10⁻¹¹ becomes a factor of 1.6 × 10⁵ in that matrix, and every rounding error in the Hamiltonian is amplified by it. Below some point the answer is noise, and a careful calculation of the inverse square root declines rather than returning it.

The right refusal names linear dependence rather than a singular matrix, which is the same fact stated in the language of the problem instead of the language of the arithmetic.

The cliff is not where the trouble starts #

The two curves in the sweep figure cross the interesting region at different points, and the gap between those points is the whole practical content of this essay.

The shape of the two curves is the argument of this essay in one picture, and it is worth stating in words.

The refusal is a cliff: on one side the calculation runs and on the other it does not. Everything interesting happens on the near side of it, in a region several orders of magnitude wide where the calculation runs perfectly happily, produces a number, lowers the energy, and is contributing nothing.

A user who tests a basis by whether the program completes has tested for the cliff. A user who tests by whether the energy went down has tested for something that is always true. Neither test sees the region between, and the region between is four orders of magnitude wide in the smallest eigenvalue and a factor of a hundred and sixty thousand wide in what the function is worth.

The other freedom a basis has is how tight or diffuse its functions should be, and what a search over it returns has already been measured. That search is what produced the optimised sets above, and it is the search whose answers crowd.

What a contraction is for #

This gives a second reading of the essay beside it, and the two fit together.

A contraction is a decision made once is about the practice of freezing primitive Gaussians into fixed combinations, on an isolated atom, before any molecule is in sight, and about what that costs: 0.0283 hartree at the effective charge a bond gives, with one freed coefficient recovering three quarters of it.

That essay treats contraction as a compromise made for speed. It is also a device for conditioning. A contracted function is a fixed combination of primitives, and choosing the combinations well means choosing them to be well separated from one another — so contraction removes exactly the near-dependencies that a large uncontracted set of primitives would have.

So the two costs run against each other: contracting loses variational freedom and buys independence, and the published basis sets are a compromise struck between the two rather than a compromise struck for speed alone.

What freezing a contraction costs in energy is the same design decision seen from the other side: a contraction removes freedom to buy speed, and near-dependence is what happens when the freedom is kept and the functions run out of room. One is priced in energy and the other in conditioning.

The one number to watch #

The whole practical content can be put as a pair of numbers to keep beside a calculation, which is more useful than a warning.

The first is the smallest eigenvalue of the normalised overlap matrix. Above about 10⁻⁴ nothing is wrong; between 10⁻⁴ and 10⁻⁶ the basis is working hard; below 10⁻⁸ any property other than the energy should be treated with suspicion, whatever the program returns.

The second is what the last function added was worth. If it is a chemically meaningful quantity — a kilojoule a mole, say — the basis is crowded but earning its size. If it is a millionth of that, the basis has a redundant member and removing it will improve every other number in the calculation at no cost.

Neither of those requires any new calculation. Both come out of matrices the calculation has already formed, and the reason they are not routinely printed is that the energy is the thing everyone came for.

Why the guarantee is not a lie #

Nothing above contradicts the variational principle, and it is worth being explicit about the relation, because the principle is exact and this essay is not proposing an exception.

Adding a function does lower the energy, always, in exact arithmetic. What the near-dependence costs is not the energy but everything else: the coefficients, which become large and nearly cancelling; the properties computed from the wavefunction, which are sums over those coefficients; and the numerical stability of the whole calculation.

That division has already been made in another setting. A better energy is not a better answer is the essay about the variational principle making the energy the least sensitive thing a wavefunction gets wrong — second order in the error where every other property is first — so a calculation can improve its energy and worsen its dipole moment. The same asymmetry is why a Gaussian is the wrong shape can be true of every function in a basis while the energies computed from them are excellent. What a dipole cannot tell apart is the essay about how little a moment pins down in the first place, and a moment computed from ill-conditioned coefficients is worse off still. Near-dependence is a mechanism for exactly that, and a sharp one.

The bases that never run into it, and what they pay instead #

Every basis of any kind is too strong, and the exception is worth stating because it accounts for the other half of electronic-structure practice.

The crowding above happens because the functions are localised and similar: they all live in one atom’s worth of space and they all have the same shape, so a new one is necessarily close to the span of the old ones. Neither property is required of a basis. A set of functions that are orthogonal by construction has every overlap eigenvalue equal to one, at every size, and its conditioning does not degrade at all — there is no smallest eigenvalue to watch, because they are all the same.

Plane waves are such a set. A calculation in them can use a hundred thousand functions with no linear-dependence threshold, no inverse square root to be careful about and no refusal, and the number is limited by memory rather than by arithmetic. That is not a small practical difference: it is why solid-state codes are routinely run at basis sizes that would be unthinkable in an atom-centred one.

What they pay is everything the argument about the cusp is about. A plane wave is delocalised and smooth, so describing a cusp with plane waves takes a very great many of them — the expansion converges as a power of the cutoff rather than exponentially — which is why such calculations almost always replace the region near the nucleus with a pseudopotential or a patched-in atomic solution rather than expanding it honestly.

So the two families are not better and worse. They are two ways of spending the same difficulty: an atom-centred basis is compact and becomes ill-conditioned, and an orthogonal one is well-conditioned and becomes enormous. The smallest overlap eigenvalue is the number to watch in the first family and is identically one in the second, which is as clean a statement as this essay can make about what it is a diagnostic of.

Two ways to be a bad basis #

The two figures in this essay describe failures that look alike and are not, and separating them is the practical conclusion.

Crowding is what the optimised sequence shows. Every function is doing something, and there is simply not much room in one atom’s worth of space, so the set becomes ill-conditioned by getting large. There is no repair for it beyond a threshold and careful arithmetic, and it is a cost of accuracy that has to be paid.

Redundancy is what the near-duplicate shows. A function is in the wrong place, contributing nothing and costing conditioning, and the repair is to move it. A basis suffering from this is worse in every respect than a smaller one, which is not true of the crowded case.

Distinguishing them from the smallest eigenvalue alone is impossible, because both drive it down. What distinguishes them is the marginal energy: a crowded basis’s last function is still worth something, and a redundant one’s is worth a millionth of that. Watching both numbers together tells them apart and watching either alone does not.

What is left #

The system here is one electron in a hydrogen atom, which is what makes every integral available in closed form and every number exact. What that leaves out is the two-electron integrals, which is where the cost of a large basis actually lives and where conditioning problems have their most expensive consequences: the number of such integrals grows as the fourth power of the basis size, and their contribution to the energy involves the same near-cancelling coefficients. The smallest many-electron calculation is where this collection goes instead, and it buys exactness by abandoning continuous space altogether.

The other absence is the second centre. A near-dependence in a molecular basis is more likely than in an atomic one, because a diffuse function on one atom can be very nearly a combination of functions on its neighbours — which is the defect known as basis set superposition, seen as a conditioning problem rather than as an energy one. That requires two-electron integrals over Gaussians on different centres, which this site still does not have.

What is here is the mechanism in the one case where every quantity is exact: the two effects of adding a function, the rate at which each proceeds, and the wide region between where one has stopped and the other has not.