Python linalg.lu_solve方法代碼示例

# 需要導入模塊: from scipy import linalg [as 別名]
# 或者: from scipy.linalg import lu_solve [as 別名]
def solve(self, vec):
        """Solve the impedance matrix for a source vector. Caches the
        factorised matrix for efficiently solving multiple vectors"""
        if self.part_o != self.part_s:
            raise ValueError("Can only invert a self-impedance matrix")

        Z_lu = self.factored()
        if isinstance(vec, LookupArray):
            vec = vec.simple_view()

        lookup = (self.unknowns, (self.part_s, self.basis_container))

        if len(vec.shape) > 1:
            lookup = lookup+(vec.shape[1],)

        I = LookupArray(lookup, dtype=np.complex128)
        I_simp = I.simple_view()
        I_simp[:] = la.lu_solve(Z_lu, vec)
        return I 

# 需要導入模塊: from scipy import linalg [as 別名]
# 或者: from scipy.linalg import lu_solve [as 別名]
def lu_compute_mtx_obj(Tbeta_lst, num_bands, K, lu_R_GtGinv_Rt_lst):
    """
    compute the matrix (M) in the objective function:
        min   c^H M c
        s.t.  c0^H c = 1

    Parameters
    ----------
    GtG_lst: 
        list of G^H * G
    Tbeta_lst: 
        list of Teoplitz matrices for beta-s
    Rc0: 
        right dual matrix for the annihilating filter (same of each block -> not a list)
    """
    mtx = np.zeros((K + 1, K + 1), dtype=float)  # <= assume G, Tbeta and Rc0 are real-valued

    for loop in range(num_bands):
        Tbeta_loop = Tbeta_lst[loop]
        mtx += np.dot(Tbeta_loop.T,
                      linalg.lu_solve(lu_R_GtGinv_Rt_lst[loop],
                                      Tbeta_loop, check_finite=False)
                      )

    return mtx 

# 需要導入模塊: from scipy import linalg [as 別名]
# 或者: from scipy.linalg import lu_solve [as 別名]
def compute_obj_val(GtG_inv_lst, Tbeta_lst, Rc0, c_ri_half, num_bands, K):
    """
    compute the fitting error.
    CAUTION: Here we assume use_lu = True
    """
    mtx = np.zeros((K + 1, K + 1), dtype=float)  # <= assume G, Tbeta and Rc0 are real-valued
    fitting_error = 0
    for band_count in range(num_bands):
        #GtG_loop = GtG_lst[band_count]
        Tbeta_loop = Tbeta_lst[band_count]

        mtx += np.dot(Tbeta_loop.T,
                      linalg.solve(
                          np.dot(Rc0, 
                              #linalg.lu_solve(GtG_loop, Rc0.T, check_finite=False)
                              np.dot(GtG_inv_lst[band_count], Rc0.T)
                              ),
                          Tbeta_loop, check_finite=False
                      )
                      )

    fitting_error += np.sqrt(np.dot(c_ri_half.T, np.dot(mtx, c_ri_half)).real)

    return fitting_error, mtx 

# 需要導入模塊: from scipy import linalg [as 別名]
# 或者: from scipy.linalg import lu_solve [as 別名]
def __init__(self, inv_array, inv_lu_and_piv, inv_lu_transposed):
        """
        Args:
            inv_array (array): 2D array specifying inverse matrix entries.
            inv_lu_and_piv (Tuple[array, array]): Pivoted LU factorisation
                represented as a tuple with first element a 2D array containing
                the lower and upper triangular factors (with the unit diagonal
                of the lower triangular factor not stored) and the second
                element a 1D array containing the pivot indices. Corresponds
                to the output of `scipy.linalg.lu_factor` and input to
                `scipy.linalg.lu_solve`.
            inv_lu_transposed (bool): Whether LU factorisation is of inverse of
                array or transpose of inverse of array.
        """
        super().__init__(inv_array.shape)
        self._inv_array = inv_array
        self._inv_lu_and_piv = inv_lu_and_piv
        self._inv_lu_transposed = inv_lu_transposed 

# 需要導入模塊: from scipy import linalg [as 別名]
# 或者: from scipy.linalg import lu_solve [as 別名]
def _matvec(self, x):
        return lu_solve(self.M_lu, x) 

