pyLOM.math#
Module contents#
- pyLOM.math.MAE(A: ndarray, B: ndarray) float #
Compute mean absolute error between A and B
- Parameters:
A (np.ndarray)
B (np.ndarray)
- Returns:
Mean absolute error.
- Return type:
(float)
- pyLOM.math.MRE_array(A: ndarray, B: ndarray, axis: int = 1) ndarray #
Mean relative error computed along a certain axis of the array.
- Parameters:
A (np.ndarray) – original field.
B (np.ndarray) – field which we want to compute the MRE of.
axis (int, optional) – along which axis the MRE will be computed (default
1
).
- Returns:
Mean relative error.
- Return type:
(np.ndarray)
- pyLOM.math.RMSE(A: ndarray, B: ndarray, relative: bool = True) float #
Compute the root mean square error between A and B
- pyLOM.math.argsort(v: ndarray) ndarray #
Returns the indices that sort a vector
- Parameters:
v (np.ndarray) – Vector v (M,)
- Result:
np.ndarray: Indices that sort v (M,)
- pyLOM.math.cellCenters(xyz: ndarray, conec: ndarray) ndarray #
Compute the cell centers given a list of elements.
- Parameters:
xyz (np.ndarray) – node positions
conec (np.ndarray) – connectivity array
- Returns:
center positions
- Return type:
np.ndarray
- pyLOM.math.cell_adjacency(edge_dict) dict #
Build a dictionary that maps each cell to its neighbors.
- pyLOM.math.cholesky(A: ndarray) ndarray #
Returns the Cholesky decompositon of A
- Parameters:
A (np.ndarray) – Matrix A (M,N)
- Result:
np.ndarray: Cholesky factorization of A (M,N)
- pyLOM.math.compute_truncation_residual(S, r)#
Compute the truncation residual. r must be a float precision (r<1) where:
r > 0: target residual
r < 0: fraction of cumulative energy to retain
- pyLOM.math.conj(A: ndarray) ndarray #
Conjugates complex number A
- Parameters:
A (np.ndarray) – Vector, matrix or number A
- Result:
np.ndarray: Conjugate of A
- pyLOM.math.data_splitting(Nt: int, mode: str, seed: int = -1)#
Generate random training, validation and test masks for a dataset of Nt samples.
- Parameters:
- Returns:
List of arrays containing the identifiers of the training, validation and test samples.
- Return type:
[(np.ndarray), (np.ndarray), (np.ndarray)]
- pyLOM.math.diag(A: ndarray) ndarray #
If A is a matrix it returns its diagonal, if its a vector it returns a diagonal matrix with A in its diagonal
- Parameters:
A (np.ndarray) – Matrix A (M,N)
- Result:
np.ndarray: Diagonal of A (M,)
- pyLOM.math.edge_to_cells(conec: ndarray) dict #
Build a dictionary that maps each edge to the cells that share it.
- Parameters:
conec (np.ndarray) – connectivity array
- Returns:
edges to cells connectivity dictionary
- Return type:
- pyLOM.math.eigen(A: ndarray) ndarray #
Eigenvalues and eigenvectors.
GPU implementation of this algoritm is still slow and thus will be executed purely on CPU level until cupy implements linalg.eig
- Parameters:
A (np.ndarray) – Matrix A (M,N)
- Result:
np.ndarray: the real eigenvalues, real(M) np.ndarray: the imaginary eigenvalues, imag(M) np.ndarray: the right eigenvectors, vecs(M,M)
- pyLOM.math.energy(original, rec)#
Compute reconstruction energy as in: Eivazi, H., Le Clainche, S., Hoyas, S., & Vinuesa, R. (2022). Towards extraction of orthogonal and parsimonious non-linear modes from turbulent flows. Expert Systems with Applications, 202, 117038. https://doi.org/10.1016
- pyLOM.math.euclidean_d(X)#
Compute Euclidean distances between simulations.
- In:
X: NxM Data matrix with N points in the mesh for M simulations
- Returns:
MxM distance matrix
- Return type:
D
- pyLOM.math.fft(t: ndarray, y: ndarray, equispaced: bool = True) ndarray #
Compute the PSD of a signal y. For non equispaced time samples the nfft package is required.
- Parameters:
t (numpy.ndarray) – time vector.
y (numpy.ndarray) – signal vector.
equispaced (bool) – whether the samples in the time vector are equispaced or not.
- Returns:
frequency. numpy.ndarray: power density spectra.
- Return type:
- pyLOM.math.find_random_sensors(bounds: ndarray, xyz: ndarray, nsensors: int)#
Generate a set of random points inside a bounding box and find the closest grid points to them
- Parameters:
bounds (np.ndarray) – bounds of the box in the following format: np.array([xmin, xmax, ymin, ymax, zmin, zmax])
xyz (np.ndarray) – coordinates of the grid
nsensors (int) – number of sensors to generate
- Returns:
array with the indices of the points
- Return type:
np.ndarray
- pyLOM.math.fix_normals_coherence(normals, edge_dict, adjacency, num_cells) ndarray #
Ensure that the normals of the cells are coherent. (i.e. they point all in the same direction). In:
normals: Array of normals of the cells
edge_dict: Dictionary mapping edges to cells sharing that edge.
adjacency: Dictionary mapping cells to their neighbors.
num_cells: Number of cells in the mesh
- Returns:
Array of normals of the cells
- Return type:
normals
- pyLOM.math.flip(A: ndarray) ndarray #
Returns the pointwise conjugate of A
This function is not implemented in the compiled layer and will raise an error if used
- Parameters:
A (np.ndarray) – Matrix A (M,N)
- Result:
np.ndarray: Flipped version of A (M,N)
- pyLOM.math.hammwin(N: int) ndarray #
Hamming windowing
- Parameters:
N (int) – Number of steps.
- Returns:
Hamming windowing.
- Return type:
- pyLOM.math.init_qr_streaming(Ai, r, q, seed=None)#
Ai(m,n) data matrix dispersed on each processor. r target number of modes
Qi(m,r) B (r,n)
- pyLOM.math.inv(A: ndarray) ndarray #
Computes the inverse matrix of A
- Parameters:
A (np.ndarray) – Matrix A (M,N)
- Result:
np.ndarray: Inverse of A (M,N)
- pyLOM.math.least_squares(A, b)#
Least squares regression (A^T * A)^-1 * A^T * b
- pyLOM.math.matmul(A: ndarray, B: ndarray) ndarray #
Matrix multiplication C = A x B
- Parameters:
A (np.ndarray) – Matrix A (M,Q)
B (np.ndarray) – Matrix B (Q,N)
- Result:
np.ndarray: Resulting matrix C (M,N)
- pyLOM.math.matmulp(A: ndarray, B: ndarray) ndarray #
Matrix multiplication in parallel C = A x B
A and B are distributed along the processors and C is the same for all of them
- Parameters:
A (np.ndarray) – Matrix A (M,Q)
B (np.ndarray) – Matrix B (Q,N)
- Result:
np.ndarray: Resulting matrix C (M,N)
- pyLOM.math.norm_variance(X: ndarray, X_mean: ndarray, X_var: ndarray) ndarray #
Normalizes the snapshot matrix X(m,n) with its variance
- Parameters:
X (numpy.ndarray) – Snapshot matrix (m,n).
X_mean (numpy.ndarray) – Mean of the snapshot matrix (m,n).
X_var (numpy.ndarray) – Variance of the snapshot matrix (m,n).
- Returns:
Snapshot matrix normalized by its variance(m,n).
- Return type:
- pyLOM.math.normals(xyz: ndarray, conec: ndarray) ndarray #
Compute the cell normals given a list of elements.
- Parameters:
xyz (np.ndarray) – node positions
conec (np.ndarray) – connectivity array
- Returns:
cell normals
- Return type:
np.ndarray
- pyLOM.math.polar(real: ndarray, imag: ndarray) ndarray #
Present a complex number in its polar form given its real and imaginary part
- Parameters:
real (np.ndarray) – the real component
imag (np.ndarray) – the imaginary component
- Result:
np.ndarray: the modulus np.ndarray: the argument
- pyLOM.math.qr(A)#
- QR factorization using Lapack
Q(m,n) is the Q matrix R(n,n) is the R matrix
- pyLOM.math.r2(A: ndarray, B: ndarray) float #
Compute r2 score between A and B
- Parameters:
A (np.ndarray)
B (np.ndarray)
- Returns:
r2 score .
- Return type:
(float)
- pyLOM.math.randomized_qr(Ai, r, q, seed=-1)#
Ai(m,n) data matrix dispersed on each processor. r target number of modes
Qi(m,r) B (r,n)
- pyLOM.math.randomized_svd(Ai, r, q, seed=-1)#
Ai(m,n) data matrix dispersed on each processor. r target number of modes
Ui(m,n) POD modes dispersed on each processor (must come preallocated). S(n) singular values. VT(n,n) right singular vectors (transposed).
- pyLOM.math.ridge_regresion(A, b, lam)#
Ridge regression
- pyLOM.math.subtract_mean(X: ndarray, X_mean: ndarray) ndarray #
Computes out(m,n) = X(m,n) - X_mean(m) where m is the spatial coordinates and n is the number of snapshots.
- Parameters:
X (numpy.ndarray) – Snapshot matrix (m,n).
X_mean (numpy.ndarray) – Averaged snapshot matrix (m,)
- Returns:
Snapshot matrix without the average(m,n).
- Return type:
- pyLOM.math.svd(A, method='gesdd')#
- Single value decomposition (SVD) using numpy.
U(m,n) are the POD modes. S(n) are the singular values. V(n,n) are the right singular vectors.
- pyLOM.math.temporal_mean(X: ndarray) ndarray #
Temporal mean of matrix X(m,n) where m is the spatial coordinates and n is the number of snapshots.
- Parameters:
X (numpy.ndarray) – Snapshot matrix (m,n).
- Returns:
Averaged snapshot matrix (m,).
- Return type:
- pyLOM.math.temporal_variance(X: ndarray, X_mean: ndarray) ndarray #
Variance of matrix X(m,n) on time where m is the spatial coordinates and n is the number of snapshots.
- Parameters:
X (numpy.ndarray) – Snapshot matrix (m,n).
X_mean (numpy.ndarray) – Mean of the snapshot matrix (m,n).
- Returns:
Variance of the snapshot matrix (m,).
- Return type:
- pyLOM.math.time_delay_embedding(X, dimension=50)#
Extract time-series of lenght equal to lag from longer time series in data, whose dimension is (number of time series, sequence length, data shape) Inputs: [Points, Time] Output: [Points, Time, Delays]
- pyLOM.math.transpose(A: ndarray) ndarray #
Transposed of matrix A
- Parameters:
A (np.ndarray) – Matrix to be transposed
- Results
np.ndarray: Transposed matrix
- pyLOM.math.tsqr(Ai)#
- Parallel QR factorization of a real array using Lapack
Q(m,n) is the Q matrix R(n,n) is the R matrix
- pyLOM.math.tsqr_svd(Ai)#
Single value decomposition (SVD) using TSQR algorithm from J. Demmel, L. Grigori, M. Hoemmen, and J. Langou, ‘Communication-optimal Parallel and Sequential QR and LU Factorizations’, SIAM J. Sci. Comput., vol. 34, no. 1, pp. A206–A239, Jan. 2012,
doi: 10.1137/080731992.
Ai(m,n) data matrix dispersed on each processor.
Ui(m,n) POD modes dispersed on each processor (must come preallocated). S(n) singular values. VT(n,n) right singular vectors (transposed).
- pyLOM.math.update_qr_streaming(Ai, Q1, B1, Yo, r, q)#
Ai(m,n) data matrix dispersed on each processor. r target number of modes
Qi(m,r) B (r,n)
- pyLOM.math.vandermonde(real: ndarray, imag: ndarray, m: int, n: int) ndarray #
Builds a Vandermonde matrix of (m x n) with the real and imaginary parts of the eigenvalues
- Parameters:
- Result:
np.ndarray: the Vandermonde matrix
- pyLOM.math.vandermondeTime(real: ndarray, imag: ndarray, m: int, time: ndarray) ndarray #
Builds a Vandermonde matrix of (m x n) with the real and imaginary parts of the eigenvalues
- Parameters:
real (np.ndarray) – the real component
imag (np.ndarray) – the imaginary component
m (int) – number of rows for the Vandermode matrix
time (np.ndarray) – the time vector
- Result:
np.ndarray: the Vandermonde matrix
- pyLOM.math.vecmat(v: ndarray, A: ndarray) ndarray #
Vector times a matrix C = v x A
- Parameters:
v (np.ndarray) – Vector v (M,)
A (np.ndarray) – Matrix A (M,N)
- Result:
np.ndarray: Resulting matrix C (M,N)
- pyLOM.math.vector_mean(v: ndarray, start: int = 0) float #
Mean of a vector
- Parameters:
v (np.ndarray) – a vector
start (int) – position of the vector where to start the mean
- Result:
float: mean of the vector
- pyLOM.math.vector_norm(v: ndarray, start: int = 0) float #
L2 norm of a vector
- Parameters:
v (np.ndarray) – a vector
start (int) – position of the vector where to start the norm
- Result:
float: norm of the vector
- pyLOM.math.vector_sum(v: ndarray, start: int = 0) float #
Sum of a vector
- Parameters:
v (np.ndarray) – a vector
start (int) – position of the vector where to start the sum
- Result:
float: sum of the vector
- pyLOM.math.wall_normals(nodes_idx, nodes_xyz, surf_normal) list #
Compute the unitary normals to the cell walls (only for 2D cells). Example: For a triangle, the wall normals are the three vectors normal to each of the sides. The wall normals are always contained in the cell plane, thus are orthogonal themselves to the cell surface normal. As a convention, wall normals are always pointing outwards the cell.
- In:
nodes_idx: List or array of the node indices of the cell
nodes_xyz: List or array of the node coordinates of the cell
surf_normal: Normal to the plane of the cell
- Returns:
List of graph edges representing the element walls (node indices) - wall_normals: List of the unitary wall normals
- Return type:
cell_edges