Python linalg.solve方法代碼示例

# 需要導入模塊: from numpy import linalg [as 別名]
# 或者: from numpy.linalg import solve [as 別名]
def si_c2(self,ww):
    """
    This computes the correlation part of the screened interaction using LinearOpt and lgmres
    lgmres method is much slower than np.linalg.solve !!
    """
    import numpy as np
    from scipy.sparse.linalg import lgmres
    from scipy.sparse.linalg import LinearOperator
    rf0 = si0 = self.rf0(ww)    
    for iw,w in enumerate(ww):                                
      k_c = np.dot(self.kernel_sq, rf0[iw,:,:])                                         
      b = np.dot(k_c, self.kernel_sq)               
      self.comega_current = w
      k_c_opt = LinearOperator((self.nprod,self.nprod), matvec=self.gw_vext2veffmatvec, dtype=self.dtypeComplex)  
      for m in range(self.nprod): 
         si0[iw,m,:],exitCode = lgmres(k_c_opt, b[m,:], atol=self.gw_iter_tol, maxiter=self.maxiter)   
      if exitCode != 0: print("LGMRES has not achieved convergence: exitCode = {}".format(exitCode))
      #np.allclose(np.dot(k_c, si0), b, atol=1e-05) == True  #Test   
    return si0 

# 需要導入模塊: from numpy import linalg [as 別名]
# 或者: from numpy.linalg import solve [as 別名]
def si_c_check (self, tol = 1e-5):
    """
    This compares np.solve and LinearOpt-lgmres methods for solving linear equation (1-v\chi_{0}) * W_c = v\chi_{0}v
    """
    import time
    import numpy as np
    ww = 1j*self.ww_ia
    t = time.time()
    si0_1 = self.si_c(ww)      #method 1:  numpy.linalg.solve
    t1 = time.time() - t
    print('numpy: {} sec'.format(t1))
    t2 = time.time()
    si0_2 = self.si_c2(ww)     #method 2:  scipy.sparse.linalg.lgmres
    t3 = time.time() - t2
    print('lgmres: {} sec'.format(t3))
    summ = abs(si0_1 + si0_2).sum()
    diff = abs(si0_1 - si0_2).sum() 
    if diff/summ < tol and diff/si0_1.size < tol:
       print('OK! scipy.lgmres methods and np.linalg.solve have identical results')
    else:
       print('Results (W_c) are NOT similar!')     
    return [[diff/summ] , [np.amax(abs(diff))] ,[tol]]

  #@profile 

# 需要導入模塊: from numpy import linalg [as 別名]
# 或者: from numpy.linalg import solve [as 別名]
def getBarycentricCoords(A, B, C, X, checkValidity = True):
    T = np.array( [ [A.x - C.x, B.x - C.x ], [A.y - C.y, B.y - C.y] ] )
    y = np.array( [ [X.x - C.x], [X.y - C.y] ] )
    lambdas = linalg.solve(T, y)
    lambdas = lambdas.flatten()
    lambdas = np.append(lambdas, 1 - (lambdas[0] + lambdas[1]))
    if checkValidity:
        if (lambdas[0] < 0 or lambdas[1] < 0 or lambdas[2] < 0):
            print "ERROR: Not a convex combination; lambda = %s"%lambdas
            print "pointInsideConvexPolygon2D = %s"%pointInsideConvexPolygon2D([A, B, C], X, 0)
            plt.plot([A.x, B.x, C.x, A.x], [A.y, B.y, C.y, A.y], 'r')
            plt.hold(True)
            plt.plot([X.x], [X.y], 'b.')
            plt.show()
        assert (lambdas[0] >= 0 and lambdas[1] >= 0 and lambdas[2] >= 0)
    else:
        lambdas[0] = max(lambdas[0], 0)
        lambdas[1] = max(lambdas[1], 0)
        lambdas[2] = max(lambdas[2], 0)
    return lambdas 

# 需要導入模塊: from numpy import linalg [as 別名]
# 或者: from numpy.linalg import solve [as 別名]
def test_0_size_k(self):
        # test zero multiple equation (K=0) case.
        class ArraySubclass(np.ndarray):
            pass
        a = np.arange(4).reshape(1, 2, 2)
        b = np.arange(6).reshape(3, 2, 1).view(ArraySubclass)

        expected = linalg.solve(a, b)[:, :, 0:0]
        result = linalg.solve(a, b[:, :, 0:0])
        assert_array_equal(result, expected)
        assert_(isinstance(result, ArraySubclass))

        # test both zero.
        expected = linalg.solve(a, b)[:, 0:0, 0:0]
        result = linalg.solve(a[:, 0:0, 0:0], b[:, 0:0, 0:0])
        assert_array_equal(result, expected)
        assert_(isinstance(result, ArraySubclass)) 

