Linear Variational Method

Linear Variational Method#

What You Need to Know

The variational method is a powerful tool for estimating upper bounds and approximations for ground-state energies in a wide range of quantum mechanical problems.
It plays a key role in electronic structure theories, such as Hartree-Fock.
The linear variational method seeks solutions to the Schrödinger equation by representing trial wavefunctions as linear combinations of simple, computationally efficient functions, such as Gaussians or exponentials.
By applying the linear variational method, the Schrödinger equation is transformed into a linear algebra problem, where the goal is to find eigenvalues (representing energies) and eigenvectors (providing coefficients for the linear combination).

Linearizing the problem#

How does linearization help us in practice? In a typical QM problem, we are attempting to solve for the wavefunction(s) and energy(ies) for a given Hamiltonian. If the Hamiltonian is such that solving the problem exactly is too challenging (e.g. any two or more electron problem) we can expand the wavefunction in a basis and attempt to solve the problem that way. We start with

ψ \approx ϕ = \sum_{n}^{N} c_{n} f_{n}

Truncating this expansion at any finite $n$ leads to an approximate solution. Also, what the variational principle tells us is that we can minimize the energy with respect to variational paramters $c_{n}$ and still have $E_{ϕ} \geq E_{0}$ .

Smallest example#

We will illustrate this idea and the general matrix construction with a simple example of two basis functions ( $N = 2$ )

ϕ = c_{1} f_{1} + c_{2} f_{2}

There is currently no need to define these functions explicitly so we will leave them as generic functions $f_{1}$ and $f_{2}$ . We now solve for the energy

E_{ϕ} = \frac{⟨ ϕ | \hat{H} | ϕ ⟩}{⟨ ϕ | ϕ ⟩} = \frac{⟨ c_{1} f_{1} + c_{2} f_{2} | \hat{H} | c_{1} f_{1} + c_{2} f_{2} ⟩}{⟨ c_{1} f_{1} + c_{2} f_{2} | c_{1} f_{1} + c_{2} f_{2} ⟩}

E_{ϕ} = \frac{c_{1}^{2} H_{11} + 2 c_{1} c_{2} H_{12} + c_{2}^{2} H_{22}}{c_{1}^{2} S_{11} + 2 c_{1} c_{2} S_{12} + c_{2}^{2} S_{22}}

$H_{i j} = ⟨ f_{i} | \hat{H} | f_{j} ⟩$ Hamiltonian matrix element expressed in basis of $f$
$S_{i j} = ⟨ f_{i} | f_{j} ⟩$ S matrix element or overlap integral expressed in basis of $f$ . This one measures how similiar $f_{i}$ is to $f_{j}$ .

Eigenvalue Problem#

Since $E_{ϕ} \geq E_{0}$ for any trial function $ϕ$ , we can minimize the energy $E_{ϕ}$ by varying the parameters $c_{1}$ and $c_{2}$ .
To minimize with respect to $c_{1}$ , we differentiate $E_{ϕ}$ with respect to $c_{1}$ and set the derivative equal to zero:

$\frac{\partial E_{ϕ}}{\partial c_{1}} = 0 = c_{1} (H_{11} - E S_{11}) + c_{2} (H_{12} - E S_{12})$
Similarly, minimizing with respect to $c_{2}$ gives:

$\frac{\partial E_{ϕ}}{\partial c_{2}} = 0 = c_{1} (H_{12} - E S_{12}) + c_{2} (H_{22} - E S_{22})$
These two coupled linear equations can be expressed compactly as a matrix equation:

$\begin{array}{r} [\begin{array}{c} H_{11} - E S_{11} & H_{12} - E S_{12} \\ H_{12} - E S_{12} & H_{22} - E S_{22} \end{array}] [\begin{array}{c} c_{1} \\ c_{2} \end{array}] = 0 \end{array}$
The matrix on the left can be rewritten as the difference between two matrices:

$\begin{array}{r} ([\begin{array}{c} H_{11} & H_{12} \\ H_{12} & H_{22} \end{array}] - E [\begin{array}{c} S_{11} & S_{12} \\ S_{12} & S_{22} \end{array}]) [\begin{array}{c} c_{1} \\ c_{2} \end{array}] = 0 \end{array}$
Rearranging, we can write:

$\begin{array}{r} [\begin{array}{c} H_{11} & H_{12} \\ H_{12} & H_{22} \end{array}] [\begin{array}{c} c_{1} \\ c_{2} \end{array}] = E [\begin{array}{c} S_{11} & S_{12} \\ S_{12} & S_{22} \end{array}] [\begin{array}{c} c_{1} \\ c_{2} \end{array}] \end{array}$
In more compact matrix notation, this becomes:

