1 /* 2 * Licensed to the Apache Software Foundation (ASF) under one or more 3 * contributor license agreements. See the NOTICE file distributed with 4 * this work for additional information regarding copyright ownership. 5 * The ASF licenses this file to You under the Apache License, Version 2.0 6 * (the "License"); you may not use this file except in compliance with 7 * the License. You may obtain a copy of the License at 8 * 9 * http://www.apache.org/licenses/LICENSE-2.0 10 * 11 * Unless required by applicable law or agreed to in writing, software 12 * distributed under the License is distributed on an "AS IS" BASIS, 13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 14 * See the License for the specific language governing permissions and 15 * limitations under the License. 16 */ 17 package org.apache.commons.rng.sampling.distribution; 18 19 import org.apache.commons.rng.UniformRandomProvider; 20 21 /** 22 * Sampling from a uniform distribution. 23 * 24 * <p>Sampling uses {@link UniformRandomProvider#nextDouble()}.</p> 25 * 26 * @since 1.0 27 */ 28 public class ContinuousUniformSampler 29 extends SamplerBase 30 implements SharedStateContinuousSampler { 31 /** The minimum ULP gap for the open interval when the doubles have the same sign. */ 32 private static final int MIN_ULP_SAME_SIGN = 2; 33 /** The minimum ULP gap for the open interval when the doubles have the opposite sign. */ 34 private static final int MIN_ULP_OPPOSITE_SIGN = 3; 35 36 /** Lower bound. */ 37 private final double lo; 38 /** Higher bound. */ 39 private final double hi; 40 /** Underlying source of randomness. */ 41 private final UniformRandomProvider rng; 42 43 /** 44 * Specialization to sample from an open interval {@code (lo, hi)}. 45 */ 46 private static class OpenIntervalContinuousUniformSampler extends ContinuousUniformSampler { 47 /** 48 * @param rng Generator of uniformly distributed random numbers. 49 * @param lo Lower bound. 50 * @param hi Higher bound. 51 */ 52 OpenIntervalContinuousUniformSampler(UniformRandomProvider rng, double lo, double hi) { 53 super(rng, lo, hi); 54 } 55 56 @Override 57 public double sample() { 58 final double x = super.sample(); 59 // Due to rounding using a variate u in the open interval (0,1) with the original 60 // algorithm may generate a value at the bound. Thus the bound is explicitly tested 61 // and the sample repeated if necessary. 62 if (x == getHi() || x == getLo()) { 63 return sample(); 64 } 65 return x; 66 } 67 68 @Override 69 public SharedStateContinuousSampler withUniformRandomProvider(UniformRandomProvider rng) { 70 return new OpenIntervalContinuousUniformSampler(rng, getLo(), getHi()); 71 } 72 } 73 74 /** 75 * @param rng Generator of uniformly distributed random numbers. 76 * @param lo Lower bound. 77 * @param hi Higher bound. 78 */ 79 public ContinuousUniformSampler(UniformRandomProvider rng, 80 double lo, 81 double hi) { 82 super(null); 83 this.rng = rng; 84 this.lo = lo; 85 this.hi = hi; 86 } 87 88 /** {@inheritDoc} */ 89 @Override 90 public double sample() { 91 final double u = rng.nextDouble(); 92 return u * hi + (1 - u) * lo; 93 } 94 95 /** 96 * Gets the lower bound. This is deliberately scoped as package private. 97 * 98 * @return the lower bound 99 */ 100 double getLo() { 101 return lo; 102 } 103 104 /** 105 * Gets the higher bound. This is deliberately scoped as package private. 106 * 107 * @return the higher bound 108 */ 109 double getHi() { 110 return hi; 111 } 112 113 /** {@inheritDoc} */ 114 @Override 115 public String toString() { 116 return "Uniform deviate [" + rng.toString() + "]"; 117 } 118 119 /** 120 * {@inheritDoc} 121 * 122 * @since 1.3 123 */ 124 @Override 125 public SharedStateContinuousSampler withUniformRandomProvider(UniformRandomProvider rng) { 126 return new ContinuousUniformSampler(rng, lo, hi); 127 } 128 129 /** 130 * Creates a new continuous uniform distribution sampler. 131 * 132 * @param rng Generator of uniformly distributed random numbers. 133 * @param lo Lower bound. 134 * @param hi Higher bound. 135 * @return the sampler 136 * @since 1.3 137 */ 138 public static SharedStateContinuousSampler of(UniformRandomProvider rng, 139 double lo, 140 double hi) { 141 return new ContinuousUniformSampler(rng, lo, hi); 142 } 143 144 /** 145 * Creates a new continuous uniform distribution sampler. 146 * 147 * <p>The bounds can be optionally excluded to sample from the open interval 148 * {@code (lower, upper)}. In this case if the bounds have the same sign the 149 * open interval must contain at least 1 double value between the limits; if the 150 * bounds have opposite signs the open interval must contain at least 2 double values 151 * between the limits excluding {@code -0.0}. Thus the interval {@code (-x,x)} will 152 * raise an exception when {@code x} is {@link Double#MIN_VALUE}. 153 * 154 * @param rng Generator of uniformly distributed random numbers. 155 * @param lo Lower bound. 156 * @param hi Higher bound. 157 * @param excludeBounds Set to {@code true} to use the open interval 158 * {@code (lower, upper)}. 159 * @return the sampler 160 * @throws IllegalArgumentException If the open interval is invalid. 161 * @since 1.4 162 */ 163 public static SharedStateContinuousSampler of(UniformRandomProvider rng, 164 double lo, 165 double hi, 166 boolean excludeBounds) { 167 if (excludeBounds) { 168 if (!validateOpenInterval(lo, hi)) { 169 throw new IllegalArgumentException("Invalid open interval (" + 170 lo + "," + hi + ")"); 171 } 172 return new OpenIntervalContinuousUniformSampler(rng, lo, hi); 173 } 174 return new ContinuousUniformSampler(rng, lo, hi); 175 } 176 177 /** 178 * Check that the open interval is valid. It must contain at least one double value 179 * between the limits if the signs are the same, or two double values between the limits 180 * if the signs are different (excluding {@code -0.0}). 181 * 182 * @param lo Lower bound. 183 * @param hi Higher bound. 184 * @return false is the interval is invalid 185 */ 186 private static boolean validateOpenInterval(double lo, double hi) { 187 // Get the raw bit representation. 188 long bitsx = Double.doubleToRawLongBits(lo); 189 long bitsy = Double.doubleToRawLongBits(hi); 190 // xor will set the sign bit if the signs are different 191 if ((bitsx ^ bitsy) < 0) { 192 // Opposite signs. Drop the sign bit to represent the count of 193 // bits above +0.0. When combined this should be above 2. 194 // This prevents the range (-Double.MIN_VALUE, Double.MIN_VALUE) 195 // which cannot be sampled unless the uniform deviate u=0.5. 196 // (MAX_VALUE has all bits set except the most significant sign bit.) 197 bitsx &= Long.MAX_VALUE; 198 bitsy &= Long.MAX_VALUE; 199 return !lessThanUnsigned(bitsx + bitsy, MIN_ULP_OPPOSITE_SIGN); 200 } 201 // Same signs, subtraction will count the ULP difference. 202 // This should be above 1. 203 return Math.abs(bitsx - bitsy) >= MIN_ULP_SAME_SIGN; 204 } 205 206 /** 207 * Compares two {@code long} values numerically treating the values as unsigned 208 * to test if the first value is less than the second value. 209 * 210 * <p>See Long.compareUnsigned(long, long) in JDK 1.8. 211 * 212 * @param x the first value 213 * @param y the second value 214 * @return true if {@code x < y} 215 */ 216 private static boolean lessThanUnsigned(long x, long y) { 217 return x + Long.MIN_VALUE < y + Long.MIN_VALUE; 218 } 219 }