SciPy - expm_frechet() Function

---

---

The scipy.linalg.expm_frechet method determines how sensitive the matrix exponential expm(A) is to changes in matrix 'A', given a direction matrix 'E'. The matrix exponential expm(A) is a mathematical transformation that is commonly used in systems such as physics and engineering to simulate changes over time. The Frechet derivative depicts how much the outcome of expm(A) varies when A is slightly modified by 'E'.

This method is very beneficial when evaluating systems that need understanding the impact of minor alterations or faults, such as stability analysis or control systems. It can also calculate the matrix exponential expm(A) and the Frechet derivative at the same time, which saves time and improves efficiency.

Errors may occur if the input matrices 'A' and 'E' are not square or if their dimensions differ. Ill-conditioned matrices can cause numerical errors in the Frechet derivative or exponential computation. Validation of input ensures accurate findings.

The expm_frechet() approach can be combined with matrix decomposition techniques such as SVD or eigenvalue analysis to study the sensitivity of specific components in dynamic systems. It is also used alongside with logm() and sqrtm() to do complex computations regarding system stability and behavior over time.

Syntax

The syntax for the SciPy expm_frechet() method is as follows −

.expm_frechet(A, E, method=None, compute_expm=True, check_finite=True)

Parameters

This method accepts the following parameters −

• A (N, N) array_like − Input square matrix (nn), real or complex.

• E (N, N) array_like − Perturbation matrix of the same dimensions as A.

• method (str, optional) − Specifies the computation method. Options include 'SPS' (Scaling-Power-Series) or 'blockEnlarge'.

• compute_expm (bool, optional) − If True, the method also computes the matrix exponential expm(A).

• check_finite (bool, optional) − If True, the input matrices A and E are checked for NaNs or infinite values. If these checks are unnecessary for performance reasons, set this to False.

Return Value

• expm_A (ndarray) − The matrix exponential of A, computed as expm(A). This represents the exponential transformation applied to the matrix A.

• expm_frechet_AE (ndarray) − The Frechet derivative of the matrix exponential expm(A) in the direction of E. This indicates how sensitive expm(A) is to perturbations along.

Example 1

The expm_frechet method calculates both the matrix exponential expm(A) and its Frechet derivative in the direction 'E'.

In the following code matrix, 'A' represents the system while 'E' denotes perturbations. The procedure returns expm(A), which depicts the system's exponential behavior, as well as the Frechet derivative, which assesses the influence of changes in A on expm.

Larger values in the output highlight areas where the system is more sensitive to changes, while smaller values indicate regions of stability.

import numpy as np
from scipy.linalg import expm_frechet

# Input: Matrices A and E
A = np.array([[1, 2], 
              [3, 4]])
E = np.array([[0, 1], 
              [1, 0]])
expm_A, frechet_AE = expm_frechet(A, E, compute_expm=True)

print("Matrix Exponential (expm(A)):\n", expm_A)
print("Frechet Derivative (Sensitivity to E):\n", frechet_AE)

When we run above program, it produces following result

Matrix Exponential (expm(A)):
 [[ 51.9689562   74.73656457]
 [112.10484685 164.07380305]]
Frechet Derivative (Sensitivity to E):
 [[ 61.30566146  80.18834128]
 [101.59837077 125.53574995]]

Example 2

In the following code, two random 33 matrices, A and E, are defined. Then, a 66 block matrix M is constructed by placing A in the diagonal blocks and E in the upper-right off-diagonal block. This structure captures both the behavior of the system (A) and its sensitivity to perturbation (E).

The matrix exponential of M is calculated and separated into diagonal and off-diagonal components. The diagonal block corresponds to expm(A), whereas the off-diagonal block represents the Frechet derivative expm_frechet_AE. This tests the correctness of expm_frechet by comparing it to the block matrix technique.

import numpy as np
from scipy import linalg

rng = np.random.default_rng()
A = rng.standard_normal((3, 3))
E = rng.standard_normal((3, 3))

expm_A, expm_frechet_AE = linalg.expm_frechet(A, E)

# Construct a 6x6 block matrix containing A and E
M = np.zeros((6, 6))
M[:3, :3] = A
M[:3, 3:] = E
M[3:, 3:] = A

# Compute the exponential of the block matrix
expm_M = linalg.expm(M)

# Verify results
print("Matrix Exponential (expm(A)) matches block diagonal:", np.allclose(expm_A, expm_M[:3, :3]))
print("Frechet Derivative matches block off-diagonal:", np.allclose(expm_frechet_AE, expm_M[:3, 3:]))

Following is an output of the above code

Matrix Exponential (expm(A)) matches block diagonal: True
Frechet Derivative matches block off-diagonal: True

Example 3

The method parameter in scipy.linalg.expm_frechet allows the usage of algorithms such as 'blockEnlarge' to efficiently compute the Frechet derivative in structured matrices, notably diagonal or block matrices.

In the example below, we created a 3x3 diagonal matrix and a perturbation matrix E and used them as inputs. The 'blockEnlarge' approach optimizes processing, by using A's structure to reduce computational cost. The resulting Frechet derivative reflects E's exponential influence on A's matrix.

import numpy as np
from scipy.linalg import expm_frechet

A = np.array([[1, 0, 0], 
              [0, 2, 0], 
              [0, 0, 3]])
E = np.array([[0, 1, 0], 
              [1, 0, 1], 
              [0, 1, 0]])

# Compute the matrix exponential and Frechet derivative using 'blockEnlarge' method
expm_A, frechet_AE = expm_frechet(A, E, method='blockEnlarge')
print("Matrix Exponential (expm(A)):\n", expm_A)
print("\nFrechet Derivative (expm_frechet_AE):\n", frechet_AE)

Output of the above code is as follows

Matrix Exponential (expm(A)):
 [[ 2.71828183  0.          0.        ]
 [ 0.          7.3890561   0.        ]
 [ 0.          0.         20.08553692]]

Frechet Derivative (expm_frechet_AE):
 [[ 0.          4.67077427  0.        ]
 [ 4.67077427  0.         12.69648082]
 [ 0.         12.69648082  0.        ]]