Merge branch 'master' into dsystem_solver_api

bin/strip_whitespace

@@ -11,8 +11,8 @@ def strip_file(filename, write, report):
     # to re-add \n to a few right-stripped lines. The hit flag will keep us
     # from unecessarily re-writing files with no changes.
 
-    # without "b" all lines appear as \n terminated
-    with open(filename, 'rb') as f:
+    # newline="" keeps current line endings
+    with open(filename, newline="") as f:
         lines = f.readlines()
     hit = False
     cr = False
@@ -51,8 +51,8 @@ def strip_file(filename, write, report):
         print("%s, %d extra newlines at eof" % (filename, extra))
 
     if write and hit:
-        # without "b" the lines may be written in sys-dep format
-        with open(filename, "wb") as f:
+        # newline="" leaves line endings unchanged
+        with open(filename, "w", newline="") as f:
             f.writelines(lines)
 
 

doc/src/modules/physics/mechanics/masses.rst

@@ -227,6 +227,8 @@ the reference frame is created and the kinematics are done. ::
   >>> from sympy.physics.mechanics import dynamicsymbols, ReferenceFrame
   >>> from sympy.physics.mechanics import RigidBody, Particle, Point, outer
   >>> from sympy.physics.mechanics import linear_momentum, angular_momentum
+  >>> from sympy.physics.vector import init_vprinting
+  >>> init_vprinting(pretty_print=False)
   >>> m, M, l1 = symbols('m M l1')
   >>> q1d = dynamicsymbols('q1d')
   >>> N = ReferenceFrame('N')
@@ -331,6 +333,8 @@ by going through an identical procedure. ::
   >>> from sympy.physics.mechanics import dynamicsymbols, ReferenceFrame, outer
   >>> from sympy.physics.mechanics import RigidBody, Particle
   >>> from sympy.physics.mechanics import kinetic_energy, potential_energy, Point
+  >>> from sympy.physics.vector import init_vprinting
+  >>> init_vprinting(pretty_print=False)
   >>> m, M, l1, g, h, H = symbols('m M l1 g h H')
   >>> omega = dynamicsymbols('omega')
   >>> N = ReferenceFrame('N')
@@ -352,9 +356,9 @@ The user can then determine the kinetic energy of any number of entities of the
 system: ::
 
   >>> kinetic_energy(N, Pa)
-  2*l1**2*m*omega(t)**2
+  2*l1**2*m*omega**2
   >>> kinetic_energy(N, Pa, A)
-  M*l1**2*omega(t)**2/2 + 2*l1**2*m*omega(t)**2 + omega(t)**2/2
+  M*l1**2*omega**2/2 + 2*l1**2*m*omega**2 + omega**2/2
 
 It should be noted that the user can determine either kinetic energy relative
 to any frame in :mod:`sympy.physics.mechanics` as the user is allowed to specify the
@@ -375,7 +379,9 @@ be determined: ::
 One can also determine the Lagrangian for this system: ::
 
   >>> from sympy.physics.mechanics import Lagrangian
+  >>> from sympy.physics.vector import init_vprinting
+  >>> init_vprinting(pretty_print=False)
   >>> Lagrangian(N, Pa, A)
-  -H*M*g + M*l1**2*omega(t)**2/2 - g*h*m + 2*l1**2*m*omega(t)**2 + omega(t)**2/2
+  -H*M*g + M*l1**2*omega**2/2 - g*h*m + 2*l1**2*m*omega**2 + omega**2/2
 
 Please refer to the docstrings to learn more about each function.

doc/src/modules/physics/mechanics/symsystem.rst

@@ -30,6 +30,8 @@ frame is rotated 90 degrees.
     >>> from sympy.physics.mechanics import (dynamicsymbols, ReferenceFrame,
     ...                                      Particle, Point)
     >>> import sympy.physics.mechanics.system as system
+    >>> from sympy.physics.vector import init_vprinting
+    >>> init_vprinting(pretty_print=False)
 
 The first step will be to initialize all of the dynamic and constant symbols. ::
 
