diff --git a/lbmpy/fluctuatinglb.py b/lbmpy/fluctuatinglb.py
index 86591879ca18bdf4682b9d6b9e02171d3ca5363c..833b93fe9699d5cd379b5312a441bfef47e9cd68 100644
--- a/lbmpy/fluctuatinglb.py
+++ b/lbmpy/fluctuatinglb.py
@@ -63,7 +63,8 @@ def fluctuation_correction(method, rng_generator, variances=SymbolGen("variance"
     conserved_moments = {sp.sympify(1), *MOMENT_SYMBOLS}
 
     # A diagonal matrix containing the random fluctuations
-    random_matrix = sp.Matrix([0 if m in conserved_moments else next(rng_generator) for m in method.moments])
+    random_matrix = sp.Matrix([0 if m in conserved_moments else (next(rng_generator) - 0.5) * sp.sqrt(12)
+                               for m in method.moments])
     random_variance = sp.diag(*[v for v, _ in zip(iter(variances), method.moments)])
 
     # corrections are applied in real space hence we need to convert to real space here