minimal set of changes for returning failure to find non-zero pivot in FFGJ

.github/workflows/ci.yml

@@ -132,11 +132,12 @@ jobs:
             CC: gcc
 
           - BUILD_TYPE: Release
-            PYTHON_VERSION: '3.8'
+            PYTHON_VERSION: '3.11'
             OS: ubuntu-20.04
             WITH_MPC: yes
             WITH_MPFR: yes
             WITH_FLINT: yes
+            WITH_FLINT_PY: yes
             WITH_SCIPY: yes
             WITH_DOCS: yes
             INTEGER_CLASS: flint

bin/install_travis.sh

@@ -16,6 +16,10 @@ if [[ "${WITH_DOCS}" == "yes" ]]; then
     export conda_pkgs="${conda_pkgs} sphinx recommonmark";
 fi
 
+if [[ "${WITH_FLINT_PY}" == "yes" ]]; then
+    export conda_pkgs="${conda_pkgs} python-flint";  # python-flint affects sympy, see e.g. sympy/sympy#26645
+fi
+
 if [[ "${WITH_SAGE}" == "yes" ]]; then
     # This is split to avoid the 10 minute limit
     conda install -q sagelib=8.1

symengine/__init__.py

@@ -52,14 +52,15 @@
     from .lib.symengine_wrapper import ComplexMPC
 
 if have_numpy:
-    from .lib.symengine_wrapper import (Lambdify, LambdifyCSE)
-
+    from .lib.symengine_wrapper import Lambdify, LambdifyCSE
     def lambdify(args, exprs, **kwargs):
+        if isinstance(exprs, Basic):
+            exprs = [exprs]
         return Lambdify(args, *exprs, **kwargs)
 else:
     def __getattr__(name):
-        if name == 'lambdify':
-            raise AttributeError("Cannot import numpy, which is required for `lambdify` to work")
+        if name in ('lambdify', 'Lambdify', 'LambdifyCSE'):
+            raise AttributeError(f"Cannot import numpy, which is required for `{name}` to work")
         raise AttributeError(f"module 'symengine' has no attribute '{name}'")
 
 

symengine/lib/symengine.pxd

@@ -687,7 +687,7 @@ cdef extern from "<symengine/matrix.h>" namespace "SymEngine":
         DenseMatrix &B) except+ nogil
     void FFLU_solve "SymEngine::fraction_free_LU_solve"(const DenseMatrix &A,
         const DenseMatrix &b, DenseMatrix &x) except+ nogil
-    void FFGJ_solve "SymEngine::fraction_free_gauss_jordan_solve"(const DenseMatrix &A,
+    int FFGJ_solve "SymEngine::fraction_free_gauss_jordan_solve"(const DenseMatrix &A,
         const DenseMatrix &b, DenseMatrix &x) except+ nogil
     void LDL_solve "SymEngine::LDL_solve"(const DenseMatrix &A, const DenseMatrix &b,
         DenseMatrix &x) except+ nogil

symengine/lib/symengine_wrapper.in.pxd

@@ -2,6 +2,7 @@
 
 cimport symengine
 from symengine cimport RCP, map_basic_basic, rcp_const_basic
+from libcpp.memory cimport unique_ptr
 from libcpp.vector cimport vector
 from libcpp.string cimport string
 from libcpp cimport bool as cppbool
@@ -44,7 +45,7 @@ cdef class _Lambdify(object):
     cpdef unsafe_eval(sef, inp, out, unsigned nbroadcast=*)
 
 cdef class LambdaDouble(_Lambdify):
-    cdef vector[symengine.LambdaRealDoubleVisitor] lambda_double
+    cdef unique_ptr[symengine.LambdaRealDoubleVisitor] lambda_visitor
     cdef _init(self, symengine.vec_basic& args_, symengine.vec_basic& outs_, cppbool cse)
     cpdef unsafe_real(self, double[::1] inp, double[::1] out, int inp_offset=*, int out_offset=*)
     cpdef as_scipy_low_level_callable(self)
@@ -54,7 +55,7 @@ cdef class LambdaDouble(_Lambdify):
                       int inp_offset=*, int out_offset=*)
 
 cdef class LambdaComplexDouble(_Lambdify):
