def test_lie_derivative(self):
Lfh = pm.lie_derivatives(self.h, self.f, self.x, 0)
self.assertEqual(Lfh, [self.h])
Lfh = pm.lie_derivatives(self.h, self.f, self.x, 1)
self.assertEqual(Lfh, [self.h,
sp.Matrix([-2*self._x1*self._x2**2
- sp.sin(self._x1)*sp.cos(self._x2)])
])
Lfh = pm.lie_derivatives(self.h, self.f, self.x, 2)
self.assertEqual(Lfh, [self.h,
sp.Matrix([-2*self._x1*self._x2**2
- sp.sin(self._x1)*sp.cos(self._x2)]),
sp.Matrix([
-self._x2**2*(
-2*self._x2**2
- sp.cos(self._x1)*sp.cos(self._x2)
)
+ sp.sin(self._x1)*(
- 4*self._x1*self._x2
+ sp.sin(self._x1)*sp.sin(self._x2)
)])
])
python类Matrix()的实例源码
def test_quantum_fourier():
assert QFT(0, 3).decompose() == \
SwapGate(0, 2)*HadamardGate(0)*CGate((0,), PhaseGate(1)) * \
HadamardGate(1)*CGate((0,), TGate(2))*CGate((1,), PhaseGate(2)) * \
HadamardGate(2)
assert IQFT(0, 3).decompose() == \
HadamardGate(2)*CGate((1,), RkGate(2, -2))*CGate((0,), RkGate(2, -3)) * \
HadamardGate(1)*CGate((0,), RkGate(1, -2))*HadamardGate(0)*SwapGate(0, 2)
assert represent(QFT(0, 3), nqubits=3) == \
Matrix([[exp(2*pi*I/8)**(i*j % 8)/sqrt(8) for i in range(8)] for j in range(8)])
assert QFT(0, 4).decompose() # non-trivial decomposition
assert qapply(QFT(0, 3).decompose()*Qubit(0, 0, 0)).expand() == qapply(
HadamardGate(0)*HadamardGate(1)*HadamardGate(2)*Qubit(0, 0, 0)
).expand()
def test_matrix_to_density():
mat = Matrix([[0, 0], [0, 1]])
assert matrix_to_density(mat) == Density([Qubit('1'), 1])
mat = Matrix([[1, 0], [0, 0]])
assert matrix_to_density(mat) == Density([Qubit('0'), 1])
mat = Matrix([[0, 0], [0, 0]])
assert matrix_to_density(mat) == 0
mat = Matrix([[0, 0, 0, 0],
[0, 0, 0, 0],
[0, 0, 1, 0],
[0, 0, 0, 0]])
assert matrix_to_density(mat) == Density([Qubit('10'), 1])
mat = Matrix([[1, 0, 0, 0],
[0, 0, 0, 0],
[0, 0, 0, 0],
[0, 0, 0, 0]])
assert matrix_to_density(mat) == Density([Qubit('00'), 1])
def cache_matrix(self, name, m):
"""Cache a matrix by its name.
Parameters
----------
name : str
A descriptive name for the matrix, like "identity2".
m : list of lists
The raw matrix data as a sympy Matrix.
"""
try:
self._sympy_matrix(name, m)
except ImportError:
pass
try:
self._numpy_matrix(name, m)
except ImportError:
pass
try:
self._scipy_sparse_matrix(name, m)
except ImportError:
pass
def get_matrix(self, name, format):
"""Get a cached matrix by name and format.
Parameters
----------
name : str
A descriptive name for the matrix, like "identity2".
format : str
The format desired ('sympy', 'numpy', 'scipy.sparse')
"""
m = self._cache.get((name, format))
if m is not None:
return m
raise NotImplementedError(
'Matrix with name %s and format %s is not available.' %
(name, format)
)
def test_parallel_axis_theorem():
# This tests the parallel axis theorem matrix by comparing to test
# matrices.
# First case, 1 in all directions.
mat1 = Matrix(((2, -1, -1), (-1, 2, -1), (-1, -1, 2)))
assert pat_matrix(1, 1, 1, 1) == mat1
assert pat_matrix(2, 1, 1, 1) == 2*mat1
# Second case, 1 in x, 0 in all others
mat2 = Matrix(((0, 0, 0), (0, 1, 0), (0, 0, 1)))
assert pat_matrix(1, 1, 0, 0) == mat2
assert pat_matrix(2, 1, 0, 0) == 2*mat2
# Third case, 1 in y, 0 in all others
mat3 = Matrix(((1, 0, 0), (0, 0, 0), (0, 0, 1)))
assert pat_matrix(1, 0, 1, 0) == mat3
assert pat_matrix(2, 0, 1, 0) == 2*mat3
# Fourth case, 1 in z, 0 in all others
mat4 = Matrix(((1, 0, 0), (0, 1, 0), (0, 0, 0)))
assert pat_matrix(1, 0, 0, 1) == mat4
assert pat_matrix(2, 0, 0, 1) == 2*mat4
def test_Pauli():
#this and the following test are testing both Pauli and Dirac matrices
#and also that the general Matrix class works correctly in a real world
#situation
sigma1 = msigma(1)
sigma2 = msigma(2)
sigma3 = msigma(3)
assert sigma1 == sigma1
assert sigma1 != sigma2
# sigma*I -> I*sigma (see #354)
assert sigma1*sigma2 == sigma3*I
assert sigma3*sigma1 == sigma2*I
assert sigma2*sigma3 == sigma1*I
assert sigma1*sigma1 == eye(2)
assert sigma2*sigma2 == eye(2)
assert sigma3*sigma3 == eye(2)
assert sigma1*2*sigma1 == 2*eye(2)
assert sigma1*sigma3*sigma1 == -sigma3
def test_one_dof():
# This is for a 1 dof spring-mass-damper case.
# It is described in more detail in the KanesMethod docstring.
q, u = dynamicsymbols('q u')
qd, ud = dynamicsymbols('q u', 1)
m, c, k = symbols('m c k')
N = ReferenceFrame('N')
P = Point('P')
P.set_vel(N, u * N.x)
kd = [qd - u]
FL = [(P, (-k * q - c * u) * N.x)]
pa = Particle('pa', P, m)
BL = [pa]
KM = KanesMethod(N, [q], [u], kd)
KM.kanes_equations(FL, BL)
MM = KM.mass_matrix
forcing = KM.forcing
rhs = MM.inv() * forcing
assert expand(rhs[0]) == expand(-(q * k + u * c) / m)
assert KM.linearize() == \
(Matrix([[0, 1], [-k, -c]]), Matrix([]), Matrix([]))
def _mat_inv_mul(A, B):
"""
Computes A^-1 * B symbolically w/ substitution, where B is not
necessarily a vector, but can be a matrix.
"""
r1, c1 = A.shape
r2, c2 = B.shape
temp1 = Matrix(r1, c1, lambda i, j: Symbol('x' + str(j) + str(r1 * i)))
temp2 = Matrix(r2, c2, lambda i, j: Symbol('y' + str(j) + str(r2 * i)))
for i in range(len(temp1)):
if A[i] == 0:
temp1[i] = 0
for i in range(len(temp2)):
if B[i] == 0:
temp2[i] = 0
temp3 = []
for i in range(c2):
temp3.append(temp1.LDLsolve(temp2[:, i]))
temp3 = Matrix([i.T for i in temp3]).T
return temp3.subs(dict(list(zip(temp1, A)))).subs(dict(list(zip(temp2, B))))
def transpose(self):
"""Return transpose of matrix.
Examples
========
>>> from sympy import MatrixSymbol, BlockMatrix, ZeroMatrix
>>> from sympy.abc import l, m, n
>>> X = MatrixSymbol('X', n, n)
>>> Y = MatrixSymbol('Y', m ,m)
>>> Z = MatrixSymbol('Z', n, m)
>>> B = BlockMatrix([[X, Z], [ZeroMatrix(m,n), Y]])
>>> B.transpose()
Matrix([
[X', 0],
[Z', Y']])
>>> _.transpose()
Matrix([
[X, Z],
[0, Y]])
"""
return self._eval_transpose()
def deblock(B):
""" Flatten a BlockMatrix of BlockMatrices """
if not isinstance(B, BlockMatrix) or not B.blocks.has(BlockMatrix):
return B
wrap = lambda x: x if isinstance(x, BlockMatrix) else BlockMatrix([[x]])
bb = B.blocks.applyfunc(wrap) # everything is a block
from sympy import Matrix
try:
MM = Matrix(0, sum(bb[0, i].blocks.shape[1] for i in range(bb.shape[1])), [])
for row in range(0, bb.shape[0]):
M = Matrix(bb[row, 0].blocks)
for col in range(1, bb.shape[1]):
M = M.row_join(bb[row, col].blocks)
MM = MM.col_join(M)
return BlockMatrix(MM)
except ShapeError:
return B
def blockcut(expr, rowsizes, colsizes):
""" Cut a matrix expression into Blocks
>>> from sympy import ImmutableMatrix, blockcut
>>> M = ImmutableMatrix(4, 4, range(16))
>>> B = blockcut(M, (1, 3), (1, 3))
>>> type(B).__name__
'BlockMatrix'
>>> ImmutableMatrix(B.blocks[0, 1])
Matrix([[1, 2, 3]])
"""
rowbounds = bounds(rowsizes)
colbounds = bounds(colsizes)
return BlockMatrix([[MatrixSlice(expr, rowbound, colbound)
for colbound in colbounds]
for rowbound in rowbounds])
def test_valued_tensor_get_matrix():
numpy = import_module("numpy")
if numpy is None:
skip("numpy not installed.")
(A, B, AB, BA, C, Lorentz, E, px, py, pz, LorentzD, mu0, mu1, mu2, ndm, n0, n1,
n2, NA, NB, NC, minkowski, ba_matrix, ndm_matrix, i0, i1, i2, i3, i4) = _get_valued_base_test_variables()
matab = AB(i0, i1).get_matrix()
assert matab == Matrix([
[1, 0, 0, 0],
[0, -1, 0, 0],
[0, 0, -1, 0],
[0, 0, 0, -1],
])
# when alternating contravariant/covariant with [1, -1, -1, -1] metric
# it becomes the identity matrix:
assert AB(i0, -i1).get_matrix() == eye(4)
# covariant and contravariant forms:
assert A(i0).get_matrix() == Matrix([E, px, py, pz])
assert A(-i0).get_matrix() == Matrix([E, -px, -py, -pz])
def as_matrix(self):
"""Returns the data of the table in Matrix form.
Examples
========
>>> from sympy import TableForm
>>> t = TableForm([[5, 7], [4, 2], [10, 3]], headings='automatic')
>>> t
| 1 2
--------
1 | 5 7
2 | 4 2
3 | 10 3
>>> t.as_matrix()
Matrix([
[ 5, 7],
[ 4, 2],
[10, 3]])
"""
from sympy import Matrix
return Matrix(self._lines)
def test_Matrix_mul():
M = Matrix([[1, 2, 3], [x, y, x]])
m = matrix([[2, 4], [x, 6], [x, z**2]])
assert M*m == Matrix([
[ 2 + 5*x, 16 + 3*z**2],
[2*x + x*y + x**2, 4*x + 6*y + x*z**2],
])
assert m*M == Matrix([
[ 2 + 4*x, 4 + 4*y, 6 + 4*x],
[ 7*x, 2*x + 6*y, 9*x],
[x + x*z**2, 2*x + y*z**2, 3*x + x*z**2],
])
a = array([2])
assert a[0] * M == 2 * M
assert M * a[0] == 2 * M
def test_threaded():
@threaded
def function(expr, *args):
return 2*expr + sum(args)
assert function(Matrix([[x, y], [1, x]]), 1, 2) == \
Matrix([[2*x + 3, 2*y + 3], [5, 2*x + 3]])
assert function(Eq(x, y), 1, 2) == Eq(2*x + 3, 2*y + 3)
assert function([x, y], 1, 2) == [2*x + 3, 2*y + 3]
assert function((x, y), 1, 2) == (2*x + 3, 2*y + 3)
assert function(set([x, y]), 1, 2) == set([2*x + 3, 2*y + 3])
@threaded
def function(expr, n):
return expr**n
assert function(x + y, 2) == x**2 + y**2
assert function(x, 2) == x**2
def main():
Format()
a = Matrix(2, 2, (1, 2, 3, 4))
b = Matrix(2, 1, (5, 6))
c = a*b
print(a, b, '=', c)
x, y = symbols('x, y')
d = Matrix(1, 2, (x**3, y**3))
e = Matrix(2, 2, (x**2, 2*x*y, 2*x*y, y**2))
f = d*e
print('%', d, e, '=', f)
xdvi()
return
def block_matrix(blocks):
"""
Construct a matrix where the elements are specified by the block structure
by joining the blocks appropriately.
Parameters
----------
blocks : two level deep iterable of sympy Matrix objects
The block specification of the matrices used to construct the block
matrix.
Returns
-------
matrix : sympy Matrix
A matrix whose elements are the elements of the blocks with the
specified block structure.
"""
return sp.Matrix.col_join(
*tuple(
sp.Matrix.row_join(
*tuple(mat for mat in row)) for row in blocks
)
)
def two_wheel_3d_deriv():
""" """
x1, x2, x3, x4, x5, x6, x7 = sympy.symbols("x1,x2,x3,x4,x5,x6,x7")
dt = sympy.symbols("dt")
# x1 - x
# x2 - y
# x3 - z
# x4 - theta
# x5 - v
# x6 - omega
# x7 - vz
# x, y, z, theta, v, omega, vz
f1 = x1 + x5 * sympy.cos(x4) * dt
f2 = x2 + x5 * sympy.sin(x4) * dt
f3 = x3 + x7 * dt
f4 = x4 + x6 * dt
f5 = x5
f6 = x6
f7 = x7
F = sympy.Matrix([f1, f2, f3, f4, f5, f6, f7])
pprint(F.jacobian([x1, x2, x3, x4, x5, x6, x7]))
def generate_jac_lambda(self,do_cse=False):
"""
translates the symbolic Jacobain to a function using SymEngines `Lambdify` or `LambdifyCSE`. If the symbolic Jacobian has not been generated, it is generated by calling `generate_jac_sym`.
Parameters
----------
do_cse : boolean
Whether a common-subexpression detection, namely `LambdifyCSE`, should be used.
"""
self._prepare_lambdas()
self._generate_jac_sym()
jac_matrix = symengine.Matrix([
[ordered_subs(entry,self._lambda_subs) for entry in line]
for line in self.jac_sym
])
lambdify = symengine.LambdifyCSE if do_cse else symengine.Lambdify
core_jac = lambdify(self._lambda_args,jac_matrix)
# self.jac = lambda t,Y,*c_pars: core_jac(np.hstack([t,Y,c_pars]))
self.jac = _lambda_wrapper(core_jac,self.control_pars)
def test_quantum_fourier():
assert QFT(0, 3).decompose() == \
SwapGate(0, 2)*HadamardGate(0)*CGate((0,), PhaseGate(1)) * \
HadamardGate(1)*CGate((0,), TGate(2))*CGate((1,), PhaseGate(2)) * \
HadamardGate(2)
assert IQFT(0, 3).decompose() == \
HadamardGate(2)*CGate((1,), RkGate(2, -2))*CGate((0,), RkGate(2, -3)) * \
HadamardGate(1)*CGate((0,), RkGate(1, -2))*HadamardGate(0)*SwapGate(0, 2)
assert represent(QFT(0, 3), nqubits=3) == \
Matrix([[exp(2*pi*I/8)**(i*j % 8)/sqrt(8) for i in range(8)] for j in range(8)])
assert QFT(0, 4).decompose() # non-trivial decomposition
assert qapply(QFT(0, 3).decompose()*Qubit(0, 0, 0)).expand() == qapply(
HadamardGate(0)*HadamardGate(1)*HadamardGate(2)*Qubit(0, 0, 0)
).expand()
def test_matrix_to_density():
mat = Matrix([[0, 0], [0, 1]])
assert matrix_to_density(mat) == Density([Qubit('1'), 1])
mat = Matrix([[1, 0], [0, 0]])
assert matrix_to_density(mat) == Density([Qubit('0'), 1])
mat = Matrix([[0, 0], [0, 0]])
assert matrix_to_density(mat) == 0
mat = Matrix([[0, 0, 0, 0],
[0, 0, 0, 0],
[0, 0, 1, 0],
[0, 0, 0, 0]])
assert matrix_to_density(mat) == Density([Qubit('10'), 1])
mat = Matrix([[1, 0, 0, 0],
[0, 0, 0, 0],
[0, 0, 0, 0],
[0, 0, 0, 0]])
assert matrix_to_density(mat) == Density([Qubit('00'), 1])
def cache_matrix(self, name, m):
"""Cache a matrix by its name.
Parameters
----------
name : str
A descriptive name for the matrix, like "identity2".
m : list of lists
The raw matrix data as a sympy Matrix.
"""
try:
self._sympy_matrix(name, m)
except ImportError:
pass
try:
self._numpy_matrix(name, m)
except ImportError:
pass
try:
self._scipy_sparse_matrix(name, m)
except ImportError:
pass
def get_matrix(self, name, format):
"""Get a cached matrix by name and format.
Parameters
----------
name : str
A descriptive name for the matrix, like "identity2".
format : str
The format desired ('sympy', 'numpy', 'scipy.sparse')
"""
m = self._cache.get((name, format))
if m is not None:
return m
raise NotImplementedError(
'Matrix with name %s and format %s is not available.' %
(name, format)
)
def _qsympify_sequence(seq):
"""Convert elements of a sequence to standard form.
This is like sympify, but it performs special logic for arguments passed
to QExpr. The following conversions are done:
* (list, tuple, Tuple) => _qsympify_sequence each element and convert
sequence to a Tuple.
* basestring => Symbol
* Matrix => Matrix
* other => sympify
Strings are passed to Symbol, not sympify to make sure that variables like
'pi' are kept as Symbols, not the SymPy built-in number subclasses.
Examples
========
>>> from sympy.physics.quantum.qexpr import _qsympify_sequence
>>> _qsympify_sequence((1,2,[3,4,[1,]]))
(1, 2, (3, 4, (1,)))
"""
return tuple(__qsympify_sequence_helper(seq))
def test_parallel_axis_theorem():
# This tests the parallel axis theorem matrix by comparing to test
# matrices.
# First case, 1 in all directions.
mat1 = Matrix(((2, -1, -1), (-1, 2, -1), (-1, -1, 2)))
assert pat_matrix(1, 1, 1, 1) == mat1
assert pat_matrix(2, 1, 1, 1) == 2*mat1
# Second case, 1 in x, 0 in all others
mat2 = Matrix(((0, 0, 0), (0, 1, 0), (0, 0, 1)))
assert pat_matrix(1, 1, 0, 0) == mat2
assert pat_matrix(2, 1, 0, 0) == 2*mat2
# Third case, 1 in y, 0 in all others
mat3 = Matrix(((1, 0, 0), (0, 0, 0), (0, 0, 1)))
assert pat_matrix(1, 0, 1, 0) == mat3
assert pat_matrix(2, 0, 1, 0) == 2*mat3
# Fourth case, 1 in z, 0 in all others
mat4 = Matrix(((1, 0, 0), (0, 1, 0), (0, 0, 0)))
assert pat_matrix(1, 0, 0, 1) == mat4
assert pat_matrix(2, 0, 0, 1) == 2*mat4
def test_Pauli():
#this and the following test are testing both Pauli and Dirac matrices
#and also that the general Matrix class works correctly in a real world
#situation
sigma1 = msigma(1)
sigma2 = msigma(2)
sigma3 = msigma(3)
assert sigma1 == sigma1
assert sigma1 != sigma2
# sigma*I -> I*sigma (see #354)
assert sigma1*sigma2 == sigma3*I
assert sigma3*sigma1 == sigma2*I
assert sigma2*sigma3 == sigma1*I
assert sigma1*sigma1 == eye(2)
assert sigma2*sigma2 == eye(2)
assert sigma3*sigma3 == eye(2)
assert sigma1*2*sigma1 == 2*eye(2)
assert sigma1*sigma3*sigma1 == -sigma3
def transpose(self):
"""Return transpose of matrix.
Examples
========
>>> from sympy import MatrixSymbol, BlockMatrix, ZeroMatrix
>>> from sympy.abc import l, m, n
>>> X = MatrixSymbol('X', n, n)
>>> Y = MatrixSymbol('Y', m ,m)
>>> Z = MatrixSymbol('Z', n, m)
>>> B = BlockMatrix([[X, Z], [ZeroMatrix(m,n), Y]])
>>> B.transpose()
Matrix([
[X.T, 0],
[Z.T, Y.T]])
>>> _.transpose()
Matrix([
[X, Z],
[0, Y]])
"""
return self._eval_transpose()
def deblock(B):
""" Flatten a BlockMatrix of BlockMatrices """
if not isinstance(B, BlockMatrix) or not B.blocks.has(BlockMatrix):
return B
wrap = lambda x: x if isinstance(x, BlockMatrix) else BlockMatrix([[x]])
bb = B.blocks.applyfunc(wrap) # everything is a block
from sympy import Matrix
try:
MM = Matrix(0, sum(bb[0, i].blocks.shape[1] for i in range(bb.shape[1])), [])
for row in range(0, bb.shape[0]):
M = Matrix(bb[row, 0].blocks)
for col in range(1, bb.shape[1]):
M = M.row_join(bb[row, col].blocks)
MM = MM.col_join(M)
return BlockMatrix(MM)
except ShapeError:
return B
def blockcut(expr, rowsizes, colsizes):
""" Cut a matrix expression into Blocks
>>> from sympy import ImmutableMatrix, blockcut
>>> M = ImmutableMatrix(4, 4, range(16))
>>> B = blockcut(M, (1, 3), (1, 3))
>>> type(B).__name__
'BlockMatrix'
>>> ImmutableMatrix(B.blocks[0, 1])
Matrix([[1, 2, 3]])
"""
rowbounds = bounds(rowsizes)
colbounds = bounds(colsizes)
return BlockMatrix([[MatrixSlice(expr, rowbound, colbound)
for colbound in colbounds]
for rowbound in rowbounds])
