Branch data Line data Source code
# 1 : : // Copyright (c) 2017-2021 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 : 136120069 : {
# 26 : 136120069 : n = c0;
# 27 : 136120069 : c0 = c1;
# 28 : 136120069 : c1 = c2;
# 29 : 136120069 : c2 = 0;
# 30 : 136120069 : }
# 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 : 1227823 : {
# 35 : 1227823 : double_limb_t t = (double_limb_t)a * b;
# 36 : 1227823 : c1 = t >> LIMB_SIZE;
# 37 : 1227823 : c0 = t;
# 38 : 1227823 : }
# 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 : 133286629 : {
# 43 : 133286629 : double_limb_t t = (double_limb_t)d0 * n + c0;
# 44 : 133286629 : c0 = t;
# 45 : 133286629 : t >>= LIMB_SIZE;
# 46 : 133286629 : t += (double_limb_t)d1 * n + c1;
# 47 : 133286629 : c1 = t;
# 48 : 133286629 : t >>= LIMB_SIZE;
# 49 : 133286629 : c2 = t + d2 * n;
# 50 : 133286629 : }
# 51 : :
# 52 : : /* [c0,c1] *= n */
# 53 : : inline void muln2(limb_t& c0, limb_t& c1, const limb_t& n)
# 54 : 2836070 : {
# 55 : 2836070 : double_limb_t t = (double_limb_t)c0 * n;
# 56 : 2836070 : c0 = t;
# 57 : 2836070 : t >>= LIMB_SIZE;
# 58 : 2836070 : t += (double_limb_t)c1 * n;
# 59 : 2836070 : c1 = t;
# 60 : 2836070 : }
# 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 : 193819568 : {
# 65 : 193819568 : double_limb_t t = (double_limb_t)a * b;
# 66 : 193819568 : limb_t th = t >> LIMB_SIZE;
# 67 : 193819568 : limb_t tl = t;
# 68 : :
# 69 : 193819568 : c0 += tl;
# 70 [ + + ]: 193819568 : th += (c0 < tl) ? 1 : 0;
# 71 : 193819568 : c1 += th;
# 72 [ + + ]: 193819568 : c2 += (c1 < th) ? 1 : 0;
# 73 : 193819568 : }
# 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 : 3166472178 : {
# 78 : 3166472178 : double_limb_t t = (double_limb_t)a * b;
# 79 : 3166472178 : limb_t th = t >> LIMB_SIZE;
# 80 : 3166472178 : limb_t tl = t;
# 81 : :
# 82 : 3166472178 : c0 += tl;
# 83 [ + + ]: 3166472178 : limb_t tt = th + ((c0 < tl) ? 1 : 0);
# 84 : 3166472178 : c1 += tt;
# 85 [ + + ]: 3166472178 : c2 += (c1 < tt) ? 1 : 0;
# 86 : 3166472178 : c0 += tl;
# 87 [ + + ]: 3166472178 : th += (c0 < tl) ? 1 : 0;
# 88 : 3166472178 : c1 += th;
# 89 [ + + ]: 3166472178 : c2 += (c1 < th) ? 1 : 0;
# 90 : 3166472178 : }
# 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 : 136121239 : {
# 98 : 136121239 : limb_t c2 = 0;
# 99 : :
# 100 : : // add
# 101 : 136121239 : c0 += a;
# 102 [ + + ]: 136121239 : if (c0 < a) {
# 103 : 195524 : c1 += 1;
# 104 : :
# 105 : : // Handle case when c1 has overflown
# 106 [ - + ]: 195524 : if (c1 == 0)
# 107 : 0 : c2 = 1;
# 108 : 195524 : }
# 109 : :
# 110 : : // extract
# 111 : 136121239 : n = c0;
# 112 : 136121239 : c0 = c1;
# 113 : 136121239 : c1 = c2;
# 114 : 136121239 : }
# 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 : 12810 : {
# 119 [ + + ]: 949769 : for (int j = 0; j < sq; ++j) in_out.Square();
# 120 : 12810 : in_out.Multiply(mul);
# 121 : 12810 : }
# 122 : :
# 123 : : } // namespace
# 124 : :
# 125 : : /** Indicates whether d is larger than the modulus. */
# 126 : : bool Num3072::IsOverflow() const
# 127 : 2838805 : {
# 128 [ + + ]: 2838805 : if (this->limbs[0] <= std::numeric_limits<limb_t>::max() - MAX_PRIME_DIFF) return false;
# 129 [ + + ]: 1825 : for (int i = 1; i < LIMBS; ++i) {
# 130 [ - + ]: 1786 : if (this->limbs[i] != std::numeric_limits<limb_t>::max()) return false;
# 131 : 1786 : }
# 132 : 39 : return true;
# 133 : 39 : }
# 134 : :
# 135 : : void Num3072::FullReduce()
# 136 : 38 : {
# 137 : 38 : limb_t c0 = MAX_PRIME_DIFF;
# 138 : 38 : limb_t c1 = 0;
# 139 [ + + ]: 1862 : for (int i = 0; i < LIMBS; ++i) {
# 140 : 1824 : addnextract2(c0, c1, this->limbs[i], this->limbs[i]);
# 141 : 1824 : }
# 142 : 38 : }
# 143 : :
# 144 : : Num3072 Num3072::GetInverse() const
# 145 : 915 : {
# 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 : 915 : Num3072 p[12]; // p[i] = a^(2^(2^i)-1)
# 151 : 915 : Num3072 out;
# 152 : :
# 153 : 915 : p[0] = *this;
# 154 : :
# 155 [ + + ]: 10980 : for (int i = 0; i < 11; ++i) {
# 156 : 10065 : p[i + 1] = p[i];
# 157 [ + + ]: 1883049 : for (int j = 0; j < (1 << i); ++j) p[i + 1].Square();
# 158 : 10065 : p[i + 1].Multiply(p[i]);
# 159 : 10065 : }
# 160 : :
# 161 : 915 : out = p[11];
# 162 : :
# 163 : 915 : square_n_mul(out, 512, p[9]);
# 164 : 915 : square_n_mul(out, 256, p[8]);
# 165 : 915 : square_n_mul(out, 128, p[7]);
# 166 : 915 : square_n_mul(out, 64, p[6]);
# 167 : 915 : square_n_mul(out, 32, p[5]);
# 168 : 915 : square_n_mul(out, 8, p[3]);
# 169 : 915 : square_n_mul(out, 2, p[1]);
# 170 : 915 : square_n_mul(out, 1, p[0]);
# 171 : 915 : square_n_mul(out, 5, p[2]);
# 172 : 915 : square_n_mul(out, 3, p[0]);
# 173 : 915 : square_n_mul(out, 2, p[0]);
# 174 : 915 : square_n_mul(out, 4, p[0]);
# 175 : 915 : square_n_mul(out, 4, p[1]);
# 176 : 915 : square_n_mul(out, 3, p[0]);
# 177 : :
# 178 : 915 : return out;
# 179 : 915 : }
# 180 : :
# 181 : : void Num3072::Multiply(const Num3072& a)
# 182 : 26124 : {
# 183 : 26124 : limb_t c0 = 0, c1 = 0, c2 = 0;
# 184 : 26124 : Num3072 tmp;
# 185 : :
# 186 : : /* Compute limbs 0..N-2 of this*a into tmp, including one reduction. */
# 187 [ + + ]: 1253945 : for (int j = 0; j < LIMBS - 1; ++j) {
# 188 : 1227821 : limb_t d0 = 0, d1 = 0, d2 = 0;
# 189 : 1227821 : mul(d0, d1, this->limbs[1 + j], a.limbs[LIMBS + j - (1 + j)]);
# 190 [ + + ]: 29467729 : for (int i = 2 + j; i < LIMBS; ++i) muladd3(d0, d1, d2, this->limbs[i], a.limbs[LIMBS + j - i]);
# 191 : 1227821 : mulnadd3(c0, c1, c2, d0, d1, d2, MAX_PRIME_DIFF);
# 192 [ + + ]: 30695519 : for (int i = 0; i < j + 1; ++i) muladd3(c0, c1, c2, this->limbs[i], a.limbs[j - i]);
# 193 : 1227821 : extract3(c0, c1, c2, tmp.limbs[j]);
# 194 : 1227821 : }
# 195 : :
# 196 : : /* Compute limb N-1 of a*b into tmp. */
# 197 : 26124 : assert(c2 == 0);
# 198 [ + + ]: 1280076 : for (int i = 0; i < LIMBS; ++i) muladd3(c0, c1, c2, this->limbs[i], a.limbs[LIMBS - 1 - i]);
# 199 : 26124 : extract3(c0, c1, c2, tmp.limbs[LIMBS - 1]);
# 200 : :
# 201 : : /* Perform a second reduction. */
# 202 : 26124 : muln2(c0, c1, MAX_PRIME_DIFF);
# 203 [ + + ]: 1280076 : for (int j = 0; j < LIMBS; ++j) {
# 204 : 1253952 : addnextract2(c0, c1, tmp.limbs[j], this->limbs[j]);
# 205 : 1253952 : }
# 206 : :
# 207 : 26124 : assert(c1 == 0);
# 208 : 0 : 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 [ + + ]: 26124 : if (this->IsOverflow()) this->FullReduce();
# 215 [ - + ]: 26124 : if (c0) this->FullReduce();
# 216 : 26124 : }
# 217 : :
# 218 : : void Num3072::Square()
# 219 : 2809933 : {
# 220 : 2809933 : limb_t c0 = 0, c1 = 0, c2 = 0;
# 221 : 2809933 : Num3072 tmp;
# 222 : :
# 223 : : /* Compute limbs 0..N-2 of this*this into tmp, including one reduction. */
# 224 [ + + ]: 134867854 : for (int j = 0; j < LIMBS - 1; ++j) {
# 225 : 132057921 : limb_t d0 = 0, d1 = 0, d2 = 0;
# 226 [ + + ]: 1682329182 : 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 [ + + ]: 132057921 : 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 : 132057921 : mulnadd3(c0, c1, c2, d0, d1, d2, MAX_PRIME_DIFF);
# 229 [ + + ]: 1682340487 : for (int i = 0; i < (j + 1) / 2; ++i) muldbladd3(c0, c1, c2, this->limbs[i], this->limbs[j - i]);
# 230 [ + + ]: 132057921 : if ((j + 1) & 1) muladd3(c0, c1, c2, this->limbs[(j + 1) / 2], this->limbs[j - (j + 1) / 2]);
# 231 : 132057921 : extract3(c0, c1, c2, tmp.limbs[j]);
# 232 : 132057921 : }
# 233 : :
# 234 : 2809933 : assert(c2 == 0);
# 235 [ + + ]: 70243883 : for (int i = 0; i < LIMBS / 2; ++i) muldbladd3(c0, c1, c2, this->limbs[i], this->limbs[LIMBS - 1 - i]);
# 236 : 2809933 : extract3(c0, c1, c2, tmp.limbs[LIMBS - 1]);
# 237 : :
# 238 : : /* Perform a second reduction. */
# 239 : 2809933 : muln2(c0, c1, MAX_PRIME_DIFF);
# 240 [ + + ]: 137675637 : for (int j = 0; j < LIMBS; ++j) {
# 241 : 134865704 : addnextract2(c0, c1, tmp.limbs[j], this->limbs[j]);
# 242 : 134865704 : }
# 243 : :
# 244 : 2809933 : assert(c1 == 0);
# 245 : 0 : 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 [ - + ]: 2809933 : if (this->IsOverflow()) this->FullReduce();
# 252 [ - + ]: 2809933 : if (c0) this->FullReduce();
# 253 : 2809933 : }
# 254 : :
# 255 : : void Num3072::SetToOne()
# 256 : 2850880 : {
# 257 : 2850880 : this->limbs[0] = 1;
# 258 [ + + ]: 136838457 : for (int i = 1; i < LIMBS; ++i) this->limbs[i] = 0;
# 259 : 2850880 : }
# 260 : :
# 261 : : void Num3072::Divide(const Num3072& a)
# 262 : 915 : {
# 263 [ + + ]: 915 : if (this->IsOverflow()) this->FullReduce();
# 264 : :
# 265 : 915 : Num3072 inv{};
# 266 [ - + ]: 915 : if (a.IsOverflow()) {
# 267 : 0 : Num3072 b = a;
# 268 : 0 : b.FullReduce();
# 269 : 0 : inv = b.GetInverse();
# 270 : 915 : } else {
# 271 : 915 : inv = a.GetInverse();
# 272 : 915 : }
# 273 : :
# 274 : 915 : this->Multiply(inv);
# 275 [ - + ]: 915 : if (this->IsOverflow()) this->FullReduce();
# 276 : 915 : }
# 277 : :
# 278 : 1894 : Num3072::Num3072(const unsigned char (&data)[BYTE_SIZE]) {
# 279 [ + + ]: 92806 : for (int i = 0; i < LIMBS; ++i) {
# 280 : 90912 : if (sizeof(limb_t) == 4) {
# 281 : 0 : this->limbs[i] = ReadLE32(data + 4 * i);
# 282 : 90912 : } else if (sizeof(limb_t) == 8) {
# 283 : 90912 : this->limbs[i] = ReadLE64(data + 8 * i);
# 284 : 90912 : }
# 285 : 90912 : }
# 286 : 1894 : }
# 287 : :
# 288 : 915 : void Num3072::ToBytes(unsigned char (&out)[BYTE_SIZE]) {
# 289 [ + + ]: 44835 : for (int i = 0; i < LIMBS; ++i) {
# 290 : 43920 : if (sizeof(limb_t) == 4) {
# 291 : 0 : WriteLE32(out + i * 4, this->limbs[i]);
# 292 : 43920 : } else if (sizeof(limb_t) == 8) {
# 293 : 43920 : WriteLE64(out + i * 8, this->limbs[i]);
# 294 : 43920 : }
# 295 : 43920 : }
# 296 : 915 : }
# 297 : :
# 298 : 1894 : Num3072 MuHash3072::ToNum3072(Span<const unsigned char> in) {
# 299 : 1894 : unsigned char tmp[Num3072::BYTE_SIZE];
# 300 : :
# 301 : 1894 : uint256 hashed_in = (CHashWriter(SER_DISK, 0) << in).GetSHA256();
# 302 : 1894 : ChaCha20(hashed_in.data(), hashed_in.size()).Keystream(tmp, Num3072::BYTE_SIZE);
# 303 : 1894 : Num3072 out{tmp};
# 304 : :
# 305 : 1894 : return out;
# 306 : 1894 : }
# 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 : 915 : {
# 315 : 915 : m_numerator.Divide(m_denominator);
# 316 : 915 : m_denominator.SetToOne(); // Needed to keep the MuHash object valid
# 317 : :
# 318 : 915 : unsigned char data[Num3072::BYTE_SIZE];
# 319 : 915 : m_numerator.ToBytes(data);
# 320 : :
# 321 : 915 : out = (CHashWriter(SER_DISK, 0) << data).GetSHA256();
# 322 : 915 : }
# 323 : :
# 324 : : MuHash3072& MuHash3072::operator*=(const MuHash3072& mul) noexcept
# 325 : 232 : {
# 326 : 232 : m_numerator.Multiply(mul.m_numerator);
# 327 : 232 : m_denominator.Multiply(mul.m_denominator);
# 328 : 232 : return *this;
# 329 : 232 : }
# 330 : :
# 331 : : MuHash3072& MuHash3072::operator/=(const MuHash3072& div) noexcept
# 332 : 174 : {
# 333 : 174 : m_numerator.Multiply(div.m_denominator);
# 334 : 174 : m_denominator.Multiply(div.m_numerator);
# 335 : 174 : return *this;
# 336 : 174 : }
# 337 : :
# 338 : 1493 : MuHash3072& MuHash3072::Insert(Span<const unsigned char> in) noexcept {
# 339 : 1493 : m_numerator.Multiply(ToNum3072(in));
# 340 : 1493 : return *this;
# 341 : 1493 : }
# 342 : :
# 343 : 29 : MuHash3072& MuHash3072::Remove(Span<const unsigned char> in) noexcept {
# 344 : 29 : m_denominator.Multiply(ToNum3072(in));
# 345 : 29 : return *this;
# 346 : 29 : }
|
