How To Simulate the Universe with DFT (Script) (Patreon)

If you used every particle in the observable universe to do a full quantum simulation, how big would that simulation be? At best a large molecule. That’s how insanely information dense the quantum wavefunction really is. And yet we routinely simulate systems with thousands, even millions of particles. How? By cheating. Let’s learn how to cheat the universe.

—-

Quantum mechanics allows us to predict the behavior of the subatomic world with truly incredible precision. For example, let’s say you want to know the intricate behavior of an electron, whether it’s bouncing around inside a box, or part of a hydrogen atom, or moving through a double slit experiment. You just scratch the good ol’ Schrodinger equation on the blackboard and solve it to learn the probability distribution of the electron locations, or the energy levels of the system, or whatever.

That’s neat if you’re really, really interested in only that single electron. But for almost every practical use you’d need to do that math for multiple quantum particles interacting - and then the blackboard doesn’t cut it. You need exponentially more computing power and storage the more particles you have. It gets pretty bad. In fact you need more particles that exist in the solar system to store the wavefunction of the electrons in a single iron atom. And yet we can do this for thousands of atoms.

First up, let’s take a look at the Schrodinger equation. This is the time-independent version. It describes how the wavefunction of a quantum particle - that’s this psi thing - changes over space, assuming the particle is in some environment described by a sum of potential energies - that’s the V. V could result from the electromagnetic field inside a hydrogen atom, or the EM fields defining the walls of a box, etc. And this E is just the total energy of the system. For an electron in a box or in a hydrogen atom, the different possible values of E that form solutions to this equation define the energy levels of the system. You can also solve the Schrodinger equation to find the wavefunction, and the square of that wavefunction gives the the probability distribution for where you’ll find the electron if you try to measure it.

The time independent Schrodinger equation isn’t the be all and end all of quantum mechanics - it’s an approximation that works for slower moving particles that don’t change over time. But it’s where we start learning quantum mechanics, and it works for a lot of simple cases. There’s one other complexity we avoided. Our Schrodinger equation just dealt with one dimension of space, x. That’s nice because then the wavefunction is just a 1-D array of values. And remember, that array stores the wavefunction - the distribution of possible locations - it doesn’t store the actual location, which isn’t even defined until you look. Of course the real universe has 3 spatial dimensions, so for most real applications we’d want Psi(x,y,z). That’s a big increase in the information - literally the cube of the 1-D case.

Let’s say we have a fairly crude grid of 10 data points in 1-D. The 3-D equivalent would need 10^3 - 1,000 data points. OK, now let’s say we want to do the 26 electrons in an iron atom instead of the 1 electron in hydrogen.

That’s 26 times the data right? Nope. Those electrons are all interacting with each other, and that increases the information content just a tad. Every new electron doesn’t just add a new set of values to the same 3-D array. Every new electron adds 3 new entire dimensions. So 26 electrons means 78 dimensions, which for our 10-point grid is 10^78 numbers. This little estimation came from the physicist Douglas Hartree. He’s the one who pointed out that to describe the iron atom on this course grid, you’d need to store more numbers than there are particles in the solar system.

Let’s dig into why we really need this much information. We’ll do it with an analogy. Roll a single six-sided die. We know we can get outcomes 1-6. That’s a 1-dimensional space of possibilities. Now roll a second die, again giving you 1-6. Our total outcome space includes 6 possible outcomes of the second die for each of the 6 possible outcomes of the first - so 36 outcomes. We can represent this by first drawing a line to represent the 1-D outcome space of the first die. The second die expands this space to 2-D.

Adding particles to the Schrodinger equation is like adding dice to this system- every time we add a particle we increase the dimensionality of the wavefunction. The x, y, and z in the Psi of the first electron isn’t the same x, y, and z in the Psi of the second. For every coordinate point for electron one, we need to consider separately every coordinate point for electron two. And for a given pair of coordinate points for electrons one and two, we need to consider every possible point for electron 3 … and so on to electron 26.

So it seems like even for a single atom of iron, a fairly run of the mill element, we can’t even store the wavefunction let alone calculate it. So, we pack up and go home? Not at all. Like I said, if you can’t win, cheat.

Cheating in quantum mechanics often means solving a problem that’s much simpler than the one you’re actually interested in, and then finding ways to tweak that solution to get what you really want. That’s the heart of perturbation theory, but it’s also what’s done in density functional theory.

So to solve the impossible case of many interacting quantum particles, we should start by thinking about the completely solvable case of many non-quantum or classical particles. In other words, Newtonian mechanics. When we use Newton to solve, say, the three-body problem, we can solve the equations for each of the 3 particles separately. We don’t need a 3-cubed or 9 dimensional equation. In fact astrophysicists do huge galaxy simulations of millions of particles without doing millions-of-dimensions calculations. Why is that possible?

