Python linalg.expm方法代码示例

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import expm [as 别名]
def test_cubic_fermionic_simulation_gate_consistency_docstring(
        weights, exponent):
    generator = np.zeros((8, 8), dtype=np.complex128)
    # w0 |110><101| + h.c.
    generator[6, 5] = weights[0]
    generator[5, 6] = weights[0].conjugate()
    # w1 |110><011| + h.c.
    generator[6, 3] = weights[1]
    generator[3, 6] = weights[1].conjugate()
    # w2 |101><011| + h.c.
    generator[5, 3] = weights[2]
    generator[3, 5] = weights[2].conjugate()
    expected_unitary = la.expm(-1j * exponent * generator)

    gate = ofc.CubicFermionicSimulationGate(weights, exponent=exponent)
    actual_unitary = cirq.unitary(gate)

    assert np.allclose(expected_unitary, actual_unitary) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import expm [as 别名]
def test_quadratic_fermionic_simulation_gate_unitary(weights, exponent):
    generator = np.zeros((4, 4), dtype=np.complex128)
    # w0 |10><01| + h.c.
    generator[2, 1] = weights[0]
    generator[1, 2] = weights[0].conjugate()
    # w1 |11><11|
    generator[3, 3] = weights[1]
    expected_unitary = la.expm(-1j * exponent * generator)

    gate = ofc.QuadraticFermionicSimulationGate(weights, exponent=exponent)
    actual_unitary = cirq.unitary(gate)

    assert np.allclose(expected_unitary, actual_unitary)

    symbolic_gate = (ofc.QuadraticFermionicSimulationGate(
        (sympy.Symbol('w0'), sympy.Symbol('w1')), exponent=sympy.Symbol('t')))
    qubits = cirq.LineQubit.range(2)
    circuit = cirq.Circuit(symbolic_gate._decompose_(qubits))
    resolver = {'w0': weights[0], 'w1': weights[1], 't': exponent}
    resolved_circuit = cirq.resolve_parameters(circuit, resolver)
    decomp_unitary = resolved_circuit.unitary(qubit_order=qubits)

    assert np.allclose(expected_unitary, decomp_unitary) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import expm [as 别名]
def test_quartic_fermionic_simulation_unitary(weights, exponent):
    generator = np.zeros((1 << 4,) * 2, dtype=np.complex128)

    # w0 |1001><0110| + h.c.
    generator[9, 6] = weights[0]
    generator[6, 9] = weights[0].conjugate()
    # w1 |1010><0101| + h.c.
    generator[10, 5] = weights[1]
    generator[5, 10] = weights[1].conjugate()
    # w2 |1100><0011| + h.c.
    generator[12, 3] = weights[2]
    generator[3, 12] = weights[2].conjugate()
    expected_unitary = la.expm(-1j * exponent * generator)

    gate = ofc.QuarticFermionicSimulationGate(weights, exponent=exponent)
    actual_unitary = cirq.unitary(gate)

    assert np.allclose(expected_unitary, actual_unitary) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import expm [as 别名]
def __expm__(self, matrix, symmetric):
        r"""Calculates matrix exponential.

        Args:
            matrix (Tensor): Matrix to take exponential of.
            symmetric (bool): Specifies whether the matrix is symmetric.

        :rtype: (:class:`Tensor`)
        """
        if symmetric:
            e, V = torch.symeig(matrix, eigenvectors=True)
            diff_mat = V @ torch.diag(e.exp()) @ V.t()
        else:
            diff_mat_np = expm(matrix.cpu().numpy())
            diff_mat = torch.Tensor(diff_mat_np).to(matrix.device)
        return diff_mat 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import expm [as 别名]
def test_lanczos_evolve(n, N_cache, tol=5.e-14):
    # generate Hermitian test array
    leg = gen_random_legcharge(ch, n)
    H = npc.Array.from_func_square(rmat.GUE, leg) - npc.diag(1., leg)
    H_flat = H.to_ndarray()
    H_Op = H  # use `matvec` of the array
    qtotal = leg.to_qflat()[0]
    psi_init = npc.Array.from_func(np.random.random, [leg], qtotal=qtotal)
    psi_init /= npc.norm(psi_init)
    psi_init_flat = psi_init.to_ndarray()
    lanc = lanczos.LanczosEvolution(H_Op, psi_init, {'verbose': 1, 'N_cache': N_cache})
    for delta in [-0.1j, 0.1j, 1.j]:  #, 0.1, 1.]:
        psi_final_flat = expm(H_flat * delta).dot(psi_init_flat)
        psi_final, N = lanc.run(delta)
        ov = np.inner(psi_final.to_ndarray().conj(), psi_final_flat)
        ov /= np.linalg.norm(psi_final_flat)
        print("<psi1|psi1_flat>/norm=", ov)
        assert (abs(1. - abs(ov)) < tol) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import expm [as 别名]
