scalar optimizations: added some scalar patterns, Canonizer, group_powers

_test_elemwise.py

@@ -164,27 +164,34 @@ class _test_CAReduce(unittest.TestCase):
 
 if __name__ == '__main__':
     unittest.main()
-#     x = modes.build(Tensor('float64', [0, 0], name = 'x'))
-#     y = modes.build(Tensor('float64', [0, 0], name = 'y'))
-#     e = Broadcast(SquareDiff, (x, y), {0:0}).out
+#     x = modes.build(Tensor('int32', [0, 0], name = 'x'))
+#     y = modes.build(Tensor('int32', [0, 0], name = 'y'))
+# #     x = modes.build(Tensor('float64', [0, 0], name = 'x'))
+# #     y = modes.build(Tensor('float64', [0, 0], name = 'y'))
+#     e = Broadcast(Pow, (x, y)).out
 #     f = gof.CLinker(env([x, y], [e])).make_function(inplace = False)
-#     xv = numpy.random.rand(1000, 1000)
-#     yv = numpy.random.rand(1000, 1000)
-#     zv = numpy.random.rand(1000, 1000)
+# #     xv = numpy.random.rand(1000, 1000)
+# #     yv = numpy.random.rand(1000, 1000)
+# #     zv = numpy.random.rand(1000, 1000)
+#     xv = numpy.random.randint(1, 5, (1000, 1000))
+#     yv = numpy.random.randint(1, 5, (1000, 1000))
 #     add = numpy.frompyfunc(lambda x, y: x + y, 2, 1)
 
-#     t0 = time.time()
-#     for i in xrange(100):
-#         xv -= yv
-#         xv *= xv
-# #        xv += yv
-#     print time.time() - t0
+# #     t0 = time.time()
+# #     for i in xrange(100):
+# #         xv / yv
+# #     print time.time() - t0
 
 #     t0 = time.time()
 #     for i in xrange(100):
 #         f(xv, yv)
 #     print time.time() - t0
-    
+
+    # speed ratios:
+    # add : 1
+    # mul : 1
+    # div : 2
+    # pow : 20
 
 
 

_test_scalar_opt.py

@@ -13,8 +13,16 @@ def inputs():
     x = Scalar('float64', name = 'x')
     y = Scalar('float64', name = 'y')
     z = Scalar('float64', name = 'z')
+    a = Scalar('float64', name = 'a')
     return x, y, z
 
+def more_inputs():
+    a = Scalar('float64', name = 'a')
+    b = Scalar('float64', name = 'b')
+    c = Scalar('float64', name = 'c')
+    d = Scalar('float64', name = 'd')
+    return a, b, c, d
+
 
 class _test_opts(unittest.TestCase):
 
@@ -24,9 +32,106 @@ class _test_opts(unittest.TestCase):
         g = Env([x], [e])
         assert str(g) == "[Pow(x, 2.0)]"
         gof.ConstantFinder().optimize(g)
-        opt2.optimize(g)
+        pow2sqr_float.optimize(g)
         assert str(g) == "[Sqr(x)]"
 
 
