Skip to content

Repository navigation

optimize coefficent combinations #176

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
odraccir1 merged 1 commit into from
Oct 17, 2026
Merged

optimize coefficent combinations #176

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
28 changes: 27 additions & 1 deletion cmake/CheckAVX.cmake
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -37,6 +37,19 @@ set(AVX512_CODE "
}
")

set(AVX512IFMA_CODE "
#include <immintrin.h>
int main() {
__m256i a = _mm256_set1_epi64x(1);
a = _mm256_madd52lo_epu64(a, a, a);
a = _mm256_madd52hi_epu64(a, a, a);
__m512i b = _mm512_set1_epi64(1);
b = _mm512_madd52lo_epu64(b, b, b);
b = _mm512_madd52hi_epu64(b, b, b);
return _mm256_extract_epi64(a, 0) + _mm_cvtsi128_si32(_mm512_castsi512_si128(b));
}
")

macro(check_avx_lang lang type flags)
set(__FLAG_I 1)
set(CMAKE_REQUIRED_FLAGS_SAVE ${CMAKE_REQUIRED_FLAGS})
Expand All @@ -45,7 +58,7 @@ macro(check_avx_lang lang type flags)
if(NOT ${lang}_${type}_FOUND)
set(CMAKE_REQUIRED_QUIET TRUE)
set(CMAKE_REQUIRED_FLAGS ${__FLAG})
if(lang STREQUAL "CXX")
if("${lang}" STREQUAL "CXX")
check_cxx_source_compiles("${${type}_CODE}" ${lang}_HAS_${type}_${__FLAG_I})
else()
check_c_source_compiles("${${type}_CODE}" ${lang}_HAS_${type}_${__FLAG_I})
Expand All @@ -69,6 +82,19 @@ macro(check_avx_lang lang type flags)

endmacro()

# Both widths use the same host contract; VL is required by the four-lane kernel.
macro(check_avx512ifma)
check_host_cpu_feature(AVX512F "avx512f")
check_host_cpu_feature(AVX512IFMA "avx512ifma")
check_host_cpu_feature(AVX512VL "avx512vl")
if(HOST_CPU_AVX512F_FOUND AND HOST_CPU_AVX512IFMA_FOUND AND HOST_CPU_AVX512VL_FOUND)
check_avx_lang(CXX "AVX512IFMA" " ;-mavx512f -mavx512vl -mavx512ifma;/arch:AVX512")
else()
set(CXX_AVX512IFMA_FOUND FALSE CACHE BOOL "CXX AVX512IFMA support" FORCE)
set(CXX_AVX512IFMA_FLAGS "" CACHE STRING "CXX AVX512IFMA flags" FORCE)
endif()
endmacro()

macro(check_avx)
check_host_cpu_feature(AVX "avx")
check_host_cpu_feature(AVX2 "avx2")
Expand Down
4 changes: 4 additions & 0 deletions cmake/CheckHostCPUFeature.cmake
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -30,6 +30,10 @@ function(check_host_cpu_feature feature builtin_name)
set(msvc_test "max_leaf >= 7 && (leaf1[2] & (1 << 27)) && (leaf1[2] & (1 << 28)) && ((_xgetbv(0) & 0x6) == 0x6) && (leaf7[1] & (1 << 5))")
elseif(feature_upper STREQUAL "AVX512F")
set(msvc_test "max_leaf >= 7 && (leaf1[2] & (1 << 27)) && ((_xgetbv(0) & 0xE6) == 0xE6) && (leaf7[1] & (1 << 16))")
elseif(feature_upper STREQUAL "AVX512IFMA")
set(msvc_test "max_leaf >= 7 && (leaf1[2] & (1 << 27)) && ((_xgetbv(0) & 0xE6) == 0xE6) && (leaf7[1] & (1 << 16)) && (leaf7[1] & (1 << 21))")
elseif(feature_upper STREQUAL "AVX512VL")
set(msvc_test "max_leaf >= 7 && (leaf1[2] & (1 << 27)) && ((_xgetbv(0) & 0xE6) == 0xE6) && (leaf7[1] & (1 << 16)) && (leaf7[1] & (1u << 31))")
elseif(feature_upper STREQUAL "BMI2")
set(msvc_test "max_leaf >= 7 && (leaf7[1] & (1 << 8))")
elseif(feature_upper STREQUAL "ADX")
Expand Down
8 changes: 8 additions & 0 deletions libs/algebra/CMakeLists.txt
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,7 @@ cm_project(algebra WORKSPACE_NAME ${CMAKE_WORKSPACE_NAME} LANGUAGES ASM C CXX)
include(CMDeploy)
include(CheckADX)
include(CheckBMI2)
include(CheckAVX)

# Blueprint components are using point {0, 0} as a point in infinity, while the standard is {0, 1}
# for most curves. We will use compatibility mode with blueprint as default. Setting this flag to
Expand Down Expand Up @@ -46,6 +47,7 @@ cm_deploy(TARGETS ${CMAKE_WORKSPACE_NAME}_${CURRENT_PROJECT_NAME}
if(CMAKE_TARGET_ARCHITECTURE STREQUAL "x86_64" OR CMAKE_TARGET_ARCHITECTURE STREQUAL "x86")
check_bmi2()
check_adx()
check_avx512ifma()

if(CXX_BMI2_FOUND)
target_compile_definitions(${CMAKE_WORKSPACE_NAME}_${CURRENT_PROJECT_NAME} INTERFACE
Expand All @@ -57,6 +59,12 @@ if(CMAKE_TARGET_ARCHITECTURE STREQUAL "x86_64" OR CMAKE_TARGET_ARCHITECTURE STRE
"${CMAKE_UPPER_WORKSPACE_NAME}_HAS_ADX")
target_compile_options(${CMAKE_WORKSPACE_NAME}_${CURRENT_PROJECT_NAME} INTERFACE ${CXX_ADX_FLAGS})
endif()
if(CXX_AVX512IFMA_FOUND)
target_compile_definitions(${CMAKE_WORKSPACE_NAME}_${CURRENT_PROJECT_NAME} INTERFACE
"${CMAKE_UPPER_WORKSPACE_NAME}_HAS_AVX512IFMA")
separate_arguments(CXX_AVX512IFMA_OPTIONS NATIVE_COMMAND "${CXX_AVX512IFMA_FLAGS}")
target_compile_options(${CMAKE_WORKSPACE_NAME}_${CURRENT_PROJECT_NAME} INTERFACE ${CXX_AVX512IFMA_OPTIONS})
endif()
endif()

include(CMTest)
Expand Down
103 changes: 103 additions & 0 deletions libs/algebra/include/nil/crypto3/algebra/fields/detail/element/fp12_fast/fp12_fast.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,12 @@

#include <array>
#include <cstddef>
#include <iterator>
#include <tuple>
#include <type_traits>
#include <utility>

#include <boost/iterator/transform_iterator.hpp>

#include <nil/crypto3/algebra/fields/detail/element/fp12_fast/types.hpp>
#include <nil/crypto3/algebra/fields/detail/element/fp12_fast/limb_ops.hpp>
Expand All @@ -12,6 +17,15 @@ namespace nil::crypto3::algebra::fields::detail::fp12_fast {
template<Fp12FastParams Params>
struct fp12_fast {
using limb_array = Params::limb_array;
using base_value_type = typename Params::base_field_type::value_type;
using base_backend_type =
std::remove_cvref_t<decltype(std::declval<base_value_type &>().data.backend().base_data())>;

// Three Fp6 products, six Fp2 products each, and three base-field products per Fp2 product.
// Every entry is a canonical Montgomery residue. These are reusable linear forms of one operand;
// multiplying corresponding entries and applying recombine_products gives the Fp12 product.
static constexpr std::size_t product_coordinate_count = 54;
using prepared_product_operand = std::array<base_backend_type, product_coordinate_count>;

// Packed input-side Fp2 value: two canonical Fp coefficients.
struct fp2_base {
Expand Down Expand Up @@ -63,6 +77,15 @@ namespace nil::crypto3::algebra::fields::detail::fp12_fast {
fp2_base::add_mod(*this, *this, other);
return *this;
}

void to_non_residue(Params::non_residue_type &result) const {
constexpr std::size_t N = Params::base_value_limb_count;
for (std::size_t coefficient = 0; coefficient < 2; ++coefficient) {
for (std::size_t i = 0; i < N; ++i) {
result.data[coefficient].data.backend().base_data().limbs()[i] = data[coefficient * N + i];
}
}
}
};

struct fp2_dbl {
Expand Down Expand Up @@ -271,6 +294,86 @@ namespace nil::crypto3::algebra::fields::detail::fp12_fast {
}
};

template<typename Fp12Value>
static prepared_product_operand prepare_product_operand(const Fp12Value &x) {
using fp2_value_type = typename Params::non_residue_type;
using fp6_value_type = typename Params::underlying_type;
const std::array<fp6_value_type, 3> outer = {x.data[0], x.data[1], x.data[0] + x.data[1]};
prepared_product_operand result;
for (std::size_t w = 0; w < outer.size(); ++w) {
const auto &a = outer[w].data;
const std::array<fp2_value_type, 6> forms = {a[0], a[1], a[2], a[0] + a[1], a[0] + a[2], a[1] + a[2]};
for (std::size_t v = 0; v < forms.size(); ++v) {
const auto &b = forms[v].data;
const base_value_type sum = b[0] + b[1];
const std::size_t offset = 18 * w + 3 * v;
result[offset] = b[0].data.backend().base_data();
result[offset + 1] = b[1].data.backend().base_data();
result[offset + 2] = sum.data.backend().base_data();
}
}
return result;
}

// This map is linear over Fp, so it also recombines sums of the 54 coordinate products. All corrections
// take place after the dot products, with canonical field arithmetic rather than unbounded lazy values.
template<typename Fp12Value>
static Fp12Value recombine_products(const std::array<base_value_type, product_coordinate_count> &products) {
using fp2_value_type = typename Params::non_residue_type;
using fp6_value_type = typename Params::underlying_type;
const auto mul_xi = [](const fp2_value_type &x) {
fp2_base packed;
fp2_base::mul_xi(packed, fp2_base(x));
fp2_value_type result;
packed.to_non_residue(result);
return result;
};
std::array<fp6_value_type, 3> outer;
for (std::size_t w = 0; w < outer.size(); ++w) {
std::array<fp2_value_type, 6> p;
for (std::size_t v = 0; v < p.size(); ++v) {
const std::size_t offset = 18 * w + 3 * v;
const auto &ac = products[offset];
const auto &bd = products[offset + 1];
auto real = ac - bd;
if constexpr (Params::u_squared == -5) {
real -= bd.doubled().doubled();
} else {
static_assert(Params::u_squared == -1, "unsupported quadratic non-residue");
}
p[v] = fp2_value_type(real, products[offset + 2] - ac - bd);
}
// p[0..2] are the diagonal Fp2 products; p[3..5] are the pairwise sum products.
outer[w] = fp6_value_type(p[0] + mul_xi(p[5] - p[1] - p[2]), p[3] - p[0] - p[1] + mul_xi(p[2]),
p[4] - p[0] - p[2] + p[1]);
}
const auto &ac = outer[0];
const auto &bd = outer[1];
const fp6_value_type v_bd(mul_xi(bd.data[2]), bd.data[0], bd.data[1]);
return Fp12Value(ac + v_bd, outer[2] - ac - bd);
}

// Preconditions: both ranges contain prepared canonical operands; the right range is at least as long
// as the left. Multipass iterators are required because each coordinate traverses the same operands.
template<typename Fp12Value, std::forward_iterator LeftIterator, std::forward_iterator RightIterator>
static Fp12Value inner_product(LeftIterator first, LeftIterator last, RightIterator right) {
if (first == last) {
return Fp12Value::zero();
}
std::array<base_value_type, product_coordinate_count> products;
for (std::size_t coordinate = 0; coordinate < product_coordinate_count; ++coordinate) {
const auto select = [coordinate](const prepared_product_operand &operand) -> const base_backend_type & {
return operand[coordinate];
};
// Reuse the existing 2*N+1-limb accumulator: canonical inputs are < p < R, so any size_t-length
// sum fits. Its final canonicalization handles sums larger than p*R, including 123-term dots.
montgomery_inner_product<Params>(
products[coordinate].data.backend().base_data(), boost::make_transform_iterator(first, select),
boost::make_transform_iterator(last, select), boost::make_transform_iterator(right, select));
}
return recombine_products<Fp12Value>(products);
}

template<typename Fp12Value>
static Fp12Value square(const Fp12Value &x) {
// Write x = a + b*w in Fp6[w]/(w^2 - v). The quadratic-extension
Expand Down
Loading
Loading