Python for Quantum Mechanics: Numpy Linear Algebra

Python for Quantum Mechanics: Numpy Linear Algebra#

from IPython.display import YouTubeVideo
YouTubeVideo('Y4kJA38vwe8',width=700, height=400)

import numpy as np

Scalar Array/Matrix Operations#

We can perform operations on an array, element-wise, using a scalar.

A = np.array([[1,2],[3,4]], dtype='float')
print(A, '\n')

print(2*A,'\n')

print(A**2,'\n')

print(A/2,'\n')

print(A//2,'\n')

print(A+2,'\n')

print(A-2,'\n')

[[1. 2.]
 [3. 4.]] 

[[2. 4.]
 [6. 8.]] 

[[ 1.  4.]
 [ 9. 16.]] 

[[0.5 1. ]
 [1.5 2. ]] 

[[0. 1.]
 [1. 2.]] 

[[3. 4.]
 [5. 6.]] 

[[-1.  0.]
 [ 1.  2.]]

Array-Array Operations#

Recall the section on basic operators in 4.1 Numpy Basics covering the use of the basic mathematical operators ( +, -, *, /, ** ) between Numpy arrays

arr1 = np.array([1,2,3,4,5])
arr2 = np.array([6,7,8,9,10], dtype='float')

print('arr1 + arr2 = ', arr1+arr2)
print('arr1 - arr2 = ', arr1-arr2)
print('arr1 * arr2 = ', arr1*arr2)
print('arr1 / arr2 = ', arr1/arr2)
print('arr1 ** arr2 = ', arr1**arr2)

arr1 + arr2 =  [ 7.  9. 11. 13. 15.]
arr1 - arr2 =  [-5. -5. -5. -5. -5.]
arr1 * arr2 =  [ 6. 14. 24. 36. 50.]
arr1 / arr2 =  [0.16666667 0.28571429 0.375      0.44444444 0.5       ]
arr1 ** arr2 =  [1.000000e+00 1.280000e+02 6.561000e+03 2.621440e+05 9.765625e+06]

While the above arrays are vectors, we can perform the same operations on higher-dimesnional matrices in the same way:

A = np.array([[1,2],[3,4]])
print(A, '\n')

A = A*A
print(A)

[[1 2]
 [3 4]] 

[[ 1  4]
 [ 9 16]]

Matrix Algebra#

The fundamental matrix operations of matrix multiplication and vector products can be computed in two ways. The first uses the \(np.dot()\) function. In the second, we can cast the arrays as matrices using \(np.matrix()\), which changes the operators above ( +, -, *) to those of matrix algebra.

np.dot(array1,array2)

m = np.array([[1,2],[3,4]])
v = np.array([[1],[2]])

print(m, '\n')
print(v, '\n')

mv = np.dot(m,v)
print(mv)

[[1 2]
 [3 4]] 

[[1]
 [2]] 

[[ 5]
 [11]]

Alternatively, we can use the @ symbol, which matrix multiplies numpy arrays

print(m@v)

[[ 5]
 [11]]

np.vdot(row_vector1,col_vector2)

v1 = np.array([[1,2]])
v2 = np.array([[1],[3]])

print(v1, '\n')
print(v2, '\n')

vv = np.dot(v1,v2)
print(vv)

[[1 2]] 

[[1]
 [3]] 

[[7]]

np.matrix(array)

m = np.array([[1,2],[3,4]])
v = np.array([[1],[2]])

m = np.matrix(m)
v = np.matrix(v)

mv = m*v

print(mv)

[[ 5]
 [11]]

Or, written more compactly

m = np.matrix([[1,2],[3,4]])
v = np.matrix([[1],[2]])

print(m*v)

[[ 5]
 [11]]

The np.linalg Routine#

We can also calculate various matrix properties using \(np.linalg\). All of the \(linalg\) methods can be found at https://numpy.org/doc/stable/reference/routines.linalg.html

np.linalg.inv(matrix)

m = np.array([[1,2],[3,4]])

print(m, '\n')


minv = np.linalg.inv(m)
print(minv, '\n')

print(np.dot(m,minv))

[[1 2]
 [3 4]] 

[[-2.   1. ]
 [ 1.5 -0.5]] 

[[1.0000000e+00 0.0000000e+00]
 [8.8817842e-16 1.0000000e+00]]

np.linalg.det(matrix)

m = np.array([[1,2],[3,4]])

print(m, '\n')


mdet = np.linalg.det(m)
print(mdet, '\n')

[[1 2]
 [3 4]] 

-2.0000000000000004

np.linalg.eig(matrix)

This gives us the eigenvalues and eigenvectors of the matrix in a list. The first element is an array of eigenvalues, the second is a matrix array of the eigenvectors.

m = np.array([[1,2],[3,4]])

print(m, '\n')


meig = np.linalg.eig(m)
print(meig, '\n')

print('eigenvalues = ', meig[0])
print('eigenvectors = ', meig[1])

[[1 2]
 [3 4]] 

EigResult(eigenvalues=array([-0.37228132,  5.37228132]), eigenvectors=array([[-0.82456484, -0.41597356],
       [ 0.56576746, -0.90937671]])) 

eigenvalues =  [-0.37228132  5.37228132]
eigenvectors =  [[-0.82456484 -0.41597356]
 [ 0.56576746 -0.90937671]]

We can thus diagonalise the matrix, giving the eigenvalues along the diagonal:

m = np.matrix(m)
meigv = np.matrix(meig[1])
meigvinv = np.matrix(np.linalg.inv(meig[1]))

print(np.round(meigvinv*m*meigv,3))

[[-0.372  0.   ]
 [ 0.     5.372]]

np.linalg.matrix_power

We can find matrix powers:

m = np.array([[1,2],[3,4]])

print('m*m = ',np.dot(m,m),'\n')

print('m**2 = ',np.linalg.matrix_power(m,2),'\n')

print('m**10 = ',np.linalg.matrix_power(m,10),'\n')

m*m =  [[ 7 10]
 [15 22]] 

m**2 =  [[ 7 10]
 [15 22]] 

m**10 =  [[ 4783807  6972050]
 [10458075 15241882]]

Matrix Transformations#

We can transform matrices/arrays in various ways, finidng transposes and complex conjugates.

A = np.array([[1+1j,2+ 2j], [3+3j, 4+4j]])

print('A = ',A, '\n')


print('transpose(A) = ',np.transpose(A),'\n')
print('conjugate(A) = ',np.conjugate(A),'\n')

A =  [[1.+1.j 2.+2.j]
 [3.+3.j 4.+4.j]] 

transpose(A) =  [[1.+1.j 3.+3.j]
 [2.+2.j 4.+4.j]] 

conjugate(A) =  [[1.-1.j 2.-2.j]
 [3.-3.j 4.-4.j]] 

A = np.array([[1+1j,2+ 2j], [3+3j, 4+4j]])

print('transpose(A) = ',A.T,'\n')

transpose(A) =  [[1.+1.j 3.+3.j]
 [2.+2.j 4.+4.j]]

We can also use the attributes \(.T\) for \(np.array\), and \(.T\) and \(.H\) for \(np.matrix\)

A = np.matrix([[1+1j,2+ 2j], [3+3j, 4+4j]])

print('transpose(A) = ',A.T,'\n')
print('Hermitian(A) = ',A.H,'\n')

transpose(A) =  [[1.+1.j 3.+3.j]
 [2.+2.j 4.+4.j]] 

Hermitian(A) =  [[1.-1.j 3.-3.j]
 [2.-2.j 4.-4.j]]
Contents