# 需要導入模塊: from numpy import linalg [as 別名]
# 或者: from numpy.linalg import solve [as 別名]
def test_0_size_k(self):
        # test zero multiple equation (K=0) case.
        class ArraySubclass(np.ndarray):
            pass
        a = np.arange(4).reshape(1, 2, 2)
        b = np.arange(6).reshape(3, 2, 1).view(ArraySubclass)

        expected = linalg.solve(a, b)[:,:, 0:0]
        result = linalg.solve(a, b[:,:, 0:0])
        assert_array_equal(result, expected)
        assert_(isinstance(result, ArraySubclass))

        # test both zero.
        expected = linalg.solve(a, b)[:, 0:0, 0:0]
        result = linalg.solve(a[:, 0:0, 0:0], b[:,0:0, 0:0])
        assert_array_equal(result, expected)
        assert_(isinstance(result, ArraySubclass)) 

# 需要導入模塊: from numpy import linalg [as 別名]
# 或者: from numpy.linalg import solve [as 別名]
def solve(self):
        # known and unknown values
        self.VU = [node[key] for node in self.U.values() for key in ("ux",)]
        self.VF = [node[key] for node in self.F.values() for key in ("fx",)]
        knw = [pos for pos,value in enumerate(self.VU) if not value is np.nan]
        unknw = [pos for pos,value in enumerate(self.VU) if value is np.nan]
        # Matrices to solve
        self.K2S = np.delete(np.delete(self.KG,knw,0),knw,1)
        self.F2S = np.delete(self.VF,knw,0)
        # For displacements
        self.solved_u = la.solve(self.K2S,self.F2S)
        # Updating U (displacements vector)
        for k,ic in enumerate(unknw):
            nd, var = self.index2key(ic)
            self.U[nd][var] = self.solved_u[k]
            self.nodes[ic].ux = self.solved_u[k]
        # For nodal forces/reactions
        self.NF = self.F.copy()
        self.VU = [node[key] for node in self.U.values() for key in ("ux",)]
        nf_calc = np.dot(self.KG, self.VU)
        for k,ic in enumerate(range(self.get_number_of_nodes())):
            nd, var = self.index2key(ic, ("fx",))
            self.NF[nd][var] = nf_calc[k]
            self.nodes[ic].fx = nf_calc[k] 

# 需要導入模塊: from numpy import linalg [as 別名]
# 或者: from numpy.linalg import solve [as 別名]
def _single_sample_update(self, X, i, w):
        n_features = X.shape[1]
        if X.indptr[i + 1] - X.indptr[i] != 0:
            X_subset = X.data[X.indptr[i]:X.indptr[i + 1]]
            subset = X.indices[X.indptr[i]:X.indptr[i + 1]]
            len_subset = subset.shape[0]
            reduction = n_features / len_subset
            self.feature_n_iter_[subset] += 1
            components_subset = self.components_[:, subset]
            Dx = components_subset.dot(X_subset)
            G = components_subset.dot(components_subset.T)
            G.flat[::self.n_components + 1] += self.alpha / reduction
            self.code_[i] = linalg.solve(G, Dx)
            code = self.code_[i]
            w_B = np.minimum(1,
                             w * self.n_iter_ / self.feature_n_iter_[subset])
            self.B_[:, subset] *= 1 - w_B
            self.B_[:, subset] += np.outer(code, X_subset * w_B) 

# 需要導入模塊: from numpy import linalg [as 別名]
# 或者: from numpy.linalg import solve [as 別名]
def enet_regression_multi_gram_(G, Dx, X, code, l1_ratio, alpha,
                                positive):
    batch_size = code.shape[0]
    if l1_ratio == 0:
        n_components = G.shape[1]
        for i in range(batch_size):
            G.flat[::n_components + 1] += alpha
            code[i] = linalg.solve(G[i], Dx[i])
            G.flat[::n_components + 1] -= alpha
    else:
        # Unused but unfortunate API
        random_state = check_random_state(0)
        for i in range(batch_size):
            cd_fast.enet_coordinate_descent_gram(
                code[i],
                alpha * l1_ratio,
                alpha * (
                    1 - l1_ratio),
                G[i], Dx[i], X[i], 100, 1e-2,
                random_state,
                False, positive)
    return code 