# 需要導入模塊: from scipy import linalg [as 別名]
# 或者: from scipy.linalg import lu_solve [as 別名]
def solve_nonlinear(self, inputs, outputs):
        self.lu = lu_factor(inputs['mtx'])

        outputs['circulations'] = lu_solve(self.lu, inputs['rhs']) 

# 需要導入模塊: from scipy import linalg [as 別名]
# 或者: from scipy.linalg import lu_solve [as 別名]
def solve_linear(self, d_outputs, d_residuals, mode):
        if mode == 'fwd':
            d_outputs['circulations'] = lu_solve(self.lu, d_residuals['circulations'], trans=0)
        else:
            d_residuals['circulations'] = lu_solve(self.lu, d_outputs['circulations'], trans=1) 

# 需要導入模塊: from scipy import linalg [as 別名]
# 或者: from scipy.linalg import lu_solve [as 別名]
def solve_system(self, rhs):
        # To do preconditioning below, we will need to scale each column of A elementwise by the entries of some vector
        col_scale = lambda A, scale: (A.T * scale).T  # Uses the trick that A*x scales the 0th column of A by x[0], etc.

        if self.factorisation_current:
            # A(preconditioned) = diag(left_scaling) * A(original) * diag(right_scaling)
            # Solve A(original)\rhs
            return col_scale(LA.lu_solve((self.lu, self.piv), col_scale(rhs, self.left_scaling)), self.right_scaling)
        else:
            if self.do_logging:
                logging.warning("model.solve_system not using factorisation")
            A, left_scaling, right_scaling = self.interpolation_matrix()
            return col_scale(LA.solve(A, col_scale(rhs, left_scaling)), right_scaling) 

# 需要導入模塊: from scipy import linalg [as 別名]
# 或者: from scipy.linalg import lu_solve [as 別名]
def _left_matrix_multiply(self, other):
        return sla.lu_solve(
            self._inv_lu_and_piv, other, self._inv_lu_transposed, 
            check_finite=False) 

# 需要導入模塊: from scipy import linalg [as 別名]
# 或者: from scipy.linalg import lu_solve [as 別名]
def _right_matrix_multiply(self, other):
        return sla.lu_solve(
            self._inv_lu_and_piv, other.T, not self._inv_lu_transposed,
            check_finite=False).T 

# 需要導入模塊: from scipy import linalg [as 別名]
# 或者: from scipy.linalg import lu_solve [as 別名]
def lu_solve_ATAI(A, rho, b, lu, piv):
    r"""
    Solve the linear system :math:`(A^T A + \rho I)\mathbf{x} = \mathbf{b}`
    or :math:`(A^T A + \rho I)X = B` using :func:`scipy.linalg.lu_solve`.

    Parameters
    ----------
    A : array_like
      Matrix :math:`A`
    rho : float
      Scalar :math:`\rho`
    b : array_like
      Vector :math:`\mathbf{b}` or matrix :math:`B`
    lu : array_like
      Matrix containing U in its upper triangle, and L in its lower triangle,
      as returned by :func:`scipy.linalg.lu_factor`
    piv : array_like
      Pivot indices representing the permutation matrix P, as returned by
      :func:`scipy.linalg.lu_factor`

    Returns
    -------
    x : ndarray
      Solution to the linear system.
    """

    N, M = A.shape
    if N >= M:
        x = linalg.lu_solve((lu, piv), b)
    else:
        x = (b - A.T.dot(linalg.lu_solve((lu, piv), A.dot(b), 1))) / rho
    return x 

# 需要導入模塊: from scipy import linalg [as 別名]
# 或者: from scipy.linalg import lu_solve [as 別名]
def lu_solve_AATI(A, rho, b, lu, piv):
    r"""
    Solve the linear system :math:`(A A^T + \rho I)\mathbf{x} = \mathbf{b}`
    or :math:`(A A^T + \rho I)X = B` using :func:`scipy.linalg.lu_solve`.

    Parameters
    ----------
    A : array_like
      Matrix :math:`A`
    rho : float
      Scalar :math:`\rho`
    b : array_like
      Vector :math:`\mathbf{b}` or matrix :math:`B`
    lu : array_like
      Matrix containing U in its upper triangle, and L in its lower triangle,
      as returned by :func:`scipy.linalg.lu_factor`
    piv : array_like
      Pivot indices representing the permutation matrix P, as returned by
      :func:`scipy.linalg.lu_factor`

    Returns
    -------
    x : ndarray
      Solution to the linear system.
    """

    N, M = A.shape
    if N >= M:
        x = (b - linalg.lu_solve((lu, piv), b.dot(A).T).T.dot(A.T)) / rho
    else:
        x = linalg.lu_solve((lu, piv), b.T).T
    return x 