@@ -146,54 +148,54 @@ equations of motion formats. ::
 
     >>> symsystem1.states
     Matrix([
-    [     x(t)],
-    [     y(t)],
-    [     u(t)],
-    [     v(t)],
-    [lambda(t)]])
+    [     x],
+    [     y],
+    [     u],
+    [     v],
+    [lambda]])
     >>> symsystem2.coordinates
     Matrix([
-    [x(t)],
-    [y(t)]])
+    [x],
+    [y]])
     >>> symsystem3.speeds
     Matrix([
-    [u(t)],
-    [v(t)]])
+    [u],
+    [v]])
     >>> symsystem1.comb_explicit_rhs
     Matrix([
-    [                                   u(t)],
-    [                                   v(t)],
-    [(-g*y(t) + u(t)**2 + v(t)**2)*x(t)/l**2],
-    [(-g*y(t) + u(t)**2 + v(t)**2)*y(t)/l**2],
-    [   m*(-g*y(t) + u(t)**2 + v(t)**2)/l**2]])
+    [                          u],
+    [                          v],
+    [(-g*y + u**2 + v**2)*x/l**2],
+    [(-g*y + u**2 + v**2)*y/l**2],
+    [m*(-g*y + u**2 + v**2)/l**2]])
     >>> symsystem2.comb_implicit_rhs
     Matrix([
-    [                       u(t)],
-    [                       v(t)],
-    [                          0],
-    [                          0],
-    [-g*y(t) + u(t)**2 + v(t)**2]])
+    [                 u],
+    [                 v],
+    [                 0],
+    [                 0],
+    [-g*y + u**2 + v**2]])
     >>> symsystem2.comb_implicit_mat
     Matrix([
-    [1, 0, 0, 0,       0],
-    [0, 1, 0, 0,       0],
-    [0, 0, 1, 0, -x(t)/m],
-    [0, 0, 0, 1, -y(t)/m],
-    [0, 0, 0, 0,  l**2/m]])
+    [1, 0, 0, 0,      0],
+    [0, 1, 0, 0,      0],
+    [0, 0, 1, 0,   -x/m],
+    [0, 0, 0, 1,   -y/m],
+    [0, 0, 0, 0, l**2/m]])
     >>> symsystem3.dyn_implicit_rhs
     Matrix([
-    [                          0],
-    [                          0],
-    [-g*y(t) + u(t)**2 + v(t)**2]])
+    [                 0],
+    [                 0],
+    [-g*y + u**2 + v**2]])
     >>> symsystem3.dyn_implicit_mat
     Matrix([
-    [1, 0, -x(t)/m],
-    [0, 1, -y(t)/m],
-    [0, 0,  l**2/m]])
+    [1, 0,   -x/m],
+    [0, 1,   -y/m],
+    [0, 0, l**2/m]])
     >>> symsystem3.kin_explicit_rhs
     Matrix([
-    [u(t)],
-    [v(t)]])
+    [u],
+    [v]])
     >>> symsystem1.alg_con
     [4]
     >>> symsystem1.bodies

doc/src/modules/physics/vector/kinematics.rst

@@ -430,6 +430,7 @@ velocity definition. ::
 
   >>> from sympy import Symbol, sin, cos
   >>> from sympy.physics.vector import *
+  >>> init_vprinting(pretty_print=False)
   >>> N = ReferenceFrame('N')
   >>> q1 = dynamicsymbols('q1')
   >>> A = N.orientnew('A', 'Axis', [q1, N.x])

doc/src/modules/polys/index.rst

@@ -27,3 +27,4 @@ Contents
     internals.rst
     ringseries.rst
     literature.rst
+    solvers.rst

doc/src/modules/polys/solvers.rst

@@ -0,0 +1,31 @@
+.. _poly_solvers-docs:
+
+Poly solvers
+============
+
+This module provides functions for solving systems of linear equations that
+are used internally in sympy.
+
+.. module::sympy.polys.solvers
+
+.. automodule:: sympy.polys.solvers
+
+solve_lin_sys
+^^^^^^^^^^^^^
+.. autofunction:: sympy.polys.solvers.solve_lin_sys
+
+eqs_to_matrix
+^^^^^^^^^^^^^
+.. autofunction:: sympy.polys.solvers.eqs_to_matrix
+
+sympy_eqs_to_ring
+^^^^^^^^^^^^^^^^^
+.. autofunction:: sympy.polys.solvers.sympy_eqs_to_ring
+
+_solve_lin_sys
+^^^^^^^^^^^^^^
+.. autofunction:: sympy.polys.solvers._solve_lin_sys
+
+_solve_lin_sys_component
+^^^^^^^^^^^^^^^^^^^^^^^^
+.. autofunction:: sympy.polys.solvers._solve_lin_sys_component

