# HG changeset patch
# User Jordi Gutiérrez Hermoso <jordigh@octave.org>
# Date 1349109487 14400
# Node ID 930117c97760fe72097e7cba1a64cc4f29ffc8ab
# Parent  6a05cad5b694a9c380bb1a173befb066ee6d1ee2
mean.m: improve geometric mean algorithm for small numbers. Add test.

diff --git a/scripts/statistics/base/mean.m b/scripts/statistics/base/mean.m
--- a/scripts/statistics/base/mean.m
+++ b/scripts/statistics/base/mean.m
@@ -113,7 +113,11 @@
   if (strcmp (opt, "a"))
     y = sum (x, dim) / n;
   elseif (strcmp (opt, "g"))
-    y = prod (x, dim) .^ (1/n);
+    if (all (x(:) >= 0))
+      y = exp (sum (log (x), dim) ./ n);
+    else
+      error ("mean: X must not contain any negative values");
+    endif
   elseif (strcmp (opt, "h"))
     y = n ./ sum (1 ./ x, dim);
   else
@@ -131,6 +135,9 @@
 %! assert (mean (y), 0);
 %! assert (mean (z), [0, 10]);
 
+## Test small numbers
+%!assert (mean (repmat (0.1,1,1000), "g"), 0.1, 20*eps)
+
 %!assert (mean (magic (3), 1), [5, 5, 5])
 %!assert (mean (magic (3), 2), [5; 5; 5])
 %!assert (mean ([2 8], "g"), 4)