qdyn.linalg module¶

Linear algebra helper routines

Summary¶

Functions:

`banded_to_full`	Convert a rectangular matrix in the Lapack “banded” format (http://www.netlib.org/lapack/lug/node124.html) into a (square) full matrix.
`choose_sparsity_model`	Estimate which might be the best storage format for the given matrix.
`inner`	Calculate the inner product of the two vectors or matrices v1, v2.
`is_hermitian`	Check, if a matrix is Hermitian.
`iscomplexobj`	Check whether the (multidimensional x object) has a type that allows for complex entries.
`norm`	Calculate the norm of a vector or matrix v, matching the inner product defined in the inner routine.
`trace`	Return the trace of the given matrix.
`tril`	Return the lower triangle of the given matrix.
`triu`	Return the upper triangle of the given matrix.
`vectorize`	Return vectorization of multi-dimensional array a.

Reference¶

qdyn.linalg.inner(v1, v2)[source]¶

Calculate the inner product of the two vectors or matrices v1, v2.

For vectors, the inner product is the standard Euclidian inner product.
For matrices, the innner product is the Hilbert-Schmidt overlap.

Note that the inner product of the vectorization of two matrices is the same as the inner product of the original matrices, and that two m x 1 matrices have the same inner product as the corresponding m-dimensional vectors.

The inner routine corresponds to overlap in QDYN.

Parameters

v1 (numpy.ndarray) – First argument. The conjugate transposed of v1 is taken before the calculation of the inner product.
v2 (numpy.ndarray) – Second argument.

Returns

Inner product of v1 and v2.

Return type

float

Examples

>>> v1 = np.array([1.0, 1.0j, 1.0, 1.0j])
>>> v2 = np.array([1.0j, 1.0j, 1.0j, 1.0j])
>>> inner(v1, v2)
(2+2j)
>>> m1 = v1.reshape((2,2))
>>> m2 = v2.reshape((2,2))
>>> inner(m1, m2)
(2+2j)
>>> m1 = v1.reshape((4,1))
>>> m2 = v2.reshape((4,1))
>>> inner(m1, m2)
(2+2j)
>>> m1 = v1.reshape((4,1), order='F')
>>> m2 = v2.reshape((4,1), order='F')
>>> inner(m1, m2)
(2+2j)
>>> m1 = v1.reshape((1,4), order='F')
>>> m2 = v2.reshape((1,4), order='F')
>>> inner(m1, m2)
(2+2j)

qdyn.linalg.trace(m)[source]¶

Return the trace of the given matrix.

Parameters: m (list, numpy.ndarray) – Input array from which the diagonals are taken.
Returns: If m is a matrix, the sum along the diagonal is returned. If m has larger dimensions, then an array of sums along diagonals is returned.
Return type: float

qdyn.linalg.norm(v)[source]¶

Calculate the norm of a vector or matrix v, matching the inner product defined in the inner routine. An algorithm like Gram-Schmidt-Orthonormalization will only work if the choice of norm and inner product are compatible.

If v is a vector, the norm is the 2-norm (i.e. the standard Euclidian vector norm).
If v is a matrix, the norm is the Hilbert-Schmidt (aka Frobenius) norm. Note that the HS norm of a matrix is identical to the 2-norm of any vectorization of that matrix (e.g. writing the columns of the matrix underneat each other). Also, the HS norm of the m x 1 matrix is the same as the 2-norm of the equivalent m-dimensional vector.

Parameters: v (numpy.ndarray, scipy.sparse.spmatrix, qutip.Qobj) – Input vector or matrix.
Returns: Norm of v.
Return type: float

qdyn.linalg.vectorize(a, order='F')[source]¶

Return vectorization of multi-dimensional array a.

Parameters

a (list, numpy.ndarray, qutip.Qobj) – Array to be vectorized.
order (str) – One of ‘C’, ‘F’, ‘A’. Read the elements of a using this index order, and place the elements into the reshaped array using this index order. ‘C’ means to read / write the elements using C-like index order, with the last axis index changing fastest, back to the first axis index changing slowest. ‘F’ means to read / write the elements using Fortran-like index order, with the first index changing fastest, and the last index changing slowest. Note that the ‘C’ and ‘F’ options take no account of the memory layout of the underlying array, and only refer to the order of indexing. ‘A’ means to read / write the elements in Fortran-like index order if a is Fortran contiguous in memory, C-like order otherwise.

Returns

The input array, but with all or a subset of the dimensions of length 1 removed. This is always a itself or a view into a.

Return type

numpy.ndarray

Examples

>>> a = np.array([1,2,3,4])
>>> vectorize(a)
array([1, 2, 3, 4])

>>> a = np.array([[1,2],[3,4]])
>>> vectorize(a)
array([1, 3, 2, 4])
>>> vectorize(a, order='C')
array([1, 2, 3, 4])

qdyn.linalg.is_hermitian(matrix)[source]¶

Check, if a matrix is Hermitian.

matrix can be a numpy array, a scipy sparse matrix, or a qutip.Qobj instance.

Parameters: matrix (list, numpy.ndarray, qutip.Qobj) – Input array.
Returns: Returns True if matrix is Hermitian, False otherwise.
Return type: bool

Examples

>>> m = np.array([[0, 1j], [-1j, 1]])
>>> is_hermitian(m)
True

>>> m = np.array([[0, 1j], [-1j, 1j]])
>>> is_hermitian(m)
False

>>> m = np.array([[0, -1j], [-1j, 1]])
>>> is_hermitian(m)
False

>>> from scipy.sparse import coo_matrix
>>> m  = coo_matrix(np.array([[0, 1j], [-1j, 0]]))
>>> is_hermitian(m)
True

qdyn.linalg.iscomplexobj(x)[source]¶

Check whether the (multidimensional x object) has a type that allows for complex entries.

Parameters: x – Multidimensional array like object.
Returns: Returns True, if x allows for complex entries, False otherwise.
Return type: bool

Examples

>>> iscomplexobj(1)
False
>>> iscomplexobj(1+0j)
True
>>> iscomplexobj([3, 1+0j, True])
True
>>> iscomplexobj(np.array([3, 1j]))
True
>>> iscomplexobj(scipy.sparse.csr_matrix([[1, 2], [4, 5]]))
False
>>> iscomplexobj(scipy.sparse.csr_matrix([[1, 2], [4, 5j]]))
True

qdyn.linalg.choose_sparsity_model(matrix)[source]¶

Estimate which might be the best storage format for the given matrix.

Parameters: matrix (numpy.matrix, scipy.sparse.spmatrix, qutip.Qobj) – Input matrix, which must be quadratic.
Returns: Returns one of ‘full’, ‘banded’, ‘dia’, or ‘indexed’.
Return type: str

Examples

>>> m = scipy.sparse.random(100, 100, 0.01)
>>> choose_sparsity_model(m)
'indexed'

>>> m = scipy.sparse.random(100, 100, 0.5)
>>> choose_sparsity_model(m)
'full'

>>> m = np.diag(np.ones(20))
>>> choose_sparsity_model(m)
'banded'

>>> m = np.diag(np.zeros(20))
>>> m[2,2] = m[10,10] = m[11,11] = 1
>>> choose_sparsity_model(m)
'indexed'

>>> td_data = np.array([[1, 2, 3, 4]]).repeat(3, axis=0)
>>> off = np.array([0, -1, 1])
>>> m = scipy.sparse.dia_matrix((td_data,off), shape=(5,5)).toarray()
>>> choose_sparsity_model(m) # should eventually be 'dia'
'indexed'

>>> m = scipy.sparse.dia_matrix((td_data,off), shape=(20,20)).toarray()
>>> m[19,19] = 1
>>> choose_sparsity_model(m)
'indexed'

qdyn.linalg.triu(matrix)[source]¶

Return the upper triangle of the given matrix.

If matrix’ is a numpy object or scipy sparse matrix, the returned matrix will have the same type as the input `matrix.
If matrix is a QuTiP operator, the type is not preserved: the result is equivalent to triu(matrix.data).

Parameters: matrix (numpy.ndarray, scipy.sparse.spmatrix, qutip.Qobj) – Input array.
Returns: Return a copy of the matrix with the elements below the diagonal zeroed.
Return type: numpy.ndarray, scipy.sparse.spmatrix
Raises: TypeError – If matrix has an invalid type.

qdyn.linalg.tril(matrix)[source]¶

Return the lower triangle of the given matrix.

If matrix’ is a numpy object or scipy sparse matrix, the returned matrix will have the same type as the input `matrix.
If matrix is a QuTiP operator, the type is not preserved: the result is equivalent to tril(matrix.data).

Parameters: matrix (numpy.ndarray, scipy.sparse.spmatrix, qutip.Qobj) – Input array.
Returns: Return a copy of the matrix with the elements above the diagonal zeroed.
Return type: numpy.ndarray, scipy.sparse.spmatrix
Raises: TypeError – If matrix has an invalid type.

qdyn.linalg.banded_to_full(banded, n, kl, ku, mode)[source]¶

Convert a rectangular matrix in the Lapack “banded” format (http://www.netlib.org/lapack/lug/node124.html) into a (square) full matrix.

Parameters

banded (numpy.ndarray) – Rectangular matrix in banded format
n (int) – The dimension of the (full) matrix
kl (int) – The number of lower diagonals (kl=0 for diagonal/upper-triangular matrix)
ku (int) – The number of upper diagonals (ku=0 for diagonal/lower-triangular matrix)
mode (str) – On of ‘g’, ‘h’, ‘s’, ‘t’ corresponding to “general”, “Hermitian”, “symmetric”, and ‘triangular’. The values ‘g’ and ‘s’ are dequivalent, except that for ‘s’ iether kl or ku must be zero. For Hermitian or symmetric storage, exactly one of kl, ku must be zero. Which one determines whether banded is assumed to contain the data for the upper or lower triangle

Returns

Numpy array of same type as banded.

Return type

numpy.ndarray