contrib/ntp/tests/libntp/g_lfpfunc.cpp

   1 #include "g_libntptest.h"
   2 #include "g_timestructs.h"
   3
   4 extern "C" {
   5 #include "ntp_fp.h"
   6 }
   7
   8 #include <float.h>
   9 #include <math.h>
  10
  11 #include <string>
  12 #include <sstream>
  13
  14 class lfpTest : public libntptest
  15 {
  16         // nothing new right now
  17 };
  18
  19 struct lfp_hl {
  20         uint32_t h, l;
  21 };
  22
  23 //----------------------------------------------------------------------
  24 // OO-wrapper for 'l_fp'
  25 //----------------------------------------------------------------------
  26
  27 class LFP
  28 {
  29 public:
  30         ~LFP();
  31         LFP();
  32         LFP(const LFP& rhs);
  33         LFP(int32 i, u_int32 f);
  34
  35         LFP  operator+ (const LFP &rhs) const;
  36         LFP& operator+=(const LFP &rhs);
  37
  38         LFP  operator- (const LFP &rhs) const;
  39         LFP& operator-=(const LFP &rhs);
  40
  41         LFP& operator=(const LFP &rhs);
  42         LFP  operator-() const;
  43
  44         bool operator==(const LFP &rhs) const;
  45
  46         LFP  neg() const;
  47         LFP  abs() const;
  48         int  signum() const;
  49
  50         bool l_isgt (const LFP &rhs) const
  51                 { return L_ISGT(&_v, &rhs._v); }
  52         bool l_isgtu(const LFP &rhs) const
  53                 { return L_ISGTU(&_v, &rhs._v); }
  54         bool l_ishis(const LFP &rhs) const
  55                 { return L_ISHIS(&_v, &rhs._v); }
  56         bool l_isgeq(const LFP &rhs) const
  57                 { return L_ISGEQ(&_v, &rhs._v); }
  58         bool l_isequ(const LFP &rhs) const
  59                 { return L_ISEQU(&_v, &rhs._v); }
  60
  61         int  ucmp(const LFP & rhs) const;
  62         int  scmp(const LFP & rhs) const;
  63
  64         std::string   toString() const;
  65         std::ostream& toStream(std::ostream &oo) const;
  66
  67         operator double() const;
  68         explicit LFP(double);
  69
  70 protected:
  71         LFP(const l_fp &rhs);
  72
  73         static int cmp_work(u_int32 a[3], u_int32 b[3]);
  74
  75         l_fp _v;
  76 };
  77
  78 static std::ostream& operator<<(std::ostream &oo, const LFP& rhs)
  79 {
  80         return rhs.toStream(oo);
  81 }
  82
  83 //----------------------------------------------------------------------
  84 // reference comparision
  85 // This is implementad as a full signed MP-subtract in 3 limbs, where
  86 // the operands are zero or sign extended before the subtraction is
  87 // executed.
  88 //----------------------------------------------------------------------
  89 int  LFP::scmp(const LFP & rhs) const
  90 {
  91         u_int32 a[3], b[3];
  92         const l_fp &op1(_v), &op2(rhs._v);
  93
  94         a[0] = op1.l_uf; a[1] = op1.l_ui; a[2] = 0;
  95         b[0] = op2.l_uf; b[1] = op2.l_ui; b[2] = 0;
  96
  97         a[2] -= (op1.l_i < 0);
  98         b[2] -= (op2.l_i < 0);
  99
 100         return cmp_work(a,b);
 101 }
 102
 103 int  LFP::ucmp(const LFP & rhs) const
 104 {
 105         u_int32 a[3], b[3];
 106         const l_fp &op1(_v), &op2(rhs._v);
 107
 108         a[0] = op1.l_uf; a[1] = op1.l_ui; a[2] = 0;
 109         b[0] = op2.l_uf; b[1] = op2.l_ui; b[2] = 0;
 110
 111         return cmp_work(a,b);
 112 }
 113
 114 int LFP::cmp_work(u_int32 a[3], u_int32 b[3])
 115 {
 116         u_int32 cy, idx, tmp;
 117         for (cy = idx = 0; idx < 3; ++idx) {
 118                 tmp = a[idx]; cy  = (a[idx] -=   cy  ) > tmp;
 119                 tmp = a[idx]; cy |= (a[idx] -= b[idx]) > tmp;
 120         }
 121         if (a[2])
 122                 return -1;
 123         return a[0] || a[1];
 124 }
 125
 126 //----------------------------------------------------------------------
 127 // imlementation of the LFP stuff
 128 // This should be easy enough...
 129 //----------------------------------------------------------------------
 130
 131 LFP::~LFP()
 132 {
 133         // NOP
 134 }
 135
 136 LFP::LFP()
 137 {
 138         _v.l_ui = 0;
 139         _v.l_uf = 0;
 140 }
 141
 142 LFP::LFP(int32 i, u_int32 f)
 143 {
 144         _v.l_i  = i;
 145         _v.l_uf = f;
 146 }
 147
 148 LFP::LFP(const LFP &rhs)
 149 {
 150         _v = rhs._v;
 151 }
 152
 153 LFP::LFP(const l_fp & rhs)
 154 {
 155         _v = rhs;
 156 }
 157
 158 LFP& LFP::operator=(const LFP & rhs)
 159 {
 160         _v = rhs._v;
 161         return *this;
 162 }
 163
 164 LFP& LFP::operator+=(const LFP & rhs)
 165 {
 166         L_ADD(&_v, &rhs._v);
 167         return *this;
 168 }
 169
 170 LFP& LFP::operator-=(const LFP & rhs)
 171 {
 172         L_SUB(&_v, &rhs._v);
 173         return *this;
 174 }
 175
 176 LFP LFP::operator+(const LFP &rhs) const
 177 {
 178         LFP tmp(*this);
 179         return tmp += rhs;
 180 }
 181
 182 LFP LFP::operator-(const LFP &rhs) const
 183 {
 184         LFP tmp(*this);
 185         return tmp -= rhs;
 186 }
 187
 188 LFP LFP::operator-() const
 189 {
 190         LFP tmp(*this);
 191         L_NEG(&tmp._v);
 192         return tmp;
 193 }
 194
 195 LFP
 196 LFP::neg() const
 197 {
 198         LFP tmp(*this);
 199         L_NEG(&tmp._v);
 200         return tmp;
 201 }
 202
 203 LFP
 204 LFP::abs() const
 205 {
 206         LFP tmp(*this);
 207         if (L_ISNEG(&tmp._v))
 208                 L_NEG(&tmp._v);
 209         return tmp;
 210 }
 211
 212 int
 213 LFP::signum() const
 214 {
 215         if (_v.l_ui & 0x80000000u)
 216                 return -1;
 217         return (_v.l_ui || _v.l_uf);
 218 }
 219
 220 std::string
 221 LFP::toString() const
 222 {
 223         std::ostringstream oss;
 224         toStream(oss);
 225         return oss.str();
 226 }
 227
 228 std::ostream&
 229 LFP::toStream(std::ostream &os) const
 230 {
 231         return os
 232             << mfptoa(_v.l_ui, _v.l_uf, 9)
 233             << " [$" << std::setw(8) << std::setfill('0') << std::hex << _v.l_ui
 234             <<  ':'  << std::setw(8) << std::setfill('0') << std::hex << _v.l_uf
 235             << ']';
 236 }
 237
 238 bool LFP::operator==(const LFP &rhs) const
 239 {
 240         return L_ISEQU(&_v, &rhs._v);
 241 }
 242
 243
 244 LFP::operator double() const
 245 {
 246         double res;
 247         LFPTOD(&_v, res);
 248         return res;
 249 }
 250
 251 LFP::LFP(double rhs)
 252 {
 253         DTOLFP(rhs, &_v);
 254 }
 255
 256
 257 //----------------------------------------------------------------------
 258 // testing the relational macros works better with proper predicate
 259 // formatting functions; it slows down the tests a bit, but makes for
 260 // readable failure messages.
 261 //----------------------------------------------------------------------
 262
 263 testing::AssertionResult isgt_p(
 264         const LFP &op1, const LFP &op2)
 265 {
 266         if (op1.l_isgt(op2))
 267                 return testing::AssertionSuccess()
 268                     << "L_ISGT(" << op1 << "," << op2 << ") is true";
 269         else
 270                 return testing::AssertionFailure()
 271                     << "L_ISGT(" << op1 << "," << op2 << ") is false";
 272 }
 273
 274 testing::AssertionResult isgeq_p(
 275         const LFP &op1, const LFP &op2)
 276 {
 277         if (op1.l_isgeq(op2))
 278                 return testing::AssertionSuccess()
 279                     << "L_ISGEQ(" << op1 << "," << op2 << ") is true";
 280         else
 281                 return testing::AssertionFailure()
 282                     << "L_ISGEQ(" << op1 << "," << op2 << ") is false";
 283 }
 284
 285 testing::AssertionResult isgtu_p(
 286         const LFP &op1, const LFP &op2)
 287 {
 288         if (op1.l_isgtu(op2))
 289                 return testing::AssertionSuccess()
 290                     << "L_ISGTU(" << op1 << "," << op2 << ") is true";
 291         else
 292                 return testing::AssertionFailure()
 293                     << "L_ISGTU(" << op1 << "," << op2 << ") is false";
 294 }
 295
 296 testing::AssertionResult ishis_p(
 297         const LFP &op1, const LFP &op2)
 298 {
 299         if (op1.l_ishis(op2))
 300                 return testing::AssertionSuccess()
 301                     << "L_ISHIS(" << op1 << "," << op2 << ") is true";
 302         else
 303                 return testing::AssertionFailure()
 304                     << "L_ISHIS(" << op1 << "," << op2 << ") is false";
 305 }
 306
 307 testing::AssertionResult isequ_p(
 308         const LFP &op1, const LFP &op2)
 309 {
 310         if (op1.l_isequ(op2))
 311                 return testing::AssertionSuccess()
 312                     << "L_ISEQU(" << op1 << "," << op2 << ") is true";
 313         else
 314                 return testing::AssertionFailure()
 315                     << "L_ISEQU(" << op1 << "," << op2 << ") is false";
 316 }
 317
 318 //----------------------------------------------------------------------
 319 // test data table for add/sub and compare
 320 //----------------------------------------------------------------------
 321
 322 static const lfp_hl addsub_tab[][3] = {
 323         // trivial idendity:
 324         {{0 ,0         }, { 0,0         }, { 0,0}},
 325         // with carry from fraction and sign change:
 326         {{-1,0x80000000}, { 0,0x80000000}, { 0,0}},
 327         // without carry from fraction
 328         {{ 1,0x40000000}, { 1,0x40000000}, { 2,0x80000000}},
 329         // with carry from fraction:
 330         {{ 1,0xC0000000}, { 1,0xC0000000}, { 3,0x80000000}},
 331         // with carry from fraction and sign change:
 332         {{0x7FFFFFFF, 0x7FFFFFFF}, {0x7FFFFFFF,0x7FFFFFFF}, {0xFFFFFFFE,0xFFFFFFFE}},
 333         // two tests w/o carry (used for l_fp<-->double):
 334         {{0x55555555,0xAAAAAAAA}, {0x11111111,0x11111111}, {0x66666666,0xBBBBBBBB}},
 335         {{0x55555555,0x55555555}, {0x11111111,0x11111111}, {0x66666666,0x66666666}},
 336         // wide-range test, triggers compare trouble
 337         {{0x80000000,0x00000001}, {0xFFFFFFFF,0xFFFFFFFE}, {0x7FFFFFFF,0xFFFFFFFF}}
 338 };
 339 static const size_t addsub_cnt(sizeof(addsub_tab)/sizeof(addsub_tab[0]));
 340 static const size_t addsub_tot(sizeof(addsub_tab)/sizeof(addsub_tab[0][0]));
 341
 342
 343 //----------------------------------------------------------------------
 344 // epsilon estimation for the precision of a conversion double --> l_fp
 345 //
 346 // The error estimation limit is as follows:
 347 //  * The 'l_fp' fixed point fraction has 32 bits precision, so we allow
 348 //    for the LSB to toggle by clamping the epsilon to be at least 2^(-31)
 349 //
 350 //  * The double mantissa has a precsion 54 bits, so the other minimum is
 351 //    dval * (2^(-53))
 352 //
 353 //  The maximum of those two boundaries is used for the check.
 354 //
 355 // Note: once there are more than 54 bits between the highest and lowest
 356 // '1'-bit of the l_fp value, the roundtrip *will* create truncation
 357 // errors. This is an inherent property caused by the 54-bit mantissa of
 358 // the 'double' type.
 359 double eps(double d)
 360 {
 361         return std::max<double>(ldexp(1.0, -31), ldexp(fabs(d), -53));
 362 }
 363
 364 //----------------------------------------------------------------------
 365 // test addition
 366 //----------------------------------------------------------------------
 367 TEST_F(lfpTest, AdditionLR) {
 368         for (size_t idx=0; idx < addsub_cnt; ++idx) {
 369                 LFP op1(addsub_tab[idx][0].h, addsub_tab[idx][0].l);
 370                 LFP op2(addsub_tab[idx][1].h, addsub_tab[idx][1].l);
 371                 LFP exp(addsub_tab[idx][2].h, addsub_tab[idx][2].l);
 372                 LFP res(op1 + op2);
 373
 374                 ASSERT_EQ(exp, res);
 375         }
 376 }
 377
 378 TEST_F(lfpTest, AdditionRL) {
 379         for (size_t idx=0; idx < addsub_cnt; ++idx) {
 380                 LFP op2(addsub_tab[idx][0].h, addsub_tab[idx][0].l);
 381                 LFP op1(addsub_tab[idx][1].h, addsub_tab[idx][1].l);
 382                 LFP exp(addsub_tab[idx][2].h, addsub_tab[idx][2].l);
 383                 LFP res(op1 + op2);
 384
 385                 ASSERT_EQ(exp, res);
 386         }
 387 }
 388
 389 //----------------------------------------------------------------------
 390 // test subtraction
 391 //----------------------------------------------------------------------
 392 TEST_F(lfpTest, SubtractionLR) {
 393         for (size_t idx=0; idx < addsub_cnt; ++idx) {
 394                 LFP op2(addsub_tab[idx][0].h, addsub_tab[idx][0].l);
 395                 LFP exp(addsub_tab[idx][1].h, addsub_tab[idx][1].l);
 396                 LFP op1(addsub_tab[idx][2].h, addsub_tab[idx][2].l);
 397                 LFP res(op1 - op2);
 398
 399                 ASSERT_EQ(exp, res);
 400         }
 401 }
 402
 403 TEST_F(lfpTest, SubtractionRL) {
 404         for (size_t idx=0; idx < addsub_cnt; ++idx) {
 405                 LFP exp(addsub_tab[idx][0].h, addsub_tab[idx][0].l);
 406                 LFP op2(addsub_tab[idx][1].h, addsub_tab[idx][1].l);
 407                 LFP op1(addsub_tab[idx][2].h, addsub_tab[idx][2].l);
 408                 LFP res(op1 - op2);
 409
 410                 ASSERT_EQ(exp, res);
 411         }
 412 }
 413
 414 //----------------------------------------------------------------------
 415 // test negation
 416 //----------------------------------------------------------------------
 417 TEST_F(lfpTest, Negation) {
 418         for (size_t idx=0; idx < addsub_cnt; ++idx) {
 419                 LFP op1(addsub_tab[idx][0].h, addsub_tab[idx][0].l);
 420                 LFP op2(-op1);
 421                 LFP sum(op1 + op2);
 422
 423                 ASSERT_EQ(LFP(0,0), sum);
 424         }
 425 }
 426
 427 //----------------------------------------------------------------------
 428 // test absolute value
 429 //----------------------------------------------------------------------
 430 TEST_F(lfpTest, Absolute) {
 431         for (size_t idx=0; idx < addsub_cnt; ++idx) {
 432                 LFP op1(addsub_tab[idx][0].h, addsub_tab[idx][0].l);
 433                 LFP op2(op1.abs());
 434
 435                 ASSERT_TRUE(op2.signum() >= 0);
 436
 437                 if (op1.signum() >= 0)
 438                         op1 -= op2;
 439                 else
 440                         op1 += op2;
 441                 ASSERT_EQ(LFP(0,0), op1);
 442         }
 443
 444         // There is one special case we have to check: the minimum
 445         // value cannot be negated, or, to be more precise, the
 446         // negation reproduces the original pattern.
 447         LFP minVal(0x80000000, 0x00000000);
 448         LFP minAbs(minVal.abs());
 449         ASSERT_EQ(-1, minVal.signum());
 450         ASSERT_EQ(minVal, minAbs);
 451 }
 452
 453 //----------------------------------------------------------------------
 454 // fp -> double -> fp rountrip test
 455 //----------------------------------------------------------------------
 456 TEST_F(lfpTest, FDF_RoundTrip) {
 457         // since a l_fp has 64 bits in it's mantissa and a double has
 458         // only 54 bits available (including the hidden '1') we have to
 459         // make a few concessions on the roundtrip precision. The 'eps()'
 460         // function makes an educated guess about the avilable precision
 461         // and checks the difference in the two 'l_fp' values against
 462         // that limit.
 463         for (size_t idx=0; idx < addsub_cnt; ++idx) {
 464                 LFP    op1(addsub_tab[idx][0].h, addsub_tab[idx][0].l);
 465                 double op2(op1);
 466                 LFP    op3(op2);
 467                 // for manual checks only:
 468                 // std::cout << std::setprecision(16) << op2 << std::endl;
 469                 ASSERT_LE(fabs(op1-op3), eps(op2));
 470         }
 471 }
 472
 473 //----------------------------------------------------------------------
 474 // test the compare stuff
 475 //
 476 // This uses the local compare and checks if the operations using the
 477 // macros in 'ntp_fp.h' produce mathing results.
 478 // ----------------------------------------------------------------------
 479 TEST_F(lfpTest, SignedRelOps) {
 480         const lfp_hl * tv(&addsub_tab[0][0]);
 481         for (size_t lc=addsub_tot-1; lc; --lc,++tv) {
 482                 LFP op1(tv[0].h,tv[0].l);
 483                 LFP op2(tv[1].h,tv[1].l);
 484                 int cmp(op1.scmp(op2));
 485
 486                 switch (cmp) {
 487                 case -1:
 488                         std::swap(op1, op2);
 489                 case 1:
 490                         EXPECT_TRUE (isgt_p(op1,op2));
 491                         EXPECT_FALSE(isgt_p(op2,op1));
 492
 493                         EXPECT_TRUE (isgeq_p(op1,op2));
 494                         EXPECT_FALSE(isgeq_p(op2,op1));
 495
 496                         EXPECT_FALSE(isequ_p(op1,op2));
 497                         EXPECT_FALSE(isequ_p(op2,op1));
 498                         break;
 499                 case 0:
 500                         EXPECT_FALSE(isgt_p(op1,op2));
 501                         EXPECT_FALSE(isgt_p(op2,op1));
 502
 503                         EXPECT_TRUE (isgeq_p(op1,op2));
 504                         EXPECT_TRUE (isgeq_p(op2,op1));
 505
 506                         EXPECT_TRUE (isequ_p(op1,op2));
 507                         EXPECT_TRUE (isequ_p(op2,op1));
 508                         break;
 509                 default:
 510                         FAIL() << "unexpected SCMP result: " << cmp;
 511                 }
 512         }
 513 }
 514
 515 TEST_F(lfpTest, UnsignedRelOps) {
 516         const lfp_hl * tv(&addsub_tab[0][0]);
 517         for (size_t lc=addsub_tot-1; lc; --lc,++tv) {
 518                 LFP op1(tv[0].h,tv[0].l);
 519                 LFP op2(tv[1].h,tv[1].l);
 520                 int cmp(op1.ucmp(op2));
 521
 522                 switch (cmp) {
 523                 case -1:
 524                         std::swap(op1, op2);
 525                 case 1:
 526                         EXPECT_TRUE (isgtu_p(op1,op2));
 527                         EXPECT_FALSE(isgtu_p(op2,op1));
 528
 529                         EXPECT_TRUE (ishis_p(op1,op2));
 530                         EXPECT_FALSE(ishis_p(op2,op1));
 531                         break;
 532                 case 0:
 533                         EXPECT_FALSE(isgtu_p(op1,op2));
 534                         EXPECT_FALSE(isgtu_p(op2,op1));
 535
 536                         EXPECT_TRUE (ishis_p(op1,op2));
 537                         EXPECT_TRUE (ishis_p(op2,op1));
 538                         break;
 539                 default:
 540                         FAIL() << "unexpected UCMP result: " << cmp;
 541                 }
 542         }
 543 }
 544
 545 //----------------------------------------------------------------------
 546 // that's all folks... but feel free to add things!
 547 //----------------------------------------------------------------------