Merge branch 'master' into release

Author: isuruf

symengine/__init__.py

@@ -55,6 +55,11 @@
 
     def lambdify(args, exprs, **kwargs):
         return Lambdify(args, *exprs, **kwargs)
+else:
+    def __getattr__(name):
+        if name == 'lambdify':
+            raise AttributeError("Cannot import numpy, which is required for `lambdify` to work")
+        raise AttributeError(f"module 'symengine' has no attribute '{name}'")
 
 
 __version__ = "0.10.0"

symengine/lib/pywrapper.cpp

@@ -281,47 +281,61 @@ inline PyObject* get_pickle_module() {
     return module;
 }
 
+PyObject* pickle_loads(const std::string &pickle_str) {
+    PyObject *module = get_pickle_module();
+    PyObject *pickle_bytes = PyBytes_FromStringAndSize(pickle_str.data(), pickle_str.size());
+    PyObject *obj = PyObject_CallMethod(module, "loads", "O", pickle_bytes);
+    Py_XDECREF(pickle_bytes);
+    if (obj == NULL) {
+        throw SerializationError("error when loading pickled symbol subclass object");
+    }
+    return obj;
+}
+
 RCP<const Basic> load_basic(cereal::PortableBinaryInputArchive &ar, RCP<const Symbol> &)
 {
     bool is_pysymbol;
+    bool store_pickle;
     std::string name;
     ar(is_pysymbol);
     ar(name);
     if (is_pysymbol) {
         std::string pickle_str;
         ar(pickle_str);
-        PyObject *module = get_pickle_module();
-        PyObject *pickle_bytes = PyBytes_FromStringAndSize(pickle_str.data(), pickle_str.size());
-        PyObject *obj = PyObject_CallMethod(module, "loads", "O", pickle_bytes);
-        if (obj == NULL) {
-            throw SerializationError("error when loading pickled symbol subclass object");
-        }
-        RCP<const Basic> result = make_rcp<PySymbol>(name, obj);
-        Py_XDECREF(pickle_bytes);
+        ar(store_pickle);
+        PyObject *obj = pickle_loads(pickle_str);
+        RCP<const Basic> result = make_rcp<PySymbol>(name, obj, store_pickle);
+        Py_XDECREF(obj);
         return result;
     } else {
         return symbol(name);
     }
 }
 
+std::string pickle_dumps(const PyObject * obj) {
+    PyObject *module = get_pickle_module();
+    PyObject *pickle_bytes = PyObject_CallMethod(module, "dumps", "O", obj);
+    if (pickle_bytes == NULL) {
+        throw SerializationError("error when pickling symbol subclass object");
+    }
+    Py_ssize_t size;
+    char* buffer;
+    PyBytes_AsStringAndSize(pickle_bytes, &buffer, &size);
+    return std::string(buffer, size);
+}
+
 void save_basic(cereal::PortableBinaryOutputArchive &ar, const Symbol &b)
 {
     bool is_pysymbol = is_a_sub<PySymbol>(b);
     ar(is_pysymbol);
     ar(b.__str__());
     if (is_pysymbol) {
         RCP<const PySymbol> p = rcp_static_cast<const PySymbol>(b.rcp_from_this());
-        PyObject *module = get_pickle_module();
-        PyObject *pickle_bytes = PyObject_CallMethod(module, "dumps", "O", p->get_py_object());
-        if (pickle_bytes == NULL) {
-            throw SerializationError("error when pickling symbol subclass object");
-        }
-        Py_ssize_t size;
-        char* buffer;
-        PyBytes_AsStringAndSize(pickle_bytes, &buffer, &size);
-        std::string pickle_str(buffer, size);
+        PyObject *obj = p->get_py_object();
+        std::string pickle_str = pickle_dumps(obj);
         ar(pickle_str);
-        Py_XDECREF(pickle_bytes);
+        ar(p->store_pickle);
+        Py_XDECREF(obj);
     }
 }
 

symengine/lib/pywrapper.h

