tools/source/generic/fract.cxx

   1 /* -*- Mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
   2 /*
   3  * This file is part of the LibreOffice project.
   4  *
   5  * This Source Code Form is subject to the terms of the Mozilla Public
   6  * License, v. 2.0. If a copy of the MPL was not distributed with this
   7  * file, You can obtain one at http://mozilla.org/MPL/2.0/.
   8  *
   9  * This file incorporates work covered by the following license notice:
  10  *
  11  *   Licensed to the Apache Software Foundation (ASF) under one or more
  12  *   contributor license agreements. See the NOTICE file distributed
  13  *   with this work for additional information regarding copyright
  14  *   ownership. The ASF licenses this file to you under the Apache
  15  *   License, Version 2.0 (the "License"); you may not use this file
  16  *   except in compliance with the License. You may obtain a copy of
  17  *   the License at http://www.apache.org/licenses/LICENSE-2.0 .
  18  */
  19
  20 #include <tools/fract.hxx>
  21 #include <tools/debug.hxx>
  22 #include <o3tl/hash_combine.hxx>
  23 #include <o3tl/safeint.hxx>
  24 #include <sal/log.hxx>
  25 #include <osl/diagnose.h>
  26
  27 #include <algorithm>
  28 #include <cmath>
  29 #include <numeric>
  30
  31 #include <boost/rational.hpp>
  32
  33 #ifdef _MSC_VER
  34 #include <intrin.h>
  35 #endif
  36
  37 static boost::rational<sal_Int32> rational_FromDouble(double dVal);
  38 static void rational_ReduceInaccurate(boost::rational<sal_Int32>& rRational, unsigned nSignificantBits);
  39 static int impl_NumberOfBits( sal_uInt32 nNum );
  40
  41 static boost::rational<sal_Int32> toRational(sal_Int32 n, sal_Int32 d)
  42 {
  43     // https://github.com/boostorg/boost/issues/335 when these are std::numeric_limits<sal_Int32>::min
  44     if (n == d)
  45         return 1;
  46     // tdf#144319 avoid boost::bad_rational e.g. if numerator=-476741369, denominator=-2147483648
  47     if (d < -std::numeric_limits<sal_Int32>::max())
  48         return 0;
  49     return boost::rational<sal_Int32>(n, d);
  50 }
  51
  52 static constexpr bool isOutOfRange(sal_Int64 nNum)
  53 {
  54     return nNum < std::numeric_limits<sal_Int32>::min()
  55             || nNum > std::numeric_limits<sal_Int32>::max();
  56 }
  57
  58 // Initialized by setting nNum as nominator and nDen as denominator
  59 // Negative values in the denominator are invalid and cause the
  60 // inversion of both nominator and denominator signs
  61 // in order to return the correct value.
  62 Fraction::Fraction( sal_Int64 nNum, sal_Int64 nDen ) : mnNumerator(nNum), mnDenominator(nDen)
  63 {
  64     if ( isOutOfRange(nNum) || isOutOfRange(nDen) )
  65     {
  66         // tdf#143200
  67         SAL_WARN("tools.fraction", "values outside of range we can represent, doing reduction, which will reduce precision");
  68         do
  69         {
  70             mnNumerator /= 2;
  71             mnDenominator /= 2;
  72         } while (isOutOfRange(mnNumerator) || isOutOfRange(mnDenominator));
  73     }
  74     if ( mnDenominator == 0 )
  75     {
  76         mbValid = false;
  77         SAL_WARN( "tools.fraction", "'Fraction(" << nNum << ",0)' invalid fraction created" );
  78         return;
  79     }
  80     else if ((nDen == -1 && nNum == std::numeric_limits<sal_Int32>::min()) ||
  81             (nNum == -1 && nDen == std::numeric_limits<sal_Int32>::min()))
  82     {
  83         mbValid = false;
  84         SAL_WARN("tools.fraction", "'Fraction(" << nNum << "," << nDen << ")' invalid fraction created");
  85         return;
  86     }
  87 }
  88
  89 /**
  90  * only here to prevent passing of NaN
  91  */
  92 Fraction::Fraction( double nNum, double nDen )
  93     : Fraction(sal_Int64(nNum), sal_Int64(nDen))
  94 {}
  95
  96 Fraction::Fraction( double dVal )
  97 {
  98     try
  99     {
 100         boost::rational<sal_Int32> v = rational_FromDouble( dVal );
 101         mnNumerator = v.numerator();
 102         mnDenominator = v.denominator();
 103     }
 104     catch (const boost::bad_rational&)
 105     {
 106         mbValid = false;
 107         SAL_WARN( "tools.fraction", "'Fraction(" << dVal << ")' invalid fraction created" );
 108     }
 109 }
 110
 111 Fraction::operator double() const
 112 {
 113     if (!mbValid)
 114     {
 115         SAL_WARN( "tools.fraction", "'double()' on invalid fraction" );
 116         return 0.0;
 117     }
 118
 119     return boost::rational_cast<double>(toRational(mnNumerator, mnDenominator));
 120 }
 121
 122 // This methods first validates both values.
 123 // If one of the arguments is invalid, the whole operation is invalid.
 124 // After computation detect if result overflows a sal_Int32 value
 125 // which cause the operation to be marked as invalid
 126 Fraction& Fraction::operator += ( const Fraction& rVal )
 127 {
 128     if ( !rVal.mbValid )
 129         mbValid = false;
 130
 131     if ( !mbValid )
 132     {
 133         SAL_WARN( "tools.fraction", "'operator +=' with invalid fraction" );
 134         return *this;
 135     }
 136
 137     boost::rational<sal_Int32> a = toRational(mnNumerator, mnDenominator);
 138     a += toRational(rVal.mnNumerator, rVal.mnDenominator);
 139     mnNumerator = a.numerator();
 140     mnDenominator = a.denominator();
 141
 142     return *this;
 143 }
 144
 145 Fraction& Fraction::operator -= ( const Fraction& rVal )
 146 {
 147     if ( !rVal.mbValid )
 148         mbValid = false;
 149
 150     if ( !mbValid )
 151     {
 152         SAL_WARN( "tools.fraction", "'operator -=' with invalid fraction" );
 153         return *this;
 154     }
 155
 156     boost::rational<sal_Int32> a = toRational(mnNumerator, mnDenominator);
 157     a -= toRational(rVal.mnNumerator, rVal.mnDenominator);
 158     mnNumerator = a.numerator();
 159     mnDenominator = a.denominator();
 160
 161     return *this;
 162 }
 163
 164 namespace
 165 {
 166     bool checked_multiply_by(boost::rational<sal_Int32>& i, const boost::rational<sal_Int32>& r)
 167     {
 168         // Protect against self-modification
 169         sal_Int32 num = r.numerator();
 170         sal_Int32 den = r.denominator();
 171
 172         // Fast-path if the number of bits in input is < the number of bits in the output, overflow cannot happen
 173         // This is considerably faster than repeated std::gcd() operations
 174         if ((impl_NumberOfBits(std::abs(i.numerator())) + impl_NumberOfBits(std::abs(r.numerator()))) < 32 &&
 175             (impl_NumberOfBits(std::abs(i.denominator())) + impl_NumberOfBits(std::abs(r.denominator()))) < 32)
 176         {
 177             i *= r;
 178             return false;
 179         }
 180
 181         // Avoid overflow and preserve normalization
 182         sal_Int32 gcd1 = std::gcd(i.numerator(), den);
 183         sal_Int32 gcd2 = std::gcd(num, i.denominator());
 184
 185         if (!gcd1 || !gcd2)
 186             return true;
 187
 188         bool fail = false;
 189         fail |= o3tl::checked_multiply(i.numerator() / gcd1, num / gcd2, num);
 190         fail |= o3tl::checked_multiply(i.denominator() / gcd2, den / gcd1, den);
 191
 192         if (!fail)
 193             i.assign(num, den);
 194
 195         return fail;
 196     }
 197 }
 198
 199 Fraction& Fraction::operator *= ( const Fraction& rVal )
 200 {
 201     if ( !rVal.mbValid )
 202         mbValid = false;
 203
 204     if ( !mbValid )
 205     {
 206         SAL_WARN( "tools.fraction", "'operator *=' with invalid fraction" );
 207         return *this;
 208     }
 209
 210     boost::rational<sal_Int32> a = toRational(mnNumerator, mnDenominator);
 211     boost::rational<sal_Int32> b = toRational(rVal.mnNumerator, rVal.mnDenominator);
 212     bool bFail = checked_multiply_by(a, b);
 213     mnNumerator = a.numerator();
 214     mnDenominator = a.denominator();
 215
 216     if (bFail)
 217     {
 218         mbValid = false;
 219     }
 220
 221     return *this;
 222 }
 223
 224 Fraction& Fraction::operator /= ( const Fraction& rVal )
 225 {
 226     if ( !rVal.mbValid )
 227         mbValid = false;
 228
 229     if ( !mbValid )
 230     {
 231         SAL_WARN( "tools.fraction", "'operator /=' with invalid fraction" );
 232         return *this;
 233     }
 234
 235     boost::rational<sal_Int32> a = toRational(mnNumerator, mnDenominator);
 236     a /= toRational(rVal.mnNumerator, rVal.mnDenominator);
 237     mnNumerator = a.numerator();
 238     mnDenominator = a.denominator();
 239
 240     return *this;
 241 }
 242
 243 /** Inaccurate cancellation for a fraction.
 244
 245     Clip both nominator and denominator to said number of bits. If
 246     either of those already have equal or less number of bits used,
 247     this method does nothing.
 248
 249     @param nSignificantBits denotes, how many significant binary
 250     digits to maintain, in both nominator and denominator.
 251
 252     @example ReduceInaccurate(8) has an error <1% [1/2^(8-1)] - the
 253     largest error occurs with the following pair of values:
 254
 255     binary    1000000011111111111111111111111b/1000000000000000000000000000000b
 256     =         1082130431/1073741824
 257     = approx. 1.007812499
 258
 259     A ReduceInaccurate(8) yields 1/1.
 260 */
 261 void Fraction::ReduceInaccurate( unsigned nSignificantBits )
 262 {
 263     if ( !mbValid )
 264     {
 265         SAL_WARN( "tools.fraction", "'ReduceInaccurate' on invalid fraction" );
 266         return;
 267     }
 268
 269     if ( !mnNumerator )
 270         return;
 271
 272     auto a = toRational(mnNumerator, mnDenominator);
 273     rational_ReduceInaccurate(a, nSignificantBits);
 274     mnNumerator = a.numerator();
 275     mnDenominator = a.denominator();
 276 }
 277
 278 sal_Int32 Fraction::GetNumerator() const
 279 {
 280     if ( !mbValid )
 281     {
 282         SAL_WARN( "tools.fraction", "'GetNumerator()' on invalid fraction" );
 283         return 0;
 284     }
 285     return mnNumerator;
 286 }
 287
 288 sal_Int32 Fraction::GetDenominator() const
 289 {
 290     if ( !mbValid )
 291     {
 292         SAL_WARN( "tools.fraction", "'GetDenominator()' on invalid fraction" );
 293         return -1;
 294     }
 295     return mnDenominator;
 296 }
 297
 298 Fraction::operator sal_Int32() const
 299 {
 300     if ( !mbValid )
 301     {
 302         SAL_WARN( "tools.fraction", "'operator sal_Int32()' on invalid fraction" );
 303         return 0;
 304     }
 305     return boost::rational_cast<sal_Int32>(toRational(mnNumerator, mnDenominator));
 306 }
 307
 308 Fraction operator+( const Fraction& rVal1, const Fraction& rVal2 )
 309 {
 310     Fraction aErg( rVal1 );
 311     aErg += rVal2;
 312     return aErg;
 313 }
 314
 315 Fraction operator-( const Fraction& rVal1, const Fraction& rVal2 )
 316 {
 317     Fraction aErg( rVal1 );
 318     aErg -= rVal2;
 319     return aErg;
 320 }
 321
 322 Fraction operator*( const Fraction& rVal1, const Fraction& rVal2 )
 323 {
 324     Fraction aErg( rVal1 );
 325     aErg *= rVal2;
 326     return aErg;
 327 }
 328
 329 Fraction operator/( const Fraction& rVal1, const Fraction& rVal2 )
 330 {
 331     Fraction aErg( rVal1 );
 332     aErg /= rVal2;
 333     return aErg;
 334 }
 335
 336 bool operator !=( const Fraction& rVal1, const Fraction& rVal2 )
 337 {
 338     return !(rVal1 == rVal2);
 339 }
 340
 341 bool operator <=( const Fraction& rVal1, const Fraction& rVal2 )
 342 {
 343     return !(rVal1 > rVal2);
 344 }
 345
 346 bool operator >=( const Fraction& rVal1, const Fraction& rVal2 )
 347 {
 348     return !(rVal1 < rVal2);
 349 }
 350
 351 bool operator == ( const Fraction& rVal1, const Fraction& rVal2 )
 352 {
 353     if ( !rVal1.mbValid || !rVal2.mbValid )
 354     {
 355         SAL_WARN( "tools.fraction", "'operator ==' with an invalid fraction" );
 356         return false;
 357     }
 358
 359     return toRational(rVal1.mnNumerator, rVal1.mnDenominator) == toRational(rVal2.mnNumerator, rVal2.mnDenominator);
 360 }
 361
 362 bool operator < ( const Fraction& rVal1, const Fraction& rVal2 )
 363 {
 364     if ( !rVal1.mbValid || !rVal2.mbValid )
 365     {
 366         SAL_WARN( "tools.fraction", "'operator <' with an invalid fraction" );
 367         return false;
 368     }
 369
 370     return toRational(rVal1.mnNumerator, rVal1.mnDenominator) < toRational(rVal2.mnNumerator, rVal2.mnDenominator);
 371 }
 372
 373 bool operator > ( const Fraction& rVal1, const Fraction& rVal2 )
 374 {
 375     if ( !rVal1.mbValid || !rVal2.mbValid )
 376     {
 377         SAL_WARN( "tools.fraction", "'operator >' with an invalid fraction" );
 378         return false;
 379     }
 380
 381     return toRational(rVal1.mnNumerator, rVal1.mnDenominator) > toRational(rVal2.mnNumerator, rVal2.mnDenominator);
 382 }
 383
 384 // If dVal > LONG_MAX or dVal < LONG_MIN, the rational throws a boost::bad_rational.
 385 // Otherwise, dVal and denominator are multiplied by 10, until one of them
 386 // is larger than (LONG_MAX / 10).
 387 //
 388 // NOTE: here we use 'sal_Int32' due that only values in sal_Int32 range are valid.
 389 static boost::rational<sal_Int32> rational_FromDouble(double dVal)
 390 {
 391     if ( dVal > std::numeric_limits<sal_Int32>::max() ||
 392          dVal < std::numeric_limits<sal_Int32>::min() ||
 393          std::isnan(dVal) )
 394         throw boost::bad_rational();
 395
 396     const sal_Int32 nMAX = std::numeric_limits<sal_Int32>::max() / 10;
 397     sal_Int32 nDen = 1;
 398     while ( std::abs( dVal ) < nMAX && nDen < nMAX ) {
 399         dVal *= 10;
 400         nDen *= 10;
 401     }
 402     return boost::rational<sal_Int32>( sal_Int32(dVal), nDen );
 403 }
 404
 405 /**
 406  * Find the number of bits required to represent this number, using the CLZ intrinsic
 407  */
 408 static int impl_NumberOfBits( sal_uInt32 nNum )
 409 {
 410     if (nNum == 0)
 411         return 0;
 412 #ifdef _MSC_VER
 413     unsigned long r = 0;
 414     _BitScanReverse(&r, nNum);
 415     return r + 1;
 416 #else
 417     return 32 - __builtin_clz(nNum);
 418 #endif
 419 }
 420
 421 /** Inaccurate cancellation for a fraction.
 422
 423     Clip both nominator and denominator to said number of bits. If
 424     either of those already have equal or less number of bits used,
 425     this method does nothing.
 426
 427     @param nSignificantBits denotes, how many significant binary
 428     digits to maintain, in both nominator and denominator.
 429
 430     @example ReduceInaccurate(8) has an error <1% [1/2^(8-1)] - the
 431     largest error occurs with the following pair of values:
 432
 433     binary    1000000011111111111111111111111b/1000000000000000000000000000000b
 434     =         1082130431/1073741824
 435     = approx. 1.007812499
 436
 437     A ReduceInaccurate(8) yields 1/1.
 438 */
 439 static void rational_ReduceInaccurate(boost::rational<sal_Int32>& rRational, unsigned nSignificantBits)
 440 {
 441     if ( !rRational )
 442         return;
 443
 444     // http://www.boost.org/doc/libs/release/libs/rational/rational.html#Internal%20representation
 445     sal_Int32 nMul = rRational.numerator();
 446     if (nMul == std::numeric_limits<sal_Int32>::min())
 447     {
 448         // ofz#32973 Integer-overflow
 449         return;
 450     }
 451     const bool bNeg = nMul < 0;
 452     if (bNeg)
 453         nMul = -nMul;
 454     sal_Int32 nDiv = rRational.denominator();
 455
 456     DBG_ASSERT(nSignificantBits<65, "More than 64 bit of significance is overkill!");
 457
 458     // How much bits can we lose?
 459     const int nMulBitsToLose = std::max( ( impl_NumberOfBits( nMul ) - int( nSignificantBits ) ), 0 );
 460     const int nDivBitsToLose = std::max( ( impl_NumberOfBits( nDiv ) - int( nSignificantBits ) ), 0 );
 461
 462     const int nToLose = std::min( nMulBitsToLose, nDivBitsToLose );
 463
 464     // Remove the bits
 465     nMul >>= nToLose;
 466     nDiv >>= nToLose;
 467
 468     if ( !nMul || !nDiv ) {
 469         // Return without reduction
 470         OSL_FAIL( "Oops, we reduced too much..." );
 471         return;
 472     }
 473
 474     rRational.assign( bNeg ? -nMul : nMul, nDiv );
 475 }
 476
 477 size_t Fraction::GetHashValue() const
 478 {
 479     size_t hash = 0;
 480     o3tl::hash_combine( hash, mnNumerator );
 481     o3tl::hash_combine( hash, mnDenominator );
 482     o3tl::hash_combine( hash, mbValid );
 483     return hash;
 484 }
 485
 486 Fraction Fraction::MakeFraction( tools::Long nN1, tools::Long nN2, tools::Long nD1, tools::Long nD2 )
 487 {
 488     if( nD1 == 0 || nD2 == 0 ) //under these bad circumstances the following while loop will be endless
 489     {
 490         SAL_WARN("tools.fraction", "Invalid parameter for ImplMakeFraction");
 491         return Fraction( 1, 1 );
 492     }
 493
 494     tools::Long i = 1;
 495
 496     if ( nN1 < 0 ) { i = -i; nN1 = -nN1; }
 497     if ( nN2 < 0 ) { i = -i; nN2 = -nN2; }
 498     if ( nD1 < 0 ) { i = -i; nD1 = -nD1; }
 499     if ( nD2 < 0 ) { i = -i; nD2 = -nD2; }
 500     // all positive; i sign
 501
 502     assert( nN1 >= std::numeric_limits<sal_Int32>::min() );
 503     assert( nN1 <= std::numeric_limits<sal_Int32>::max( ));
 504     assert( nD1 >= std::numeric_limits<sal_Int32>::min() );
 505     assert( nD1 <= std::numeric_limits<sal_Int32>::max( ));
 506     assert( nN2 >= std::numeric_limits<sal_Int32>::min() );
 507     assert( nN2 <= std::numeric_limits<sal_Int32>::max( ));
 508     assert( nD2 >= std::numeric_limits<sal_Int32>::min() );
 509     assert( nD2 <= std::numeric_limits<sal_Int32>::max( ));
 510
 511     boost::rational<sal_Int32> a = toRational(i*nN1, nD1);
 512     boost::rational<sal_Int32> b = toRational(nN2, nD2);
 513     bool bFail = checked_multiply_by(a, b);
 514
 515     while ( bFail ) {
 516         if ( nN1 > nN2 )
 517             nN1 = (nN1 + 1) / 2;
 518         else
 519             nN2 = (nN2 + 1) / 2;
 520         if ( nD1 > nD2 )
 521             nD1 = (nD1 + 1) / 2;
 522         else
 523             nD2 = (nD2 + 1) / 2;
 524
 525         a = toRational(i*nN1, nD1);
 526         b = toRational(nN2, nD2);
 527         bFail = checked_multiply_by(a, b);
 528     }
 529
 530     return Fraction(a.numerator(), a.denominator());
 531 }
 532
 533 /* vim:set shiftwidth=4 softtabstop=4 expandtab: */