By left-multiplying both sides by $S^{- 1}$ , we transform this into a standard eigenvalue problem:

$S^{- 1} H c = E I c$
Therefore, the minimum energies correspond to the eigenvalues of $S^{- 1} H$ , and the variational parameters that minimize the energies are the eigenvectors of $S^{- 1} H$ .

Breaking problem down to matrix eigenvalue eigenvector problem

In the equation

S^{- 1} H c = E I c,

$I$ represents the identity matrix. Its role in this context is essential to express the equation as a standard eigenvalue problem.

Eigenvalue Problem Form:

In linear algebra, a standard eigenvalue problem is written as:
$ $A v = λ I v,$ $
where:
- $A$ is a square matrix,
- $λ$ is a scalar eigenvalue,
- $I$ is the identity matrix, and
- $v$ is the corresponding eigenvector.
The identity matrix $I$ ensures that $λ$ scales the eigenvector $v$ without altering its direction. The eigenvalue problem is about finding the values of $λ$ and their associated $v$ .
Connecting to $S^{- 1} H c = E I c$ :
Here, $S^{- 1} H$ acts as the operator $A$ in the standard eigenvalue problem.
- $S^{- 1} H$ is a matrix resulting from left-multiplying $H$ by the inverse of $S$ .
- $c$ represents the eigenvector.
- $E$ is the eigenvalue (corresponding to the energy in the quantum mechanical system).
The identity matrix $I$ is explicitly included to highlight that $E$ is a scalar multiplying the vector $c$ . This ensures that the left-hand side (a matrix operation) matches the right-hand side (a scaled vector).
Why $S^{- 1}$ Appears:
Initially, we had:
$ $H c = E S c,$ $w h i c h c a n n o t d i r e c t l y b e i n t e r p r e t e d a s a n e i g e n v a l u e p r o b l e m b e c a u s e o f t h e p r e s e n c e o f$ \mathbf{S} $(t h e o v e r l a p m a t r i x) . T o t r a n s f o r m t h i s i n t o a s t a n d a r d f o r m, w e p r e - m u l t i p l y b o t h s i d e s b y$ \mathbf{S}^{-1} $:$ $S^{- 1} H c = E S^{- 1} S c .$ $S i n c e$ \mathbf{S}^{-1}\mathbf{S} = \mathbf{I} $, t h i s s i m p l i f i e s t o :$ $S^{- 1} H c = E I c .$ $
How to Interpret This as an Eigenvalue Problem:
The equation now has the form of a standard eigenvalue problem:
$ $A v = λ I v,$ $
where:
- $A = S^{- 1} H$ is the effective matrix to diagonalize,
- $λ = E$ are the eigenvalues, corresponding to the energy levels,
- $v = c$ are the eigenvectors, containing the coefficients of the trial wavefunctions.
Physical Interpretation:
Solving the eigenvalue problem $S^{- 1} H c = E I c$ gives the approximate energies ( $E$ ) of the quantum system as eigenvalues and the corresponding variational parameters ( $c$ ) as eigenvectors. The identity matrix $I$ is crucial for preserving the standard form of the eigenvalue problem, ensuring proper mathematical and physical interpretation.

Example: Particle in a Box#

Let’s consider a free particle in 1D bounded to $0 \leq x \leq a$ . The Hamiltonian for such a system is simply the kinetic energy operator

\hat{H} = - \frac{ℏ^{2}}{2 m} \frac{d^{2}}{d x^{2}}

While we can have solved this problem analytically, it will be instructive to see how the variational solution works. We start by approximating $ψ (x)$ as an expansion in two basis functions

ψ (x) \approx c_{1} x (a - x) + c_{2} x^{2} (a - x)^{2}

where the basis functions are $f_{1} (x) = x (a - x)$ and $f_{2} (x) = x^{2} (a - x)^{2}$ and $c_{1}$ and $c_{2}$ are the variational parameters (/linear coefficients of basis functions). In order to solve for the variational energies and wavefunctions given this expansion we must construct the Hamiltonian matrix, $H$ , and the basis function overlap matrix, $S$ . We will then compute and diagonlize the $S^{- 1} H$ matrix. Recall,

H_{i j} = ⟨ f_{i} | \hat{H} | f_{j} ⟩

and

S_{i j} = ⟨ f_{i} | f_{j} ⟩

In this problem, we have two basis functions