+# class _test_canonize(unittest.TestCase):
+
+#     def test_muldiv(self):
+#         x, y, z = inputs()
+#         a, b, c, d = more_inputs()
+# #        e = (2.0 * x) / (2.0 * y)
+# #        e = (2.0 * x) / (4.0 * y)
+# #        e = x / (y / z)
+# #        e = (x * y) / x
+# #        e = (x / y) * (y / z) * (z / x)
+# #        e = (a / b) * (b / c) * (c / d)
+# #        e = (a * b) / (b * c) / (c * d)
+# #        e = 2 * x / 2
+# #        e = x / y / x
+#         g = Env([x, y, z, a, b, c, d], [e])
+#         print g
+#         gof.ConstantFinder().optimize(g)
+#         mulfn = lambda *inputs: reduce(lambda x, y: x * y, (1,) + inputs)
+#         divfn = lambda x, y: x / y
+#         invfn = lambda x: 1 / x
+#         Canonizer(Mul, Div, Inv, mulfn, divfn, invfn).optimize(g)
+#         print g
+        
+#     def test_plusmin(self):
+#         x, y, z = inputs()
+#         a, b, c, d = more_inputs()
+# #        e = x - x
+# #        e = (2.0 + x) - (2.0 + y)
+# #        e = (2.0 + x) - (4.0 + y)
+# #        e = x - (y - z)
+# #        e = (x + y) - x
+# #        e = (x - y) + (y - z) + (z - x)
+# #        e = (a - b) + (b - c) + (c - d)
+# #        e = x + -y
+# #        e = a - b - b + a + b + c + b - c
+#         e = x + log(y) - x + y
+#         g = Env([x, y, z, a, b, c, d], [e])
+#         print g
+#         gof.ConstantFinder().optimize(g)
+#         addfn = lambda *inputs: reduce(lambda x, y: x + y, (0,) + inputs)
+#         subfn = lambda x, y: x - y
+#         negfn = lambda x: -x
+#         Canonizer(Add, Sub, Neg, addfn, subfn, negfn).optimize(g)
+#         print g
+
+#     def test_both(self):
+#         x, y, z = inputs()
+#         a, b, c, d = more_inputs()
+#         e0 = (x * y / x)
+#         e = e0 + e0 - e0
+#         g = Env([x, y, z, a, b, c, d], [e])
+#         print g
+#         gof.ConstantFinder().optimize(g)
+#         mulfn = lambda *inputs: reduce(lambda x, y: x * y, (1,) + inputs)
+#         divfn = lambda x, y: x / y
+#         invfn = lambda x: 1 / x
+#         Canonizer(Mul, Div, Inv, mulfn, divfn, invfn).optimize(g)
+#         addfn = lambda *inputs: reduce(lambda x, y: x + y, (0,) + inputs)
+#         subfn = lambda x, y: x - y
+#         negfn = lambda x: -x
+#         Canonizer(Add, Sub, Neg, addfn, subfn, negfn).optimize(g)
+#         print g
+
+#     def test_group_powers(self):
+#         x, y, z = inputs()
+#         a, b, c, d = more_inputs()
+# #        e = x * exp(y) * exp(z)
+# #        e = x * pow(x, y) * pow(x, z)
+# #        e = pow(x, y) / pow(x, z)
+# #        e = pow(x, 2.0) * pow(x, y) / pow(x, 7.0)
+# #        e = pow(x - x, y)
+# #        e = pow(x, 2.0 + y - 7.0)
+# #        e = pow(x, 2.0) * pow(x, y) / pow(x, 7.0) / pow(x, z)
+# #        e = pow(x, 2.0 + y - 7.0 - z)
+# #        e = x ** y / x ** y
+# #        e = x ** y / x ** (y - 1.0)
+#         e = exp(x) * a * exp(y) / exp(z)
+#         g = Env([x, y, z, a, b, c, d], [e])
+#         print g
+#         gof.ConstantFinder().optimize(g)
+#         mulfn = lambda *inputs: reduce(lambda x, y: x * y, (1,) + inputs)
+#         divfn = lambda x, y: x / y
+#         invfn = lambda x: 1 / x
+#         Canonizer(Mul, Div, Inv, mulfn, divfn, invfn, group_powers).optimize(g)
+#         print g
+#         addfn = lambda *inputs: reduce(lambda x, y: x + y, (0,) + inputs)
+#         subfn = lambda x, y: x - y
+#         negfn = lambda x: -x
+#         Canonizer(Add, Sub, Neg, addfn, subfn, negfn).optimize(g)
+#         print g
+#         pow2one_float.optimize(g)
+#         pow2x_float.optimize(g)
+#         print g
+        
+
+
 if __name__ == '__main__':
     unittest.main()
+

_test_tensor.py

@@ -219,6 +219,47 @@ def make_broadcast_tester(op_class, expected, checks = {}, **kwargs):
     return make_tester(name, op_class, expected, checks, **kwargs)
 
 