@@ -8,6 +8,9 @@
 
 namespace SymEngine {
 
+std::string pickle_dumps(const PyObject *);
+PyObject* pickle_loads(const std::string &);
+
 /*
  * PySymbol is a subclass of Symbol that keeps a reference to a Python object.
  * When subclassing a Symbol from Python, the information stored in subclassed
@@ -27,16 +30,30 @@ namespace SymEngine {
 class PySymbol : public Symbol {
 private:
     PyObject* obj;
+    std::string bytes;
 public:
-    PySymbol(const std::string& name, PyObject* obj) : Symbol(name), obj(obj) {
-        Py_INCREF(obj);
+    const bool store_pickle;
+    PySymbol(const std::string& name, PyObject* obj, bool store_pickle) :
+            Symbol(name), obj(obj), store_pickle(store_pickle) {
+        if (store_pickle) {
+            bytes = pickle_dumps(obj);
+        } else {
+            Py_INCREF(obj);
+        }
     }
     PyObject* get_py_object() const {
-        return obj;
+        if (store_pickle) {
+            return pickle_loads(bytes);
+        } else {
+            Py_INCREF(obj);
+            return obj;
+        }
     }
     virtual ~PySymbol() {
-        // TODO: This is never called because of the cyclic reference.
-        Py_DECREF(obj);
+        if (not store_pickle) {
+            // TODO: This is never called because of the cyclic reference.
+            Py_DECREF(obj);
+        }
     }
 };
 

symengine/lib/symengine.pxd

@@ -197,7 +197,7 @@ cdef extern from "<symengine/basic.h>" namespace "SymEngine":
     bool neq(const Basic &a, const Basic &b) nogil except +
 
     RCP[const Symbol] rcp_static_cast_Symbol "SymEngine::rcp_static_cast<const SymEngine::Symbol>"(rcp_const_basic &b) nogil
-    RCP[const PySymbol] rcp_static_cast_PySymbol "SymEngine::rcp_static_cast<const SymEngine::PySymbol>"(rcp_const_basic &b) nogil
+    RCP[const PySymbol] rcp_static_cast_PySymbol "SymEngine::rcp_static_cast<const SymEngine::PySymbol>"(rcp_const_basic &b) nogil except +
     RCP[const Integer] rcp_static_cast_Integer "SymEngine::rcp_static_cast<const SymEngine::Integer>"(rcp_const_basic &b) nogil
     RCP[const Rational] rcp_static_cast_Rational "SymEngine::rcp_static_cast<const SymEngine::Rational>"(rcp_const_basic &b) nogil
     RCP[const Complex] rcp_static_cast_Complex "SymEngine::rcp_static_cast<const SymEngine::Complex>"(rcp_const_basic &b) nogil
@@ -369,8 +369,8 @@ cdef extern from "pywrapper.h" namespace "SymEngine":
 
 cdef extern from "pywrapper.h" namespace "SymEngine":
     cdef cppclass PySymbol(Symbol):
-        PySymbol(string name, PyObject* pyobj)
-        PyObject* get_py_object()
+        PySymbol(string name, PyObject* pyobj, bool use_pickle) except +
+        PyObject* get_py_object() except +
 
     string wrapper_dumps(const Basic &x) nogil except +
     rcp_const_basic wrapper_loads(const string &s) nogil except +
@@ -479,7 +479,7 @@ cdef extern from "<symengine/basic.h>" namespace "SymEngine":
     rcp_const_basic make_rcp_Symbol "SymEngine::make_rcp<const SymEngine::Symbol>"(string name) nogil
     rcp_const_basic make_rcp_Dummy "SymEngine::make_rcp<const SymEngine::Dummy>"() nogil
     rcp_const_basic make_rcp_Dummy "SymEngine::make_rcp<const SymEngine::Dummy>"(string name) nogil
-    rcp_const_basic make_rcp_PySymbol "SymEngine::make_rcp<const SymEngine::PySymbol>"(string name, PyObject * pyobj) nogil
+    rcp_const_basic make_rcp_PySymbol "SymEngine::make_rcp<const SymEngine::PySymbol>"(string name, PyObject * pyobj, bool use_pickle) except +
     rcp_const_basic make_rcp_Constant "SymEngine::make_rcp<const SymEngine::Constant>"(string name) nogil
     rcp_const_basic make_rcp_Infty "SymEngine::make_rcp<const SymEngine::Infty>"(RCP[const Number] i) nogil
     rcp_const_basic make_rcp_NaN "SymEngine::make_rcp<const SymEngine::NaN>"() nogil

symengine/lib/symengine_wrapper.pyx

@@ -46,6 +46,7 @@ cpdef void assign_to_capsule(object capsule, object value):
 
 cdef object c2py(rcp_const_basic o):
     cdef Basic r
+    cdef PyObject *obj
     if (symengine.is_a_Add(deref(o))):
         r = Expr.__new__(Add)
     elif (symengine.is_a_Mul(deref(o))):
@@ -74,7 +75,10 @@ cdef object c2py(rcp_const_basic o):
         r = Dummy.__new__(Dummy)
     elif (symengine.is_a_Symbol(deref(o))):
         if (symengine.is_a_PySymbol(deref(o))):
-            return <object>(deref(symengine.rcp_static_cast_PySymbol(o)).get_py_object())
+            obj = deref(symengine.rcp_static_cast_PySymbol(o)).get_py_object()
+            result = <object>(obj)
+            Py_XDECREF(obj);
+            return result
         r = Symbol.__new__(Symbol)
     elif (symengine.is_a_Constant(deref(o))):
         r = S.Pi
@@ -1216,16 +1220,26 @@ cdef class Expr(Basic):
 
 
 cdef class Symbol(Expr):
-
     """
     Symbol is a class to store a symbolic variable with a given name.
+    Subclassing Symbol leads to a memory leak due to a cycle in reference counting.
+    To avoid this with a performance penalty, set the kwarg store_pickle=True
+    in the constructor and support the pickle protocol in the subclass by
+    implmenting __reduce__.
     """
 
     def __init__(Basic self, name, *args, **kwargs):
+        cdef cppbool store_pickle;
         if type(self) == Symbol:
             self.thisptr = symengine.make_rcp_Symbol(name.encode("utf-8"))
         else:
-            self.thisptr = symengine.make_rcp_PySymbol(name.encode("utf-8"), <PyObject*>self)
+            store_pickle = kwargs.pop("store_pickle", False)
+            if store_pickle:
+                # First set the pointer to a regular symbol so that when pickle.dumps
+                # is called when the PySymbol is created, methods like name works.
+                self.thisptr = symengine.make_rcp_Symbol(name.encode("utf-8"))
+            self.thisptr = symengine.make_rcp_PySymbol(name.encode("utf-8"), <PyObject*>self,
+                store_pickle)
 
     def _sympy_(self):
         import sympy
@@ -2635,6 +2649,14 @@ class atan2(Function):
         cdef Basic Y = sympify(y)
         return c2py(symengine.atan2(X.thisptr, Y.thisptr))
 
+    def _sympy_(self):
+        import sympy
+        return sympy.atan2(*self.args_as_sympy())
+
+    def _sage_(self):
+        import sage.all as sage
+        return sage.atan2(*self.args_as_sage())
+
 # For backwards compatibility
 
 Sin = sin
@@ -3895,7 +3917,7 @@ cdef class DenseMatrixBase(MatrixBase):
                 l.append(c2py(A.get(i, j))._sympy_())
             s.append(l)
         import sympy
-        return sympy.Matrix(s)
+        return sympy.ImmutableMatrix(s)
 
     def _sage_(self):
         s = []
@@ -3906,7 +3928,7 @@ cdef class DenseMatrixBase(MatrixBase):
                 l.append(c2py(A.get(i, j))._sage_())
             s.append(l)
         import sage.all as sage
-        return sage.Matrix(s)
+        return sage.Matrix(s, immutable=True)
 
     def dump_real(self, double[::1] out):
         cdef size_t ri, ci, nr, nc
@@ -4046,6 +4068,12 @@ cdef class ImmutableDenseMatrix(DenseMatrixBase):
     def __setitem__(self, key, value):
         raise TypeError("Cannot set values of {}".format(self.__class__))
 
+    def _applyfunc(self, f):
+        res = DenseMatrix(self)
+        res._applyfunc(f)
+        return ImmutableDenseMatrix(res)
+
+
 ImmutableMatrix = ImmutableDenseMatrix
 
 
@@ -5203,7 +5231,7 @@ def Lambdify(args, *exprs, cppbool real=True, backend=None, order='C',
         Whether datatype is ``double`` (``double complex`` otherwise).
     backend : str
         'llvm' or 'lambda'. When ``None`` the environment variable
-        'SYMENGINE_LAMBDIFY_BACKEND' is used (taken as 'lambda' if unset).
+        'SYMENGINE_LAMBDIFY_BACKEND' is used (taken as 'llvm' if available, otherwise 'lambda').
     order : 'C' or 'F'
         C- or Fortran-contiguous memory layout. Note that this affects
         broadcasting: e.g. a (m, n) matrix taking 3 arguments and given a
@@ -5235,7 +5263,11 @@ def Lambdify(args, *exprs, cppbool real=True, backend=None, order='C',
 
     """
     if backend is None:
-        backend = os.getenv('SYMENGINE_LAMBDIFY_BACKEND', "lambda")
+        IF HAVE_SYMENGINE_LLVM:
+            backend_default = 'llvm' if real else 'lambda'
+        ELSE:
+            backend_default = 'lambda'
+        backend = os.getenv('SYMENGINE_LAMBDIFY_BACKEND', backend_default)
     if backend == "llvm":
         IF HAVE_SYMENGINE_LLVM:
             if dtype == None:

symengine/tests/test_matrices.py

@@ -3,13 +3,23 @@
     Rational, function_symbol, I, NonSquareMatrixError, ShapeError, zeros,
     ones, eye, ImmutableMatrix)
 from symengine.test_utilities import raises
+import unittest
 
 
 try:
     import numpy as np
-    HAVE_NUMPY = True
+    have_numpy = True
 except ImportError:
-    HAVE_NUMPY = False
+    have_numpy = False
+    
+try:
+    import sympy
+    from sympy.core.cache import clear_cache
+    import atexit
+    atexit.register(clear_cache)
+    have_sympy = True
+except ImportError:
+    have_sympy = False
 
 
 def test_init():
@@ -520,21 +530,18 @@ def test_reshape():
     assert C != A
 
 
-# @pytest.mark.skipif(not HAVE_NUMPY, reason='requires numpy')
+@unittest.skipIf(not have_numpy, 'requires numpy')
 def test_dump_real():
-    if not HAVE_NUMPY:  # nosetests work-around
-        return
     ref = [1, 2, 3, 4]
     A = DenseMatrix(2, 2, ref)
     out = np.empty(4)
     A.dump_real(out)
     assert np.allclose(out, ref)
 
 
-# @pytest.mark.skipif(not HAVE_NUMPY, reason='requires numpy')
+
+@unittest.skipIf(not have_numpy, 'requires numpy')
 def test_dump_complex():
-    if not HAVE_NUMPY:  # nosetests work-around
-        return
     ref = [1j, 2j, 3j, 4j]
     A = DenseMatrix(2, 2, ref)
     out = np.empty(4, dtype=np.complex128)
@@ -741,3 +748,8 @@ def test_repr_latex():
     latex_string = testmat._repr_latex_()
     assert isinstance(latex_string, str)
     init_printing(False)
+
+@unittest.skipIf(not have_sympy, "SymPy not installed")
+def test_simplify():
+    A = ImmutableMatrix([1])
+    assert type(A.simplify()) == type(A)

symengine/tests/test_sympy_conv.py

@@ -2,7 +2,7 @@
         function_symbol, I, E, pi, oo, zoo, nan, true, false,
         exp, gamma, have_mpfr, have_mpc, DenseMatrix, sin, cos, tan, cot,
         csc, sec, asin, acos, atan, acot, acsc, asec, sinh, cosh, tanh, coth,
-        asinh, acosh, atanh, acoth, Add, Mul, Pow, diff, GoldenRatio,
+        asinh, acosh, atanh, acoth, atan2, Add, Mul, Pow, diff, GoldenRatio,
         Catalan, EulerGamma, UnevaluatedExpr, RealDouble)
 from symengine.lib.symengine_wrapper import (Subs, Derivative, RealMPFR,
         ComplexMPC, PyNumber, Function, LambertW, zeta, dirichlet_eta,
@@ -171,6 +171,7 @@ def test_conv7():
     assert acot(x/3) == acot(sympy.Symbol("x") / 3)
     assert acsc(x/3) == acsc(sympy.Symbol("x") / 3)
     assert asec(x/3) == asec(sympy.Symbol("x") / 3)
+    assert atan2(x/3, y) == atan2(sympy.Symbol("x") / 3, sympy.Symbol("y"))
 
     assert sin(x/3)._sympy_() == sympy.sin(sympy.Symbol("x") / 3)
     assert sin(x/3)._sympy_() != sympy.cos(sympy.Symbol("x") / 3)
@@ -185,6 +186,22 @@ def test_conv7():
     assert acot(x/3)._sympy_() == sympy.acot(sympy.Symbol("x") / 3)
     assert acsc(x/3)._sympy_() == sympy.acsc(sympy.Symbol("x") / 3)
     assert asec(x/3)._sympy_() == sympy.asec(sympy.Symbol("x") / 3)
+    assert atan2(x/3, y)._sympy_() == sympy.atan2(sympy.Symbol("x") / 3, sympy.Symbol("y"))
+
+    assert sympy.sympify(sin(x/3)) == sympy.sin(sympy.Symbol("x") / 3)
+    assert sympy.sympify(sin(x/3)) != sympy.cos(sympy.Symbol("x") / 3)
+    assert sympy.sympify(cos(x/3)) == sympy.cos(sympy.Symbol("x") / 3)
+    assert sympy.sympify(tan(x/3)) == sympy.tan(sympy.Symbol("x") / 3)
+    assert sympy.sympify(cot(x/3)) == sympy.cot(sympy.Symbol("x") / 3)
+    assert sympy.sympify(csc(x/3)) == sympy.csc(sympy.Symbol("x") / 3)
+    assert sympy.sympify(sec(x/3)) == sympy.sec(sympy.Symbol("x") / 3)
+    assert sympy.sympify(asin(x/3)) == sympy.asin(sympy.Symbol("x") / 3)
+    assert sympy.sympify(acos(x/3)) == sympy.acos(sympy.Symbol("x") / 3)
+    assert sympy.sympify(atan(x/3)) == sympy.atan(sympy.Symbol("x") / 3)
+    assert sympy.sympify(acot(x/3)) == sympy.acot(sympy.Symbol("x") / 3)
+    assert sympy.sympify(acsc(x/3)) == sympy.acsc(sympy.Symbol("x") / 3)
+    assert sympy.sympify(asec(x/3)) == sympy.asec(sympy.Symbol("x") / 3)
+    assert sympy.sympify(atan2(x/3, y)) == sympy.atan2(sympy.Symbol("x") / 3, sympy.Symbol("y"))
 
 
 @unittest.skipIf(not have_sympy, "SymPy not installed")
@@ -204,6 +221,7 @@ def test_conv7b():
     assert sympify(sympy.acot(x/3)) == acot(Symbol("x") / 3)
     assert sympify(sympy.acsc(x/3)) == acsc(Symbol("x") / 3)
     assert sympify(sympy.asec(x/3)) == asec(Symbol("x") / 3)
+    assert sympify(sympy.atan2(x/3, y)) == atan2(Symbol("x") / 3, Symbol("y"))
 
 
 @unittest.skipIf(not have_sympy, "SymPy not installed")