Branch data Line data Source code
# 1 : : // Copyright (c) 2017-2020 The Bitcoin Core developers
# 2 : : // Distributed under the MIT software license, see the accompanying
# 3 : : // file COPYING or http://www.opensource.org/licenses/mit-license.php.
# 4 : :
# 5 : : #include <crypto/muhash.h>
# 6 : :
# 7 : : #include <crypto/chacha20.h>
# 8 : : #include <crypto/common.h>
# 9 : : #include <hash.h>
# 10 : :
# 11 : : #include <cassert>
# 12 : : #include <cstdio>
# 13 : : #include <limits>
# 14 : :
# 15 : : namespace {
# 16 : :
# 17 : : using limb_t = Num3072::limb_t;
# 18 : : using double_limb_t = Num3072::double_limb_t;
# 19 : : constexpr int LIMB_SIZE = Num3072::LIMB_SIZE;
# 20 : : /** 2^3072 - 1103717, the largest 3072-bit safe prime number, is used as the modulus. */
# 21 : : constexpr limb_t MAX_PRIME_DIFF = 1103717;
# 22 : :
# 23 : : /** Extract the lowest limb of [c0,c1,c2] into n, and left shift the number by 1 limb. */
# 24 : : inline void extract3(limb_t& c0, limb_t& c1, limb_t& c2, limb_t& n)
# 25 : 119479152 : {
# 26 : 119479152 : n = c0;
# 27 : 119479152 : c0 = c1;
# 28 : 119479152 : c1 = c2;
# 29 : 119479152 : c2 = 0;
# 30 : 119479152 : }
# 31 : :
# 32 : : /** [c0,c1] = a * b */
# 33 : : inline void mul(limb_t& c0, limb_t& c1, const limb_t& a, const limb_t& b)
# 34 : 1087392 : {
# 35 : 1087392 : double_limb_t t = (double_limb_t)a * b;
# 36 : 1087392 : c1 = t >> LIMB_SIZE;
# 37 : 1087392 : c0 = t;
# 38 : 1087392 : }
# 39 : :
# 40 : : /* [c0,c1,c2] += n * [d0,d1,d2]. c2 is 0 initially */
# 41 : : inline void mulnadd3(limb_t& c0, limb_t& c1, limb_t& c2, limb_t& d0, limb_t& d1, limb_t& d2, const limb_t& n)
# 42 : 116990003 : {
# 43 : 116990003 : double_limb_t t = (double_limb_t)d0 * n + c0;
# 44 : 116990003 : c0 = t;
# 45 : 116990003 : t >>= LIMB_SIZE;
# 46 : 116990003 : t += (double_limb_t)d1 * n + c1;
# 47 : 116990003 : c1 = t;
# 48 : 116990003 : t >>= LIMB_SIZE;
# 49 : 116990003 : c2 = t + d2 * n;
# 50 : 116990003 : }
# 51 : :
# 52 : : /* [c0,c1] *= n */
# 53 : : inline void muln2(limb_t& c0, limb_t& c1, const limb_t& n)
# 54 : 2489149 : {
# 55 : 2489149 : double_limb_t t = (double_limb_t)c0 * n;
# 56 : 2489149 : c0 = t;
# 57 : 2489149 : t >>= LIMB_SIZE;
# 58 : 2489149 : t += (double_limb_t)c1 * n;
# 59 : 2489149 : c1 = t;
# 60 : 2489149 : }
# 61 : :
# 62 : : /** [c0,c1,c2] += a * b */
# 63 : : inline void muladd3(limb_t& c0, limb_t& c1, limb_t& c2, const limb_t& a, const limb_t& b)
# 64 : 170586576 : {
# 65 : 170586576 : double_limb_t t = (double_limb_t)a * b;
# 66 : 170586576 : limb_t th = t >> LIMB_SIZE;
# 67 : 170586576 : limb_t tl = t;
# 68 : :
# 69 : 170586576 : c0 += tl;
# 70 [ + + ]: 170586576 : th += (c0 < tl) ? 1 : 0;
# 71 : 170586576 : c1 += th;
# 72 [ + + ]: 170586576 : c2 += (c1 < th) ? 1 : 0;
# 73 : 170586576 : }
# 74 : :
# 75 : : /** [c0,c1,c2] += 2 * a * b */
# 76 : : inline void muldbladd3(limb_t& c0, limb_t& c1, limb_t& c2, const limb_t& a, const limb_t& b)
# 77 : 2781662664 : {
# 78 : 2781662664 : double_limb_t t = (double_limb_t)a * b;
# 79 : 2781662664 : limb_t th = t >> LIMB_SIZE;
# 80 : 2781662664 : limb_t tl = t;
# 81 : :
# 82 : 2781662664 : c0 += tl;
# 83 [ + + ]: 2781662664 : limb_t tt = th + ((c0 < tl) ? 1 : 0);
# 84 : 2781662664 : c1 += tt;
# 85 [ + + ]: 2781662664 : c2 += (c1 < tt) ? 1 : 0;
# 86 : 2781662664 : c0 += tl;
# 87 [ + + ]: 2781662664 : th += (c0 < tl) ? 1 : 0;
# 88 : 2781662664 : c1 += th;
# 89 [ + + ]: 2781662664 : c2 += (c1 < th) ? 1 : 0;
# 90 : 2781662664 : }
# 91 : :
# 92 : : /**
# 93 : : * Add limb a to [c0,c1]: [c0,c1] += a. Then extract the lowest
# 94 : : * limb of [c0,c1] into n, and left shift the number by 1 limb.
# 95 : : * */
# 96 : : inline void addnextract2(limb_t& c0, limb_t& c1, const limb_t& a, limb_t& n)
# 97 : 119480208 : {
# 98 : 119480208 : limb_t c2 = 0;
# 99 : :
# 100 : : // add
# 101 : 119480208 : c0 += a;
# 102 [ + + ]: 119480208 : if (c0 < a) {
# 103 : 200412 : c1 += 1;
# 104 : :
# 105 : : // Handle case when c1 has overflown
# 106 [ - + ]: 200412 : if (c1 == 0)
# 107 : 0 : c2 = 1;
# 108 : 200412 : }
# 109 : :
# 110 : : // extract
# 111 : 119480208 : n = c0;
# 112 : 119480208 : c0 = c1;
# 113 : 119480208 : c1 = c2;
# 114 : 119480208 : }
# 115 : :
# 116 : : /** in_out = in_out^(2^sq) * mul */
# 117 : : inline void square_n_mul(Num3072& in_out, const int sq, const Num3072& mul)
# 118 : 11242 : {
# 119 [ + + ]: 833514 : for (int j = 0; j < sq; ++j) in_out.Square();
# 120 : 11242 : in_out.Multiply(mul);
# 121 : 11242 : }
# 122 : :
# 123 : : } // namespace
# 124 : :
# 125 : : /** Indicates whether d is larger than the modulus. */
# 126 : : bool Num3072::IsOverflow() const
# 127 : 2491558 : {
# 128 [ + + ]: 2491558 : if (this->limbs[0] <= std::numeric_limits<limb_t>::max() - MAX_PRIME_DIFF) return false;
# 129 [ + + ]: 1056 : for (int i = 1; i < LIMBS; ++i) {
# 130 [ - + ]: 1034 : if (this->limbs[i] != std::numeric_limits<limb_t>::max()) return false;
# 131 : 1034 : }
# 132 : 22 : return true;
# 133 : 22 : }
# 134 : :
# 135 : : void Num3072::FullReduce()
# 136 : 22 : {
# 137 : 22 : limb_t c0 = MAX_PRIME_DIFF;
# 138 : 22 : limb_t c1 = 0;
# 139 [ + + ]: 1078 : for (int i = 0; i < LIMBS; ++i) {
# 140 : 1056 : addnextract2(c0, c1, this->limbs[i], this->limbs[i]);
# 141 : 1056 : }
# 142 : 22 : }
# 143 : :
# 144 : : Num3072 Num3072::GetInverse() const
# 145 : 803 : {
# 146 : : // For fast exponentiation a sliding window exponentiation with repunit
# 147 : : // precomputation is utilized. See "Fast Point Decompression for Standard
# 148 : : // Elliptic Curves" (Brumley, Järvinen, 2008).
# 149 : :
# 150 : 803 : Num3072 p[12]; // p[i] = a^(2^(2^i)-1)
# 151 : 803 : Num3072 out;
# 152 : :
# 153 : 803 : p[0] = *this;
# 154 : :
# 155 [ + + ]: 9636 : for (int i = 0; i < 11; ++i) {
# 156 : 8833 : p[i + 1] = p[i];
# 157 [ + + ]: 1652574 : for (int j = 0; j < (1 << i); ++j) p[i + 1].Square();
# 158 : 8833 : p[i + 1].Multiply(p[i]);
# 159 : 8833 : }
# 160 : :
# 161 : 803 : out = p[11];
# 162 : :
# 163 : 803 : square_n_mul(out, 512, p[9]);
# 164 : 803 : square_n_mul(out, 256, p[8]);
# 165 : 803 : square_n_mul(out, 128, p[7]);
# 166 : 803 : square_n_mul(out, 64, p[6]);
# 167 : 803 : square_n_mul(out, 32, p[5]);
# 168 : 803 : square_n_mul(out, 8, p[3]);
# 169 : 803 : square_n_mul(out, 2, p[1]);
# 170 : 803 : square_n_mul(out, 1, p[0]);
# 171 : 803 : square_n_mul(out, 5, p[2]);
# 172 : 803 : square_n_mul(out, 3, p[0]);
# 173 : 803 : square_n_mul(out, 2, p[0]);
# 174 : 803 : square_n_mul(out, 4, p[0]);
# 175 : 803 : square_n_mul(out, 4, p[1]);
# 176 : 803 : square_n_mul(out, 3, p[0]);
# 177 : :
# 178 : 803 : return out;
# 179 : 803 : }
# 180 : :
# 181 : : void Num3072::Multiply(const Num3072& a)
# 182 : 23136 : {
# 183 : 23136 : limb_t c0 = 0, c1 = 0, c2 = 0;
# 184 : 23136 : Num3072 tmp;
# 185 : :
# 186 : : /* Compute limbs 0..N-2 of this*a into tmp, including one reduction. */
# 187 [ + + ]: 1110528 : for (int j = 0; j < LIMBS - 1; ++j) {
# 188 : 1087392 : limb_t d0 = 0, d1 = 0, d2 = 0;
# 189 : 1087392 : mul(d0, d1, this->limbs[1 + j], a.limbs[LIMBS + j - (1 + j)]);
# 190 [ + + ]: 26097408 : for (int i = 2 + j; i < LIMBS; ++i) muladd3(d0, d1, d2, this->limbs[i], a.limbs[LIMBS + j - i]);
# 191 : 1087392 : mulnadd3(c0, c1, c2, d0, d1, d2, MAX_PRIME_DIFF);
# 192 [ + + ]: 27184800 : for (int i = 0; i < j + 1; ++i) muladd3(c0, c1, c2, this->limbs[i], a.limbs[j - i]);
# 193 : 1087392 : extract3(c0, c1, c2, tmp.limbs[j]);
# 194 : 1087392 : }
# 195 : :
# 196 : : /* Compute limb N-1 of a*b into tmp. */
# 197 : 23136 : assert(c2 == 0);
# 198 [ + + ]: 1133664 : for (int i = 0; i < LIMBS; ++i) muladd3(c0, c1, c2, this->limbs[i], a.limbs[LIMBS - 1 - i]);
# 199 : 23136 : extract3(c0, c1, c2, tmp.limbs[LIMBS - 1]);
# 200 : :
# 201 : : /* Perform a second reduction. */
# 202 : 23136 : muln2(c0, c1, MAX_PRIME_DIFF);
# 203 [ + + ]: 1133664 : for (int j = 0; j < LIMBS; ++j) {
# 204 : 1110528 : addnextract2(c0, c1, tmp.limbs[j], this->limbs[j]);
# 205 : 1110528 : }
# 206 : :
# 207 : 23136 : assert(c1 == 0);
# 208 : 23136 : assert(c0 == 0 || c0 == 1);
# 209 : :
# 210 : : /* Perform up to two more reductions if the internal state has already
# 211 : : * overflown the MAX of Num3072 or if it is larger than the modulus or
# 212 : : * if both are the case.
# 213 : : * */
# 214 [ + + ]: 23136 : if (this->IsOverflow()) this->FullReduce();
# 215 [ - + ]: 23136 : if (c0) this->FullReduce();
# 216 : 23136 : }
# 217 : :
# 218 : : void Num3072::Square()
# 219 : 2466013 : {
# 220 : 2466013 : limb_t c0 = 0, c1 = 0, c2 = 0;
# 221 : 2466013 : Num3072 tmp;
# 222 : :
# 223 : : /* Compute limbs 0..N-2 of this*this into tmp, including one reduction. */
# 224 [ + + ]: 118368624 : for (int j = 0; j < LIMBS - 1; ++j) {
# 225 : 115902611 : limb_t d0 = 0, d1 = 0, d2 = 0;
# 226 [ + + ]: 1477141787 : for (int i = 0; i < (LIMBS - 1 - j) / 2; ++i) muldbladd3(d0, d1, d2, this->limbs[i + j + 1], this->limbs[LIMBS - 1 - i]);
# 227 [ + + ]: 115902611 : if ((j + 1) & 1) muladd3(d0, d1, d2, this->limbs[(LIMBS - 1 - j) / 2 + j + 1], this->limbs[LIMBS - 1 - (LIMBS - 1 - j) / 2]);
# 228 : 115902611 : mulnadd3(c0, c1, c2, d0, d1, d2, MAX_PRIME_DIFF);
# 229 [ + + ]: 1477141787 : for (int i = 0; i < (j + 1) / 2; ++i) muldbladd3(c0, c1, c2, this->limbs[i], this->limbs[j - i]);
# 230 [ + + ]: 115902611 : if ((j + 1) & 1) muladd3(c0, c1, c2, this->limbs[(j + 1) / 2], this->limbs[j - (j + 1) / 2]);
# 231 : 115902611 : extract3(c0, c1, c2, tmp.limbs[j]);
# 232 : 115902611 : }
# 233 : :
# 234 : 2466013 : assert(c2 == 0);
# 235 [ + + ]: 61650325 : for (int i = 0; i < LIMBS / 2; ++i) muldbladd3(c0, c1, c2, this->limbs[i], this->limbs[LIMBS - 1 - i]);
# 236 : 2466013 : extract3(c0, c1, c2, tmp.limbs[LIMBS - 1]);
# 237 : :
# 238 : : /* Perform a second reduction. */
# 239 : 2466013 : muln2(c0, c1, MAX_PRIME_DIFF);
# 240 [ + + ]: 120834637 : for (int j = 0; j < LIMBS; ++j) {
# 241 : 118368624 : addnextract2(c0, c1, tmp.limbs[j], this->limbs[j]);
# 242 : 118368624 : }
# 243 : :
# 244 : 2466013 : assert(c1 == 0);
# 245 : 2466013 : assert(c0 == 0 || c0 == 1);
# 246 : :
# 247 : : /* Perform up to two more reductions if the internal state has already
# 248 : : * overflown the MAX of Num3072 or if it is larger than the modulus or
# 249 : : * if both are the case.
# 250 : : * */
# 251 [ - + ]: 2466013 : if (this->IsOverflow()) this->FullReduce();
# 252 [ - + ]: 2466013 : if (c0) this->FullReduce();
# 253 : 2466013 : }
# 254 : :
# 255 : : void Num3072::SetToOne()
# 256 : 2502725 : {
# 257 : 2502725 : this->limbs[0] = 1;
# 258 [ + + ]: 120130775 : for (int i = 1; i < LIMBS; ++i) this->limbs[i] = 0;
# 259 : 2502725 : }
# 260 : :
# 261 : : void Num3072::Divide(const Num3072& a)
# 262 : 803 : {
# 263 [ + + ]: 803 : if (this->IsOverflow()) this->FullReduce();
# 264 : :
# 265 : 803 : Num3072 inv{};
# 266 [ - + ]: 803 : if (a.IsOverflow()) {
# 267 : 0 : Num3072 b = a;
# 268 : 0 : b.FullReduce();
# 269 : 0 : inv = b.GetInverse();
# 270 : 803 : } else {
# 271 : 803 : inv = a.GetInverse();
# 272 : 803 : }
# 273 : :
# 274 : 803 : this->Multiply(inv);
# 275 [ - + ]: 803 : if (this->IsOverflow()) this->FullReduce();
# 276 : 803 : }
# 277 : :
# 278 : 1818 : Num3072::Num3072(const unsigned char (&data)[BYTE_SIZE]) {
# 279 [ + + ]: 89082 : for (int i = 0; i < LIMBS; ++i) {
# 280 : 87264 : if (sizeof(limb_t) == 4) {
# 281 : 0 : this->limbs[i] = ReadLE32(data + 4 * i);
# 282 : 87264 : } else if (sizeof(limb_t) == 8) {
# 283 : 87264 : this->limbs[i] = ReadLE64(data + 8 * i);
# 284 : 87264 : }
# 285 : 87264 : }
# 286 : 1818 : }
# 287 : :
# 288 : 803 : void Num3072::ToBytes(unsigned char (&out)[BYTE_SIZE]) {
# 289 [ + + ]: 39347 : for (int i = 0; i < LIMBS; ++i) {
# 290 : 38544 : if (sizeof(limb_t) == 4) {
# 291 : 0 : WriteLE32(out + i * 4, this->limbs[i]);
# 292 : 38544 : } else if (sizeof(limb_t) == 8) {
# 293 : 38544 : WriteLE64(out + i * 8, this->limbs[i]);
# 294 : 38544 : }
# 295 : 38544 : }
# 296 : 803 : }
# 297 : :
# 298 : 1818 : Num3072 MuHash3072::ToNum3072(Span<const unsigned char> in) {
# 299 : 1818 : unsigned char tmp[Num3072::BYTE_SIZE];
# 300 : :
# 301 : 1818 : uint256 hashed_in = (CHashWriter(SER_DISK, 0) << in).GetSHA256();
# 302 : 1818 : ChaCha20(hashed_in.data(), hashed_in.size()).Keystream(tmp, Num3072::BYTE_SIZE);
# 303 : 1818 : Num3072 out{tmp};
# 304 : :
# 305 : 1818 : return out;
# 306 : 1818 : }
# 307 : :
# 308 : : MuHash3072::MuHash3072(Span<const unsigned char> in) noexcept
# 309 : 372 : {
# 310 : 372 : m_numerator = ToNum3072(in);
# 311 : 372 : }
# 312 : :
# 313 : : void MuHash3072::Finalize(uint256& out) noexcept
# 314 : 803 : {
# 315 : 803 : m_numerator.Divide(m_denominator);
# 316 : 803 : m_denominator.SetToOne(); // Needed to keep the MuHash object valid
# 317 : :
# 318 : 803 : unsigned char data[Num3072::BYTE_SIZE];
# 319 : 803 : m_numerator.ToBytes(data);
# 320 : :
# 321 : 803 : out = (CHashWriter(SER_DISK, 0) << data).GetSHA256();
# 322 : 803 : }
# 323 : :
# 324 : : MuHash3072& MuHash3072::operator*=(const MuHash3072& mul) noexcept
# 325 : 248 : {
# 326 : 248 : m_numerator.Multiply(mul.m_numerator);
# 327 : 248 : m_denominator.Multiply(mul.m_denominator);
# 328 : 248 : return *this;
# 329 : 248 : }
# 330 : :
# 331 : : MuHash3072& MuHash3072::operator/=(const MuHash3072& div) noexcept
# 332 : 158 : {
# 333 : 158 : m_numerator.Multiply(div.m_denominator);
# 334 : 158 : m_denominator.Multiply(div.m_numerator);
# 335 : 158 : return *this;
# 336 : 158 : }
# 337 : :
# 338 : 1409 : MuHash3072& MuHash3072::Insert(Span<const unsigned char> in) noexcept {
# 339 : 1409 : m_numerator.Multiply(ToNum3072(in));
# 340 : 1409 : return *this;
# 341 : 1409 : }
# 342 : :
# 343 : 37 : MuHash3072& MuHash3072::Remove(Span<const unsigned char> in) noexcept {
# 344 : 37 : m_denominator.Multiply(ToNum3072(in));
# 345 : 37 : return *this;
# 346 : 37 : }
|