Backport of fixes from master branch.

Fix: NextDoubleNonZero() could return zero (#5)

Redzen/Random/IRandomSource.cs

@@ -34,7 +34,7 @@ public interface IRandomSource
 
         /// <summary>
         /// Generate a random Int32 over the interval [minValue, maxValue), i.e. excluding maxValue.
-        /// maxValue must be >= minValue. minValue may be negative.
+        /// maxValue must be > minValue. minValue may be negative.
         /// </summary>
         int Next(int minValue, int maxValue);
 
@@ -44,11 +44,11 @@ public interface IRandomSource
         double NextDouble();
 
         /// <summary>
-        /// Generate a random double over the interval [0, 1), i.e. inclusive of 0.0 and exclusive of 1.0.
+        /// Generate a random double over the interval [0, 1], i.e. inclusive of both 0.0 and 1.0.
         /// </summary>
         /// <remarks>
         /// Uses an alternative sampling method that is capable of generating all possible values in the
-        /// interval [0,1) that can be represented by a double precision float. Note however that this method 
+        /// interval [0,1] that can be represented by a double precision float. Note however that this method 
         /// is significantly slower than NextDouble().
         /// </remarks>
         double NextDoubleHighRes();
@@ -82,7 +82,7 @@ public interface IRandomSource
         ulong NextULong();
 
         /// <summary>
-        /// Generates a random double over the interval (0, 1), i.e. exclusive of both 0.0 and 1.0
+        /// Generate a random double over the interval (0, 1], i.e. exclusive 0.0 and inclusive of 1.0.
         /// </summary>
         double NextDoubleNonZero();
 

Redzen/Random/RandomSourceBase.cs

@@ -59,10 +59,6 @@ public int Next(int maxValue)
                 throw new ArgumentOutOfRangeException(nameof(maxValue), maxValue, "maxValue must be > 0");
             }
 
-            if(1 == maxValue) {
-                return 0;
-            }
-
             return NextInner(maxValue);
         }
 
@@ -234,19 +230,13 @@ public ulong NextULong()
         }
 
         /// <summary>
-        /// Generate a random double over the interval (0, 1), i.e. exclusive of both 0.0 and 1.0
+        /// Generate a random double over the interval (0, 1], i.e. exclusive 0.0 and inclusive of 1.0.
         /// </summary>
         public double NextDoubleNonZero()
         {
-            // Here we generate a random value in the interval [0, 0x1f_ffff_ffff_fffe], and add one
-            // to generate a random value in the interval [1, 0x1f_ffff_ffff_ffff].
-            //
-            // We then multiply by the fractional unit 1.0 / 2^53 to obtain a floating point value 
-            // in the interval [ 1/(2^53-1) , 1.0].
-            //
-            // Note. the bit shift right here may appear redundant, however, the high significance
-            // bits have better randomness than the low bits, thus this approach is preferred.
-            return ((NextULongInner() >> 11) & 0x1f_ffff_ffff_fffe) * INCR_DOUBLE;
+            // Here we generate a random double in the interval [0, 1 - (1 / 2^53)], and add INCR_DOUBLE
+            // to produce a value in the interval [1 / 2^53, 1]
+            return NextDoubleInner() + INCR_DOUBLE;
         }
 
         /// <summary>

Redzen/Random/XorShiftRandom.cs

@@ -138,10 +138,6 @@ public int Next(int maxValue)
                 throw new ArgumentOutOfRangeException(nameof(maxValue), maxValue, "maxValue must be > 0");
             }
 
-            if(1 == maxValue) {
-                return 0;
-            }
-
             return NextInner(maxValue);
         }
 
@@ -375,16 +371,13 @@ public ulong NextULong()
         }
 
         /// <summary>
-        /// Generate a random double over the interval (0, 1), i.e. exclusive of both 0.0 and 1.0
+        /// Generate a random double over the interval (0, 1], i.e. exclusive of both 0.0 and 1.0
         /// </summary>
         public double NextDoubleNonZero()
         {
-            // Here we generate a random value in the interval [0, 0xffff_fffe], and add one
-            // to generate a random value in the interval [1, 0xffff_ffff].
-            //
-            // We then multiply by the fractional unit 1.0 / 2^32 to obtain a floating point value 
-            // in the interval [ 1/(2^32-1) , 1.0].
-            return ((NextInner() & 0xffff_fffe) + 1) * INCR_DOUBLE;
+            // Here we generate a random double in the interval [0, 1 - (1 / 2^32)], and add INCR_DOUBLE
+            // to produce a value in the interval [1 / 2^32, 1]
+            return NextDoubleInner() + INCR_DOUBLE;
         }
 
         /// <summary>
@@ -478,15 +471,15 @@ private long NextInner(long maxValue)
         private double NextDoubleInner()
         {
             // Notes. 
-            // Here we generate a random integer in the interval [0, 2^53-1]  (i.e. the max value is 53 binary 1s),
-            // and multiply by the fractional value 1.0 / 2^53, thus the result has a min value of 0.0 and a max value of 
-            // 1.0 - (1.0 / 2^53), or 0.99999999999999989 in decimal.
+            // Here we generate a random integer in the interval [0, 2^32-1]  (i.e. the max value is 32 binary 1s),
+            // and multiply by the fractional value 1.0 / 2^32, thus the result has a min value of 0.0 and a max value of 
+            // 1.0 - (1.0 / 2^32), or 0.99999999976716936 in decimal.
             //
-            // I.e. we break the interval [0,1) into 2^53 uniformly distributed discrete values, and thus the interval between
-            // two adjacent values is 1.0 / 2^53. This increment is chosen because it is the smallest value at which each 
-            // distinct value in the full range (from 0.0 to 1.0 exclusive) can be represented directly by a double precision
-            // float, and thus no rounding occurs in the representation of these values, which in turn ensures no bias in the 
-            // random samples.                   
+            // I.e. we break the interval [0,1) into 2^32 uniformly distributed discrete values, and thus the interval between
+            // two adjacent values is 1.0 / 2^32.
+            // Use of 32 bits was a historical choice based on that fact that the underlying XorShift PRNG produces 32 bits of randomness
+            // per invocation/cycle. This approach is maintained here for backwards compatibility, however, this class is deprecated in 
+            // favour of RandomSourceBase, which uses 53 bits of entropy per double precision float instead of 32.
             return NextInner() * INCR_DOUBLE;
         }
 

Redzen/Redzen.csproj

@@ -2,13 +2,14 @@
 
   <PropertyGroup>
     <TargetFramework>netstandard2.0</TargetFramework>
-    <Version>9.1.0.0</Version>
+    <Version>9.1.1.0</Version>
     <Copyright>Copyright Colin D. Green 2015-2020</Copyright>
     <Description>General purpose C# code library.</Description>
     <PackageProjectUrl>https://github.com/colgreen/Redzen</PackageProjectUrl>
     <PackageLicenseExpression>MIT</PackageLicenseExpression>
     <PackageTags>C# .NET rng prng ziggurat gaussian io sorting timsort xoshiro</PackageTags>
-    <PackageReleaseNotes>Added MathArrayUtils.Clip(). nuget ref updates.</PackageReleaseNotes>
+    <PackageReleaseNotes>Fix: NextDoubleNonZero() could return zero (very infrequently). (#5)
+</PackageReleaseNotes>
     <SignAssembly>true</SignAssembly>
     <AssemblyOriginatorKeyFile>colindgreen.pfx</AssemblyOriginatorKeyFile>
     <Authors>Colin D. Green</Authors>