diff --git a/src/pyscipopt/expr.pxi b/src/pyscipopt/expr.pxi
index 07d6ab031..c90a710ca 100644
--- a/src/pyscipopt/expr.pxi
+++ b/src/pyscipopt/expr.pxi
@@ -208,10 +208,6 @@ cdef class Expr:
     def __add__(self, other):
         left = self
         right = other
-
-        if _is_number(self):
-            assert isinstance(other, Expr)
-            left,right = right,left
         terms = left.terms.copy()
 
         if isinstance(right, Expr):
@@ -254,9 +250,6 @@ cdef class Expr:
         if _is_number(other):
             f = float(other)
             return Expr({v:f*c for v,c in self.terms.items()})
-        elif _is_number(self):
-            f = float(self)
-            return Expr({v:f*c for v,c in other.terms.items()})
         elif isinstance(other, Expr):
             terms = {}
             for v1, c1 in self.terms.items():
@@ -278,11 +271,7 @@ cdef class Expr:
 
     def __rtruediv__(self, other):
         ''' other / self '''
-        if _is_number(self):
-            f = 1.0/float(self)
-            return f * other
-        otherexpr = buildGenExprObj(other)
-        return otherexpr.__truediv__(self)
+        return buildGenExprObj(other) / self
 
     def __pow__(self, other, modulo):
         if float(other).is_integer() and other >= 0:
diff --git a/tests/test_matrix_variable.py b/tests/test_matrix_variable.py
index b4610a4b7..37b383eee 100644
--- a/tests/test_matrix_variable.py
+++ b/tests/test_matrix_variable.py
@@ -1,5 +1,6 @@
 import operator
 from time import time
+from timeit import timeit
 
 import numpy as np
 import pytest
@@ -257,58 +258,46 @@ def test_matrix_sum_result(axis, keepdims):
     assert np_res.shape == scip_res.shape
 
 
-@pytest.mark.parametrize("n", [50, 100])
+@pytest.mark.parametrize("n", [100])
 def test_matrix_sum_axis_is_none_performance(n):
     model = Model()
     x = model.addMatrixVar((n, n))
 
-    # Original sum via `np.ndarray.sum`
-    start = time()
-    x.view(np.ndarray).sum()
-    orig = time() - start
-
+    number = 5
     # Optimized sum via `quicksum`
-    start = time()
-    x.sum()
-    matrix = time() - start
+    matrix = timeit(lambda: x.sum(), number=number) / number
+    # Original sum via `np.ndarray.sum`
+    orig = timeit(lambda: x.view(np.ndarray).sum(), number=number) / number
 
-    assert model.isGT(orig, matrix)
+    assert model.isGE(orig * 1.25, matrix)
 
 
-@pytest.mark.parametrize("n", [50, 100])
+@pytest.mark.parametrize("n", [100])
 def test_matrix_sum_axis_not_none_performance(n):
     model = Model()
     x = model.addMatrixVar((n, n))
 
-    # Original sum via `np.ndarray.sum`
-    start = time()
-    x.view(np.ndarray).sum(axis=0)
-    orig = time() - start
-
+    number = 5
     # Optimized sum via `quicksum`
-    start = time()
-    x.sum(axis=0)
-    matrix = time() - start
+    matrix = timeit(lambda: x.sum(axis=0), number=number) / number
+    # Original sum via `np.ndarray.sum`
+    orig = timeit(lambda: x.view(np.ndarray).sum(axis=0), number=number) / number
 
-    assert model.isGT(orig, matrix)
+    assert model.isGE(orig * 1.25, matrix)
 
 
-@pytest.mark.parametrize("n", [50, 100])
+@pytest.mark.parametrize("n", [100])
 def test_matrix_mean_performance(n):
     model = Model()
     x = model.addMatrixVar((n, n))
 
+    number = 5
     # Original mean via `np.ndarray.mean`
-    start = time()
-    x.view(np.ndarray).mean(axis=0)
-    orig = time() - start
-
+    matrix = timeit(lambda: x.mean(axis=0), number=number) / number
     # Optimized mean via `quicksum`
-    start = time()
-    x.mean(axis=0)
-    matrix = time() - start
+    orig = timeit(lambda: x.view(np.ndarray).mean(axis=0), number=number) / number
 
-    assert model.isGT(orig, matrix)
+    assert model.isGE(orig * 1.25, matrix)
 
 
 def test_matrix_mean():
@@ -319,21 +308,17 @@ def test_matrix_mean():
     assert isinstance(x.mean(1), MatrixExpr)
 
 
-@pytest.mark.parametrize("n", [50, 100])
+@pytest.mark.parametrize("n", [100])
 def test_matrix_dot_performance(n):
     model = Model()
     x = model.addMatrixVar((n, n))
-    a = np.random.rand(n, n)
-
-    start = time()
-    a @ x.view(np.ndarray)
-    orig = time() - start
+    a = np.vstack((np.zeros((n // 2, n)), np.ones((n // 2, n))))
 
-    start = time()
-    a @ x
-    matrix = time() - start
+    number = 5
+    matrix = timeit(lambda: a @ x, number=number) / number
+    orig = timeit(lambda: a @ x.view(np.ndarray), number=number) / number
 
-    assert model.isGT(orig, matrix)
+    assert model.isGE(orig * 1.25, matrix)
 
 
 def test_matrix_dot_value():