Repository: pyccel/pyccel
Pull request #235 (2020-02-17 13:41:06)

Added matmul and some tests
Here, `numpy.prod` and `numpy.matmul` are added together with according unit tests. Also, `.gitignore` is modified to ignore `__pyccel__` and `__epyccel__` entries and some linter stuff.


Patch:
diff --git a/.gitignore b/.gitignore
index 81623916..afeaa7a8 100644
--- a/.gitignore
+++ b/.gitignore
@@ -12,7 +12,10 @@
 *.pyo
 .ipynb_checkpoints/
 .cache
+.mypy_cache
 __pycache__
+__pyccel__
+__epyccel__
 
 build
 *egg*
diff --git a/pyccel/ast/datatypes.py b/pyccel/ast/datatypes.py
index 907adc3c..e2a990a0 100644
--- a/pyccel/ast/datatypes.py
+++ b/pyccel/ast/datatypes.py
@@ -43,6 +43,7 @@ class NativeBool(DataType):
 
 class NativeInteger(DataType):
     _name = 'Int'
+    pass
 
 class NativeReal(DataType):
     _name = 'Real'
diff --git a/pyccel/ast/numpyext.py b/pyccel/ast/numpyext.py
index 44d5eb57..f7b4a345 100644
--- a/pyccel/ast/numpyext.py
+++ b/pyccel/ast/numpyext.py
@@ -126,7 +126,7 @@ class Sum(Function):
     """
 
     def __new__(cls, arg):
-        if not isinstance(arg, (list, tuple, Tuple, List, Variable)):
+        if not isinstance(arg, (list, tuple, Tuple, List, Variable, Mul, Add, Pow, sp_Rational)):
             raise TypeError('Uknown type of  %s.' % type(arg))
         return Basic.__new__(cls, arg)
 
@@ -151,6 +151,92 @@ class Sum(Function):
             return '{0} = sum({1})'.format(lhs_code, rhs_code)
         return 'sum({0})'.format(rhs_code)
 
+#=======================================================================================
+
+class Product(Function):
+    """Represents a call to  numpy.prod for code generation.
+
+    arg : list , tuple , Tuple, List, Variable
+    """
+
+    def __new__(cls, arg):
+        if not isinstance(arg, (list, tuple, Tuple, List, Variable, Mul, Add, Pow, sp_Rational)):
+            raise TypeError('Uknown type of  %s.' % type(arg))
+        return Basic.__new__(cls, arg)
+
+    @property
+    def arg(self):
+        return self._args[0]
+
+    @property
+    def dtype(self):
+        return self._args[0].dtype
+
+    @property
+    def rank(self):
+        return 0
+
+    def fprint(self, printer, lhs=None):
+        """Fortran print."""
+
+        rhs_code = printer(self.arg)
+        if lhs:
+            lhs_code = printer(lhs)
+            return '{0} = product({1})'.format(lhs_code, rhs_code)
+        return 'product({0})'.format(rhs_code)
+
+#=======================================================================================
+
+class Matmul(Application):
+    """Represents a call to numpy.matmul for code generation.
+    arg : list , tuple , Tuple, List, Variable
+    """
+
+    def __new__(cls, a, b):
+        if not isinstance(a, (list, tuple, Tuple, List, Variable, Mul, Add, Pow, sp_Rational)):
+            raise TypeError('Uknown type of  %s.' % type(a))
+        if not isinstance(b, (list, tuple, Tuple, List, Variable, Mul, Add, Pow, sp_Rational)):
+            raise TypeError('Uknown type of  %s.' % type(a))
+        return Basic.__new__(cls, a, b)
+
+    @property
+    def a(self):
+        return self._args[0]
+
+    @property
+    def b(self):
+        return self._args[1]
+
+    @property
+    def dtype(self):
+        return self._args[0].dtype
+
+    @property
+    def rank(self):
+        return 1 # TODO: make this general
+
+    def fprint(self, printer, lhs=None):
+        """Fortran print."""
+        a_code = printer(self.a)
+        b_code = printer(self.b)
+
+        if lhs:
+            lhs_code = printer(lhs)
+
+        if self.a.order and self.b.order:
+            if self.a.order != self.b.order:
+                raise NotImplementedError("Mixed order matmul not supported.")
+
+        # Fortran ordering
+        if self.a.order == 'F':
+            if lhs:
+                return '{0} = matmul({1},{2})'.format(lhs_code, a_code, b_code)
+            return 'matmul({0},{1})'.format(a_code, b_code)
+
+        # C ordering
+        if lhs:
+            return '{0} = matmul({2},{1})'.format(lhs_code, a_code, b_code)
+        return 'matmul({1},{0})'.format(a_code, b_code)
 
 #=======================================================================================
 
diff --git a/pyccel/ast/utilities.py b/pyccel/ast/utilities.py
index 9e10f896..cbba1719 100644
--- a/pyccel/ast/utilities.py
+++ b/pyccel/ast/utilities.py
@@ -10,10 +10,11 @@ from .core import FunctionDef, Return, Assign
 from .core import Constant, Variable, IndexedVariable
 from .numpyext import Zeros, Ones, Empty, ZerosLike, FullLike, Diag, Cross
 from .numpyext import Min, Max, Abs, Norm, EmptyLike, Where
-from .numpyext import Array, Shape, Int, Rand, Sum, Real, Complex, Imag, Mod
+from .numpyext import Array, Shape, Int, Rand, Sum, Matmul, Real, Complex, Imag, Mod
 from .numpyext import Int64, Int32, Float32, Float64, Complex64, Complex128
 from .numpyext import Sqrt, Asin, Acsc, Acos, Asec, Atan, Acot, Sinh, Cosh, Tanh, Log
 from .numpyext import numpy_constants, Linspace
+from .numpyext import Product as Prod
 from pyccel.symbolic import lambdify
 from sympy import Symbol, Lambda, floor
 from sympy import Not, Float
@@ -73,7 +74,10 @@ numpy_functions = {
     'int64'     : Int64,
     'complex128': Complex128,
     'complex64' : Complex64,
+    'matmul'    : Matmul,
     'sum'       : Sum,
+    'prod'      : Prod,
+    'product'   : Prod,
     'rand'      : Rand,
     'random'    : Rand,
     'linspace'  : Linspace,
diff --git a/pyccel/codegen/printing/fcode.py b/pyccel/codegen/printing/fcode.py
index 37e949d5..36ae4c6c 100644
--- a/pyccel/codegen/printing/fcode.py
+++ b/pyccel/codegen/printing/fcode.py
@@ -623,6 +623,12 @@ class FCodePrinter(CodePrinter):
     def _print_Sum(self, expr):
         return expr.fprint(self._print)
 
+    def _print_Product(self, expr):
+        return expr.fprint(self._print)
+
+    def _print_Matmul(self, expr):
+        return expr.fprint(self._print)
+
     def _print_Cross(self, expr):
         return expr.fprint(self._print)
 
diff --git a/pyccel/parser/semantic.py b/pyccel/parser/semantic.py
index b2f62f37..f0db9686 100644
--- a/pyccel/parser/semantic.py
+++ b/pyccel/parser/semantic.py
@@ -825,12 +825,37 @@ class SemanticParser(BasicParser):
 
                 d_var['is_target'] = True # ISSUE 177: TODO this should be done using update_variable
                 
-            elif name in ['Len', 'Sum', 'Rand', 'Min', 'Max']:
+            elif name in ['Len', 'Sum', 'Product', 'Rand', 'Min', 'Max']:
                 d_var['datatype'   ] = sp_dtype(expr)
                 d_var['rank'       ] = 0
                 d_var['allocatable'] = False
                 d_var['is_pointer' ] = False
 
+            elif name in ['Matmul']:
+
+                d_vars = [self._infere_type(arg,**settings) for arg in expr.args]
+
+                var0_is_vector = d_vars[0]['rank'] < 2
+                var1_is_vector = d_vars[1]['rank'] < 2
+
+                if(d_vars[0]['shape'] is None or d_vars[1]['shape'] is None):
+                    d_var['shape'] = None
+                else:
+
+                    m = 1 if var0_is_vector else d_vars[0]['shape'][0]
+                    n = 1 if var1_is_vector else d_vars[1]['shape'][1]
+                    d_var['shape'] = [m, n]
+
+                d_var['datatype'   ] = d_vars[0]['datatype']
+                if var0_is_vector or var1_is_vector:
+                    d_var['rank'   ] = 1
+                else:
+                    d_var['rank'   ] = 2
+                d_var['allocatable'] = False
+                d_var['is_pointer' ] = False
+                d_var['precision'  ] = max(d_vars[0]['precision'],
+                                           d_vars[1]['precision'])
+
             elif name in ['Int','Int32','Int64','Real','Imag',
                           'Float32','Float64','Complex',
                           'Complex128','Complex64']:
@@ -1900,7 +1925,6 @@ class SemanticParser(BasicParser):
                 if str(lhs.dtype) != str(var.dtype):
                     txt = '|{name}| {old} <-> {new}'
                     txt = txt.format(name=name, old=var.dtype, new=lhs.dtype)
-                 
                     errors.report(INCOMPATIBLE_TYPES_IN_ASSIGNMENT,
                     symbol=txt,bounding_box=self._current_fst_node.absolute_bounding_box,
                     severity='error', blocker=False)


Test patch:
diff --git a/tests/epyccel/modules/arrays.py b/tests/epyccel/modules/arrays.py
index b345626d..05a9da3f 100644
--- a/tests/epyccel/modules/arrays.py
+++ b/tests/epyccel/modules/arrays.py
@@ -241,3 +241,49 @@ def array_real_1d_div_stack_array():
     for i in range(10):
         s += 1.0 / a[i]
     return s
+
+#==============================================================================
+# TEST: Product and matrix multiplication
+#==============================================================================
+
+@types('real[:], real[:], real[:]')
+def array_real_1d_1d_prod(x, out):
+    from numpy import prod
+    out[:] = prod(x)
+
+@types('real[:,:], real[:], real[:]')
+def array_real_2d_1d_matmul(A, x, out):
+    from numpy import matmul
+    out[:] = matmul(A, x)
+
+@types('real[:,:](order=F), real[:], real[:]')
+def array_real_2d_1d_matmul_order_F(A, x, out):
+    from numpy import matmul
+    out[:] = matmul(A, x)
+
+@types('real[:], real[:,:], real[:]')
+def array_real_1d_2d_matmul(x, A, out):
+    from numpy import matmul
+    out[:] = matmul(x, A)
+
+@types('real[:,:], real[:,:], real[:,:]')
+def array_real_2d_2d_matmul(A, B, out):
+    from numpy import matmul
+    out[:,:] = matmul(A, B)
+
+@types('real[:,:](order=F), real[:,:](order=F), real[:,:](order=F)')
+def array_real_2d_2d_matmul_F_F(A, B, out):
+    from numpy import matmul
+    out[:,:] = matmul(A, B)
+
+# Mixed order, not supported currently, see #244
+@types('real[:,:], real[:,:](order=F), real[:,:]')
+def array_real_2d_2d_matmul_mixorder(A, B, out):
+    from numpy import matmul
+    out[:,:] = matmul(A, B)
+
+@types('real[:], real[:], real[:]')
+def array_real_loopdiff(x, y, out):
+    dxy = x - y
+    for k in range(len(x)):
+        out[k] = dxy[k]
diff --git a/tests/epyccel/modules/loops.py b/tests/epyccel/modules/loops.py
index 0897e79a..7ecac6a5 100644
--- a/tests/epyccel/modules/loops.py
+++ b/tests/epyccel/modules/loops.py
@@ -138,3 +138,15 @@ def zip_prod( m ):
         res += i1*i2
 
     return res
+
+# ...
+@types( 'real[:], real[:]' )
+def product_loop_on_real_array( z, out ):
+
+    from numpy     import shape
+
+    s = shape( z )
+    n = s[0]
+
+    for i in range(n):
+        out[i] = z[i]**2
diff --git a/tests/epyccel/test_arrays.py b/tests/epyccel/test_arrays.py
index cc5a220e..53031498 100644
--- a/tests/epyccel/test_arrays.py
+++ b/tests/epyccel/test_arrays.py
@@ -720,6 +720,176 @@ def test_array_real_div_stack_array():
     x2 = f2()
     assert np.equal( x1, x2 )
 
+#==============================================================================
+# TEST: Product and matrix multiplication
+#==============================================================================
+def test_array_real_1d_1d_prod():
+    f1 = arrays.array_real_1d_1d_prod
+    f2 = epyccel( f1 )
+    x1 = np.array([3.0, 2.0, 1.0])
+    x2 = np.copy(x1)
+    y1 = np.empty(3)
+    y2 = np.empty(3)
+    f1(x1, y1)
+    f2(x2, y2)
+    assert np.array_equal(y1, y2)
+
+def test_array_real_2d_1d_matmul():
+    f1 = arrays.array_real_2d_1d_matmul
+    f2 = epyccel( f1 )
+    A1 = np.ones([3, 2])
+    A1[1,0] = 2
+    A2 = np.copy(A1)
+    x1 = np.ones([2, 1])
+    x2 = np.copy(x1)
+    y1 = np.empty([3,1])
+    y2 = np.empty([3,1])
+    f1(A1, x1, y1)
+    f2(A2, x2, y2)
+    assert np.array_equal(y1, y2)
+
+def test_array_real_2d_1d_matmul_order_C_F():
+    f1 = arrays.array_real_2d_1d_matmul
+    f2 = epyccel( f1 )
+    A1 = np.ones([3, 2], order='F')
+    A1[1,0] = 2
+    A2 = np.copy(A1)
+    x1 = np.ones([2, 1])
+    x2 = np.copy(x1)
+    y1 = np.empty([3, 1])
+    y2 = np.empty([3, 1])
+    f1(A1, x1, y1)
+    f2(A2, x2, y2)
+    assert np.array_equal(y1, y2)
+
+def test_array_real_1d_2d_matmul_order_F():
+    f1 = arrays.array_real_1d_2d_matmul
+    f2 = epyccel( f1 )
+    A1 = np.ones([3, 2], order='F')
+    A1[1, 0] = 2
+    A2 = np.copy(A1)
+    x1 = np.ones([1, 3])
+    x2 = np.copy(x1)
+    y1 = np.empty(2)
+    y2 = np.empty(2)
+    f1(x1, A1, y1)
+    f2(x2, A2, y2)
+    assert np.array_equal(y1, y2)
+
+def test_array_real_2d_1d_matmul_order_F_C():
+    f1 = arrays.array_real_2d_1d_matmul_order_F
+    f2 = epyccel( f1 )
+    A1 = np.ones([3, 2])
+    A1[1,0] = 2
+    A2 = np.copy(A1)
+    x1 = np.ones([2, 1])
+    x2 = np.copy(x1)
+    y1 = np.empty([3, 1])
+    y2 = np.empty([3, 1])
+    f1(A1, x1, y1)
+    f2(A2, x2, y2)
+    assert np.array_equal(y1, y2)
+
+def test_array_real_2d_1d_matmul_order_F_F():
+    f1 = arrays.array_real_2d_1d_matmul_order_F
+    f2 = epyccel( f1 )
+    A1 = np.ones([3, 2], order='F')
+    A1[1,0] = 2
+    A2 = np.copy(A1)
+    x1 = np.ones([2, 1])
+    x2 = np.copy(x1)
+    y1 = np.empty([3, 1])
+    y2 = np.empty([3, 1])
+    f1(A1, x1, y1)
+    f2(A2, x2, y2)
+    assert np.array_equal(y1, y2)
+
+def test_array_real_2d_2d_matmul():
+    f1 = arrays.array_real_2d_2d_matmul
+    f2 = epyccel( f1 )
+    A1 = np.ones([3, 2])
+    A1[1, 0] = 2
+    A2 = np.copy(A1)
+    B1 = np.ones([2, 3])
+    B2 = np.copy(B1)
+    C1 = np.empty([3,3])
+    C2 = np.empty([3,3])
+    f1(A1, B1, C1)
+    f2(A2, B2, C2)
+    assert np.array_equal(C1, C2)
+
+def test_array_real_2d_2d_matmul_C_C_F_F():
+    f1 = arrays.array_real_2d_2d_matmul
+    f2 = epyccel( f1 )
+    A1 = np.ones([3, 2], order='F')
+    A1[1, 0] = 2
+    A2 = np.copy(A1)
+    B1 = np.ones([2, 3], order='F')
+    B2 = np.copy(B1)
+    C1 = np.empty([3,3])
+    C2 = np.empty([3,3])
+    f1(A1, B1, C1)
+    f2(A2, B2, C2)
+    assert np.array_equal(C1, C2)
+
+def test_array_real_2d_2d_matmul_C_C_C_F():
+    f1 = arrays.array_real_2d_2d_matmul
+    f2 = epyccel( f1 )
+    A1 = np.ones([3, 2])
+    A1[1, 0] = 2
+    A2 = np.copy(A1)
+    B1 = np.ones([2, 3], order='F')
+    B2 = np.copy(B1)
+    C1 = np.empty([3,3])
+    C2 = np.empty([3,3])
+    f1(A1, B1, C1)
+    f2(A2, B2, C2)
+    assert np.array_equal(C1, C2)
+
+def test_array_real_2d_2d_matmul_F_F_F_F():
+    f1 = arrays.array_real_2d_2d_matmul_F_F
+    f2 = epyccel( f1 )
+    A1 = np.ones([3, 2], order='F')
+    A1[1, 0] = 2
+    A2 = np.copy(A1)
+    B1 = np.ones([2, 3], order='F')
+    B2 = np.copy(B1)
+    C1 = np.empty([3,3], order='F')
+    C2 = np.empty([3,3], order='F')
+    f1(A1, B1, C1)
+    f2(A2, B2, C2)
+    assert np.array_equal(C1, C2)
+
+@pytest.mark.xfail(reason="Should fail as long as mixed order not supported, see #244")
+def test_array_real_2d_2d_matmul_mixorder():
+    f1 = arrays.array_real_2d_2d_matmul_mixorder
+    f2 = epyccel( f1 )
+    A1 = np.ones([3, 2])
+    A1[1, 0] = 2
+    A2 = np.copy(A1)
+    B1 = np.ones([2, 3])
+    B2 = np.copy(B1)
+    C1 = np.empty([3,3])
+    C2 = np.empty([3,3])
+    f1(A1, B1, C1)
+    f2(A2, B2, C2)
+    assert np.array_equal(C1, C2)
+
+# TODO: Fix #245
+@pytest.mark.xfail(reason="Needs fixing, see #245")
+def test_array_real_loopdiff():
+    f1 = arrays.array_real_loopdiff
+    f2 = epyccel( f1 )
+    x1 = np.ones(5)
+    y1 = np.zeros(5)
+    x2 = np.copy(x1)
+    y2 = np.copy(y1)
+    z1 = np.empty(5)
+    z2 = np.empty(5)
+    f1(x1, y1, z1)
+    f2(x2, y2, z2)
+    assert np.array_equal(z1, z2)
+
 ##==============================================================================
 ## CLEAN UP GENERATED FILES AFTER RUNNING TESTS
 ##==============================================================================
diff --git a/tests/epyccel/test_loops.py b/tests/epyccel/test_loops.py
index 6016b9bf..2268bb37 100644
--- a/tests/epyccel/test_loops.py
+++ b/tests/epyccel/test_loops.py
@@ -99,6 +99,21 @@ def test_zip_prod():
 
     assert np.array_equal( f1(10), f2(10) )
 
+def test_loop_on_real_array():
+
+    f1 = loops.product_loop_on_real_array
+    f2 = epyccel( f1 )
+
+    z1 = np.ones(11)
+    out1 = np.empty_like(z1)
+    z2 = z1.copy()
+    out2 = out1.copy()
+
+    f1(z1, out1)
+    f2(z2, out2)
+
+    assert np.array_equal( out1, out2 )
+
 ##==============================================================================
 ## CLEAN UP GENERATED FILES AFTER RUNNING TESTS
 ##==============================================================================