-    cdef vector[symengine.LambdaComplexDoubleVisitor] lambda_double
+    cdef unique_ptr[symengine.LambdaComplexDoubleVisitor] lambda_visitor
     cdef _init(self, symengine.vec_basic& args_, symengine.vec_basic& outs_, cppbool cse)
     cpdef unsafe_complex(self, double complex[::1] inp, double complex[::1] out, int inp_offset=*, int out_offset=*)
 
@@ -63,22 +64,22 @@ IF HAVE_SYMENGINE_LLVM:
         cdef int opt_level
 
     cdef class LLVMDouble(_LLVMLambdify):
-        cdef vector[symengine.LLVMDoubleVisitor] lambda_double
+        cdef unique_ptr[symengine.LLVMDoubleVisitor] lambda_visitor
         cdef _init(self, symengine.vec_basic& args_, symengine.vec_basic& outs_, cppbool cse)
         cdef _load(self, const string &s)
         cpdef unsafe_real(self, double[::1] inp, double[::1] out, int inp_offset=*, int out_offset=*)
         cpdef as_scipy_low_level_callable(self)
         cpdef as_ctypes(self)
 
     cdef class LLVMFloat(_LLVMLambdify):
-        cdef vector[symengine.LLVMFloatVisitor] lambda_double
+        cdef unique_ptr[symengine.LLVMFloatVisitor] lambda_visitor
         cdef _init(self, symengine.vec_basic& args_, symengine.vec_basic& outs_, cppbool cse)
         cdef _load(self, const string &s)
         cpdef unsafe_real(self, float[::1] inp, float[::1] out, int inp_offset=*, int out_offset=*)
 
     IF HAVE_SYMENGINE_LLVM_LONG_DOUBLE:
         cdef class LLVMLongDouble(_LLVMLambdify):
-            cdef vector[symengine.LLVMLongDoubleVisitor] lambda_double
+            cdef unique_ptr[symengine.LLVMLongDoubleVisitor] lambda_visitor
             cdef _init(self, symengine.vec_basic& args_, symengine.vec_basic& outs_, cppbool cse)
             cdef _load(self, const string &s)
             cpdef unsafe_real(self, long double[::1] inp, long double[::1] out, int inp_offset=*, int out_offset=*)

symengine/lib/symengine_wrapper.in.pyx

@@ -1211,7 +1211,7 @@ cdef class Basic(object):
         return int(float(self))
 
     def __long__(self):
-        return long(float(self))
+        return int(float(self))
 
     def __complex__(self):
         f = self.n(real=False)
@@ -1627,7 +1627,9 @@ class Equality(Relational):
     def is_Equality(self):
         return True
 
-    func = __class__
+    @property
+    def func(self):
+        return self.__class__
 
 
 Eq = Equality
@@ -1648,7 +1650,9 @@ class Unequality(Relational):
         s = self.args_as_sage()
         return sage.ne(*s)
 
-    func = __class__
+    @property
+    def func(self):
+        return self.__class__
 
 
 Ne = Unequality
@@ -3529,6 +3533,9 @@ cdef class DenseMatrixBase(MatrixBase):
     def __neg__(self):
         return self.mul_scalar(-1)
 
+    def __abs__(self):
+        return self.applyfunc(abs)
+
     def __getitem__(self, item):
         if isinstance(item, slice):
             if (self.ncols() == 0 or self.nrows() == 0):
@@ -3998,9 +4005,11 @@ cdef class DenseMatrixBase(MatrixBase):
                 deref(symengine.static_cast_DenseMatrix(b_.thisptr)),
                 deref(symengine.static_cast_DenseMatrix(x.thisptr)))
         elif method.upper() == 'FFGJ':
