1 /////////////////////////////////////////////////////////////////////////////////// 2 /// OpenGL Mathematics (glm.g-truc.net) 3 /// 4 /// Copyright (c) 2005 - 2014 G-Truc Creation (www.g-truc.net) 5 /// Permission is hereby granted, free of charge, to any person obtaining a copy 6 /// of this software and associated documentation files (the "Software"), to deal 7 /// in the Software without restriction, including without limitation the rights 8 /// to use, copy, modify, merge, publish, distribute, sublicense, and/or sell 9 /// copies of the Software, and to permit persons to whom the Software is 10 /// furnished to do so, subject to the following conditions: 11 /// 12 /// The above copyright notice and this permission notice shall be included in 13 /// all copies or substantial portions of the Software. 14 /// 15 /// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR 16 /// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, 17 /// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE 18 /// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER 19 /// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, 20 /// OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN 21 /// THE SOFTWARE. 22 /// 23 /// @ref gtc_ulp 24 /// @file glm/gtc/ulp.inl 25 /// @date 2011-03-07 / 2012-04-07 26 /// @author Christophe Riccio 27 /////////////////////////////////////////////////////////////////////////////////// 28 /// Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. 29 /// 30 /// Developed at SunPro, a Sun Microsystems, Inc. business. 31 /// Permission to use, copy, modify, and distribute this 32 /// software is freely granted, provided that this notice 33 /// is preserved. 34 /////////////////////////////////////////////////////////////////////////////////// 35 36 #include "../detail/type_int.hpp" 37 #include <cmath> 38 #include <cfloat> 39 #include <limits> 40 41 #if(GLM_COMPILER & GLM_COMPILER_VC) 42 # pragma warning(push) 43 # pragma warning(disable : 4127) 44 #endif 45 46 typedef union 47 { 48 float value; 49 /* FIXME: Assumes 32 bit int. */ 50 unsigned int word; 51 } ieee_float_shape_type; 52 53 typedef union 54 { 55 double value; 56 struct 57 { 58 glm::detail::int32 lsw; 59 glm::detail::int32 msw; 60 } parts; 61 } ieee_double_shape_type; 62 63 #define GLM_EXTRACT_WORDS(ix0,ix1,d) \ 64 do { \ 65 ieee_double_shape_type ew_u; \ 66 ew_u.value = (d); \ 67 (ix0) = ew_u.parts.msw; \ 68 (ix1) = ew_u.parts.lsw; \ 69 } while (0) 70 71 #define GLM_GET_FLOAT_WORD(i,d) \ 72 do { \ 73 ieee_float_shape_type gf_u; \ 74 gf_u.value = (d); \ 75 (i) = gf_u.word; \ 76 } while (0) 77 78 #define GLM_SET_FLOAT_WORD(d,i) \ 79 do { \ 80 ieee_float_shape_type sf_u; \ 81 sf_u.word = (i); \ 82 (d) = sf_u.value; \ 83 } while (0) 84 85 #define GLM_INSERT_WORDS(d,ix0,ix1) \ 86 do { \ 87 ieee_double_shape_type iw_u; \ 88 iw_u.parts.msw = (ix0); \ 89 iw_u.parts.lsw = (ix1); \ 90 (d) = iw_u.value; \ 91 } while (0) 92 93 namespace glm{ 94 namespace detail 95 { 96 GLM_FUNC_QUALIFIER float nextafterf(float x, float y) 97 { 98 volatile float t; 99 glm::detail::int32 hx, hy, ix, iy; 100 101 GLM_GET_FLOAT_WORD(hx, x); 102 GLM_GET_FLOAT_WORD(hy, y); 103 ix = hx&0x7fffffff; // |x| 104 iy = hy&0x7fffffff; // |y| 105 106 if((ix>0x7f800000) || // x is nan 107 (iy>0x7f800000)) // y is nan 108 return x+y; 109 if(x==y) return y; // x=y, return y 110 if(ix==0) { // x == 0 111 GLM_SET_FLOAT_WORD(x,(hy&0x80000000)|1);// return +-minsubnormal 112 t = x*x; 113 if(t==x) return t; else return x; // raise underflow flag 114 } 115 if(hx>=0) { // x > 0 116 if(hx>hy) { // x > y, x -= ulp 117 hx -= 1; 118 } else { // x < y, x += ulp 119 hx += 1; 120 } 121 } else { // x < 0 122 if(hy>=0||hx>hy){ // x < y, x -= ulp 123 hx -= 1; 124 } else { // x > y, x += ulp 125 hx += 1; 126 } 127 } 128 hy = hx&0x7f800000; 129 if(hy>=0x7f800000) return x+x; // overflow 130 if(hy<0x00800000) { // underflow 131 t = x*x; 132 if(t!=x) { // raise underflow flag 133 GLM_SET_FLOAT_WORD(y,hx); 134 return y; 135 } 136 } 137 GLM_SET_FLOAT_WORD(x,hx); 138 return x; 139 } 140 141 GLM_FUNC_QUALIFIER double nextafter(double x, double y) 142 { 143 volatile double t; 144 glm::detail::int32 hx, hy, ix, iy; 145 glm::detail::uint32 lx, ly; 146 147 GLM_EXTRACT_WORDS(hx, lx, x); 148 GLM_EXTRACT_WORDS(hy, ly, y); 149 ix = hx & 0x7fffffff; // |x| 150 iy = hy & 0x7fffffff; // |y| 151 152 if(((ix>=0x7ff00000)&&((ix-0x7ff00000)|lx)!=0) || // x is nan 153 ((iy>=0x7ff00000)&&((iy-0x7ff00000)|ly)!=0)) // y is nan 154 return x+y; 155 if(x==y) return y; // x=y, return y 156 if((ix|lx)==0) { // x == 0 157 GLM_INSERT_WORDS(x, hy & 0x80000000, 1); // return +-minsubnormal 158 t = x*x; 159 if(t==x) return t; else return x; // raise underflow flag 160 } 161 if(hx>=0) { // x > 0 162 if(hx>hy||((hx==hy)&&(lx>ly))) { // x > y, x -= ulp 163 if(lx==0) hx -= 1; 164 lx -= 1; 165 } else { // x < y, x += ulp 166 lx += 1; 167 if(lx==0) hx += 1; 168 } 169 } else { // x < 0 170 if(hy>=0||hx>hy||((hx==hy)&&(lx>ly))){// x < y, x -= ulp 171 if(lx==0) hx -= 1; 172 lx -= 1; 173 } else { // x > y, x += ulp 174 lx += 1; 175 if(lx==0) hx += 1; 176 } 177 } 178 hy = hx&0x7ff00000; 179 if(hy>=0x7ff00000) return x+x; // overflow 180 if(hy<0x00100000) { // underflow 181 t = x*x; 182 if(t!=x) { // raise underflow flag 183 GLM_INSERT_WORDS(y,hx,lx); 184 return y; 185 } 186 } 187 GLM_INSERT_WORDS(x,hx,lx); 188 return x; 189 } 190 }//namespace detail 191 }//namespace glm 192 193 #if(GLM_COMPILER & GLM_COMPILER_VC) 194 # pragma warning(pop) 195 #endif 196 197 namespace glm 198 { 199 template <> 200 GLM_FUNC_QUALIFIER float next_float(float const & x) 201 { 202 # if((GLM_LANG & GLM_LANG_CXX11_FLAG)) 203 return std::nextafter(x, std::numeric_limits<float>::max()); 204 # elif((GLM_COMPILER & GLM_COMPILER_VC) || ((GLM_COMPILER & GLM_COMPILER_INTEL) && (GLM_PLATFORM & GLM_PLATFORM_WINDOWS))) 205 return detail::nextafterf(x, FLT_MAX); 206 # else 207 return nextafterf(x, FLT_MAX); 208 # endif 209 } 210 211 template <> 212 GLM_FUNC_QUALIFIER double next_float(double const & x) 213 { 214 # if((GLM_LANG & GLM_LANG_CXX11_FLAG)) 215 return std::nextafter(x, std::numeric_limits<double>::max()); 216 # elif((GLM_COMPILER & GLM_COMPILER_VC) || ((GLM_COMPILER & GLM_COMPILER_INTEL) && (GLM_PLATFORM & GLM_PLATFORM_WINDOWS))) 217 return detail::nextafter(x, std::numeric_limits<double>::max()); 218 # else 219 return nextafter(x, DBL_MAX); 220 # endif 221 } 222 223 template<typename T, precision P, template<typename, precision> class vecType> 224 GLM_FUNC_QUALIFIER vecType<T, P> next_float(vecType<T, P> const & x) 225 { 226 vecType<T, P> Result; 227 for(length_t i = 0; i < Result.length(); ++i) 228 Result[i] = next_float(x[i]); 229 return Result; 230 } 231 232 GLM_FUNC_QUALIFIER float prev_float(float const & x) 233 { 234 # if((GLM_LANG & GLM_LANG_CXX11_FLAG)) 235 return std::nextafter(x, std::numeric_limits<float>::min()); 236 # elif((GLM_COMPILER & GLM_COMPILER_VC) || ((GLM_COMPILER & GLM_COMPILER_INTEL) && (GLM_PLATFORM & GLM_PLATFORM_WINDOWS))) 237 return detail::nextafterf(x, FLT_MIN); 238 # else 239 return nextafterf(x, FLT_MIN); 240 # endif 241 } 242 243 GLM_FUNC_QUALIFIER double prev_float(double const & x) 244 { 245 # if((GLM_LANG & GLM_LANG_CXX11_FLAG)) 246 return std::nextafter(x, std::numeric_limits<double>::min()); 247 # elif((GLM_COMPILER & GLM_COMPILER_VC) || ((GLM_COMPILER & GLM_COMPILER_INTEL) && (GLM_PLATFORM & GLM_PLATFORM_WINDOWS))) 248 return _nextafter(x, DBL_MIN); 249 # else 250 return nextafter(x, DBL_MIN); 251 # endif 252 } 253 254 template<typename T, precision P, template<typename, precision> class vecType> 255 GLM_FUNC_QUALIFIER vecType<T, P> prev_float(vecType<T, P> const & x) 256 { 257 vecType<T, P> Result; 258 for(length_t i = 0; i < Result.length(); ++i) 259 Result[i] = prev_float(x[i]); 260 return Result; 261 } 262 263 template <typename T> 264 GLM_FUNC_QUALIFIER T next_float(T const & x, uint const & ulps) 265 { 266 T temp = x; 267 for(uint i = 0; i < ulps; ++i) 268 temp = next_float(temp); 269 return temp; 270 } 271 272 template<typename T, precision P, template<typename, precision> class vecType> 273 GLM_FUNC_QUALIFIER vecType<T, P> next_float(vecType<T, P> const & x, vecType<uint, P> const & ulps) 274 { 275 vecType<T, P> Result; 276 for(length_t i = 0; i < Result.length(); ++i) 277 Result[i] = next_float(x[i], ulps[i]); 278 return Result; 279 } 280 281 template <typename T> 282 GLM_FUNC_QUALIFIER T prev_float(T const & x, uint const & ulps) 283 { 284 T temp = x; 285 for(uint i = 0; i < ulps; ++i) 286 temp = prev_float(temp); 287 return temp; 288 } 289 290 template<typename T, precision P, template<typename, precision> class vecType> 291 GLM_FUNC_QUALIFIER vecType<T, P> prev_float(vecType<T, P> const & x, vecType<uint, P> const & ulps) 292 { 293 vecType<T, P> Result; 294 for(length_t i = 0; i < Result.length(); ++i) 295 Result[i] = prev_float(x[i], ulps[i]); 296 return Result; 297 } 298 299 template <typename T> 300 GLM_FUNC_QUALIFIER uint float_distance(T const & x, T const & y) 301 { 302 uint ulp = 0; 303 304 if(x < y) 305 { 306 T temp = x; 307 while(temp != y)// && ulp < std::numeric_limits<std::size_t>::max()) 308 { 309 ++ulp; 310 temp = next_float(temp); 311 } 312 } 313 else if(y < x) 314 { 315 T temp = y; 316 while(temp != x)// && ulp < std::numeric_limits<std::size_t>::max()) 317 { 318 ++ulp; 319 temp = next_float(temp); 320 } 321 } 322 else // == 323 { 324 325 } 326 327 return ulp; 328 } 329 330 template<typename T, precision P, template<typename, precision> class vecType> 331 GLM_FUNC_QUALIFIER vecType<uint, P> float_distance(vecType<T, P> const & x, vecType<T, P> const & y) 332 { 333 vecType<uint, P> Result; 334 for(length_t i = 0; i < Result.length(); ++i) 335 Result[i] = float_distance(x[i], y[i]); 336 return Result; 337 } 338 }//namespace glm 339