# 需要導入模塊: from scipy import linalg [as 別名]
# 或者: from scipy.linalg import lu_solve [as 別名]
def lu_solve_ATAI(A, rho, b, lu, piv, check_finite=True):
    r"""Solve the linear system :math:`(A^T A + \rho I)\mathbf{x} = \mathbf{b}`
    or :math:`(A^T A + \rho I)X = B` using :func:`scipy.linalg.lu_solve`.

    Parameters
    ----------
    A : array_like
      Matrix :math:`A`
    rho : float
      Scalar :math:`\rho`
    b : array_like
      Vector :math:`\mathbf{b}` or matrix :math:`B`
    lu : array_like
      Matrix containing U in its upper triangle, and L in its lower triangle,
      as returned by :func:`scipy.linalg.lu_factor`
    piv : array_like
      Pivot indices representing the permutation matrix P, as returned by
      :func:`scipy.linalg.lu_factor`
    check_finite : bool, optional (default False)
      Flag indicating whether the input array should be checked for Inf
      and NaN values

    Returns
    -------
    x : ndarray
      Solution to the linear system
    """

    N, M = A.shape
    if N >= M:
        x = linalg.lu_solve((lu, piv), b, check_finite=check_finite)
    else:
        x = (b - A.T.dot(linalg.lu_solve((lu, piv), A.dot(b), 1,
                                         check_finite=check_finite))) / rho
    return x 

# 需要導入模塊: from scipy import linalg [as 別名]
# 或者: from scipy.linalg import lu_solve [as 別名]
def lu_solve_AATI(A, rho, b, lu, piv, check_finite=True):
    r"""Solve the linear system :math:`(A A^T + \rho I)\mathbf{x} = \mathbf{b}`
    or :math:`(A A^T + \rho I)X = B` using :func:`scipy.linalg.lu_solve`.

    Parameters
    ----------
    A : array_like
      Matrix :math:`A`
    rho : float
      Scalar :math:`\rho`
    b : array_like
      Vector :math:`\mathbf{b}` or matrix :math:`B`
    lu : array_like
      Matrix containing U in its upper triangle, and L in its lower triangle,
      as returned by :func:`scipy.linalg.lu_factor`
    piv : array_like
      Pivot indices representing the permutation matrix P, as returned by
      :func:`scipy.linalg.lu_factor`
    check_finite : bool, optional (default False)
      Flag indicating whether the input array should be checked for Inf
      and NaN values

    Returns
    -------
    x : ndarray
      Solution to the linear system
    """

    N, M = A.shape
    if N >= M:
        x = (b - linalg.lu_solve((lu, piv), b.dot(A).T,
                                 check_finite=check_finite).T.dot(A.T)) / rho
    else:
        x = linalg.lu_solve((lu, piv), b.T, check_finite=check_finite).T
    return x 

# 需要導入模塊: from scipy import linalg [as 別名]
# 或者: from scipy.linalg import lu_solve [as 別名]
def mimo_block_arnoldi(self, frequency, r):

        n = self.ss.states

        A = self.ss.A
        B = self.ss.B
        C = self.ss.C

        for i in range(self.nfreq):

            if self.frequency[i] == np.inf:
                F = A
                G = B
            else:
                lu_a = krylovutils.lu_factor(frequency[i], A)
                F = krylovutils.lu_solve(lu_a, np.eye(n))
                G = krylovutils.lu_solve(lu_a, B)

            if i == 0:
                V = krylovutils.block_arnoldi_krylov(r, F, G)
            else:
                Vi = krylovutils.block_arnoldi_krylov(r, F, G)
                V = np.block([V, Vi])

        self.V = V

        Ar = V.T.dot(A.dot(V))
        Br = V.T.dot(B)
        Cr = C.dot(V)

        return Ar, Br, Cr