# 需要導入模塊: from numpy import linalg [as 別名]
# 或者: from numpy.linalg import solve [as 別名]
def enet_regression_single_gram_(G, Dx, X, code, code_l1_ratio, code_alpha,
                                 code_pos):
    batch_size = code.shape[0]
    if code_l1_ratio == 0:
        n_components = G.shape[0]
        G = G.copy()
        G.flat[::n_components + 1] += code_alpha
        code[:] = linalg.solve(G, Dx.T).T
        G.flat[::n_components + 1] -= code_alpha
    else:
        # Unused but unfortunate API
        random_state = check_random_state(0)
        for i in range(batch_size):
            cd_fast.enet_coordinate_descent_gram(
                code[i],
                code_alpha * code_l1_ratio,
                code_alpha * (
                    1 - code_l1_ratio),
                G, Dx[i], X[i], 100, 1e-2,
                random_state,
                False, code_pos)
    return code 

# 需要導入模塊: from numpy import linalg [as 別名]
# 或者: from numpy.linalg import solve [as 別名]
def find_line_intersection(r1, angle1, r2, angle2):
    """
    Find intersection between two lines defined by point and direction.

    :param r1: Origin of the first line.
    :param angle1: Angle of the first line.
    :param r2: Origin of the second line.
    :param angle2: Angle of the second line.
    :return: Tuple of point of intersection and distances from the origins.
    :rtype: tuple
    """
    u1 = np.array([np.cos(angle1), np.sin(angle1)])
    u2 = np.array([np.cos(angle2), np.sin(angle2)])

    a = np.array([[u1[0], -u2[0]], [u1[1], -u2[1]]])
    b = r2 - r1

    # noinspection PyTypeChecker
    x = linalg.solve(a, b)
    return r1 + u1 * x[0], x 

# 需要導入模塊: from numpy import linalg [as 別名]
# 或者: from numpy.linalg import solve [as 別名]
def MakeIaosRaw(COcc, S1, S2, S12):
    # calculate the molecule-intrinsic atomic orbital (IAO) basis
    # ref: [1] Knizia, J. Chem. Theory Comput., http://dx.doi.org/10.1021/ct400687b
    # This is the "Simple/2014" version from ibo-ref at sites.psu.edu/knizia/software
    assert(S1.shape[0] == S1.shape[1] and S1.shape[0] == S12.shape[0])
    assert(S2.shape[0] == S2.shape[1] and S2.shape[0] == S12.shape[1])
    assert(COcc.shape[0] == S1.shape[0] and COcc.shape[1] <= S2.shape[0])
    P12 = la.solve(S1, S12)   # P12 = S1^{-1} S12
    COcc2 = mdot(S12.T, COcc)              # O(N m^2)
    CTil = la.solve(S2, COcc2)             # O(m^3)
    STil = mdot(COcc2.T, CTil)             # O(m^3)
    CTil2Bar = la.solve(STil, CTil.T).T    # O(m^3)
    T4 = COcc - mdot(P12, CTil2Bar)        # O(N m^2)
    CIb = P12 + mdot(T4, COcc2.T)          # O(N m^2)
    return CIb 

# 需要導入模塊: from numpy import linalg [as 別名]
# 或者: from numpy.linalg import solve [as 別名]
def horner(self, s):
        """Evaluate the systems's transfer function for a complex variable

        Returns a matrix of values evaluated at complex variable s.
        """
        resp = np.dot(self.C, solve(s * eye(self.states) - self.A,
                                    self.B)) + self.D
        return array(resp)

    # Method for generating the frequency response of the system 

# 需要導入模塊: from numpy import linalg [as 別名]
# 或者: from numpy.linalg import solve [as 別名]
def dcgain(self):
        """Return the zero-frequency gain

        The zero-frequency gain of a continuous-time state-space
        system is given by:

        .. math: G(0) = - C A^{-1} B + D

        and of a discrete-time state-space system by:

        .. math: G(1) = C (I - A)^{-1} B + D

        Returns
        -------
        gain : ndarray
            An array of shape (outputs,inputs); the array will either
            be the zero-frequency (or DC) gain, or, if the frequency
            response is singular, the array will be filled with np.nan.
        """
        try:
            if self.isctime():
                gain = np.asarray(self.D-self.C.dot(np.linalg.solve(self.A, self.B)))
            else:
                gain = self.horner(1)
        except LinAlgError:
            # eigenvalue at DC
            gain = np.tile(np.nan, (self.outputs, self.inputs))
        return np.squeeze(gain)


# TODO: add discrete time check 

# 需要導入模塊: from numpy import linalg [as 別名]
# 或者: from numpy.linalg import solve [as 別名]
def do(self, a, b, tags):
        x = linalg.solve(a, b)
        assert_almost_equal(b, dot_generalized(a, x))
        assert_(consistent_subclass(x, b)) 

# 需要導入模塊: from numpy import linalg [as 別名]
# 或者: from numpy.linalg import solve [as 別名]
def test_types(self, dtype):
        x = np.array([[1, 0.5], [0.5, 1]], dtype=dtype)
        assert_equal(linalg.solve(x, x).dtype, dtype)