Well, these systems do sort of exist in a many-dimensional space called configuration space. But for the classical case we can throw away most of configuration space and cyonsider just the few points in that space where the particles actually exist at a given point in time. In quantum mechanics, we’re dealing with the wavefunction, and the wavefunction fills all of configuration space. For perfect precision we can’t throw any of it away. Not only that, but QM contains non-local correlations which arise because the position of one quantum particle can restrict the set of possible positions for the other particles, for example through the Pauli exclusion principle and through quantum entanglement.

In the Newtonian case, particles only interact locally, and that means the Newtonian equations of motion are what we call separable. An equation is separable if the solution along one axis - in one dimension - doesn’t depend on the solution along any other axis.

If this is true then you can take an equation for N particles in 3D and reduce it from a 3^N dimensional equation to simply N coupled equations each in 3D. In fact in Newtonian mechanics we can not only write down the equations of motion for each particle separately, we can write the x, y, and z equations of motion separately. That gives us a sufficiently sane number of coupled equations that can actually be solved.

But if we want to keep the quantum behaviour of quantum mechanics we can’t throw away most of configuration space like we do in Newtonian mechanics. We need to know the wavefunction for every particle everywhere. And we can’t reduce the dimensionality by treating particles separately because the Schrodinger equation can’t be made “separable”. Doing so destroys quantum correlations. And yet despite the apparent impossibility of solving the Schrodinger equation for more than a few particles, researchers still manage to do quantum simulations of some extremely complex systems. That’s because they have the cheat code.

The code is D-F-T. Make sure you get the order right or it won’t work. It stands for density functional theory, and it’s perhaps the most successful approach to tackling the extreme dimensionality problem when solving for realistic quantum systems. As an example, here’s a quantum simulation of the million atoms comprising the capsid of a virus done using density functional theory, and here’s resonant excitation in the many electrons of a molecule .

So how does DFT do a calculation that should need to manipulate vastly more bits of information than there are particles in the entire universe? Remember what we talked about - try to solve a problem that you CAN solve, and then tweak it to the more challenging case. In the case of DFT, what you do is just pretend the electrons aren’t interacting with each other at all and solve that case. And through a very mysterious quality of the quantum world, it's possible to map that fake solution to real answers.

Let me try to give you a sense of how this is done. The “mysterious quality” is expressed in a set of theorems devised by the physicists Walter Kohn and Pierre Hohenberg. The first theorem basically states that if you have a system of electrons in their ground state, no matter how complicated, the properties of that system are uniquely determined by position density of those electrons. The position density - or more commonly the charge density - is just a tiny fragment of the information held in the total wavefunction of all of those electrons. It’s a 3-D entity - just a map of how much “electron-ness” exists through space. The Kohn-Hohenberg theorems say that you can map from this charge density to the interesting observables - like the energy of the system, without having to go through the impossibly complex many-electron wavefunction.

So we start with a bad guess at the ground state charge distribution and then iterate closer to the truth, and we do this for the totally unrealistic case of non-interacting electrons. Because they aren’t interacting, the equations of motion for these electrons are separable, just like in Newtonian mechanics. Instead of using the Schrodinger equation, you use an analog called the Kohn-Sham equations.

There’s that Kohn guy again - he got the Nobel prize in chemistry for this stuff by the way. These equations basically let you translate between a realistic interacting system and the fictitious non-interacting case. Solving those equations becomes possible, and it gets you the ground state energy for your guess at the charge density. Then iterate until everything is consistent - when the ground state energy, the potential, and the charge density converge. And according to the theorem we started with, that ground state charge distribution is unique - it corresponds both to the fake non-interacting and the true interacting case. From there you can go ahead and calculate whatever observable you’re interested in.

The secret sauce that makes this work is something called the energy functional - it’s the functional in density functional theory. When given a ground state density, it let’s us construct this fictitious non-interacting “Kohn-Sham” system and figure out its total energy. DFT is all about finding the energy functional. We never know the exact functional, just that it exists - and so we have to approximate it.

So yeah, that’s a very crude summary. The main takeaway is that physicists realised that a tiny sliver of the full wavefunction - the density distribution - could be mapped to all sorts of useful observables. Without taking, like, the entire age of the universe and a solar-system sized computer. The key is understanding that the final solution must be self-consistent - basically, it has to obey the Schrodinger equation - and the cheat is finding a way to inch your way to that solution. Density functional theory has now been used to model the intricate quantum-level behavior of chemical reactions, of complex molecules even as far as DNA or that viral capsid from earlier, and many advanced materials like semiconductors, nanostructures, etc.

Being able to calculate this stuff is nice. I mean, it’s essential for our technological advancement. But what does DFT actually tell us about the universe? What does it mean that there exists a map between the low-information slice of the wavefunction and really all the information we want to get from it?  It’s like the ultimate compression algorithm. For one thing, it’s good news if we ever want to simulate another universe. We won’t need a computer bigger than the universe. But there are no doubt deeper truths to be found by understanding how the universal wavefunction with its insane dimensionality is connected to the narrow sliver of that wavefunction that represents our observed reality. Our computationally tractable reality, due to its very few dimensions of spacetime.