f_{1} (x) = x (a - x)

and

f_{2} (x) = x^{2} (a - x)^{2}

Computing the matrix elements:#

H_{11} = ⟨ x (a - x) | \hat{H} | x (a - x) ⟩

= ⟨ a x | \hat{H} | a x ⟩ - ⟨ a x | \hat{H} | x^{2} ⟩ - ⟨ x^{2} | \hat{H} | a x ⟩ + ⟨ x^{2} | \hat{H} | x^{2} ⟩

= - a ⟨ x | \hat{H} | x^{2} ⟩ + ⟨ x^{2} | \hat{H} | x^{2} ⟩

where the last equality holds because the second derivative of $x$ with respect to x is zero. Let’s now investigate each of these integrals

\begin{array}{r} ⟨ x | \hat{H} | x^{2} ⟩ = \int_{0}^{a} x \frac{- ℏ^{2}}{2 m} \frac{d^{2}}{d x^{2}} x^{2} d x \\ = \frac{- ℏ^{2}}{m} \int_{0}^{a} x d x = \frac{- ℏ^{2} a^{2}}{2 m} \end{array}

Similarly we should get

⟨ x^{2} | \hat{H} | x^{2} ⟩ = \frac{- ℏ^{2} a^{3}}{3 m}

Thus

H_{11} = - a \frac{- ℏ^{2} a^{2}}{2 m} + \frac{- ℏ^{2} a^{3}}{3 m} = \frac{ℏ^{2} a^{3}}{6 m}

H_{22} = ⟨ x^{2} (a - x)^{2} | \hat{H} | x^{2} (a - x)^{2} ⟩

= ⟨ a^{2} x^{2} | \hat{H} | a^{2} x^{2} ⟩ - 2 ⟨ a^{2} x^{2} | \hat{H} | 2 a x^{3} ⟩ + ⟨ a^{2} x^{2} | \hat{H} | x^{4} ⟩ - 2 ⟨ a x^{3} | \hat{H} | a^{2} x^{2} ⟩ + 4 ⟨ a x^{3} | \hat{H} | a x^{3} ⟩ - 2 ⟨ a x^{3} | \hat{H} | x^{4} ⟩ + ⟨ x^{4} | \hat{H} | a^{2} x^{2} ⟩ - 2 ⟨ x^{4} | \hat{H} | a x^{3} ⟩ + ⟨ x^{4} | \hat{H} | x^{4} ⟩

= \frac{- ℏ^{2} a^{7}}{3 m} + \frac{3 ℏ^{2} a^{7}}{m} - \frac{6 ℏ^{2} a^{7}}{5 m} + \frac{ℏ^{2} a^{7}}{2 m} - \frac{12 ℏ^{2} a^{7}}{5 m} + \frac{2 ℏ^{2} a^{7}}{m} - \frac{ℏ^{2} a^{7}}{5 m} + \frac{ℏ^{2} a^{7}}{m} - \frac{6 ℏ^{2} a^{7}}{7 m}

= \frac{ℏ^{2} a^{7}}{105 m}

H_{12} = H_{21} = \frac{ℏ^{2} a^{5}}{30 m}

We can now complete the Hamiltonian matrix, $H$ ,

\begin{array}{r} H = \frac{ℏ^{2} a^{3}}{m} [\begin{array}{c} \frac{1}{6} & \frac{a^{2}}{30} \\ \frac{a^{2}}{30} & \frac{a^{4}}{105} \end{array}] \end{array}

It can also be shown that $S$ is

\begin{array}{r} S = \frac{a^{5}}{10} [\begin{array}{c} \frac{1}{3} & \frac{a^{2}}{14} \\ \frac{a^{2}}{14} & \frac{a^{4}}{63} \end{array}] \end{array}

Eigenvalues: [ 4.93487481+0.j 51.06512519+0.j]
First eigenvector: [-0.66168489 -0.74978204]

So we see that the smallest energy in this basis is

E_{ϕ} = 4.9349 \frac{ℏ^{2}}{m}

How does this compare to the analytic solution? We need to recall that for particle in a box we have:

E_{n} = \frac{n^{2} π^{2} ℏ^{2}}{2 m a^{2}}

Plugging in for the ground state, $n = 1$ , and $a = 1$ since that is what we used numerically above we get

E_{1} = \frac{π^{2} ℏ^{2}}{2} \approx 4.9348 \frac{ℏ^{2}}{m}

So we can see that our variational solution worked out well for the energy. Now how about the wavefunction?

<matplotlib.legend.Legend at 0x7f9083657a90>