sympy/calculus/finite_diff.py

@@ -138,7 +138,7 @@ def finite_diff_weights(order, x_list, x0=S.One):
     >>> [simplify(c) for c in  mycoeffs] #doctest: +NORMALIZE_WHITESPACE
     [(h**3/2 + h**2*x - 3*h*x**2 - 4*x**3)/h**4,
     (-sqrt(2)*h**3 - 4*h**2*x + 3*sqrt(2)*h*x**2 + 8*x**3)/h**4,
-    6*x/h**2 - 8*x**3/h**4,
+    (6*h**2*x - 8*x**3)/h**4,
     (sqrt(2)*h**3 - 4*h**2*x - 3*sqrt(2)*h*x**2 + 8*x**3)/h**4,
     (-h**3/2 + h**2*x + 3*h*x**2 - 4*x**3)/h**4]
 

sympy/core/add.py

@@ -916,7 +916,13 @@ def _eval_as_leading_term(self, x, cdir=0):
         if not new_expr:
             # simple leading term analysis gave us cancelled terms but we have to send
             # back a term, so compute the leading term (via series)
-            return old.compute_leading_term(x)
+            n0 = min.getn()
+            res = Order(1)
+            incr = S.One
+            while res.is_Order:
+                res = old._eval_nseries(x, n=n0+incr, logx=None, cdir=cdir).cancel().powsimp().trigsimp()
+                incr *= 2
+            return res.as_leading_term(x, cdir=cdir)
 
         elif new_expr is S.NaN:
             return old.func._from_args(infinite)

sympy/core/mul.py

@@ -1766,7 +1766,7 @@ def ndiv(a, b):
         return co_residual*self2.func(*margs)*self2.func(*nc)
 
     def _eval_nseries(self, x, n, logx, cdir=0):
-        from sympy import Mul, Order, ceiling, powsimp
+        from sympy import degree, Mul, Order, ceiling, powsimp, PolynomialError
         from itertools import product
 
         def coeff_exp(term, x):
@@ -1821,8 +1821,19 @@ def coeff_exp(term, x):
             if power < n:
                 res += Mul(*coeffs)*(x**power)
 
-        if (res - self).expand() is not S.Zero:
-            res += Order(x**n, x)
+        if self.is_polynomial(x):
+            try:
+                if degree(self, x) != degree(res, x):
+                    res += Order(x**n, x)
+            except PolynomialError:
+                pass
+            else:
+                return res
+
+        for i in (1, 2, 3):
+            if (res - self).subs(x, i) is not S.Zero:
+                res += Order(x**n, x)
+                break
         return res
 
     def _eval_as_leading_term(self, x, cdir=0):

sympy/core/power.py

@@ -1589,10 +1589,11 @@ def mul(d1, d2):
             ex = e1 + inex
             res += terms[e1]*inco*x**(ex)
 
-        if (res - self).cancel() == S.Zero:
-            return res
-
-        return res + O(x**n, x)
+        for i in (1, 2, 3):
+            if (res - self).subs(x, i) is not S.Zero:
+                res += O(x**n, x)
+                break
+        return res
 
     def _eval_as_leading_term(self, x, cdir=0):
         from sympy import exp, I, im, log

sympy/core/tests/test_args.py

@@ -1076,9 +1076,9 @@ def test_sympy__stats__joint_rv__JointRandomSymbol():
     from sympy.stats.joint_rv import JointRandomSymbol
     assert _test_args(JointRandomSymbol(x))
 
-def test_sympy__stats__joint_rv__JointDistributionHandmade():
+def test_sympy__stats__joint_rv_types__JointDistributionHandmade():
     from sympy import Indexed
-    from sympy.stats.joint_rv import JointDistributionHandmade
+    from sympy.stats.joint_rv_types import JointDistributionHandmade
     x1, x2 = (Indexed('x', i) for i in (1, 2))
     assert _test_args(JointDistributionHandmade(x1 + x2, S.Reals**2))
 

