A Discrete-Event Network Simulator Home Tutorials  ▼ Docs    ▼ Develop ▼ API
int64x64-128.cc
Go to the documentation of this file.
1 #include "int64x64-128.h"
2 #include "abort.h"
3 #include "assert.h"
4 #include "log.h"
5
6 // Note: Logging in this file is largely avoided due to the
7 // number of calls that are made to these functions and the possibility
8 // of causing recursions leading to stack overflow
9
10 NS_LOG_COMPONENT_DEFINE ("int64x64-128");
11
12 namespace ns3 {
13
14 static inline
15 bool
17  const int128_t sb,
18  uint128_t & ua,
19  uint128_t & ub)
20 {
21  bool negA = sa < 0;
22  bool negB = sb < 0;
23  ua = negA ? -sa : sa;
24  ub = negB ? -sb : sb;
25  return (negA && !negB) || (!negA && negB);
26 }
27
28 void
30 {
31  uint128_t a, b;
32  bool negative = output_sign (_v, o._v, a, b);
33  uint128_t result = Umul (a, b);
34  _v = negative ? -result : result;
35 }
36
39 {
40  uint128_t aL = a & HP_MASK_LO;
41  uint128_t bL = b & HP_MASK_LO;
42  uint128_t aH = (a >> 64) & HP_MASK_LO;
43  uint128_t bH = (b >> 64) & HP_MASK_LO;
44
45  uint128_t result;
46  uint128_t hiPart, loPart, midPart;
47  uint128_t res1, res2;
48
49  // Multiplying (a.h 2^64 + a.l) x (b.h 2^64 + b.l) =
50  // 2^128 a.h b.h + 2^64*(a.h b.l+b.h a.l) + a.l b.l
51  // get the low part a.l b.l
52  // multiply the fractional part
53  loPart = aL * bL;
54  // compute the middle part 2^64*(a.h b.l+b.h a.l)
55  midPart = aL * bH + aH * bL;
56  // compute the high part 2^128 a.h b.h
57  hiPart = aH * bH;
58  // if the high part is not zero, put a warning
59  NS_ABORT_MSG_IF ((hiPart & HP_MASK_HI) != 0,
60  "High precision 128 bits multiplication error: multiplication overflow.");
61
62  // Adding 64-bit terms to get 128-bit results, with carries
63  res1 = loPart >> 64;
64  res2 = midPart & HP_MASK_LO;
65  result = res1 + res2;
66
67  res1 = midPart >> 64;
68  res2 = hiPart & HP_MASK_LO;
69  res1 += res2;
70  res1 <<= 64;
71
72  result += res1;
73
74  return result;
75 }
76
77 void
79 {
80  uint128_t a, b;
81  bool negative = output_sign (_v, o._v, a, b);
82  int128_t result = Udiv (a, b);
83  _v = negative ? -result : result;
84 }
85
88 {
89
90  uint128_t rem = a;
91  uint128_t den = b;
92  uint128_t quo = rem / den;
93  rem = rem % den;
94  uint128_t result = quo;
95
96  // Now, manage the remainder
97  const uint64_t DIGITS = 64; // Number of fraction digits (bits) we need
98  const uint128_t ZERO = 0;
99
100  NS_ASSERT_MSG (rem < den,
101  "Remainder not less than divisor");
102
103  uint64_t digis = 0; // Number of digits we have already
104  uint64_t shift = 0; // Number we are going to get this round
105
106  // Skip trailing zeros in divisor
107  while ( (shift < DIGITS) && !(den & 0x1))
108  {
109  ++shift;
110  den >>= 1;
111  }
112
113  while ( (digis < DIGITS) && (rem != ZERO) )
114  {
115  // Skip leading zeros in remainder
116  while ( (digis + shift < DIGITS) &&
118  {
119  ++shift;
120  rem <<= 1;
121  }
122
123  // Cast off denominator bits if:
124  // Need more digits and
125  // LSB is zero or
126  // rem < den
127  while ( (digis + shift < DIGITS) &&
128  ( !(den & 0x1) || (rem < den) ) )
129  {
130  ++shift;
131  den >>= 1;
132  }
133
134  // Do the division
135  quo = rem / den;
136  rem = rem % den;
137
138  // Add in the quotient as shift bits of the fraction
139  result <<= shift;
140  result += quo;
141
142  digis += shift;
143  shift = 0;
144  }
145  // Did we run out of remainder?
146  if (digis < DIGITS)
147  {
148  shift = DIGITS - digis;
149  result <<= shift;
150  }
151
152  return result;
153 }
154
155 void
157 {
158  bool negResult = _v < 0;
159  uint128_t a = negResult ? -_v : _v;
160  uint128_t result = UmulByInvert (a, o._v);
161
162  _v = negResult ? -result : result;
163 }
164
165 uint128_t
167 {
168  uint128_t result, ah, bh, al, bl;
169  uint128_t hi, mid;
170  ah = a >> 64;
171  bh = b >> 64;
172  al = a & HP_MASK_LO;
173  bl = b & HP_MASK_LO;
174  hi = ah * bh;
175  mid = ah * bl + al * bh;
176  mid >>= 64;
177  result = hi + mid;
178  return result;
179 }
180
181 int64x64_t
182 int64x64_t::Invert (const uint64_t v)
183 {
184  NS_ASSERT (v > 1);
185  uint128_t a;
186  a = 1;
187  a <<= 64;
188  int64x64_t result;
189  result._v = Udiv (a, v);
190  int64x64_t tmp = int64x64_t (v, false);
191  tmp.MulByInvert (result);
192  if (tmp.GetHigh () != 1)
193  {
194  result._v += 1;
195  }
196  return result;
197 }
198
199 } // namespace ns3
200
High precision numerical type, implementing Q64.64 fixed precision.
Definition: int64x64-128.h:20
static uint128_t Udiv(const uint128_t a, const uint128_t b)
Unsigned division of Q64.64 values.
Definition: int64x64-128.cc:87
static int64x64_t Invert(const uint64_t v)
Compute the inverse of an integer value.
#define NS_ASSERT(condition)
At runtime, in debugging builds, if this condition is not true, the program prints the source file...
Definition: assert.h:61
#define NS_LOG_COMPONENT_DEFINE(name)
Define a Log component with a specific name.
Definition: log.h:170
__int128_t int128_t
Definition: int64x64-128.h:10
void MulByInvert(const int64x64_t &o)
Multiply this value by a Q0.128 value, presumably representing an inverse, completing a division oper...
static uint128_t UmulByInvert(const uint128_t a, const uint128_t b)
Unsigned multiplication of Q64.64 and Q0.128 values.
__uint128_t uint128_t
Definition: int64x64-128.h:9
static uint128_t Umul(const uint128_t a, const uint128_t b)
Unsigned multiplication of Q64.64 values.
Definition: int64x64-128.cc:38
int128_t _v
The Q64.64 value.
Definition: int64x64-128.h:310
void Mul(const int64x64_t &o)
Implement *=.
Definition: int64x64-128.cc:29
int64x64_t()
Default constructor.
Definition: int64x64-128.h:60
static bool output_sign(const int128_t sa, const int128_t sb, uint128_t &ua, uint128_t &ub)
Definition: int64x64-128.cc:16
uint128_t high bit (sign bit)
Definition: int64x64-128.h:23
#define NS_ASSERT_MSG(condition, message)
At runtime, in debugging builds, if this condition is not true, the program prints the message to out...
Definition: assert.h:84