def test_round_trip_random_float(self):
        np.random.seed(1234)
        for n in range(1, 6):
            M_unscaled = np.random.randn(n, n)
            for scale in np.logspace(-4, 4, 9):
                M = M_unscaled * scale

                # Eigenvalues are related to the branch cut.
                W = np.linalg.eigvals(M)
                err_msg = 'M:{0} eivals:{1}'.format(M, W)

                # Check sqrtm round trip because it is used within logm.
                M_sqrtm, info = sqrtm(M, disp=False)
                M_sqrtm_round_trip = M_sqrtm.dot(M_sqrtm)
                assert_allclose(M_sqrtm_round_trip, M)

                # Check logm round trip.
                M_logm, info = logm(M, disp=False)
                M_logm_round_trip = expm(M_logm)
                assert_allclose(M_logm_round_trip, M, err_msg=err_msg) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import expm [as 别名]
def test_random_matrices_and_powers(self):
        # Each independent iteration of this fuzz test picks random parameters.
        # It tries to hit some edge cases.
        np.random.seed(1234)
        nsamples = 20
        for i in range(nsamples):
            # Sample a matrix size and a random real power.
            n = random.randrange(1, 5)
            p = np.random.randn()

            # Sample a random real or complex matrix.
            matrix_scale = np.exp(random.randrange(-4, 5))
            A = np.random.randn(n, n)
            if random.choice((True, False)):
                A = A + 1j * np.random.randn(n, n)
            A = A * matrix_scale

            # Check a couple of analytically equivalent ways
            # to compute the fractional matrix power.
            # These can be compared because they both use the principal branch.
            A_power = fractional_matrix_power(A, p)
            A_logm, info = logm(A, disp=False)
            A_power_expm_logm = expm(A_logm * p)
            assert_allclose(A_power, A_power_expm_logm) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import expm [as 别名]
def test_expm_frechet(self):
        # a test of the basic functionality
        M = np.array([
            [1, 2, 3, 4],
            [5, 6, 7, 8],
            [0, 0, 1, 2],
            [0, 0, 5, 6],
            ], dtype=float)
        A = np.array([
            [1, 2],
            [5, 6],
            ], dtype=float)
        E = np.array([
            [3, 4],
            [7, 8],
            ], dtype=float)
        expected_expm = scipy.linalg.expm(A)
        expected_frechet = scipy.linalg.expm(M)[:2, 2:]
        for kwargs in ({}, {'method':'SPS'}, {'method':'blockEnlarge'}):
            observed_expm, observed_frechet = expm_frechet(A, E, **kwargs)
            assert_allclose(expected_expm, observed_expm)
            assert_allclose(expected_frechet, observed_frechet) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import expm [as 别名]
def _test_variable_step_method(self, method_str):
        """Some tests for a variable step method."""

        # get method and set general options
        ode_method = method_from_string(method_str)
        options = DE_Options(method=method_str, atol=10**-10, rtol=10**-10)

        # run on matrix problem
        solver = ode_method(t0=self.t0, y0=self.y0, rhs=self.rhs, options=options)
        solver.integrate(1.)
        expected = expm(-1j * np.pi * self.X)

        # set the comparison tolerance to be somewhat lenient
        self.assertAlmostEqual(solver.y, expected, tol=10**-8)

        # test with an arbitrary problem
        def rhs(t, y):
            return np.array([t**2])

        solver = ode_method(t0=0., y0=np.array(0.), rhs={'rhs': rhs}, options=options)

        solver.integrate(1.)
        expected = 1./3
        self.assertAlmostEqual(solver.y, expected, tol=10**-9) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import expm [as 别名]
def test_RK4(self):
        """Run tests on RK4 fixed-step solver."""
        ode_method = method_from_string('RK4')
        options = DE_Options(max_dt=10**-3)

        # run on matrix problem
        solver = ode_method(t0=self.t0, y0=self.y0, rhs=self.rhs, options=options)
        solver.integrate(1.)
        expected = expm(-1j * np.pi * self.X)

        # set the comparison tolerance to be somewhat lenient
        self.assertAlmostEqual(solver.y, expected, tol=10**-8)

        # test with an arbitrary problem
        def rhs(t, y):
            return np.array([t**2])

        solver = ode_method(t0=0., y0=np.array(0.), rhs={'rhs': rhs}, options=options)

        solver.integrate(1.)
        expected = 1./3
        self.assertAlmostEqual(solver.y, expected, tol=10**-8) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import expm [as 别名]