sympy/core/tests/test_subs.py

@@ -855,3 +855,10 @@ def test_issue_17823():
 def test_issue_19326():
     x, y = [i(t) for i in map(Function, 'xy')]
     assert (x*y).subs({x: 1 + x, y: x}) == (1 + x)*x
+
+def test_issue_19558():
+    e = (7*x*cos(x) - 12*log(x)**3)*(-log(x)**4 + 2*sin(x) + 1)**2/ \
+    (2*(x*cos(x) - 2*log(x)**3)*(3*log(x)**4 - 7*sin(x) + 3)**2)
+
+    assert e.subs(x, oo) == AccumBounds(-oo, oo)
+    assert (sin(x) + cos(x)).subs(x, oo) == AccumBounds(-2, 2)

sympy/functions/elementary/complexes.py

@@ -394,6 +394,9 @@ def _eval_rewrite_as_Heaviside(self, arg, **kwargs):
         if arg.is_extended_real:
             return Heaviside(arg, H0=S(1)/2) * 2 - 1
 
+    def _eval_rewrite_as_Abs(self, arg, **kwargs):
+        return Piecewise((0, Eq(arg, 0)), (arg / Abs(arg), True))
+
     def _eval_simplify(self, **kwargs):
         return self.func(factor_terms(self.args[0]))  # XXX include doit?
 

sympy/functions/elementary/exponential.py

@@ -700,6 +700,7 @@ def eval(cls, arg, base=None):
                     from sympy.simplify import ratsimp
                     # Check for arguments involving rational multiples of pi
                     t = (i_/r_).cancel()
+                    t1 = (-t).cancel()
                     atan_table = {
                         # first quadrant only
                         sqrt(3): S.Pi/3,
@@ -728,12 +729,12 @@ def eval(cls, arg, base=None):
                             return cls(modulus) + I * atan_table[t]
                         else:
                             return cls(modulus) + I * (atan_table[t] - S.Pi)
-                    elif -t in atan_table:
+                    elif t1 in atan_table:
                         modulus = ratsimp(coeff * Abs(arg_))
                         if r_.is_positive:
-                            return cls(modulus) + I * (-atan_table[-t])
+                            return cls(modulus) + I * (-atan_table[t1])
                         else:
-                            return cls(modulus) + I * (S.Pi - atan_table[-t])
+                            return cls(modulus) + I * (S.Pi - atan_table[t1])
 
     def as_base_exp(self):
         """

sympy/functions/elementary/tests/test_complexes.py

@@ -4,7 +4,7 @@
     pi, Rational, re, S, sign, sin, sqrt, Symbol, symbols, transpose,
     zoo, exp_polar, Piecewise, Interval, comp, Integral, Matrix,
     ImmutableMatrix, SparseMatrix, ImmutableSparseMatrix, MatrixSymbol,
-    FunctionMatrix, Lambda, Derivative)
+    FunctionMatrix, Lambda, Derivative, Eq)
 from sympy.core.expr import unchanged
 from sympy.core.function import ArgumentIndexError
 from sympy.testing.pytest import XFAIL, raises
@@ -296,11 +296,14 @@ def test_sign():
     assert sign(Symbol('x', real=True, zero=False)).is_nonpositive is None
 
     x, y = Symbol('x', real=True), Symbol('y')
+    f = Function('f')
     assert sign(x).rewrite(Piecewise) == \
         Piecewise((1, x > 0), (-1, x < 0), (0, True))
     assert sign(y).rewrite(Piecewise) == sign(y)
     assert sign(x).rewrite(Heaviside) == 2*Heaviside(x, H0=S(1)/2) - 1
     assert sign(y).rewrite(Heaviside) == sign(y)
+    assert sign(y).rewrite(Abs) == Piecewise((0, Eq(y, 0)), (y/Abs(y), True))
+    assert sign(f(y)).rewrite(Abs) == Piecewise((0, Eq(f(y), 0)), (f(y)/Abs(f(y)), True))
 
     # evaluate what can be evaluated
     assert sign(exp_polar(I*pi)*pi) is S.NegativeOne