-            symengine.FFGJ_solve(deref(symengine.static_cast_DenseMatrix(self.thisptr)),
+            failure = symengine.FFGJ_solve(deref(symengine.static_cast_DenseMatrix(self.thisptr)),
                 deref(symengine.static_cast_DenseMatrix(b_.thisptr)),
                 deref(symengine.static_cast_DenseMatrix(x.thisptr)))
+            if (failure != 0):
+                raise Exception("Underdetermined system. Failed to find non-zero pivot in column: %d" % failure)
         else:
             raise Exception("Unsupported method.")
 
@@ -5152,27 +5161,27 @@ cdef class LambdaDouble(_Lambdify):
         pass
 
     cdef _init(self, symengine.vec_basic& args_, symengine.vec_basic& outs_, cppbool cse):
-        self.lambda_double.resize(1)
-        self.lambda_double[0].init(args_, outs_, cse)
+        self.lambda_visitor.reset(new symengine.LambdaRealDoubleVisitor())
+        deref(self.lambda_visitor).init(args_, outs_, cse)
 
     cpdef unsafe_real(self, double[::1] inp, double[::1] out, int inp_offset=0, int out_offset=0):
-        self.lambda_double[0].call(&out[out_offset], &inp[inp_offset])
+        deref(self.lambda_visitor).call(&out[out_offset], &inp[inp_offset])
 
     cpdef unsafe_eval(self, inp, out, unsigned nbroadcast=1):
         cdef double[::1] c_inp, c_out
         cdef unsigned idx
         c_inp = np.ascontiguousarray(inp.ravel(order=self.order), dtype=self.numpy_dtype)
         c_out = out
         for idx in range(nbroadcast):
-            self.lambda_double[0].call(&c_out[idx*self.tot_out_size], &c_inp[idx*self.args_size])
+            deref(self.lambda_visitor).call(&c_out[idx*self.tot_out_size], &c_inp[idx*self.args_size])
 
     cpdef as_scipy_low_level_callable(self):
         from ctypes import c_double, c_void_p, c_int, cast, POINTER, CFUNCTYPE
         if self.tot_out_size > 1:
             raise RuntimeError("SciPy LowLevelCallable supports only functions with 1 output")
         addr1 = cast(<size_t>&_scipy_callback_lambda_real,
                         CFUNCTYPE(c_double, c_int, POINTER(c_double), c_void_p))
-        addr2 = cast(<size_t>&self.lambda_double[0], c_void_p)
+        addr2 = cast(<size_t>self.lambda_visitor.get(), c_void_p)
         return create_low_level_callable(self, addr1, addr2)
 
     cpdef as_ctypes(self):
@@ -5187,7 +5196,7 @@ cdef class LambdaDouble(_Lambdify):
         from ctypes import c_double, c_void_p, c_int, cast, POINTER, CFUNCTYPE
         addr1 = cast(<size_t>&_ctypes_callback_lambda_real,
                         CFUNCTYPE(c_void_p, POINTER(c_double), POINTER(c_double), c_void_p))
-        addr2 = cast(<size_t>&self.lambda_double[0], c_void_p)
+        addr2 = cast(<size_t>self.lambda_visitor.get(), c_void_p)
         return addr1, addr2
 
 
@@ -5197,19 +5206,19 @@ cdef class LambdaComplexDouble(_Lambdify):
         pass
 
     cdef _init(self, symengine.vec_basic& args_, symengine.vec_basic& outs_, cppbool cse):
-        self.lambda_double.resize(1)
-        self.lambda_double[0].init(args_, outs_, cse)
+        self.lambda_visitor.reset(new symengine.LambdaComplexDoubleVisitor())
+        deref(self.lambda_visitor).init(args_, outs_, cse)
 
     cpdef unsafe_complex(self, double complex[::1] inp, double complex[::1] out, int inp_offset=0, int out_offset=0):
-        self.lambda_double[0].call(&out[out_offset], &inp[inp_offset])
+        deref(self.lambda_visitor).call(&out[out_offset], &inp[inp_offset])
 
     cpdef unsafe_eval(self, inp, out, unsigned nbroadcast=1):
         cdef double complex[::1] c_inp, c_out
         cdef unsigned idx
         c_inp = np.ascontiguousarray(inp.ravel(order=self.order), dtype=self.numpy_dtype)
         c_out = out
         for idx in range(nbroadcast):
-            self.lambda_double[0].call(&c_out[idx*self.tot_out_size], &c_inp[idx*self.args_size])
+            deref(self.lambda_visitor).call(&c_out[idx*self.tot_out_size], &c_inp[idx*self.args_size])
 
 
 IF HAVE_SYMENGINE_LLVM:
@@ -5218,31 +5227,31 @@ IF HAVE_SYMENGINE_LLVM:
             self.opt_level = opt_level
 
         cdef _init(self, symengine.vec_basic& args_, symengine.vec_basic& outs_, cppbool cse):
-            self.lambda_double.resize(1)
-            self.lambda_double[0].init(args_, outs_, cse, self.opt_level)
+            self.lambda_visitor.reset(new symengine.LLVMDoubleVisitor())
+            deref(self.lambda_visitor).init(args_, outs_, cse, self.opt_level)
 
         cdef _load(self, const string &s):
-            self.lambda_double.resize(1)
-            self.lambda_double[0].loads(s)
+            self.lambda_visitor.reset(new symengine.LLVMDoubleVisitor())
+            deref(self.lambda_visitor).loads(s)
 
         def __reduce__(self):
             """
             Interface for pickle. Note that the resulting object is platform dependent.
             """
-            cdef bytes s = self.lambda_double[0].dumps()
+            cdef bytes s = deref(self.lambda_visitor).dumps()
             return llvm_loading_func, (self.args_size, self.tot_out_size, self.out_shapes, self.real, \
                 self.n_exprs, self.order, self.accum_out_sizes, self.numpy_dtype, s)
 
         cpdef unsafe_real(self, double[::1] inp, double[::1] out, int inp_offset=0, int out_offset=0):
-            self.lambda_double[0].call(&out[out_offset], &inp[inp_offset])
+            deref(self.lambda_visitor).call(&out[out_offset], &inp[inp_offset])
 
         cpdef unsafe_eval(self, inp, out, unsigned nbroadcast=1):
             cdef double[::1] c_inp, c_out
             cdef unsigned idx
             c_inp = np.ascontiguousarray(inp.ravel(order=self.order), dtype=self.numpy_dtype)
             c_out = out
             for idx in range(nbroadcast):
-                self.lambda_double[0].call(&c_out[idx*self.tot_out_size], &c_inp[idx*self.args_size])
+                deref(self.lambda_visitor).call(&c_out[idx*self.tot_out_size], &c_inp[idx*self.args_size])
 
         cpdef as_scipy_low_level_callable(self):
             from ctypes import c_double, c_void_p, c_int, cast, POINTER, CFUNCTYPE
@@ -5252,7 +5261,7 @@ IF HAVE_SYMENGINE_LLVM:
                 raise RuntimeError("SciPy LowLevelCallable supports only functions with 1 output")
             addr1 = cast(<size_t>&_scipy_callback_llvm_real,
                             CFUNCTYPE(c_double, c_int, POINTER(c_double), c_void_p))
-            addr2 = cast(<size_t>&self.lambda_double[0], c_void_p)
+            addr2 = cast(<size_t>self.lambda_visitor.get(), c_void_p)
             return create_low_level_callable(self, addr1, addr2)
 
         cpdef as_ctypes(self):
@@ -5269,71 +5278,71 @@ IF HAVE_SYMENGINE_LLVM:
                 raise RuntimeError("Lambda function has to be real")
             addr1 = cast(<size_t>&_ctypes_callback_llvm_real,
                             CFUNCTYPE(c_void_p, POINTER(c_double), POINTER(c_double), c_void_p))
-            addr2 = cast(<size_t>&self.lambda_double[0], c_void_p)
+            addr2 = cast(<size_t>self.lambda_visitor.get(), c_void_p)
             return addr1, addr2
 
     cdef class LLVMFloat(_LLVMLambdify):
         def __cinit__(self, args, *exprs, cppbool real=True, order='C', cppbool cse=False, cppbool _load=False, opt_level=3, dtype=None):
             self.opt_level = opt_level
 
         cdef _init(self, symengine.vec_basic& args_, symengine.vec_basic& outs_, cppbool cse):
-            self.lambda_double.resize(1)
-            self.lambda_double[0].init(args_, outs_, cse, self.opt_level)
+            self.lambda_visitor.reset(new symengine.LLVMFloatVisitor())
+            deref(self.lambda_visitor).init(args_, outs_, cse, self.opt_level)
 
         cdef _load(self, const string &s):
-            self.lambda_double.resize(1)
-            self.lambda_double[0].loads(s)
+            self.lambda_visitor.reset(new symengine.LLVMFloatVisitor())
+            deref(self.lambda_visitor).loads(s)
 
         def __reduce__(self):
             """
             Interface for pickle. Note that the resulting object is platform dependent.
             """
-            cdef bytes s = self.lambda_double[0].dumps()
+            cdef bytes s = deref(self.lambda_visitor).dumps()
             return llvm_float_loading_func, (self.args_size, self.tot_out_size, self.out_shapes, self.real, \
                 self.n_exprs, self.order, self.accum_out_sizes, self.numpy_dtype, s)
 
         cpdef unsafe_real(self, float[::1] inp, float[::1] out, int inp_offset=0, int out_offset=0):
-            self.lambda_double[0].call(&out[out_offset], &inp[inp_offset])
+            deref(self.lambda_visitor).call(&out[out_offset], &inp[inp_offset])
 
         cpdef unsafe_eval(self, inp, out, unsigned nbroadcast=1):
             cdef float[::1] c_inp, c_out
             cdef unsigned idx
             c_inp = np.ascontiguousarray(inp.ravel(order=self.order), dtype=self.numpy_dtype)
             c_out = out
             for idx in range(nbroadcast):
-                self.lambda_double[0].call(&c_out[idx*self.tot_out_size], &c_inp[idx*self.args_size])
+                deref(self.lambda_visitor).call(&c_out[idx*self.tot_out_size], &c_inp[idx*self.args_size])
 
     IF HAVE_SYMENGINE_LLVM_LONG_DOUBLE:
         cdef class LLVMLongDouble(_LLVMLambdify):
             def __cinit__(self, args, *exprs, cppbool real=True, order='C', cppbool cse=False, cppbool _load=False, opt_level=3, dtype=None):
                 self.opt_level = opt_level
 
             cdef _init(self, symengine.vec_basic& args_, symengine.vec_basic& outs_, cppbool cse):
-                self.lambda_double.resize(1)
-                self.lambda_double[0].init(args_, outs_, cse, self.opt_level)
+                self.lambda_visitor.reset(new symengine.LLVMLongDoubleVisitor())
+                deref(self.lambda_visitor).init(args_, outs_, cse, self.opt_level)
 
             cdef _load(self, const string &s):
-                self.lambda_double.resize(1)
-                self.lambda_double[0].loads(s)
+                self.lambda_visitor.reset(new symengine.LLVMLongDoubleVisitor())
+                deref(self.lambda_visitor).loads(s)
 
             def __reduce__(self):
                 """
                 Interface for pickle. Note that the resulting object is platform dependent.
                 """
-                cdef bytes s = self.lambda_double[0].dumps()
+                cdef bytes s = deref(self.lambda_visitor).dumps()
                 return llvm_long_double_loading_func, (self.args_size, self.tot_out_size, self.out_shapes, self.real, \
                     self.n_exprs, self.order, self.accum_out_sizes, self.numpy_dtype, s)
 
             cpdef unsafe_real(self, long double[::1] inp, long double[::1] out, int inp_offset=0, int out_offset=0):
-                self.lambda_double[0].call(&out[out_offset], &inp[inp_offset])
+                deref(self.lambda_visitor).call(&out[out_offset], &inp[inp_offset])
 
             cpdef unsafe_eval(self, inp, out, unsigned nbroadcast=1):
                 cdef long double[::1] c_inp, c_out
                 cdef unsigned idx
                 c_inp = np.ascontiguousarray(inp.ravel(order=self.order), dtype=self.numpy_dtype)
                 c_out = out
                 for idx in range(nbroadcast):
-                    self.lambda_double[0].call(&c_out[idx*self.tot_out_size], &c_inp[idx*self.args_size])
+                    deref(self.lambda_visitor).call(&c_out[idx*self.tot_out_size], &c_inp[idx*self.args_size])
 
     def llvm_loading_func(*args):
         return LLVMDouble(args, _load=True)

symengine/tests/test_logic.py

@@ -27,6 +27,10 @@ def test_relationals():
     assert Ge(1, 1) == true
     assert Eq(I, 2) == false
     assert Ne(I, 2) == true
+    eq = Eq(x, y)
+    assert eq.func(*eq.args) == eq
+    ne = Ne(x, y)
+    assert ne.func(*ne.args) == ne
 
 
 def test_rich_cmp():
@@ -118,4 +122,3 @@ def test_Contains():
     assert Contains(x, Interval(1, 1)) != false
     assert Contains(oo, Interval(-oo, oo)) == false
     assert Contains(-oo, Interval(-oo, oo)) == false
-