Numeric arrays in Python

Numeric arrays in Python#

Examples demonstrating NumPy usage in ASE documentation.

Links to NumPy’s webpage:

ASE makes heavy use of an extension to Python called NumPy. The NumPy module defines an ndarray type that can hold large arrays of uniform multidimensional numeric data. An array is similar to a list or a tuple, but it is a lot more powerful and efficient.

Some examples from everyday ASE-life here …

import numpy as np

a = np.zeros((3, 2))
a[:, 1] = 1.0
a[1] = 2.0
print(a)
print(a.shape)
print(a.ndim)
[[0. 1.]
 [2. 2.]
 [0. 1.]]
(3, 2)
2

The conventions of numpy’s linear algebra package:

# --- Basic array example ---
a = np.zeros((3, 2))
a[:, 1] = 1.0
a[1] = 2.0
print('Array a:\n', a)
print('Shape:', a.shape)
print('Number of dimensions:', a.ndim)

# --- Linear algebra example ---
# Make a random Hermitian matrix
H = np.random.rand(6, 6) + 1.0j * np.random.rand(6, 6)
H = H + H.T.conj()

# Eigenvalues and eigenvectors
eps, U = np.linalg.eigh(H)

# Sort eigenvalues and corresponding eigenvectors
sorted_indices = eps.real.argsort()
eps = eps[sorted_indices]
U = U[:, sorted_indices]

# Verify diagonalization
print(
    'Check diagonalization:\n', np.dot(np.dot(U.T.conj(), H), U) - np.diag(eps)
)
print('All close?', np.allclose(np.dot(np.dot(U.T.conj(), H), U), np.diag(eps)))

# Check eigenvectors individually
print(
    'Eigenvector check (one column):',
    np.allclose(np.dot(H, U[:, 3]), eps[3] * U[:, 3]),
)
print('Eigenvector check (all):', np.allclose(np.dot(H, U), eps * U))
Array a:
 [[0. 1.]
 [2. 2.]
 [0. 1.]]
Shape: (3, 2)
Number of dimensions: 2
Check diagonalization:
 [[-6.66133815e-16+0.00000000e+00j -1.73472348e-18+2.98372438e-16j
  -3.19189120e-16-1.38777878e-16j -3.26128013e-16-1.11022302e-16j
   3.19582427e-16-1.60603581e-16j -2.92864862e-16+1.01597601e-17j]
 [-3.64291930e-17-2.28983499e-16j  2.22044605e-16+6.93889390e-17j
  -5.55111512e-17+4.30211422e-16j  5.68989300e-16-1.11022302e-16j
  -6.38650453e-16+1.10739173e-16j  4.09580400e-16+8.31759666e-17j]
 [-2.11636264e-16+1.04083409e-16j  4.85722573e-17-3.85108612e-16j
   3.60822483e-16-1.38777878e-17j  1.42247325e-16+2.29850861e-16j
   8.97941782e-17-1.24273607e-16j  1.59256842e-16-9.77935701e-17j]
 [-3.33066907e-16+5.55111512e-17j  5.82867088e-16+1.66533454e-16j
   8.32667268e-17-1.17961196e-16j  2.22044605e-16+0.00000000e+00j
   1.04809334e-17-1.53803417e-16j -2.12515629e-16+1.83659084e-17j]
 [ 3.33066907e-16+1.66533454e-16j -6.38378239e-16-1.80411242e-16j
   1.66533454e-16+0.00000000e+00j  2.77555756e-17+6.93889390e-17j
   1.77635684e-15+4.71381579e-17j -3.63393744e-16-3.16281843e-16j]
 [ 1.11022302e-16+0.00000000e+00j  3.33066907e-16-2.77555756e-17j
   5.55111512e-17+5.55111512e-17j  0.00000000e+00-5.55111512e-17j
  -7.11093185e-16+2.77612627e-16j  4.44089210e-15+4.29859064e-17j]]
All close? True
Eigenvector check (one column): True
Eigenvector check (all): True

The rules for multiplying 1D arrays with 2D arrays:

  • 1D arrays and treated like shape (1, N) arrays (row vectors).

  • left and right multiplications are treated identically.

  • A length \(m\) row vector can be multiplied with an \(n \times m\) matrix, producing the same result as if replaced by a matrix with \(n\) copies of the vector as rows.

  • A length \(n\) column vector can be multiplied with an \(n \times m\) matrix, producing the same result as if replaced by a matrix with \(m\) copies of the vector as columns.

Thus, for the arrays below:

# --- 1D vs 2D multiplication rules ---
M = np.arange(5 * 6).reshape(5, 6)  # A matrix of shape (5, 6)
v5 = np.arange(5) + 10  # A vector of length 5
v51 = v5[:, None]  # Column vector (5, 1)
v6 = np.arange(6) - 12  # A vector of length 6
v16 = v6[None, :]  # Row vector (1, 6)

The following identities hold:

v6 * M == v16 * M == M * v6 == M * v16 == M * v16.repeat(5, 0)
v51 * M == M * v51 == M * v51.repeat(6, 1)

The same rules apply for adding and subtracting 1D arrays to from 2D arrays.

# Identities
print('v6 * M == M * v6?', np.allclose(v6 * M, M * v6))
print('v16 * M == M * v16?', np.allclose(v16 * M, M * v16))
print('v51 * M == M * v51?', np.allclose(v51 * M, M * v51))
v6 * M == M * v6? True
v16 * M == M * v16? True
v51 * M == M * v51? True

Gallery generated by Sphinx-Gallery