$$\require{cancel}$$

# 7.3: Higher Dimensions


We can work out the finite-difference equations for higher dimensions in a similar manner. In two dimensions, for example, the wavefunction $$\psi(x,y)$$ is described with two indices:

$\psi_{mn} \equiv \psi(x_m, y_n).$

The discretization of the derivatives is carried out in the same way, using the mid-point rule for first partial derivatives in each direction, and the three-point rule for the second partial derivative in each direction. Let us suppose that the discretization spacing is equal in both directions:

$h = x_{m+1} - x_m = y_{n+1} - y_n.$

Then, for the second derivative, the Laplacian operator

$\nabla^2 \psi(x,y) \equiv \frac{\partial^2\psi}{\partial x^2} + \frac{\partial^2\psi}{\partial y^2}$

can be approximated by a five-point rule, which involves the value of the function at $$(m,n)$$ and its four nearest neighbors:

$\nabla^2\psi(x_m,y_n) \approx \frac{\psi_{m+1,n} + \psi_{m,n+1} - 4\psi_{mn} + \psi_{m-1,n} + \psi_{m,n-1}}{h^2} + O(h^2).$

For instance, the finite-difference equations for the 2D Schrödinger wave equation is

$-\frac{1}{2h^2}\, \Big[\psi_{m+1,n} + \psi_{m,n+1} - 4\psi_{mn} + \psi_{m-1,n} + \psi_{m,n-1} \Big] + V_{mn} \psi_{mn} = E \psi_{mn}.$

## 7.3.1 Matrix Reshaping

Higher-dimensional differential equations introduce one annoying complication: in order to convert between the finite-difference equation and the matrix equation, the indices have to be re-organized. For instance, the matrix form of the 2D Schrödinger wave equation should have the form

$\sum_{\nu} H_{\mu\nu} \psi_\nu = E \psi_\mu,$

where the wavefunctions are organized into a 1D array labeled by a "point index" $$\mu$$. Each point index corresponds to a pair of "grid indices", $$(m,n)$$, representing spatial coordinates on a 2D grid. We have to be careful not to mix up the two types of indices.

We will adopt the following conversion scheme between point indices and grid indices:

$\mu(m,n) = m N + n,\quad \mathrm{where}\; m \in \{ 0, \dots, M-1\}, \;\; n \in \{ 0, \dots, N-1\}.$

One good thing about this conversion scheme is that Scipy provides a reshape function which can convert a 2D array with grid indices $$(m,n)$$ into a 1D array with the point index $$\mu$$:

>>> a = array([[0,1,2],[3,4,5],[6,7,8]])
>>> a
array([[0, 1, 2],
[3, 4, 5],
[6, 7, 8]])
>>> b = reshape(a, (9))     # Reshape a into a 1D array of size 9
>>> b
array([0, 1, 2, 3, 4, 5, 6, 7, 8])


The reshape function can also convert a 1D back into the 2D array, in the right order:

>>> c = reshape(b, (3,3))   # Reshape b into a 2D array of size 3x3
>>> c
array([[0, 1, 2],
[3, 4, 5],
[6, 7, 8]])


Under point indices, the discretized derivatives take the following forms:

$\frac{\partial \psi}{\partial x}(\vec{r}_\mu)\;\, \approx \frac{1}{2h} \left(\psi_{\mu+N} - \psi_{\mu-N}\right)$

$\frac{\partial \psi}{\partial y}(\vec{r}_\mu)\;\, \approx \frac{1}{2h} \left(\psi_{\mu+1} - \psi_{\mu-1}\right)$

$\nabla^2\psi(\vec{r}_\mu) \approx \frac{1}{h^2} \left(\psi_{\mu+N} + \psi_{\mu+1} - 4\psi_{\mu} + \psi_{\mu-N} + \psi_{\mu-1}\right).$

The role of boundary conditions is left as an exercise. There are now two sets of boundaries, at $$m \in \{0,M-1\}$$ and $$n \in \{0, N-1\}$$. By examining the finite-difference equations along each boundary, we can (i) assign the right discretization coordinates and (ii) modify the finite-difference matrix elements to fit the boundary conditions. The details are slightly tedious to work out, but the logic is essentially the same as in the previously-discussed 1D cases.

7.3: Higher Dimensions is shared under a CC BY-SA 4.0 license and was authored, remixed, and/or curated by Y. D. Chong via source content that was edited to conform to the style and standards of the LibreTexts platform; a detailed edit history is available upon request.