go-blockchain/crypto/upstream/msm.h

// Copyright (c) 2014-2018 Zano Project
// Copyright (c) 2014-2018 The Louisdor Project
// Copyright (c) 2012-2013 The Boolberry developers
// Copyright (c) 2017-2025 Lethean (https://lt.hn)
//
// Licensed under the European Union Public Licence (EUPL) version 1.2.
// You may obtain a copy of the licence at:
//
//     https://joinup.ec.europa.eu/software/page/eupl/licence-eupl
//
// The EUPL is a copyleft licence that is compatible with the MIT/X11
// licence used by the original projects; the MIT terms are therefore
// considered “grandfathered” under the EUPL for this code.
//
// SPDX‑License‑Identifier: EUPL-1.2
//
#pragma once

// This file contains Multi-Scalar Multiplication routines

#include "epee/include/misc_log_ex.h"
#include "crypto-sugar.h"

namespace crypto
{

  template<typename CT>
  bool msm_and_check_zero_naive(const scalar_vec_t& g_scalars, const scalar_vec_t& h_scalars, const point_t& summand)
  {
    CHECK_AND_ASSERT_MES(g_scalars.size() <= CT::c_bpp_mn_max, false, "g_scalars oversized");
    CHECK_AND_ASSERT_MES(h_scalars.size() <= CT::c_bpp_mn_max, false, "h_scalars oversized");

    point_t result = summand;

    for (size_t i = 0; i < g_scalars.size(); ++i)
      result += g_scalars[i] * CT::get_generator(false, i);

    for (size_t i = 0; i < h_scalars.size(); ++i)
      result += h_scalars[i] * CT::get_generator(true, i);

    if (!result.is_zero())
    {
      LOG_PRINT_L0("msm result is non zero: " << result);
      return false;
    }
    return true;
  }


  // https://eprint.iacr.org/2022/999.pdf
  // "Pippenger algorithm [1], and its variant that is widely used in the ZK space is called the bucket method"
  template<typename CT>
  bool msm_and_check_zero_pippenger_v3(const scalar_vec_t& g_scalars, const scalar_vec_t& h_scalars, const point_t& summand, uint8_t c)
  {
    // TODO: with c = 8 and with direct access got much worse result than with c = 7 and get_bits() for N = 128..256, consider checking again for bigger datasets (N>256) 
    // TODO: consider preparing a cached generators' points

    CHECK_AND_ASSERT_MES(g_scalars.size() <= CT::c_bpp_mn_max, false, "g_scalars oversized");
    CHECK_AND_ASSERT_MES(h_scalars.size() <= CT::c_bpp_mn_max, false, "h_scalars oversized");
    CHECK_AND_ASSERT_MES(c < 10, false, "c is too big");

    size_t C = 1ull << c;

    // k_max * c + (c-1) >= max_bit_idx
    // 
    //                 max_bit_idx - (c - 1)            max_bit_idx - (c - 1) + (c - 1)              max_bit_idx
    // k_max = ceil ( --------------------- ) = floor ( ------------------------------ )  =  floor ( ----------- )
    //                           c                                    c                                   c      
    const size_t b = 253; // the maximum number of bits in x  https://eprint.iacr.org/2022/999.pdf  TODO: we may also scan for maximum bit used in all the scalars if all the scalars are small
    const size_t max_bit_idx = b - 1;
    const size_t k_max = max_bit_idx / c;
    const size_t K = k_max + 1;

    std::vector<point_t> buckets(C * K);
    std::vector<bool> buckets_inited(C * K);
    std::vector<point_t> Sk(K);
    std::vector<bool> Sk_inited(K);
    std::vector<point_t> Gk(K);
    std::vector<bool> Gk_inited(K);

    // first loop, calculate partial bucket sums
    for (size_t n = 0; n < g_scalars.size(); ++n)
    {
      for (size_t k = 0; k < K; ++k)
      {
        uint64_t l = g_scalars[n].get_bits((uint8_t)(k * c), c); // l in [0; 2^c-1]
        if (l != 0)
        {
          size_t bucket_id = l * K + k;
          if (buckets_inited[bucket_id])
            buckets[bucket_id] += CT::get_generator(false, n);
          else
          {
            buckets[bucket_id] =  CT::get_generator(false, n);
            buckets_inited[bucket_id] = true;
          }
        }
      }
    }
    // still the first loop (continued)
    for (size_t n = 0; n < h_scalars.size(); ++n)
    {
      for (size_t k = 0; k < K; ++k)
      {
        uint64_t l = h_scalars[n].get_bits((uint8_t)(k * c), c); // l in [0; 2^c-1]
        if (l != 0)
        {
          size_t bucket_id = l * K + k;
          if (buckets_inited[bucket_id])
            buckets[bucket_id] += CT::get_generator(true, n);
          else
          {
            buckets[bucket_id] =  CT::get_generator(true, n);
            buckets_inited[bucket_id] = true;
          }
        }
      }
    }

    // the second loop
    for (size_t l = C - 1; l > 0; --l)
    {
      for (size_t k = 0; k < K; ++k)
      {
        size_t bucket_id = l * K + k;
        if (buckets_inited[bucket_id])
        {
          if (Sk_inited[k])
            Sk[k] += buckets[bucket_id];
          else
          {
            Sk[k] = buckets[bucket_id];
            Sk_inited[k] = true;
          }
        }

        if (Sk_inited[k])
        {
          if (Gk_inited[k])
            Gk[k] += Sk[k];
          else
          {
            Gk[k] = Sk[k];
            Gk_inited[k] = true;
          }
        }
      }
    }

    // the third loop: Horner’s rule
    point_t result = Gk_inited[K - 1] ? Gk[K - 1] : c_point_0;
    for (size_t k = K - 2; k != SIZE_MAX; --k)
    {
      result.modify_mul_pow_2(c);
      if (Gk_inited[k])
        result += Gk[k];
    }

    result += summand;

    if (!result.is_zero())
    {
      LOG_PRINT_L0("multiexp result is non zero: " << result);
      return false;
    }

    return true;
  }



  

  // Just switcher

  template<typename CT>
  bool msm_and_check_zero(const scalar_vec_t& g_scalars, const scalar_vec_t& h_scalars, const point_t& summand)
  {
    //return msm_and_check_zero_naive<CT>(g_scalars, h_scalars, summand);
    return msm_and_check_zero_pippenger_v3<CT>(g_scalars, h_scalars, summand, 7);
  }



} // namespace crypto