:py:mod:`FELiCS.Solvers.LinearSolver` ===================================== .. py:module:: FELiCS.Solvers.LinearSolver Module Contents --------------- Classes ~~~~~~~ .. autoapisummary:: FELiCS.Solvers.LinearSolver.LinearSolver FELiCS.Solvers.LinearSolver.ResolventOperator Attributes ~~~~~~~~~~ .. autoapisummary:: FELiCS.Solvers.LinearSolver.logger .. py:data:: logger .. py:class:: LinearSolver Linear algebra utilities using PETSc and SLEPc. Provides static methods for solving generalized eigenvalue problems (GEVP), singular value decompositions (SVD), and linear systems efficiently using PETSc and SLEPc. **Initialize the LinearSolver object** :param None: .. admonition:: Notes Consists of only static methods that don't need an instance of this class. Call the methods via "LinearSolver.method()". .. admonition:: Examples Solving a linear system: >>> from petsc4py import PETSc >>> A = PETSc.Mat().create() >>> b = PETSc.Vec().create() >>> x = LinearSolver.solveEquationSystem(A, b) Solving a generalized eigenvalue problem: >>> eigVals, eigVecs, error = LinearSolver.solveGeneralEigenproblem( ... A, B, sigma=0.0, nev=5) .. py:method:: solveGeneralEigenproblem(A, B, sigma, nev, tol=1e-12, max_it=200, adjoint=False, isForEigenProblem=True) :staticmethod: Solve the generalized eigenvalue problem (GEVP) using SLEPc and PETSc. :param A: Matrix defining the eigenvalue problem. :type A: :py:class:`PETSc.Mat` :param B: Second matrix in the generalized eigenvalue problem. :type B: :py:class:`PETSc.Mat` :param sigma: Initial guess for the eigenvalue. :type sigma: :py:class:`float` :param nev: Number of eigenvalues to compute. :type nev: :py:class:`int` :param tol: Tolerance for convergence, by default 1.e-12. :type tol: :py:class:`float`, *optional* :param max_it: Maximum number of iterations, by default 200. :type max_it: :py:class:`int`, *optional* :param adjoint: Compute adjoint eigenvalues and eigenvectors if True, by default False. :type adjoint: :py:class:`bool`, *optional* :param isForEigenProblem: Flag indicating whether this is for an eigenproblem, by default True. :type isForEigenProblem: :py:class:`bool`, *optional* :returns: * **eigVals** (:py:class:`numpy.ndarray`) -- Computed eigenvalues. * **eigVecs** (:py:class:`numpy.ndarray`) -- Computed eigenvectors. * **error** (:py:class:`numpy.ndarray`) -- Relative residuals for the eigenvalue problem. .. py:method:: solveSVDOfResolvent(resolventOperator, nev, tol=1e-16, max_it=200) :staticmethod: Solve the singular value decomposition (SVD) of a resolvent operator. :param resolventOperator: Resolvent operator object for which the SVD is computed. :type resolventOperator: :py:class:`ResolventOperator` :param nev: Number of singular values to compute. :type nev: :py:class:`int` :param tol: Tolerance for convergence, by default 1.e-16. :type tol: :py:class:`float`, *optional* :param max_it: Maximum number of iterations, by default 200. :type max_it: :py:class:`int`, *optional* :returns: * **eigVals** (:py:class:`numpy.ndarray`) -- Singular values of the operator. * **eigVecs** (:py:class:`numpy.ndarray`) -- Singular vectors corresponding to the singular values. .. py:method:: solveEquationSystem(A, b, destroy=False) :staticmethod: Solve a linear system of equations Ax = b using PETSc. :param A: Matrix representing the linear system. :type A: :py:class:`PETSc.Mat` :param b: Right-hand side of the equation. :type b: :py:class:`PETSc.Vec` :param destroy: Whether to destroy the matrix and vector after solving, by default False. :type destroy: :py:class:`bool`, *optional* :returns: **x** -- Solution vector. :rtype: :py:class:`numpy.ndarray` .. py:method:: solveTransposeEquationSystem(A, b, destroy=False) :staticmethod: Solve the transpose of a linear system A^T x = b using PETSc. :param A: Matrix representing the linear system. :type A: :py:class:`PETSc.Mat` :param b: Right-hand side of the equation. :type b: :py:class:`PETSc.Vec` :param destroy: Whether to destroy the matrix and vector after solving, by default False. :type destroy: :py:class:`bool`, *optional* :returns: **x** -- Solution vector. :rtype: :py:class:`numpy.ndarray` .. py:method:: createEquationSystemSolver(A) :staticmethod: Create a KSP PETSc solver to solve a linear equation system. This is useful if several linear equation systems with the same matrix are solved, since it stores the preconditioner and the calculation time is significantly reduced. To solve the equation system use the method "solveEquationSystemWithPredefinedSolver". :param A: Matrix of the linear system. :type A: :py:class:`PETSc.Mat` :returns: **solver** -- Preconfigured solver for the given matrix. :rtype: :py:class:`PETSc.KSP` .. py:method:: solveEquationSystemWithPredefinedSolver(solver, b, destroy=False) :staticmethod: Solve a linear equation system Ax=b using PETSc and a predefined solver. :param solver: Preconfigured solver, can be created with the method "createEquationSystemSolver". :type solver: :py:class:`PETSc.KSP` :param b: Right-hand side of the equation. :type b: :py:class:`PETSc.Vec` :param destroy: Whether to destroy the vector after solving, by default False. Can be useful by repetitive computations to avoid memory leaks. :type destroy: :py:class:`bool`, *optional* :returns: **x** -- Solution vector of the linear equation system. :rtype: :py:class:`numpy.ndarray` .. py:class:: ResolventOperator(ResolventOperator, FEMWeightMatrix_fullSystem, FEMWeightMatrix_forcingNorm, FEMWeightMatrix_responseNorm, RestrictorMatrix_forcing, RestrictorMatrix_response) Bases: :py:obj:`object` Matrix-free representation of the Resolvent operator multiplied with its conjugate transpose. Used to conduct the singular value decomposition of the system for resolvent analysis. Contains a method called "mult", which is called by the eigenvalue solver and returns a matrix-vector product of the represented matrix. The resolvent operator requires additional full-size quadratic matrices to calculate the forcing and response norms, as well as (possibly rectangular) restrictor matrices that limit the spatial domain and variable dimensions. **Initialize the ResolventOperator object** :param ResolventOperator: Matrix representing the linear system. :type ResolventOperator: :py:class:`PETSc.Mat` :param FEMWeightMatrix_fullSystem: Weight matrix for the full FEM system. :type FEMWeightMatrix_fullSystem: :py:class:`PETSc.Mat` :param FEMWeightMatrix_forcingNorm: Weight matrix for the forcing norm. Has default size of full system (surplus DOFs will be ignored). :type FEMWeightMatrix_forcingNorm: :py:class:`PETSc.Mat` :param FEMWeightMatrix_responseNorm: Weight matrix for the response norm. Has default size of full system (surplus DOFs will be ignored). :type FEMWeightMatrix_responseNorm: :py:class:`PETSc.Mat` :param RestrictorMatrix_forcing: Restrictor matrix for forcing. Rectangular matrix of appropriate size without FEM weights. Spatial restrictor values can be between 0 and 1. :type RestrictorMatrix_forcing: :py:class:`PETSc.Mat` :param RestrictorMatrix_response: Restrictor matrix for response. Rectangular matrix of appropriate size without FEM weights. Spatial restrictor values can be between 0 and 1. :type RestrictorMatrix_response: :py:class:`PETSc.Mat` .. py:method:: getSize() Get the size of the operator. :returns: The size of the operator (rows, columns). :rtype: :py:class:`tuple` .. py:method:: getVecs() Get PETSc vectors for the operator. :returns: A tuple of PETSc.Vec objects used internally by the operator. :rtype: :py:class:`tuple` .. py:method:: mult(mat, X, Y) Compute the matrix-vector product Y = mat * X. :param mat: The matrix represented by the operator. :type mat: :py:class:`PETSc.Mat` :param X: Input vector. :type X: :py:class:`PETSc.Vec` :param Y: Output vector. :type Y: :py:class:`PETSc.Vec` :returns: The result of the matrix-vector product. :rtype: :py:class:`PETSc.Vec` .. py:method:: getKSP() Get the KSP solver for the operator. :returns: KSP solver instance that solves the resolvent equation (without its conjugate transpose). Can be used to get the response to a given forcing. :rtype: :py:class:`PETSc.KSP` .. py:method:: destroySelf() Clean up resources to avoid memory leaks. This method should be called if multiple resolvent SVDs are performed consecutively.