+def make_broadcast_tester_unary(op_class, expected, checks = {}, **kwargs):
+
+    _randint = randint
+    _rand = rand
+    if kwargs.has_key('nonzero'):
+        if kwargs['nonzero']:
+            _randint = banzero(_randint)
+            _rand = banzero(_rand)
+        del kwargs['nonzero']
+
+    if kwargs.has_key('positive'):
+        if kwargs['positive']:
+            _randint = banneg(_randint)
+            _rand = banneg(_rand)
+        del kwargs['positive']
+
+    _good_broadcast = dict(normal = (_rand(2, 3), ),
+                           int = (_rand(2, 3), ))
+
+    _bad_build_broadcast = dict()
+
+    _bad_runtime_broadcast = dict()
+
+    _grad_broadcast = dict(normal = (_rand(2, 3), ),
+                           int = (_rand(2, 3), ))
+    
+    kwargs.setdefault('good', _good_broadcast)
+    kwargs.setdefault('bad_build', _bad_build_broadcast)
+    kwargs.setdefault('bad_runtime', _bad_runtime_broadcast)
+    kwargs.setdefault('grad', _grad_broadcast)
+    name = op_class.__name__ + "Tester"
+    if kwargs.has_key('inplace'):
+        if kwargs['inplace']:
+            _expected = expected
+            expected = lambda *inputs: numpy.array(_expected(*inputs), dtype = inputs[0].dtype)
+            checks = dict(checks,
+                          inplace_check = lambda inputs, outputs: inputs[0] is outputs[0])
+        del kwargs['inplace']
+    return make_tester(name, op_class, expected, checks, **kwargs)
+
+
 
 
 
@@ -264,11 +305,11 @@ def make_broadcast_tester(op_class, expected, checks = {}, **kwargs):
 #                                          good = _pow_good)
 
 
-# AbsTester = make_broadcast_tester(op_class = Abs,
-#                                   expected = lambda x: abs(x))
-# AbsInplaceTester = make_broadcast_tester(op_class = AbsInplace,
-#                                          expected = lambda x: abs(x),
-#                                          inplace = True)
+AbsTester = make_broadcast_tester_unary(op_class = Abs,
+                                        expected = lambda x: abs(x))
+AbsInplaceTester = make_broadcast_tester_unary(op_class = AbsInplace,
+                                               expected = lambda x: abs(x),
+                                               inplace = True)
 
 # ExpTester = make_broadcast_tester(op_class = Exp,
 #                                   expected = lambda x: numpy.exp(x))

gof/op.py

@@ -273,7 +273,7 @@ class GuardedOp(Op):
         try:
             if not old.same_properties(new):
                 raise TypeError("The new input must have the same properties as the previous one.")
-        except AbstractFunction:
+        except AbstractFunctionError:
             pass
         Op.set_input(self, i, new)
 

scalar_opt.py

 
 from scalar import *
-from gof import PatternOptimizer
+from gof import PatternOptimizer as Pattern
+from gof import utils
 
-c2 = constant(2.0)
+C = constant
 
-opt1 = PatternOptimizer((Mul, 'x', 'x'), (Sqr, 'x'))
-opt2 = PatternOptimizer((Pow, 'x', c2), (Sqr, 'x'))
+# x**2 -> x*x
+pow2sqr_float = Pattern((Pow, 'x', C(2.0)), (Sqr, 'x'))
+pow2sqr_int = Pattern((Pow, 'x', C(2)), (Sqr, 'x'))
 
+# x**0 -> 1
+pow2one_float = Pattern((Pow, 'x', C(0.0)), C(1.0))
+pow2one_int = Pattern((Pow, 'x', C(0)), C(1))
 
