cuPQC-NTT Usage#
The cuPQC-NTT library provides GPU-optimized Number Theoretic Transform (NTT) operations that integrate directly into CUDA kernels. An NTT descriptor is defined by combining cuPQC-NTT operators to specify the algorithm, direction, precision, transform length, and block configuration at compile time.
Defining an NTT Descriptor#
An NTT descriptor specifies all transform parameters as a C++ type using operator composition:
#include <ntt.hpp>
using namespace cupqc;
using ForwardNTT = decltype(Algorithm<algorithm::NTT>()
+ Direction<nttDirection::FORWARD>()
+ Precision<uint32_t>()
+ Size<1024>()
+ Block()
+ BlockDim<128>());
A descriptor for the inverse transform uses nttDirection::INVERSE:
using InverseNTT = decltype(Algorithm<algorithm::NTT>()
+ Direction<nttDirection::INVERSE>()
+ Precision<uint32_t>()
+ Size<1024>()
+ Block()
+ BlockDim<128>());
Prime Fields and Primitive Roots#
The cuPQC-NTT library accepts any suitable prime modulus. When using a 64-bit integer precision type, the prime must be at most 62 bits. For convenience, the library provides built-in compile-time constants for the BabyBear (231 - 227 + 1) and KoalaBear (231 - 224 + 1) prime fields, including pre-computed primitive roots of unity for sizes 210 through 224:
// KoalaBear prime with pre-computed primitive root for N=1024=2^10
constexpr uint32_t p = cupqc::KoalaBear;
constexpr uint32_t g = cupqc::KoalaBear_primitive_root_10;
// BabyBear prime with pre-computed primitive root for N=1024=2^10
constexpr uint32_t p_bb = cupqc::BabyBear;
constexpr uint32_t g_bb = cupqc::BabyBear_primitive_root_10;
Generating Twiddle Factors#
Twiddle factors are precomputed once per prime/generator pair and reused across all
length-N buffers in a batch. A single-threaded kernel calls make_twiddles():
template<class NTT>
__global__ void make_twiddles_kernel(uint32_t* twiddles, const uint32_t p, const uint32_t g) {
NTT().make_twiddles(twiddles, p, g);
}
// Allocate and generate forward and inverse twiddle tables on the device
uint32_t* d_twiddles;
uint32_t* d_inv_twiddles;
cudaMalloc(reinterpret_cast<void**>(&d_twiddles), 1024 * sizeof(uint32_t));
cudaMalloc(reinterpret_cast<void**>(&d_inv_twiddles), 1024 * sizeof(uint32_t));
make_twiddles_kernel<ForwardNTT><<<1, 1>>>(d_twiddles, p, g);
// Compute the inverse generator: g_inv = g^{p-2} mod p (Fermat's little theorem)
make_twiddles_kernel<InverseNTT><<<1, 1>>>(d_inv_twiddles, p, g_inv);
// N_inv = N^{-1} mod p is also required for the inverse NTT normalisation step.
For Montgomery-domain workflows, convert twiddle tables to Montgomery form after generation.
Unlike make_twiddles, transform_twiddles_to_mont is a blockwise function — launch it
with a full block of NTT::BlockDim threads rather than a single thread:
template<class NTT>
__global__ void to_mont_kernel(uint32_t* twiddles, const uint32_t p) {
NTT().transform_twiddles_to_mont(twiddles, p);
}
to_mont_kernel<ForwardNTT><<<1, ForwardNTT::BlockDim>>>(d_twiddles, p);
to_mont_kernel<InverseNTT><<<1, InverseNTT::BlockDim>>>(d_inv_twiddles, p);
Standard NTT Execution#
Each CUDA block runs one length-N transform. The execution pattern is: load data into shared
memory (converting to Montgomery form), run the transform, then store results back to global
memory (converting out of Montgomery form). execute operates entirely on the shared memory
buffer that was populated by load/load_to_mont — always pass that same shared memory
buffer (sdata below), never a global memory pointer, as the first argument to execute,
and it must already be in Montgomery form, which load_to_mont guarantees.
Use NTT::Size and NTT::BlockDim to derive the launch configuration, and
ntt_shared_workspace_size<N, Precision>() for the required shared memory, which must hold
N elements of Precision. Construct an nttConst from the prime on the host and
pass it to the kernel; twiddle tables must be in Montgomery form beforehand (see
Generating Twiddle Factors above).
constexpr size_t smem = ntt_shared_workspace_size<1024, uint32_t>();
// Build Montgomery constants on the host (usable on device too)
cupqc::nttConst<uint32_t> ntt_const(p);
__global__ void forward_ntt_kernel(uint32_t* polys, const uint32_t* twiddles,
const cupqc::nttConst<uint32_t> ntt_const) {
extern __shared__ uint32_t sdata[];
uint32_t* poly = polys + blockIdx.x * ForwardNTT::Size;
ForwardNTT().load_to_mont(sdata, poly, ntt_const);
ForwardNTT().execute(sdata, twiddles, ntt_const.p);
ForwardNTT().store_from_mont(sdata, poly, ntt_const);
}
__global__ void inverse_ntt_kernel(uint32_t* polys, const uint32_t* inv_twiddles,
const cupqc::nttConst<uint32_t> ntt_const,
const uint32_t N_inv) {
extern __shared__ uint32_t sdata[];
uint32_t* poly = polys + blockIdx.x * InverseNTT::Size;
InverseNTT().load_to_mont(sdata, poly, ntt_const);
InverseNTT().execute(sdata, inv_twiddles, ntt_const.p, N_inv);
InverseNTT().store_from_mont(sdata, poly, ntt_const);
}
// One block per length-N buffer
forward_ntt_kernel<<<batch, ForwardNTT::BlockDim, smem>>>(d_data, d_twiddles, ntt_const);
inverse_ntt_kernel<<<batch, InverseNTT::BlockDim, smem>>>(d_data, d_inv_twiddles, ntt_const, N_inv);
Polynomial Arithmetic in the NTT Domain#
After a forward NTT on Montgomery-domain data, use mul, add, or sub to perform
pointwise modular arithmetic on the transformed coefficients in shared memory before
the inverse NTT. Both operands must be in Montgomery form.
The following example shows polynomial multiplication via NTT-based cyclic convolution
(modulo xN - 1):
constexpr size_t smem = ntt_shared_workspace_size<1024, uint32_t>();
__global__ void poly_mul_kernel(const uint32_t* poly_a, const uint32_t* poly_b,
uint32_t* result,
const uint32_t* fwd_twiddles,
const uint32_t* inv_twiddles,
const cupqc::nttConst<uint32_t> ntt_const,
const uint32_t N_inv) {
extern __shared__ uint32_t sdata_a[];
// Allocate a second shared buffer for poly_b immediately after sdata_a
uint32_t* sdata_b = sdata_a + ForwardNTT::Size;
const uint32_t* a = poly_a + blockIdx.x * ForwardNTT::Size;
const uint32_t* b = poly_b + blockIdx.x * ForwardNTT::Size;
uint32_t* r = result + blockIdx.x * ForwardNTT::Size;
ForwardNTT fwd;
InverseNTT inv;
// Forward NTT on both inputs (into Montgomery form)
fwd.load_to_mont(sdata_a, a, ntt_const);
fwd.execute(sdata_a, fwd_twiddles, ntt_const.p);
fwd.load_to_mont(sdata_b, b, ntt_const);
fwd.execute(sdata_b, fwd_twiddles, ntt_const.p);
// Pointwise multiply in the NTT domain (both buffers in Montgomery form)
fwd.mul(sdata_a, sdata_b, ntt_const.p);
// Inverse NTT and write result back to global memory
inv.execute(sdata_a, inv_twiddles, ntt_const.p, N_inv);
inv.store_from_mont(sdata_a, r, ntt_const);
}
// Launch with twice the shared memory to hold both sdata_a and sdata_b
poly_mul_kernel<<<batch, ForwardNTT::BlockDim, 2 * smem>>>(
d_poly_a, d_poly_b, d_result,
d_fwd_twiddles, d_inv_twiddles,
ntt_const, N_inv);
add and sub follow the same pattern as mul.
Staged NTT for Large Transform Sizes#
When the transform length N exceeds available shared memory, use SubSize<M> to split the
transform into two kernel passes. Given M, K = N / M is derived automatically:
constexpr uint32_t NTT_N = 16384;
constexpr uint32_t NTT_M = 256;
constexpr uint32_t NTT_K = NTT_N / NTT_M;
using ForwardNTT = decltype(Algorithm<algorithm::NTT>()
+ Direction<nttDirection::FORWARD>()
+ Precision<uint32_t>()
+ Size<NTT_N>()
+ SubSize<NTT_M>()
+ Block()
+ BlockDim<128>());
using InverseNTT = decltype(Algorithm<algorithm::NTT>()
+ Direction<nttDirection::INVERSE>()
+ Precision<uint32_t>()
+ Size<NTT_N>()
+ SubSize<NTT_M>()
+ Block()
+ BlockDim<128>());
Each stage has its own shared memory requirement and launch grid. As with the standard NTT,
stage_1_execute/stage_2_execute operate on the shared memory buffer filled by the
matching stage_1_load_to_mont/stage_2_load — always pass that buffer, sized as noted below, and
never a global memory pointer, and it must already be in Montgomery form. The intermediate
buffer between stage 1 and stage 2 is always in Montgomery form, so stage_2_load never
needs a conversion; only the initial stage_1_load has a _to_mont variant for
standard-form input (see below). Effectively we convert the input data to Montgomery form in the first stage, and it remains in that form
until the final store in the second stage. Pass blockIdx.x as the id argument to the staged load
and store functions:
// fwd_s1_smem/inv_s2_smem hold K = N / M elements; fwd_s2_smem/inv_s1_smem hold M elements.
constexpr size_t fwd_s1_smem = fwd_stage_1_ntt_shared_workspace_size<NTT_N, NTT_M, uint32_t>();
constexpr size_t fwd_s2_smem = fwd_stage_2_ntt_shared_workspace_size<NTT_N, NTT_M, uint32_t>();
constexpr size_t inv_s1_smem = inv_stage_1_ntt_shared_workspace_size<NTT_N, NTT_M, uint32_t>();
constexpr size_t inv_s2_smem = inv_stage_2_ntt_shared_workspace_size<NTT_N, NTT_M, uint32_t>();
__global__ void fwd_stage_1_kernel(uint32_t* poly, const uint32_t* twiddles,
const cupqc::nttConst<uint32_t> ntt_const) {
// sdata is a shared memory buffer of size K (= NTT_N / NTT_M)
extern __shared__ uint32_t sdata[];
ForwardNTT ntt;
ntt.stage_1_load_to_mont(sdata, poly, blockIdx.x, ntt_const);
ntt.stage_1_execute(sdata, twiddles, ntt_const.p);
ntt.stage_1_store(sdata, poly, blockIdx.x);
}
__global__ void fwd_stage_2_kernel(uint32_t* poly, const uint32_t* twiddles,
const cupqc::nttConst<uint32_t> ntt_const) {
// sdata is a shared memory buffer of size M (= NTT_M)
extern __shared__ uint32_t sdata[];
ForwardNTT ntt;
ntt.stage_2_load(sdata, poly, blockIdx.x);
ntt.stage_2_execute(sdata, twiddles, ntt_const.p);
ntt.stage_2_store_from_mont(sdata, poly, blockIdx.x, ntt_const);
}
__global__ void inv_stage_1_kernel(uint32_t* poly, const uint32_t* inv_twiddles,
const cupqc::nttConst<uint32_t> ntt_const) {
// sdata is a shared memory buffer of size M (= NTT_M)
extern __shared__ uint32_t sdata[];
InverseNTT intt;
intt.stage_1_load_to_mont(sdata, poly, blockIdx.x, ntt_const);
intt.stage_1_execute(sdata, inv_twiddles, ntt_const.p);
intt.stage_1_store(sdata, poly, blockIdx.x);
}
__global__ void inv_stage_2_kernel(uint32_t* poly, const uint32_t* inv_twiddles,
const cupqc::nttConst<uint32_t> ntt_const, const uint32_t N_inv) {
// sdata is a shared memory buffer of size K (= NTT_N / NTT_M)
extern __shared__ uint32_t sdata[];
InverseNTT intt;
intt.stage_2_load(sdata, poly, blockIdx.x);
intt.stage_2_execute(sdata, inv_twiddles, ntt_const.p, N_inv);
intt.stage_2_store_from_mont(sdata, poly, blockIdx.x, ntt_const);
}
/* Another method to convert twiddle tables to Montgomery form is to use the to_mont function.
This function operates on a single piece of data, so it can be used in many contexts when converting to
Montgomery form. The function from_mont is used in an analogous manner to convert back to standard form. */
template<class NTT>
__global__ void twiddles_to_mont_kernel(uint32_t* twiddles, const cupqc::nttConst<uint32_t> ntt_const) {
const uint32_t idx = threadIdx.x + blockIdx.x * blockDim.x;
twiddles[idx] = NTT().to_mont(twiddles[idx], ntt_const);
}
// Build Montgomery constants on the host (usable on device too)
cupqc::nttConst<uint32_t> ntt_const(p);
// Convert twiddle tables to Montgomery form, one thread per element
twiddles_to_mont_kernel<ForwardNTT><<<NTT_N / ForwardNTT::BlockDim, ForwardNTT::BlockDim>>>(
d_twiddles, ntt_const);
twiddles_to_mont_kernel<InverseNTT><<<NTT_N / InverseNTT::BlockDim, InverseNTT::BlockDim>>>(
d_inv_twiddles, ntt_const);
// Forward: stage 1 (M blocks), then stage 2 (K blocks)
fwd_stage_1_kernel<<<NTT_M, ForwardNTT::BlockDim, fwd_s1_smem>>>(d_poly, d_twiddles, ntt_const);
fwd_stage_2_kernel<<<NTT_K, ForwardNTT::BlockDim, fwd_s2_smem>>>(d_poly, d_twiddles, ntt_const);
// Inverse: stage 1 (K blocks), then stage 2 (M blocks, applies 1/N normalisation)
inv_stage_1_kernel<<<NTT_K, InverseNTT::BlockDim, inv_s1_smem>>>(d_poly, d_inv_twiddles, ntt_const);
inv_stage_2_kernel<<<NTT_M, InverseNTT::BlockDim, inv_s2_smem>>>(d_poly, d_inv_twiddles, ntt_const, N_inv);
Compilation#
Include ntt.hpp and link with the cuPQC-NTT static library using LTO flags:
nvcc -std=c++17 -dlto -arch=sm_80 \
-I${PATH_TO_CUPQC_SDK_INCLUDE} \
-L${PATH_TO_STATIC_LIB} -lcupqc-ntt \
my_ntt_program.cu -o my_ntt_program
For CMake-based projects:
target_link_libraries(YourProgram PRIVATE cupqc-ntt_static)
For detailed installation and compilation instructions, see the Getting Started guide.
API Reference:
cuPQC-NTT API Reference - cuPQC-NTT API documentation