general dtype inference using solve

theano/tensor/slinalg.py

@@ -202,12 +202,12 @@ class Solve(Op):
         b = as_tensor_variable(b)
         assert A.ndim == 2
         assert b.ndim in [1, 2]
-        if ((A.dtype == 'float32' and b.dtype == 'float32')
-                or (A.dtype in ['int8', 'int16'] and b.dtype == 'float32')
-                or (b.dtype in ['int8', 'int16'] and A.dtype == 'float32')):
-            o_dtype = 'float32'
-        else:
-            o_dtype = 'float64'
+
+        # infer dtype by solving the most simple
+        # case with (1, 1) matrices
+        o_dtype = scipy.linalg.solve(
+            numpy.eye(1).astype(A.dtype),
+            numpy.eye(1).astype(b.dtype)).dtype
         x = tensor.tensor(
             broadcastable=b.broadcastable,
             dtype=o_dtype)

theano/tensor/tests/test_slinalg.py

@@ -235,12 +235,14 @@ class test_Solve(utt.InferShapeTester):
         if not imported_scipy:
             raise SkipTest("Scipy needed for the Solve op.")
 
-        dtypes = ['int8', 'int16', 'int32', 'int64',
-                  'float32', 'float64']
+        dtypes = ['uint8', 'uint16', 'uint32', 'uint64',
+                  'int8', 'int16', 'int32', 'int64',
+                  'float16', 'float32', 'float64']
 
         A_val = numpy.eye(2)
         b_val = numpy.ones((2, 1))
 
+        # try all dtype combinations
         for A_dtype, b_dtype in itertools.product(dtypes, dtypes):
             A = tensor.matrix(dtype=A_dtype)
             b = tensor.matrix(dtype=b_dtype)