001/* 002 * Licensed to the Apache Software Foundation (ASF) under one or more 003 * contributor license agreements. See the NOTICE file distributed with 004 * this work for additional information regarding copyright ownership. 005 * The ASF licenses this file to You under the Apache License, Version 2.0 006 * (the "License"); you may not use this file except in compliance with 007 * the License. You may obtain a copy of the License at 008 * 009 * http://www.apache.org/licenses/LICENSE-2.0 010 * 011 * Unless required by applicable law or agreed to in writing, software 012 * distributed under the License is distributed on an "AS IS" BASIS, 013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 014 * See the License for the specific language governing permissions and 015 * limitations under the License. 016 */ 017package org.apache.commons.rng.sampling.distribution; 018 019import org.apache.commons.rng.UniformRandomProvider; 020 021/** 022 * Sampling from a uniform distribution. 023 * 024 * <p>Sampling uses {@link UniformRandomProvider#nextDouble()}.</p> 025 * 026 * @since 1.0 027 */ 028public class ContinuousUniformSampler 029 extends SamplerBase 030 implements SharedStateContinuousSampler { 031 /** The minimum ULP gap for the open interval when the doubles have the same sign. */ 032 private static final int MIN_ULP_SAME_SIGN = 2; 033 /** The minimum ULP gap for the open interval when the doubles have the opposite sign. */ 034 private static final int MIN_ULP_OPPOSITE_SIGN = 3; 035 036 /** Lower bound. */ 037 private final double lo; 038 /** Higher bound. */ 039 private final double hi; 040 /** Underlying source of randomness. */ 041 private final UniformRandomProvider rng; 042 043 /** 044 * Specialization to sample from an open interval {@code (lo, hi)}. 045 */ 046 private static class OpenIntervalContinuousUniformSampler extends ContinuousUniformSampler { 047 /** 048 * @param rng Generator of uniformly distributed random numbers. 049 * @param lo Lower bound. 050 * @param hi Higher bound. 051 */ 052 OpenIntervalContinuousUniformSampler(UniformRandomProvider rng, double lo, double hi) { 053 super(rng, lo, hi); 054 } 055 056 @Override 057 public double sample() { 058 final double x = super.sample(); 059 // Due to rounding using a variate u in the open interval (0,1) with the original 060 // algorithm may generate a value at the bound. Thus the bound is explicitly tested 061 // and the sample repeated if necessary. 062 if (x == getHi() || x == getLo()) { 063 return sample(); 064 } 065 return x; 066 } 067 068 @Override 069 public SharedStateContinuousSampler withUniformRandomProvider(UniformRandomProvider rng) { 070 return new OpenIntervalContinuousUniformSampler(rng, getLo(), getHi()); 071 } 072 } 073 074 /** 075 * @param rng Generator of uniformly distributed random numbers. 076 * @param lo Lower bound. 077 * @param hi Higher bound. 078 */ 079 public ContinuousUniformSampler(UniformRandomProvider rng, 080 double lo, 081 double hi) { 082 super(null); 083 this.rng = rng; 084 this.lo = lo; 085 this.hi = hi; 086 } 087 088 /** {@inheritDoc} */ 089 @Override 090 public double sample() { 091 final double u = rng.nextDouble(); 092 return u * hi + (1 - u) * lo; 093 } 094 095 /** 096 * Gets the lower bound. This is deliberately scoped as package private. 097 * 098 * @return the lower bound 099 */ 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}