+# x**1 -> x
+pow2x_float = Pattern((Pow, 'x', C(1.0)), 'x')
+pow2x_int = Pattern((Pow, 'x', C(1)), 'x')
+
+# log(x**y) -> y*log(x)
+logpow = Pattern((Log, (Pow, 'x', 'y')),
+                 (Mul, 'y', (Log, 'x')))
+
+
+class Canonizer(gof.Optimizer):
+
+    def __init__(self, main, inverse, reciprocal, mainfn, invfn, recfn, transform = None):
+        self.main = main
+        self.inverse = inverse
+        self.reciprocal = reciprocal
+        self.mainfn = mainfn
+        self.invfn = invfn
+        self.recfn = recfn
+        self.neutral = mainfn()
+        self.transform = transform
+
+    def apply(self, env):
+
+        def canonize(r):
+            if r in env.inputs or r in env.orphans():
+                return
+            
+            def flatten(r, nclients_check = True):
+                op = r.owner
+                if op is None or r in env.inputs or r in env.orphans():
+                    return [r], []
+                
+                results = [r2.dtype == r.dtype and flatten(r2) or ([r2], []) for r2 in op.inputs]
+                if isinstance(op, self.main) and (not nclients_check or env.nclients(r) == 1):
+                    nums = [x[0] for x in results]
+                    denums = [x[1] for x in results]
+                elif isinstance(op, self.inverse) and (not nclients_check or env.nclients(r) == 1):
+                    nums = [results[0][0], results[1][1]]
+                    denums = [results[0][1], results[1][0]]
+                elif isinstance(op, self.reciprocal) and (not nclients_check or env.nclients(r) == 1):
+                    nums = [results[0][1]]
+                    denums = [results[0][0]]
+                else:
+                    return [r], []
+
+                return reduce(list.__add__, nums), reduce(list.__add__, denums)
+
+            num, denum = flatten(r, False)
+
+            if (num, denum) == ([r], []):
+                if r.owner is None:
+                    return
+                else:
+                    for input in r.owner.inputs:
+                        canonize(input)
+                    return
+            
+            for d in list(denum):
+                if d in list(num):
+                    num.remove(d)
+                    denum.remove(d)
+
+            numct, num = utils.partition(lambda factor: getattr(factor, 'constant', False) and factor.data is not None, num)
+            denumct, denum = utils.partition(lambda factor: getattr(factor, 'constant', False) and factor.data is not None, denum)
+            
+            v = self.invfn(self.mainfn(*[x.data for x in numct]), self.mainfn(*[x.data for x in denumct]))
+            if v != self.neutral:
+                num.insert(0, C(v))
+
+            if self.transform is not None:
+                num, denum = self.transform(env, num, denum)
+
+            def make(factors):
+                n = len(factors)
+                if n == 0:
+                    return None
+                elif n == 1:
+                    return factors[0]
+                else:
+                    return self.main(*factors).out
+
+            numr, denumr = make(num), make(denum)
+            
+            if numr is None:
+                if denumr is None:
+                    new_r = Scalar(dtype = r.dtype)
+                    new_r.constant = True
+                    new_r.data = self.neutral
+                else:
+                    new_r = self.reciprocal(denumr).out
+            else:
+                if denumr is None:
+                    new_r = numr
+                else:
+                    new_r = self.inverse(numr, denumr).out
+
+            env.replace(r, new_r)
+
+            for factor in num + denum:
+                canonize(factor)
+
+        for output in env.outputs:
+            canonize(output)
+
+
+def group_powers(env, num, denum):
+
+    num_powers = {}
+    denum_powers = {}
+
+    def populate(d, seq):
+        for factor in list(seq):
+            op = factor.owner
+            if op is None or factor in env.inputs or factor in env.orphans():
+                continue
+            if isinstance(op, Exp):
+                d.setdefault('e', []).append(op.inputs[0])
+                seq.remove(factor)
+            elif isinstance(op, Pow):
+                d.setdefault(op.inputs[0], []).append(op.inputs[1])
+                seq.remove(factor)
+
+    populate(num_powers, num)
+    populate(denum_powers, denum)
+    
+    for x in set(num_powers.keys() + denum_powers.keys()):
+        
+        try: num_ys = num_powers.pop(x)
+        except KeyError: num_ys = []
+        
+        try: denum_ys = denum_powers.pop(x)
+        except KeyError: denum_ys = []
+
+        num_r = num_ys and add(*num_ys) or C(0)
+        denum_r = denum_ys and add(*denum_ys) or C(0)
+        if x == 'e':
+            num.append(exp(num_r - denum_r))
+        else:
+            num.append(pow(x, num_r - denum_r))
+
+    return num, denum
+
+
+def simple_factorize(env, num, denum):
+
+    # a*b + a*c -> a*(b+c)
+    # a*b + a*c + b*c -> a*(b+c) + b*c
+    #                 -> a*b + (a+b)*c
+    #  => a: {b, c}, b: {a, c}, c: {a, b}
+    # a*c + a*d + b*c + b*d
+    #  => a: {c, d}, b: {c, d}, c: {a, b}, d: {a, b}
+    # (a+b*x)*(c+d) --> a*c + a*d + b*x*c + b*x*d
+    #  => a: {c, d}, b: {xc, xd}, c: {a, bx}, d: {a, bx}, x: {bc, bd}
+    
+    pass