def _cphase_legacy(angle=np.pi, leakage=0.):
    """A perfect controlled phase rotation.
    First qubit is low-frequency, second qubit is high-frequency (it leaks).

    Parameters
    ----------
    angle: float, optional
        Rotation angle in radians. Default is :math:`\\pi`.
    leakage: float, optional
        Leakage rate of a CPhase gate

    Returns
    -------
    Operation
        An operation, that corresponds to the rotation.
    """
    dcphase = np.zeros((9, 9))
    dcphase[2, 4] = 1
    dcphase[4, 2] = 1
    angle_frac = 1 - np.arcsin(np.sqrt(leakage)) / np.pi
    unitary = expm(-1j * angle * angle_frac * dcphase)
    return Operation.from_kraus(unitary, bases2_default) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import expm [as 别名]
def _line_search(Y, W, density, direction, signs, current_loss, ls_tries,
                 verbose, ortho, extended):
    '''
    Performs a backtracking line search, starting from Y and W, in the
    direction direction. I
    '''
    N = W.shape[0]
    alpha = 1.
    if current_loss is None:
        current_loss = _loss(Y, W, density, signs, ortho, extended)
    for _ in range(ls_tries):
        if ortho:
            transform = expm(alpha * direction)
        else:
            transform = np.eye(N) + alpha * direction
        Y_new = np.dot(transform, Y)
        W_new = np.dot(transform, W)
        new_loss = _loss(Y_new, W_new, density, signs, ortho, extended)
        if new_loss < current_loss:
            return True, Y_new, W_new, new_loss, alpha * direction
        alpha /= 2.
    else:
        if verbose:
            print('line search failed, falling back to gradient')
        return False, Y_new, W_new, new_loss, alpha * direction 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import expm [as 别名]
def calc_likelihood(args):
	[delta_t,r_vec,Q,rho_at_present,r_vec_indexes,sign_list,sp_OrigTimeIndex]=args
	PvDes= rho_at_present
	#print rho_at_present
	recursive = np.arange(sp_OrigTimeIndex,len(delta_t))[::-1]
	for i in recursive:
		# get time span
		t=delta_t[i] 
		# get rho vector
		r_ind= r_vec_indexes[i]
		sign=  sign_list[i]
		rho_vec= np.prod(abs(sign-r_vec[r_ind]),axis=1)
		# prob of at least 1  1-exp(-rho_vec*t)
		#print rho_vec,r_vec[r_ind]
		#print obs_area_series[i],"\t",rho_vec				
		Pt= linalg.expm(Q.T *(t))
		condLik_temp= np.dot(PvDes,Pt)
		PvDes= condLik_temp *rho_vec
		
		#print "temp,",t,PvDes,log(condLik_temp),rho_vec
	return np.log(np.sum(PvDes)) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import expm [as 别名]
def calc_likelihood_mQ(args):
	[delta_t,r_vec_list,Q_list,rho_at_present,r_vec_indexes,sign_list,sp_OrigTimeIndex,index_r,index_q]=args
	PvDes= rho_at_present
	#print rho_at_present
	recursive = np.arange(sp_OrigTimeIndex,len(delta_t))[::-1]
	for i in recursive:
		#print "here",i, Q_index[i]
		Q = Q_list[index_q[i]]
		r_vec=r_vec_list[index_r[i]]
		# get time span
		t=delta_t[i] 
		# get rho vector
		r_ind= r_vec_indexes[i]
		sign=  sign_list[i]
		rho_vec= np.prod(abs(sign-r_vec[r_ind]),axis=1)
		# prob of at least 1  1-exp(-rho_vec*t)
		Pt= linalg.expm(Q.T *(t))
		condLik_temp= np.dot(PvDes,Pt)
		PvDes= condLik_temp *rho_vec
		
		#print "temp,",t,PvDes,log(condLik_temp),rho_vec
	if np.sum(PvDes) <= 0: 
		print np.sum(PvDes), list(PvDes)
	
	return np.log(np.sum(PvDes)) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import expm [as 别名]
def v2U(v,exp=True):
    """
    Transform a vector into a unitary matrix
    """
    n = int(np.sqrt(v.shape[0]+1)) # dimension of the matrix
    M = np.zeros((n,n),dtype=np.complex) # define the matrix
    ii = 0 # counter
    d = np.zeros(n,dtype=np.complex) # diagonal
    for i in range(n-1):
        d[i] = v[ii]
        ii += 1 # increase counter
    d[n-1] = -np.sum(d) # ensure the trace is zero
    for i in range(n): M[i,i] = d[i] # put the diagonal elements
    for j in range(n):
        for i in range(j+1,n):
            M[i,j] = v[ii] + 1j*v[ii+1] # store
            M[j,i] = v[ii] - 1j*v[ii+1] # store
            ii += 2 # increase
    if exp: return lg.expm(1j*M) # return unitary matrix
    else: return M # return the Lie algebra matrix