def test_commutator():
A = Operator('A')
B = Operator('B')
c = Commutator(A, B)
c_tall = Commutator(A**2, B)
assert str(c) == '[A,B]'
assert pretty(c) == '[A,B]'
assert upretty(c) == u('[A,B]')
assert latex(c) == r'\left[A,B\right]'
sT(c, "Commutator(Operator(Symbol('A')),Operator(Symbol('B')))")
assert str(c_tall) == '[A**2,B]'
ascii_str = \
"""\
[ 2 ]\n\
[A ,B]\
"""
ucode_str = \
u("""\
? 2 ?\n\
?A ,B?\
""")
assert pretty(c_tall) == ascii_str
assert upretty(c_tall) == ucode_str
assert latex(c_tall) == r'\left[\left(A\right)^{2},B\right]'
sT(c_tall, "Commutator(Pow(Operator(Symbol('A')), Integer(2)),Operator(Symbol('B')))")
python类Pow()的实例源码
def _expand_powers(factors):
"""
Helper function for normal_ordered_form and normal_order: Expand a
power expression to a multiplication expression so that that the
expression can be handled by the normal ordering functions.
"""
new_factors = []
for factor in factors.args:
if isinstance(factor, Pow):
for n in range(factor.args[1]):
new_factors.append(factor.args[0])
else:
new_factors.append(factor)
return new_factors
def _extractFunctionsFromSympy(self, expr, top=False):
newVars = []
newArgs = []
isExpensive = False
for arg in expr.args:
(newVarsFromArg, newArg,dontcare) = self._extractFunctionsFromSympy(arg)
newVars += newVarsFromArg
newArgs.append(newArg)
if newVars:
expr = expr.func(*newArgs)
if isinstance(expr.func, type(sympy.Function)) or (
expr.func == sympy.Pow and (
not expr.exp.is_constant or
int(expr.exp) != expr.exp)):
if top:
isExpensive = True
else:
newSym = self.addSSA("_expensive_functions", expr)
expr = newSym
newVars.append(newSym)
return (newVars, expr, isExpensive)
def _collect_factor_and_dimension(expr):
if isinstance(expr, Quantity):
return expr.scale_factor, expr.dimension
elif isinstance(expr, Mul):
factor = 1
dimension = 1
for arg in expr.args:
arg_factor, arg_dim = Quantity._collect_factor_and_dimension(arg)
factor *= arg_factor
dimension *= arg_dim
return factor, dimension
elif isinstance(expr, Pow):
factor, dim = Quantity._collect_factor_and_dimension(expr.base)
return factor ** expr.exp, dim ** expr.exp
elif isinstance(expr, Add):
raise NotImplementedError
else:
return 1, 1
def _expand_powers(factors):
"""
Helper function for normal_ordered_form and normal_order: Expand a
power expression to a multiplication expression so that that the
expression can be handled by the normal ordering functions.
"""
new_factors = []
for factor in factors.args:
if (isinstance(factor, Pow)
and isinstance(factor.args[1], Integer)
and factor.args[1] > 0):
for n in range(factor.args[1]):
new_factors.append(factor.args[0])
else:
new_factors.append(factor)
return new_factors
def _print_Pow(self, power):
from sympy.simplify.simplify import fraction
b, e = power.as_base_exp()
if power.is_commutative:
if e is S.NegativeOne:
return prettyForm("1")/self._print(b)
n, d = fraction(e)
if n is S.One and d.is_Atom and not e.is_Integer:
return self._print_nth_root(b, e)
if e.is_Rational and e < 0:
return prettyForm("1")/self._print(Pow(b, -e, evaluate=False))
if b.is_Relational:
return prettyForm(*self._print(b).parens()).__pow__(self._print(e))
return self._print(b)**self._print(e)
def _real_to_rational(expr, tolerance=None):
"""
Replace all reals in expr with rationals.
>>> from sympy import nsimplify
>>> from sympy.abc import x
>>> nsimplify(.76 + .1*x**.5, rational=True)
sqrt(x)/10 + 19/25
"""
p = expr
reps = {}
reduce_num = None
if tolerance is not None and tolerance < 1:
reduce_num = ceiling(1/tolerance)
for float in p.atoms(C.Float):
key = float
if reduce_num is not None:
r = Rational(float).limit_denominator(reduce_num)
elif (tolerance is not None and tolerance >= 1 and
float.is_Integer is False):
r = Rational(tolerance*round(float/tolerance)
).limit_denominator(int(tolerance))
else:
r = nsimplify(float, rational=False)
# e.g. log(3).n() -> log(3) instead of a Rational
if not r.is_Rational:
if float < 0:
float = -float
d = Pow(10, int((mpmath.log(float)/mpmath.log(10))))
r = -Rational(str(float/d))*d
elif float > 0:
d = Pow(10, int((mpmath.log(float)/mpmath.log(10))))
r = Rational(str(float/d))*d
else:
r = Integer(0)
reps[key] = r
return p.subs(reps, simultaneous=True)
def _get_min_qubits(a_gate):
if isinstance(a_gate, Pow):
return a_gate.base.min_qubits
else:
return a_gate.min_qubits
def _sympystr(self, printer, *args):
from sympy.printing.str import sstr
length = len(self.args)
s = ''
for i in range(length):
if isinstance(self.args[i], (Add, Pow, Mul)):
s = s + '('
s = s + sstr(self.args[i])
if isinstance(self.args[i], (Add, Pow, Mul)):
s = s + ')'
if i != length - 1:
s = s + 'x'
return s
def _print_Pow(self, power):
from sympy.simplify.simplify import fraction
b, e = power.as_base_exp()
if power.is_commutative:
if e is S.NegativeOne:
return prettyForm("1")/self._print(b)
n, d = fraction(e)
if n is S.One and d.is_Atom and not e.is_Integer:
return self._print_nth_root(b, e)
if e.is_Rational and e < 0:
return prettyForm("1")/self._print(C.Pow(b, -e, evaluate=False))
return self._print(b)**self._print(e)
def _print_ProductSet(self, p):
if len(p.sets) > 1 and not has_variety(p.sets):
from sympy import Pow
return self._print(Pow(p.sets[0], len(p.sets), evaluate=False))
else:
prod_char = u('\xd7')
return self._print_seq(p.sets, None, None, ' %s ' % prod_char,
parenthesize=lambda set: set.is_Union or set.is_Intersection)
def _eval_nseries(self, x, n, logx):
# NOTE this is fairly inefficient
from sympy import log, EulerGamma, Pow
n += 1
if self.args[0].subs(x, 0) != 0:
return super(TrigonometricIntegral, self)._eval_nseries(x, n, logx)
baseseries = self._trigfunc(x)._eval_nseries(x, n, logx)
if self._trigfunc(0) != 0:
baseseries -= 1
baseseries = baseseries.replace(Pow, lambda t, n: t**n/n, simultaneous=False)
if self._trigfunc(0) != 0:
baseseries += EulerGamma + log(x)
return baseseries.subs(x, self.args[0])._eval_nseries(x, n, logx)
def main():
x = Pow(2, 50, evaluate=False)
y = Pow(10, -50, evaluate=False)
# A large, unevaluated expression
m = Mul(x, y, evaluate=False)
# Evaluating the expression
e = S(2)**50/S(10)**50
print("%s == %s" % (m, e))
def expr(self):
"""
Process the atomic expression, by returning a sympy expression
"""
if DEBUG:
print("> Power ")
expr = self.op.expr
if self.exponent is None:
return expr
else:
exponent = self.exponent.expr
return Pow(expr, exponent)
def expr(self):
"""
Process the term, by returning a sympy expression
"""
if DEBUG:
print("> Term ")
ret = self.op[0].expr
for operation, operand in zip(self.op[1::2], self.op[2::2]):
if operation == '*':
ret = Mul(ret, operand.expr)
else:
a = Pow(operand.expr, -1)
ret = Mul(ret, a)
return ret
def powerProcess(self, x, y):
if (not x.astUnit.isNull()) and y.sympy.is_constant():
newUnit = x.astUnit ** float(y.sympy)
else:
newUnit = ASTUnit.null()
return AST(sympy.Pow(x.sympy, y.sympy), newUnit)
############################################
def p_multiplicitiveExpr_impl(self, p):
'''multiplicitiveExpr : multiplicitiveExpr '*' unaryExpr
| multiplicitiveExpr '/' unaryExpr
'''
lhs = p[1]
rhs = p[3]
if p[2] == '/':
rhs = AST(sympy.Pow(rhs.sympy,sympy.Integer(-1)), rhs.astUnit ** -1)
ret = AST(sympy.Mul(lhs.sympy,rhs.sympy), lhs.astUnit*rhs.astUnit)
p[0] = ret
def test_Pow_is_algebraic():
e = Symbol('e', algebraic=True)
assert Pow(1, e, evaluate=False).is_algebraic
assert Pow(0, e, evaluate=False).is_algebraic
a = Symbol('a', algebraic=True)
na = Symbol('na', algebraic=False)
ia = Symbol('ia', algebraic=True, irrational=True)
ib = Symbol('ib', algebraic=True, irrational=True)
r = Symbol('r', rational=True)
x = Symbol('x')
assert (a**r).is_algebraic
assert (a**x).is_algebraic is None
assert (na**r).is_algebraic is None
assert (ia**r).is_algebraic
assert (ia**ib).is_algebraic is False
assert (a**e).is_algebraic is None
# Gelfond-Schneider constant:
assert Pow(2, sqrt(2), evaluate=False).is_algebraic is False
assert Pow(S.GoldenRatio, sqrt(3), evaluate=False).is_algebraic is False
# issue 8649
t = Symbol('t', real=True, transcendental=True)
n = Symbol('n', integer=True)
assert (t**n).is_algebraic is None
assert (t**n).is_integer is None
def test_Pow_is_prime():
from sympy import Pow
x = Symbol('x', positive=True, integer=True)
y = Symbol('y', positive=True, integer=True)
assert (x**y).is_prime is None
x = Symbol('x', positive=True)
assert (x**y).is_prime is None
assert Pow(6, S.One, evaluate=False).is_prime is False
assert Pow(9, S.Half, evaluate=False).is_prime is True
assert Pow(5, S.One, evaluate=False).is_prime is True
def test_dim_simplify_pow():
assert dim_simplify(Pow(L, 2)) == L**2
assert dim_simplify(L**2) == L**2
def dim_simplify(expr):
"""
NOTE: this function could be deprecated in the future.
Simplify expression by recursively evaluating the dimension arguments.
This function proceeds to a very rough dimensional analysis. It tries to
simplify expression with dimensions, and it deletes all what multiplies a
dimension without being a dimension. This is necessary to avoid strange
behavior when Add(L, L) be transformed into Mul(2, L).
"""
if isinstance(expr, Dimension):
return expr
if isinstance(expr, Pow):
return dim_simplify(expr.base)**dim_simplify(expr.exp)
elif isinstance(expr, Function):
return dim_simplify(expr.args[0])
elif isinstance(expr, Add):
if (all(isinstance(arg, Dimension) for arg in expr.args) or
all(arg.is_dimensionless for arg in expr.args if isinstance(arg, Dimension))):
return reduce(lambda x, y: x.add(y), expr.args)
else:
raise ValueError("Dimensions cannot be added: %s" % expr)
elif isinstance(expr, Mul):
return Dimension(Mul(*[dim_simplify(i).name for i in expr.args if isinstance(i, Dimension)]))
raise ValueError("Cannot be simplifed: %s", expr)
def _get_min_qubits(a_gate):
if isinstance(a_gate, Pow):
return a_gate.base.min_qubits
else:
return a_gate.min_qubits
def test_anticommutator():
A = Operator('A')
B = Operator('B')
ac = AntiCommutator(A, B)
ac_tall = AntiCommutator(A**2, B)
assert str(ac) == '{A,B}'
assert pretty(ac) == '{A,B}'
assert upretty(ac) == u'{A,B}'
assert latex(ac) == r'\left\{A,B\right\}'
sT(ac, "AntiCommutator(Operator(Symbol('A')),Operator(Symbol('B')))")
assert str(ac_tall) == '{A**2,B}'
ascii_str = \
"""\
/ 2 \\\n\
<A ,B>\n\
\\ /\
"""
ucode_str = \
u("""\
? 2 ?\n\
?A ,B?\n\
? ?\
""")
assert pretty(ac_tall) == ascii_str
assert upretty(ac_tall) == ucode_str
assert latex(ac_tall) == r'\left\{A^{2},B\right\}'
sT(ac_tall, "AntiCommutator(Pow(Operator(Symbol('A')), Integer(2)),Operator(Symbol('B')))")
def _sympystr(self, printer, *args):
from sympy.printing.str import sstr
length = len(self.args)
s = ''
for i in range(length):
if isinstance(self.args[i], (Add, Pow, Mul)):
s = s + '('
s = s + sstr(self.args[i])
if isinstance(self.args[i], (Add, Pow, Mul)):
s = s + ')'
if i != length - 1:
s = s + 'x'
return s
def _print_ProductSet(self, p):
if len(p.sets) > 1 and not has_variety(p.sets):
from sympy import Pow
return self._print(Pow(p.sets[0], len(p.sets), evaluate=False))
else:
prod_char = u"\N{MULTIPLICATION SIGN}" if self._use_unicode else 'x'
return self._print_seq(p.sets, None, None, ' %s ' % prod_char,
parenthesize=lambda set: set.is_Union or
set.is_Intersection or set.is_ProductSet)
def _eval_nseries(self, x, n, logx):
# NOTE this is fairly inefficient
from sympy import log, EulerGamma, Pow
n += 1
if self.args[0].subs(x, 0) != 0:
return super(TrigonometricIntegral, self)._eval_nseries(x, n, logx)
baseseries = self._trigfunc(x)._eval_nseries(x, n, logx)
if self._trigfunc(0) != 0:
baseseries -= 1
baseseries = baseseries.replace(Pow, lambda t, n: t**n/n, simultaneous=False)
if self._trigfunc(0) != 0:
baseseries += EulerGamma + log(x)
return baseseries.subs(x, self.args[0])._eval_nseries(x, n, logx)
def _separatevars(expr, force):
if len(expr.free_symbols) == 1:
return expr
# don't destroy a Mul since much of the work may already be done
if expr.is_Mul:
args = list(expr.args)
changed = False
for i, a in enumerate(args):
args[i] = separatevars(a, force)
changed = changed or args[i] != a
if changed:
expr = expr.func(*args)
return expr
# get a Pow ready for expansion
if expr.is_Pow:
expr = Pow(separatevars(expr.base, force=force), expr.exp)
# First try other expansion methods
expr = expr.expand(mul=False, multinomial=False, force=force)
_expr, reps = posify(expr) if force else (expr, {})
expr = factor(_expr).subs(reps)
if not expr.is_Add:
return expr
# Find any common coefficients to pull out
args = list(expr.args)
commonc = args[0].args_cnc(cset=True, warn=False)[0]
for i in args[1:]:
commonc &= i.args_cnc(cset=True, warn=False)[0]
commonc = Mul(*commonc)
commonc = commonc.as_coeff_Mul()[1] # ignore constants
commonc_set = commonc.args_cnc(cset=True, warn=False)[0]
# remove them
for i, a in enumerate(args):
c, nc = a.args_cnc(cset=True, warn=False)
c = c - commonc_set
args[i] = Mul(*c)*Mul(*nc)
nonsepar = Add(*args)
if len(nonsepar.free_symbols) > 1:
_expr = nonsepar
_expr, reps = posify(_expr) if force else (_expr, {})
_expr = (factor(_expr)).subs(reps)
if not _expr.is_Add:
nonsepar = _expr
return commonc*nonsepar
def __new__(cls, *args):
""" Construct a Trace object.
Parameters
==========
args = sympy expression
indices = tuple/list if indices, optional
"""
# expect no indices,int or a tuple/list/Tuple
if (len(args) == 2):
if not isinstance(args[1], (list, Tuple, tuple)):
indices = Tuple(args[1])
else:
indices = Tuple(*args[1])
expr = args[0]
elif (len(args) == 1):
indices = Tuple()
expr = args[0]
else:
raise ValueError("Arguments to Tr should be of form "
"(expr[, [indices]])")
if isinstance(expr, Matrix):
return expr.trace()
elif hasattr(expr, 'trace') and callable(expr.trace):
#for any objects that have trace() defined e.g numpy
return expr.trace()
elif isinstance(expr, Add):
return Add(*[Tr(arg, indices) for arg in expr.args])
elif isinstance(expr, Mul):
c_part, nc_part = expr.args_cnc()
if len(nc_part) == 0:
return Mul(*c_part)
else:
obj = Expr.__new__(cls, Mul(*nc_part), indices )
#this check is needed to prevent cached instances
#being returned even if len(c_part)==0
return Mul(*c_part)*obj if len(c_part) > 0 else obj
elif isinstance(expr, Pow):
if (_is_scalar(expr.args[0]) and
_is_scalar(expr.args[1])):
return expr
else:
return Expr.__new__(cls, expr, indices)
else:
if (_is_scalar(expr)):
return expr
return Expr.__new__(cls, expr, indices)
def gate_sort(circuit):
"""Sorts the gates while keeping track of commutation relations
This function uses a bubble sort to rearrange the order of gate
application. Keeps track of Quantum computations special commutation
relations (e.g. things that apply to the same Qubit do not commute with
each other)
circuit is the Mul of gates that are to be sorted.
"""
# Make sure we have an Add or Mul.
if isinstance(circuit, Add):
return sum(gate_sort(t) for t in circuit.args)
if isinstance(circuit, Pow):
return gate_sort(circuit.base)**circuit.exp
elif isinstance(circuit, Gate):
return circuit
if not isinstance(circuit, Mul):
return circuit
changes = True
while changes:
changes = False
circ_array = circuit.args
for i in xrange(len(circ_array) - 1):
# Go through each element and switch ones that are in wrong order
if isinstance(circ_array[i], (Gate, Pow)) and \
isinstance(circ_array[i + 1], (Gate, Pow)):
# If we have a Pow object, look at only the base
first_base, first_exp = circ_array[i].as_base_exp()
second_base, second_exp = circ_array[i + 1].as_base_exp()
# Use sympy's hash based sorting. This is not mathematical
# sorting, but is rather based on comparing hashes of objects.
# See Basic.compare for details.
if first_base.compare(second_base) > 0:
if Commutator(first_base, second_base).doit() == 0:
new_args = (circuit.args[:i] + (circuit.args[i + 1],) +
(circuit.args[i],) + circuit.args[i + 2:])
circuit = Mul(*new_args)
circ_array = circuit.args
changes = True
break
if AntiCommutator(first_base, second_base).doit() == 0:
new_args = (circuit.args[:i] + (circuit.args[i + 1],) +
(circuit.args[i],) + circuit.args[i + 2:])
sign = Integer(-1)**(first_exp*second_exp)
circuit = sign*Mul(*new_args)
circ_array = circuit.args
changes = True
break
return circuit
#-----------------------------------------------------------------------------
# Utility functions
#-----------------------------------------------------------------------------
def test_operator():
a = Operator('A')
b = Operator('B', Symbol('t'), S(1)/2)
inv = a.inv()
f = Function('f')
x = symbols('x')
d = DifferentialOperator(Derivative(f(x), x), f(x))
op = OuterProduct(Ket(), Bra())
assert str(a) == 'A'
assert pretty(a) == 'A'
assert upretty(a) == u('A')
assert latex(a) == 'A'
sT(a, "Operator(Symbol('A'))")
assert str(inv) == 'A**(-1)'
ascii_str = \
"""\
-1\n\
A \
"""
ucode_str = \
u("""\
-1\n\
A \
""")
assert pretty(inv) == ascii_str
assert upretty(inv) == ucode_str
assert latex(inv) == r'\left(A\right)^{-1}'
sT(inv, "Pow(Operator(Symbol('A')), Integer(-1))")
assert str(d) == 'DifferentialOperator(Derivative(f(x), x),f(x))'
ascii_str = \
"""\
/d \\\n\
DifferentialOperator|--(f(x)),f(x)|\n\
\dx /\
"""
ucode_str = \
u("""\
?d ?\n\
DifferentialOperator???(f(x)),f(x)?\n\
?dx ?\
""")
assert pretty(d) == ascii_str
assert upretty(d) == ucode_str
assert latex(d) == \
r'DifferentialOperator\left(\frac{d}{d x} f{\left (x \right )},f{\left (x \right )}\right)'
sT(d, "DifferentialOperator(Derivative(Function('f')(Symbol('x')), Symbol('x')),Function('f')(Symbol('x')))")
assert str(b) == 'Operator(B,t,1/2)'
assert pretty(b) == 'Operator(B,t,1/2)'
assert upretty(b) == u('Operator(B,t,1/2)')
assert latex(b) == r'Operator\left(B,t,\frac{1}{2}\right)'
sT(b, "Operator(Symbol('B'),Symbol('t'),Rational(1, 2))")
assert str(op) == '|psi><psi|'
assert pretty(op) == '|psi><psi|'
assert upretty(op) == u('??????')
assert latex(op) == r'{\left|\psi\right\rangle }{\left\langle \psi\right|}'
sT(op, "OuterProduct(Ket(Symbol('psi')),Bra(Symbol('psi')))")
