From 219721c29a2c96497bfd615aae90b328a266c90c Mon Sep 17 00:00:00 2001 From: Carlo Pascoe Date: Wed, 24 Jun 2026 11:52:34 -0400 Subject: [PATCH 1/4] Reduce CKKS bootstrapping rotation-key footprint in the homomorphic FFT Reformulate the baby-step/giant-step (BSGS) evaluation of the CoeffsToSlots and SlotsToCoeffs linear transforms so that the same slot permutations are realized with far fewer distinct rotation (automorphism) keys. The keys these transforms need dominate the memory cost of a bootstrapping context, and EvalBootstrapKeyGen now generates a smaller, more heavily shared set. The underlying math is unchanged. Core changes to EvalCoeffsToSlots / EvalSlotsToCoeffs and their precompute: - Horner single giant-step. Replace the forward outer sum sum_i Aut_{i*t}(block_i) with the nested Horner form block_0 + Aut_t(block_1 + Aut_t(...)). Both use the same rotation count, but Horner needs only one giant-step key per level (stride t = g*scale) instead of the b-1 distinct keys {t, 2t, ..., (b-1)t}. - Zero-based hoisted inner rotations. Replace the centered inner baby-step indices {(j-offset)*sigma} with zero-based {j*sigma}. The j=0 term is now always rotation 0 and is handled directly by KeySwitchExt (no key), and the per-level offset delta_s = offset*scale is pushed offline into the precomputed plaintexts as a pre-rotation. - Folded per-level corrections. The per-level zero-basing corrections (an O(levels) set of runtime EvalAtIndex calls in both transforms) are commuted forward and absorbed into the precomputed plaintexts, leaving a single accumulated rotation applied once at the end of SlotsToCoeffs. EvalMod is equivariant under slot rotations, so the CoeffsToSlots correction passes through it unchanged. - Split, not combined, accumulated correction. Apply Aut_{-(slots-1)} at the end of CoeffsToSlots so it always outputs correctly-ordered slots, rather than deferring a single combined correction to the end of SlotsToCoeffs. This fixes EvalFBTNoDecoding + EvalHomDecoding (FBT_CONSECLEV), where a user operation between the two transforms previously saw a residual rotation. - Sparse-packing index reduction. Reduce all BSGS rotation indices modulo min(2*slots, M/4). Under sparse packing the precomputed-plaintext vector has cyclic period 2*slots (the concatenated real/imaginary blocks), not M/4, so indices reduced only to [0, M/4) could otherwise be inconsistent with the period-2*slots plaintext pre-rotations. For full and half packing the modulus equals M/4 and behavior is unchanged. The single-level linear transform (EvalLinearTransform, used when the level budget is 1) is converted to the same single-giant-step Horner form, and FindLinearTransformRotationIndices no longer emits the giant-step keys {2g, 3g, ..., (h-1)g} that the forward form required. Supporting cleanups (no behavior change): - Inline the precomputed rot_in index tables in EvalCoeffsToSlots / EvalSlotsToCoeffs (compute per level, drop the 2D allocation and a redundant scale pass). - Fix the over-large reserve() in the Find*RotationIndices helpers (they reserved ~M entries for a list of a few hundred). - Hoist a redundant KeySwitchExt out of the EvalLinearTransform giant-step loop, and hoist repeated GetParams()/GetElementAtIndex() calls in ExtendCiphertext. - Take crypto parameters by const reference and compute the bootstrap scale factor with std::ldexp. - Add Doxygen for the transform functions and a note in CKKS_BOOTSTRAPPING.md. Not included: the SlotsToCoeffs decode-layer re-partition that moves the remainder group to scale 1 to mirror CoeffsToSlots. The StC remainder still sits at the last (largest-scale) position, so the additional key sharing that re-partition enables (e.g. dropping the StC-specific large-scale remainder key) is not realized here. --- src/pke/examples/CKKS_BOOTSTRAPPING.md | 2 +- src/pke/include/scheme/ckksrns/ckksrns-fhe.h | 23 + src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp | 577 ++++++++----------- 3 files changed, 255 insertions(+), 347 deletions(-) diff --git a/src/pke/examples/CKKS_BOOTSTRAPPING.md b/src/pke/examples/CKKS_BOOTSTRAPPING.md index 466a12501..23c169789 100644 --- a/src/pke/examples/CKKS_BOOTSTRAPPING.md +++ b/src/pke/examples/CKKS_BOOTSTRAPPING.md @@ -8,7 +8,7 @@ The bootstrapping procedure follows a modular framework consisting of four main The original method started with a modulus raising step, as follows: 1. **ModRaise** (Modulus Raising): The ciphertext $(c_0, c_1)\mod q$ is logically viewed as $(c_0, c_1)\mod Q$, for $Q >> q$. This process effectively changes the underlying plaintext from $m(X)$ to $m(X) + q\cdot I(X)$, where $I$ is an integer overflow. Note that the modulus $Q$ needs to be large enough to support the desired number of levels after bootstrapping, as well as the levels required for the bootstrapping procedure. Also note that the overflows are over the coefficients of the ciphertext. -2. **CoeffsToSlots** (Homomorphic Encoding): This step performs a linear transformation to move the coefficients of the ciphertext polynomials into slots, in order to enable the removal of the overflows. In OpenFHE, this is implemented using a matrix-vector multiplication approach (often using a FFT-like butterfly structure) to minimize the number of rotations. Note that this step is expensive since it happens at the largest modulus. +2. **CoeffsToSlots** (Homomorphic Encoding): This step performs a linear transformation to move the coefficients of the ciphertext polynomials into slots, in order to enable the removal of the overflows. In OpenFHE, this is implemented using a matrix-vector multiplication approach (often using a FFT-like butterfly structure) to minimize the number of rotations. The baby-step/giant-step evaluation accumulates the giant steps in Horner form with a single giant-step stride, which also minimizes the number of distinct rotation keys that must be generated and stored. Note that this step is expensive since it happens at the largest modulus. 3. **EvalMod** (Approximate Modular Reduction): This is the core functional evaluation. To recover an encryption of the message $m$ without the overflow, we have to evaluate an approximate modular reduction function homomorphically, using a sine function approximation for small angles. To handle larger ranges of $I$ (determined by the parameter $K$ depending in turn on the secret key distribution), we employ $r$ iterations of the double-angle formula, which reduces the computational complexity at the expense of a higher multiplicative depth. The underlying trigonometric function is computed using Chebyshev series interpolation. See [https://eprint.iacr.org/2019/688](https://eprint.iacr.org/2019/688), [https://eprint.iacr.org/2018/1043](https://eprint.iacr.org/2018/1043), [https://eprint.iacr.org/2020/1203](https://eprint.iacr.org/2020/1203) for details. 4. **SlotsToCoeffs** (Homomorphic Decoding): This step is the inverse of the CoeffsToSlots step and returns the underlying message to the slots encoding, resulting in a refreshed ciphertext $(c_0', c_1')(\bmod Q')$. diff --git a/src/pke/include/scheme/ckksrns/ckksrns-fhe.h b/src/pke/include/scheme/ckksrns/ckksrns-fhe.h index 48b66a9c3..ff9c1c8eb 100644 --- a/src/pke/include/scheme/ckksrns/ckksrns-fhe.h +++ b/src/pke/include/scheme/ckksrns/ckksrns-fhe.h @@ -256,13 +256,36 @@ class FHECKKSRNS : public FHERNS { //------------------------------------------------------------------------------ // EVALUATION: CoeffsToSlots and SlotsToCoeffs //------------------------------------------------------------------------------ + // The transforms below use a baby-step/giant-step (BSGS) decomposition in + // which the giant steps are accumulated in Horner form with a single giant- + // step stride. This requires only one giant-step rotation key per level + // (instead of one per giant step), minimizing the number of distinct + // rotation keys that EvalBootstrapKeyGen must generate and store. + /** + * Single-level linear transform, used for CoeffsToSlots/SlotsToCoeffs when the + * level budget is 1. Evaluated with a BSGS decomposition (Horner giant steps). + * @param A precomputed diagonal plaintexts of the linear transform + * @param ct input ciphertext + */ Ciphertext EvalLinearTransform(const std::vector& A, ConstCiphertext& ct) const; + /** + * Homomorphic encoding (CoeffsToSlots) over multiple BSGS levels, evaluated + * with Horner giant steps to minimize the number of rotation keys. + * @param A precomputed encoding plaintexts, one inner vector per BSGS level + * @param ctxt input ciphertext + */ Ciphertext EvalCoeffsToSlots(const std::vector>& A, ConstCiphertext& ctxt) const; + /** + * Homomorphic decoding (SlotsToCoeffs); the inverse of EvalCoeffsToSlots, + * evaluated with the same Horner giant-step BSGS structure. + * @param A precomputed decoding plaintexts, one inner vector per BSGS level + * @param ctxt input ciphertext + */ Ciphertext EvalSlotsToCoeffs(const std::vector>& A, ConstCiphertext& ctxt) const; diff --git a/src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp b/src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp index 98fc8d0be..09d51012b 100644 --- a/src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp +++ b/src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp @@ -34,8 +34,8 @@ #include "key/evalkeyrelin.h" #include "key/privatekey.h" #include "lattice/lat-hal.h" -#include "math/hal/basicint.h" #include "math/dftransform.h" +#include "math/hal/basicint.h" #include "scheme/ckksrns/ckksrns-cryptoparameters.h" #include "scheme/ckksrns/ckksrns-fhe.h" #include "scheme/ckksrns/ckksrns-utils.h" @@ -50,30 +50,28 @@ #include #include #include -#ifdef BOOTSTRAPTIMING - #include -#endif #include #include #include #include #ifdef BOOTSTRAPTIMING + #include #define PROFILE #endif namespace { // GetBigModulus() calculates the big modulus as the product of // the "compositeDegree" number of parameter modulus -double GetBigModulus(const std::shared_ptr cryptoParams) { - double qDouble = 1.0; +double GetBigModulus(const std::shared_ptr& cryptoParams) { + const auto& params = cryptoParams->GetElementParams()->GetParams(); uint32_t compositeDegree = cryptoParams->GetCompositeDegree(); + double qDouble = 1.0; for (uint32_t j = 0; j < compositeDegree; ++j) { - qDouble *= cryptoParams->GetElementParams()->GetParams()[j]->GetModulus().ConvertToDouble(); + qDouble *= params[j]->GetModulus().ConvertToDouble(); } return qDouble; } - } // namespace namespace lbcrypto { @@ -167,6 +165,7 @@ void FHECKKSRNS::EvalBootstrapSetup(const CryptoContextImpl& cc, std:: // computes all powers of a primitive root of unity exp(2 * M_PI/m) std::vector> ksiPows(m + 1); double ak = 2 * M_PI / m; +#pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(4)) for (uint32_t j = 0; j < m; ++j) { double angle = ak * j; ksiPows[j].real(std::cos(angle)); @@ -192,7 +191,7 @@ void FHECKKSRNS::EvalBootstrapSetup(const CryptoContextImpl& cc, std:: uint32_t compositeDegree = cryptoParams->GetCompositeDegree(); double qDouble = GetBigModulus(cryptoParams); - double factor = static_cast(1) << static_cast(std::round(std::log2(qDouble))); + double factor = std::ldexp(1.0, static_cast(std::round(std::log2(qDouble)))); double pre = (compositeDegree > 1) ? 1.0 : qDouble / factor; double scaleEnc = pre / k; // TODO: YSP Can be extended to FLEXIBLE* scaling techniques as well as the closeness of 2^p to moduli is no longer needed @@ -284,9 +283,8 @@ std::shared_ptr>> FHECKKSRNS::EvalBootstrap (*evalKeys)[M - 1] = ConjugateKeyGen(privateKey); if (cryptoParams->GetSecretKeyDist() == SPARSE_ENCAPSULATED) { - DCRTPoly::TugType tug; - // sparse key used for the modraising step + DCRTPoly::TugType tug; auto skNew = std::make_shared>(cc); skNew->SetPrivateElement(DCRTPoly(tug, cryptoParams->GetElementParams(), Format::EVALUATION, 32)); @@ -304,9 +302,9 @@ std::vector FHECKKSRNS::EvalBootstrapKeyMapIndices(const CryptoContext if (!cryptoParams) OPENFHE_THROW("Invalid crypto parameters: expected CryptoParametersCKKSRNS"); - auto M = cc->GetCyclotomicOrder(); + auto M = cc->GetCyclotomicOrder(); + slots = (slots == 0) ? M / 4 : slots; - slots = (slots == 0) ? M / 4 : slots; const auto bootstrapRotationsIndices = FindBootstrapRotationIndices(slots, M); std::set indexList; @@ -357,6 +355,7 @@ void FHECKKSRNS::EvalBootstrapPrecompute(const CryptoContextImpl& cc, // computes all powers of a primitive root of unity exp(2 * M_PI/m) std::vector> ksiPows(m + 1); double ak = 2 * M_PI / m; +#pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(4)) for (uint32_t j = 0; j < m; ++j) { double angle = ak * j; ksiPows[j].real(std::cos(angle)); @@ -382,7 +381,7 @@ void FHECKKSRNS::EvalBootstrapPrecompute(const CryptoContextImpl& cc, uint32_t compositeDegree = cryptoParams->GetCompositeDegree(); double qDouble = GetBigModulus(cryptoParams); - double factor = static_cast(1) << static_cast(std::round(std::log2(qDouble))); + double factor = std::ldexp(1.0, static_cast(std::round(std::log2(qDouble)))); double pre = (compositeDegree > 1) ? 1.0 : qDouble / factor; double scaleEnc = pre / k; // TODO: YSP Can be extended to FLEXIBLE* scaling techniques as well as the closeness of 2^p to moduli is no longer needed @@ -542,8 +541,8 @@ Ciphertext FHECKKSRNS::EvalBootstrap(ConstCiphertext& cipher if (st == FLEXIBLEAUTOEXT) elementParamsRaised.PopLastParam(); - auto paramsQ = elementParamsRaised.GetParams(); - uint32_t sizeQ = paramsQ.size(); + const auto& paramsQ = elementParamsRaised.GetParams(); + uint32_t sizeQ = paramsQ.size(); std::vector moduli(sizeQ); std::vector roots(sizeQ); for (uint32_t i = 0; i < sizeQ; ++i) { @@ -949,8 +948,8 @@ Ciphertext FHECKKSRNS::EvalBootstrapStCFirst(ConstCiphertext if (st == FLEXIBLEAUTOEXT) elementParamsRaised.PopLastParam(); - auto paramsQ = elementParamsRaised.GetParams(); - uint32_t sizeQ = paramsQ.size(); + const auto& paramsQ = elementParamsRaised.GetParams(); + uint32_t sizeQ = paramsQ.size(); std::vector moduli(sizeQ); std::vector roots(sizeQ); for (uint32_t i = 0; i < sizeQ; ++i) { @@ -1279,24 +1278,14 @@ std::vector FHECKKSRNS::FindBootstrapRotationIndices(uint32_t slots, ui // so it DOES NOT remove possible duplicates and automorphisms corresponding to 0 and M/4. // This method completely depends on FindBootstrapRotationIndices() to do that. std::vector FHECKKSRNS::FindLinearTransformRotationIndices(uint32_t slots, uint32_t M) { - // Computing the baby-step g and the giant-step h. const auto& p = GetBootPrecom(slots); const uint32_t g = (p.m_paramsEnc.g == 0) ? std::ceil(std::sqrt(slots)) : p.m_paramsEnc.g; - const uint32_t h = std::ceil(static_cast(slots) / g); - - // To avoid overflowing uint32_t variables, we do some math operations below in a specific order - // computing all indices for baby-step giant-step procedure - const int32_t indexListSz = static_cast(g) + h + M - 2; - if (indexListSz < 0) - OPENFHE_THROW("indexListSz can not be negative"); std::vector indexList; - indexList.reserve(indexListSz); + indexList.reserve(g + 32); for (uint32_t i = 1; i <= g; ++i) indexList.emplace_back(i); - for (uint32_t i = 2; i < h; ++i) - indexList.emplace_back(i * g); // additional automorphisms are needed for sparse bootstrapping if (uint32_t m = slots * 4; m != M) { @@ -1313,14 +1302,8 @@ std::vector FHECKKSRNS::FindLinearTransformRotationIndices(uint32_t sl std::vector FHECKKSRNS::FindCoeffsToSlotsRotationIndices(uint32_t slots, uint32_t M) { const auto& p = GetBootPrecom(slots).m_paramsEnc; - // To avoid overflowing uint32_t variables, we do some math operations below in a specific order - // Computing all indices for baby-step giant-step procedure for encoding and decoding - const int32_t indexListSz = static_cast(p.b) + p.g - 2 + p.bRem + p.gRem - 2 + 1 + M; - if (indexListSz < 0) - OPENFHE_THROW("indexListSz can not be negative"); - std::vector indexList; - indexList.reserve(indexListSz); + indexList.reserve(p.lvlb * (p.g + 1) + p.gRem + 32); // additional automorphisms are needed for sparse bootstrapping if (uint32_t m = slots * 4; m != M) { @@ -1329,40 +1312,35 @@ std::vector FHECKKSRNS::FindCoeffsToSlotsRotationIndices(uint32_t slot } M >>= 2; - const int32_t flagRem = (p.remCollapse == 0) ? 0 : 1; - const int32_t halfRots = 1 - (p.numRotations + 1) / 2; - const int32_t halfRotsg = halfRots + p.g; + // clamp all indices to [0, 2*slots) + const uint32_t reduceMod = std::min(M, 2 * slots); + const int32_t flagRem = (p.remCollapse == 0) ? 0 : 1; for (int32_t s = -1 + p.lvlb; s >= flagRem; --s) { const uint32_t scalingFactor = 1U << ((s - flagRem) * p.layersCollapse + p.remCollapse); - for (int32_t j = halfRots; j < halfRotsg; ++j) - indexList.emplace_back(ReduceRotation(j * scalingFactor, slots)); - for (uint32_t i = 0; i < p.b; ++i) - indexList.emplace_back(ReduceRotation(i * p.g * scalingFactor, M)); + for (uint32_t j = 1; j < p.g; ++j) + indexList.emplace_back(ReduceRotation(j * scalingFactor, reduceMod)); + if (p.b > 1) + indexList.emplace_back(ReduceRotation(p.g * scalingFactor, reduceMod)); } if (flagRem == 1) { - const int32_t halfRotsRem = (1 - (p.numRotationsRem + 1) / 2); - const int32_t halfRotsRemg = halfRotsRem + p.gRem; - for (int32_t j = halfRotsRem; j < halfRotsRemg; ++j) - indexList.emplace_back(ReduceRotation(j, slots)); - for (uint32_t i = 0; i < p.bRem; ++i) - indexList.emplace_back(ReduceRotation(i * p.gRem, M)); + for (uint32_t j = 1; j < p.gRem; ++j) + indexList.emplace_back(ReduceRotation(j, reduceMod)); + if (p.bRem > 1) + indexList.emplace_back(ReduceRotation(p.gRem, reduceMod)); } + // StC correction key Aut_{-(slots-1)} (CtS correction is now absorbed into plaintexts) + indexList.emplace_back(ReduceRotation(-static_cast(slots - 1), reduceMod)); + return indexList; } std::vector FHECKKSRNS::FindSlotsToCoeffsRotationIndices(uint32_t slots, uint32_t M) { const auto& p = GetBootPrecom(slots).m_paramsDec; - // To avoid overflowing uint32_t variables, we do some math operations below in a specific order - // Computing all indices for baby-step giant-step procedure for encoding and decoding - const int32_t indexListSz = static_cast(p.b) + p.g - 2 + p.bRem + p.gRem - 2 + 1 + M; - if (indexListSz < 0) - OPENFHE_THROW("indexListSz can not be negative"); - std::vector indexList; - indexList.reserve(indexListSz); + indexList.reserve(p.lvlb * (p.g + 1) + p.gRem + 32); // additional automorphisms are needed for sparse bootstrapping if (uint32_t m = slots * 4; m != M) { @@ -1371,26 +1349,24 @@ std::vector FHECKKSRNS::FindSlotsToCoeffsRotationIndices(uint32_t slot } M >>= 2; - const uint32_t flagRem = (p.remCollapse == 0) ? 0 : 1; - const uint32_t smax = p.lvlb - flagRem; - const int32_t halfRots = (1 - (p.numRotations + 1) / 2); - const int32_t halfRotsg = halfRots + p.g; + // clamp all indices to [0, 2*slots) + const uint32_t reduceMod = std::min(M, 2 * slots); + const uint32_t flagRem = (p.remCollapse == 0) ? 0 : 1; + const uint32_t smax = p.lvlb - flagRem; for (uint32_t s = 0; s < smax; ++s) { const uint32_t scalingFactor = 1U << (s * p.layersCollapse); - for (int32_t j = halfRots; j < halfRotsg; ++j) - indexList.emplace_back(ReduceRotation(j * scalingFactor, M)); - for (uint32_t i = 0; i < p.b; ++i) - indexList.emplace_back(ReduceRotation(i * p.g * scalingFactor, M)); + for (uint32_t j = 1; j < p.g; ++j) + indexList.emplace_back(ReduceRotation(j * scalingFactor, reduceMod)); + if (p.b > 1) + indexList.emplace_back(ReduceRotation(p.g * scalingFactor, reduceMod)); } if (flagRem == 1) { const uint32_t scalingFactor = 1U << (smax * p.layersCollapse); - const int32_t halfRotsRem = (1 - (p.numRotationsRem + 1) / 2); - const int32_t halfRotsRemg = halfRotsRem + p.gRem; - for (int32_t j = halfRotsRem; j < halfRotsRemg; ++j) - indexList.emplace_back(ReduceRotation(j * scalingFactor, M)); - for (uint32_t i = 0; i < p.bRem; ++i) - indexList.emplace_back(ReduceRotation(i * p.gRem * scalingFactor, M)); + for (uint32_t j = 1; j < p.gRem; ++j) + indexList.emplace_back(ReduceRotation(j * scalingFactor, reduceMod)); + if (p.bRem > 1) + indexList.emplace_back(ReduceRotation(p.gRem * scalingFactor, reduceMod)); } return indexList; @@ -1415,10 +1391,10 @@ std::vector FHECKKSRNS::EvalLinearTransformPrecompute( for (uint32_t i = 0; i < towersToDrop; ++i) elementParams.PopLastParam(); - auto paramsQ = elementParams.GetParams(); - uint32_t sizeQ = paramsQ.size(); - auto paramsP = cryptoParams->GetParamsP()->GetParams(); - uint32_t sizeP = paramsP.size(); + const auto& paramsQ = elementParams.GetParams(); + uint32_t sizeQ = paramsQ.size(); + const auto& paramsP = cryptoParams->GetParamsP()->GetParams(); + uint32_t sizeP = paramsP.size(); std::vector moduli(sizeQ + sizeP); std::vector roots(sizeQ + sizeP); for (uint32_t i = 0; i < sizeQ; ++i) { @@ -1461,10 +1437,10 @@ std::vector FHECKKSRNS::EvalLinearTransformPrecompute( for (uint32_t i = 0; i < towersToDrop; ++i) elementParams.PopLastParam(); - auto paramsQ = elementParams.GetParams(); - uint32_t sizeQ = paramsQ.size(); - auto paramsP = cryptoParams->GetParamsP()->GetParams(); - uint32_t sizeP = paramsP.size(); + const auto& paramsQ = elementParams.GetParams(); + uint32_t sizeQ = paramsQ.size(); + const auto& paramsP = cryptoParams->GetParamsP()->GetParams(); + uint32_t sizeP = paramsP.size(); std::vector moduli(sizeQ + sizeP); std::vector roots(sizeQ + sizeP); for (uint32_t i = 0; i < sizeQ; ++i) { @@ -1560,11 +1536,10 @@ std::vector> FHECKKSRNS::EvalCoeffsToSlotsPrecomp uint32_t level0 = towersToDrop + compositeDegree * (p.lvlb - 1); - auto paramsQ = elementParams.GetParams(); - uint32_t sizeQ = paramsQ.size(); - auto paramsP = cryptoParams->GetParamsP()->GetParams(); - uint32_t sizeP = paramsP.size(); - + const auto& paramsQ = elementParams.GetParams(); + uint32_t sizeQ = paramsQ.size(); + const auto& paramsP = cryptoParams->GetParamsP()->GetParams(); + uint32_t sizeP = paramsP.size(); std::vector moduli(sizeQ + sizeP); std::vector roots(sizeQ + sizeP); for (uint32_t i = 0; i < sizeQ; ++i) { @@ -1579,14 +1554,19 @@ std::vector> FHECKKSRNS::EvalCoeffsToSlotsPrecomp // we need to pre-compute the plaintexts in the extended basis P*Q uint32_t M = cc.GetCyclotomicOrder(); std::vector>> paramsVector(p.lvlb - stop); - for (int32_t s = -1 + p.lvlb; s >= stop; --s) { + for (int32_t s = p.lvlb - 1; s >= stop; --s) { paramsVector[s - stop] = std::make_shared>(M, moduli, roots); - for (uint32_t j = 0; j < compositeDegree; ++j, --sizeQ) { - moduli.erase(moduli.begin() + sizeQ - 1); - roots.erase(roots.begin() + sizeQ - 1); - } + sizeQ -= compositeDegree; + moduli.erase(moduli.begin() + sizeQ, moduli.begin() + sizeQ + compositeDegree); + roots.erase(roots.begin() + sizeQ, roots.begin() + sizeQ + compositeDegree); } + // zero-based inner indices require shifting the pre-rotation by +offset*scale + const int32_t offset = static_cast((p.numRotations + 1) / 2) - 1; + const int32_t offsetRem = static_cast((p.numRotationsRem + 1) / 2) - 1; + // running accumulated correction; plaintexts pre-absorb the per-level Aut_{-delta_s} + int32_t runningAcc = 0; + if (uint32_t M4 = M / 4; slots == M4 || flagStCComplex) { //------------------------------------------------------------------------------ // fully-packed mode @@ -1595,25 +1575,25 @@ std::vector> FHECKKSRNS::EvalCoeffsToSlotsPrecomp auto coeff = CoeffEncodingCollapse(A, rotGroup, p.lvlb, flag_i); for (int32_t s = -1 + p.lvlb; s > stop; --s) { - const int32_t rotScale = (1 << ((s - flagRem) * p.layersCollapse + p.remCollapse)) * p.g; - const uint32_t limit = p.b * p.g; + const int32_t shiftScale = 1 << ((s - flagRem) * p.layersCollapse + p.remCollapse); + const int32_t rotScale = shiftScale * static_cast(p.g); + const uint32_t limit = p.b * p.g; #if !defined(__MINGW32__) && !defined(__MINGW64__) #pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(limit)) #endif for (uint32_t ij = 0; ij < limit; ++ij) { if (ij != p.numRotations) { if ((flagRem == 0) && (s == stop + 1)) { - // do the scaling only at the last set of coefficients for (auto& c : coeff[s][ij]) c *= scale; } - auto rot = Rotate(coeff[s][ij], ReduceRotation(-rotScale * (ij / p.g), slots)); + auto rot = Rotate(coeff[s][ij], ReduceRotation(-rotScale * (ij / p.g) + offset * shiftScale + runningAcc, slots)); - result[s][ij] = - MakeAuxPlaintext(cc, paramsVector[s - stop], rot, 1, level0 - compositeDegree * s, rot.size()); + result[s][ij] = MakeAuxPlaintext(cc, paramsVector[s - stop], rot, 1, level0 - compositeDegree * s, rot.size()); } } + runningAcc += offset * shiftScale; } if (flagRem == 1) { @@ -1626,7 +1606,7 @@ std::vector> FHECKKSRNS::EvalCoeffsToSlotsPrecomp for (auto& c : coeff[stop][ij]) c *= scale; - auto rot = Rotate(coeff[stop][ij], ReduceRotation(-p.gRem * (ij / p.gRem), slots)); + auto rot = Rotate(coeff[stop][ij], ReduceRotation(-static_cast(p.gRem) * (ij / p.gRem) + offsetRem + runningAcc, slots)); result[stop][ij] = MakeAuxPlaintext(cc, paramsVector[0], rot, 1, level0, rot.size()); } @@ -1642,8 +1622,9 @@ std::vector> FHECKKSRNS::EvalCoeffsToSlotsPrecomp auto coeffi = CoeffEncodingCollapse(A, rotGroup, p.lvlb, true); for (int32_t s = -1 + p.lvlb; s > stop; --s) { - const int32_t rotScale = (1 << ((s - flagRem) * p.layersCollapse + p.remCollapse)) * p.g; - const uint32_t limit = p.b * p.g; + const int32_t shiftScale = 1 << ((s - flagRem) * p.layersCollapse + p.remCollapse); + const int32_t rotScale = shiftScale * static_cast(p.g); + const uint32_t limit = p.b * p.g; #if !defined(__MINGW32__) && !defined(__MINGW64__) #pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(limit)) #endif @@ -1654,17 +1635,16 @@ std::vector> FHECKKSRNS::EvalCoeffsToSlotsPrecomp auto& clearTmpi = coeffi[s][ij]; clearTmp.insert(clearTmp.end(), clearTmpi.begin(), clearTmpi.end()); if ((flagRem == 0) && (s == stop + 1)) { - // do the scaling only at the last set of coefficients for (auto& c : clearTmp) c *= scale; } - auto rot = Rotate(clearTmp, ReduceRotation(-rotScale * (ij / p.g), M4)); + auto rot = Rotate(clearTmp, ReduceRotation(-rotScale * (ij / p.g) + offset * shiftScale + runningAcc, 2 * slots)); - result[s][ij] = - MakeAuxPlaintext(cc, paramsVector[s - stop], rot, 1, level0 - compositeDegree * s, rot.size()); + result[s][ij] = MakeAuxPlaintext(cc, paramsVector[s - stop], rot, 1, level0 - compositeDegree * s, rot.size()); } } + runningAcc += offset * shiftScale; } if (flagRem == 1) { @@ -1681,7 +1661,7 @@ std::vector> FHECKKSRNS::EvalCoeffsToSlotsPrecomp for (auto& c : clearTmp) c *= scale; - auto rot = Rotate(clearTmp, ReduceRotation(-p.gRem * (ij / p.gRem), M4)); + auto rot = Rotate(clearTmp, ReduceRotation(-p.gRem * (ij / p.gRem) + offsetRem + runningAcc, 2 * slots)); result[stop][ij] = MakeAuxPlaintext(cc, paramsVector[0], rot, 1, level0, rot.size()); } @@ -1717,10 +1697,10 @@ std::vector> FHECKKSRNS::EvalSlotsToCoeffsPrecomp for (uint32_t i = 0; i < towersToDrop; ++i) elementParams.PopLastParam(); - auto paramsQ = elementParams.GetParams(); - uint32_t sizeQ = paramsQ.size(); - auto paramsP = cryptoParams->GetParamsP()->GetParams(); - uint32_t sizeP = paramsP.size(); + const auto& paramsQ = elementParams.GetParams(); + uint32_t sizeQ = paramsQ.size(); + const auto& paramsP = cryptoParams->GetParamsP()->GetParams(); + uint32_t sizeP = paramsP.size(); std::vector moduli(sizeQ + sizeP); std::vector roots(sizeQ + sizeP); for (uint32_t i = 0; i < sizeQ; ++i) { @@ -1743,35 +1723,42 @@ std::vector> FHECKKSRNS::EvalSlotsToCoeffsPrecomp } } + // zero-based inner indices require shifting the pre-rotation by +offset*shiftScale + const int32_t offset = static_cast((p.numRotations + 1) / 2) - 1; + const int32_t offsetRem = static_cast((p.numRotationsRem + 1) / 2) - 1; + // CtS output now carries no extra rotation (absorbed into CtS precompute), so start at 0 + int32_t runningAcc = 0; + if (uint32_t M4 = cc.GetCyclotomicOrder() / 4; M4 == slots || flagStCComplex) { - // fully-packed - auto coeff = CoeffDecodingCollapse(A, rotGroup, p.lvlb, flag_i); + auto coeff = CoeffDecodingCollapse(A, rotGroup, p.lvlb, flag_i); + const uint32_t smax = p.lvlb - flagRem; for (uint32_t s = 0; s < smax; ++s) { - const int32_t rotScale = (1 << (s * p.layersCollapse)) * p.g; - const uint32_t limit = p.b * p.g; + const int32_t shiftScale = 1 << (s * p.layersCollapse); + const int32_t rotScale = shiftScale * static_cast(p.g); + const uint32_t limit = p.b * p.g; #if !defined(__MINGW32__) && !defined(__MINGW64__) #pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(limit)) #endif for (uint32_t ij = 0; ij < limit; ++ij) { if (ij != p.numRotations) { if ((flagRem == 0) && (s + 1 == smax)) { - // do the scaling only at the last set of coefficients for (auto& c : coeff[s][ij]) c *= scale; } - auto rot = Rotate(coeff[s][ij], ReduceRotation(-rotScale * (ij / p.g), slots)); + auto rot = Rotate(coeff[s][ij], ReduceRotation(-rotScale * (ij / p.g) + offset * shiftScale + runningAcc, slots)); - result[s][ij] = - MakeAuxPlaintext(cc, paramsVector[s], rot, 1, towersToDrop + compositeDegree * s, rot.size()); + result[s][ij] = MakeAuxPlaintext(cc, paramsVector[s], rot, 1, towersToDrop + compositeDegree * s, rot.size()); } } + runningAcc += offset * shiftScale; } if (flagRem == 1) { - const int32_t rotScale = (1 << (smax * p.layersCollapse)) * p.gRem; - const uint32_t limit = p.bRem * p.gRem; + const int32_t shiftScaleRem = 1 << (smax * p.layersCollapse); + const int32_t rotScale = shiftScaleRem * static_cast(p.gRem); + const uint32_t limit = p.bRem * p.gRem; #if !defined(__MINGW32__) && !defined(__MINGW64__) #pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(limit)) #endif @@ -1780,10 +1767,9 @@ std::vector> FHECKKSRNS::EvalSlotsToCoeffsPrecomp for (auto& c : coeff[smax][ij]) c *= scale; - auto rot = Rotate(coeff[smax][ij], ReduceRotation(-rotScale * (ij / p.gRem), slots)); + auto rot = Rotate(coeff[smax][ij], ReduceRotation(-rotScale * (ij / p.gRem) + offsetRem * shiftScaleRem + runningAcc, slots)); - result[smax][ij] = MakeAuxPlaintext(cc, paramsVector[smax], rot, 1, - towersToDrop + compositeDegree * smax, rot.size()); + result[smax][ij] = MakeAuxPlaintext(cc, paramsVector[smax], rot, 1, towersToDrop + compositeDegree * smax, rot.size()); } } } @@ -1798,8 +1784,9 @@ std::vector> FHECKKSRNS::EvalSlotsToCoeffsPrecomp const uint32_t smax = p.lvlb - flagRem; for (uint32_t s = 0; s < smax; ++s) { - const int32_t rotScale = (1 << (s * p.layersCollapse)) * p.g; - const uint32_t limit = p.b * p.g; + const int32_t shiftScale = 1 << (s * p.layersCollapse); + const int32_t rotScale = shiftScale * static_cast(p.g); + const uint32_t limit = p.b * p.g; #if !defined(__MINGW32__) && !defined(__MINGW64__) #pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(limit)) #endif @@ -1815,17 +1802,18 @@ std::vector> FHECKKSRNS::EvalSlotsToCoeffsPrecomp c *= scale; } - auto rot = Rotate(clearTmp, ReduceRotation(-rotScale * (ij / p.g), M4)); + auto rot = Rotate(clearTmp, ReduceRotation(-rotScale * (ij / p.g) + offset * shiftScale + runningAcc, 2 * slots)); - result[s][ij] = - MakeAuxPlaintext(cc, paramsVector[s], rot, 1, towersToDrop + compositeDegree * s, rot.size()); + result[s][ij] = MakeAuxPlaintext(cc, paramsVector[s], rot, 1, towersToDrop + compositeDegree * s, rot.size()); } } + runningAcc += offset * shiftScale; } if (flagRem == 1) { - const int32_t rotScale = (1 << (smax * p.layersCollapse)) * p.g; - const uint32_t limit = p.bRem * p.gRem; + const int32_t shiftScaleRem = 1 << (smax * p.layersCollapse); + const int32_t rotScale = shiftScaleRem * static_cast(p.g); + const uint32_t limit = p.bRem * p.gRem; #if !defined(__MINGW32__) && !defined(__MINGW64__) #pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(limit)) #endif @@ -1838,10 +1826,9 @@ std::vector> FHECKKSRNS::EvalSlotsToCoeffsPrecomp for (auto& c : clearTmp) c *= scale; - auto rot = Rotate(clearTmp, ReduceRotation(-rotScale * (ij / p.gRem), M4)); + auto rot = Rotate(clearTmp, ReduceRotation(-rotScale * (ij / p.gRem) + offsetRem * shiftScaleRem + runningAcc, 2 * slots)); - result[smax][ij] = MakeAuxPlaintext(cc, paramsVector[smax], rot, 1, - towersToDrop + compositeDegree * smax, rot.size()); + result[smax][ij] = MakeAuxPlaintext(cc, paramsVector[smax], rot, 1, towersToDrop + compositeDegree * smax, rot.size()); } } } @@ -1856,53 +1843,41 @@ std::vector> FHECKKSRNS::EvalSlotsToCoeffsPrecomp Ciphertext FHECKKSRNS::EvalLinearTransform(const std::vector& A, ConstCiphertext& ct) const { // Computing the baby-step bStep and the giant-step gStep. - const uint32_t slots = A.size(); - const auto& p = GetBootPrecom(slots); - const uint32_t bStep = (p.m_paramsEnc.g == 0) ? std::ceil(std::sqrt(slots)) : p.m_paramsEnc.g; - const uint32_t gStep = std::ceil(static_cast(slots) / bStep); + const int32_t slots = A.size(); + const auto& p = GetBootPrecom(slots); + const int32_t bStep = (p.m_paramsEnc.g == 0) ? std::ceil(std::sqrt(slots)) : p.m_paramsEnc.g; + const int32_t gStep = std::ceil(static_cast(slots) / bStep); auto cc = ct->GetCryptoContext(); auto digits = cc->EvalFastRotationPrecompute(ct); - // hoisted automorphisms std::vector> fastRotation(bStep - 1); #pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(bStep - 1)) - for (uint32_t j = 1; j < bStep; ++j) + for (int32_t j = 1; j < bStep; ++j) fastRotation[j - 1] = cc->EvalFastRotationExt(ct, j, digits, true); - const uint32_t M = cc->GetCyclotomicOrder(); - const uint32_t N = cc->GetRingDimension(); - std::vector map(N); - Ciphertext result; - DCRTPoly first; - for (uint32_t j = 0; j < gStep; ++j) { - auto inner = EvalMultExt(cc->KeySwitchExt(ct, true), A[bStep * j]); - for (uint32_t i = 1; i < bStep; ++i) { - if (bStep * j + i < slots) - EvalAddExtInPlace(inner, EvalMultExt(fastRotation[i - 1], A[bStep * j + i])); - } - - if (j == 0) { - first = cc->KeySwitchDownFirstElement(inner); - auto elements = inner->GetElements(); - elements[0].SetValuesToZero(); - inner->SetElements(std::move(elements)); - result = std::move(inner); + auto ctExt = cc->KeySwitchExt(ct, true); + + // Horner backward accumulation with a single giant stride t = bStep: + Ciphertext outer; + for (int32_t j = gStep - 1; j >= 0; --j) { + const int32_t G = j * bStep; + auto inner = EvalMultExt(ctExt, A[G]); + for (int32_t i = 1; i < bStep; ++i) + if (G + i < slots) + EvalAddExtInPlace(inner, EvalMultExt(fastRotation[i - 1], A[G + i])); + + if (j == gStep - 1) { + outer = std::move(inner); } else { - inner = cc->KeySwitchDown(inner); - // Find the automorphism index that corresponds to rotation index index. - uint32_t autoIndex = FindAutomorphismIndex2nComplex(bStep * j, M); - PrecomputeAutoMap(N, autoIndex, &map); - first += inner->GetElements()[0].AutomorphismTransform(autoIndex, map); - - auto&& innerDigits = cc->EvalFastRotationPrecompute(inner); - EvalAddExtInPlace(result, cc->EvalFastRotationExt(inner, bStep * j, innerDigits, false)); + auto outerStd = cc->KeySwitchDown(outer); + auto outerDigits = cc->EvalFastRotationPrecompute(outerStd); + outer = cc->EvalFastRotationExt(outerStd, bStep, outerDigits, true); + EvalAddExtInPlace(outer, inner); } } - result = cc->KeySwitchDown(result); - result->GetElements()[0] += first; - return result; + return cc->KeySwitchDown(outer); } Ciphertext FHECKKSRNS::EvalCoeffsToSlots(const std::vector>& A, @@ -1911,46 +1886,20 @@ Ciphertext FHECKKSRNS::EvalCoeffsToSlots(const std::vector> rot_out(p.lvlb, std::vector(p.b + p.bRem)); - std::vector> rot_in(p.lvlb, std::vector(p.numRotations + 1)); - int32_t stop = -1; int32_t flagRem = 0; if (p.remCollapse != 0) { stop = 0; flagRem = 1; - - // remainder corresponds to index 0 in encoding and to last index in decoding - rot_in[0].resize(p.numRotationsRem + 1); } auto cc = ctxt->GetCryptoContext(); - const uint32_t M4 = cc->GetCyclotomicOrder() / 4; - - int32_t offset = static_cast((p.numRotations + 1) / 2) - 1; - for (int32_t s = p.lvlb - 1; s > stop; --s) { - int32_t scale = (1 << ((s - flagRem) * p.layersCollapse + p.remCollapse)); - for (uint32_t i = 0; i < p.b; ++i) - rot_out[s][i] = ReduceRotation(scale * p.g * i, M4); - for (uint32_t j = 0; j < p.g; ++j) - rot_in[s][j] = ReduceRotation(scale * (j - offset), slots); - } - - if (flagRem == 1) { - offset = static_cast((p.numRotationsRem + 1) / 2) - 1; - for (uint32_t i = 0; i < p.bRem; ++i) - rot_out[stop][i] = ReduceRotation(p.gRem * i, M4); - for (uint32_t j = 0; j < p.gRem; ++j) - rot_in[stop][j] = ReduceRotation(j - offset, slots); - } + const uint32_t M4 = cc->GetCyclotomicOrder() / 4; + const uint32_t reduceMod = std::min(M4, 2 * slots); auto result = ctxt->Clone(); - uint32_t N = cc->GetRingDimension(); - std::vector map(N); - auto algo = cc->GetScheme(); const auto cryptoParams = std::dynamic_pointer_cast(cc->GetCryptoParameters()); uint32_t compositeDegree = cryptoParams->GetCompositeDegree(); @@ -1961,104 +1910,86 @@ Ciphertext FHECKKSRNS::EvalCoeffsToSlots(const std::vectorModReduceInternalInPlace(result, compositeDegree); - // computes the NTTs for each CRT limb (for the hoisted automorphisms used later on) + const int32_t scale = 1 << ((s - flagRem) * p.layersCollapse + p.remCollapse); + auto digits = cc->EvalFastRotationPrecompute(result); std::vector> fastRotation(p.g); #pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(p.g)) for (uint32_t j = 0; j < p.g; ++j) - fastRotation[j] = (rot_in[s][j] != 0) ? cc->EvalFastRotationExt(result, rot_in[s][j], digits, true) : - cc->KeySwitchExt(result, true); + fastRotation[j] = (j != 0) ? + cc->EvalFastRotationExt(result, ReduceRotation(scale * j, reduceMod), digits, true) : + cc->KeySwitchExt(result, true); + + // Horner backward accumulation with single giant stride t = g * scale + const int32_t t = ReduceRotation(scale * p.g, reduceMod); Ciphertext outer; - DCRTPoly first; - for (uint32_t i = 0; i < p.b; ++i) { - // for the first iteration with j=0: - uint32_t G = p.g * i; + for (int32_t i = static_cast(p.b) - 1; i >= 0; --i) { + uint32_t G = p.g * static_cast(i); auto inner = EvalMultExt(fastRotation[0], A[s][G]); - // continue the loop for (uint32_t j = 1; j < p.g; ++j) { if ((G + j) != p.numRotations) EvalAddExtInPlace(inner, EvalMultExt(fastRotation[j], A[s][G + j])); } - if (i == 0) { - first = cc->KeySwitchDownFirstElement(inner); + if (i == static_cast(p.b) - 1) { outer = std::move(inner); - outer->GetElements()[0].SetValuesToZero(); } else { - if (rot_out[s][i] != 0) { - inner = cc->KeySwitchDown(inner); - // Find the automorphism index that corresponds to rotation index index. - uint32_t autoIndex = FindAutomorphismIndex2nComplex(rot_out[s][i], cc->GetCyclotomicOrder()); - PrecomputeAutoMap(N, autoIndex, &map); - first += inner->GetElements()[0].AutomorphismTransform(autoIndex, map); - auto&& innerDigits = cc->EvalFastRotationPrecompute(inner); - EvalAddExtInPlace(outer, cc->EvalFastRotationExt(inner, rot_out[s][i], innerDigits, false)); - } - else { - first += cc->KeySwitchDownFirstElement(inner); - auto& elements = inner->GetElements(); - elements[0].SetValuesToZero(); - EvalAddExtInPlace(outer, inner); + if (t != 0) { + auto outerStd = cc->KeySwitchDown(outer); + auto outerDigits = cc->EvalFastRotationPrecompute(outerStd); + outer = cc->EvalFastRotationExt(outerStd, t, outerDigits, true); } + EvalAddExtInPlace(outer, inner); } } result = cc->KeySwitchDown(outer); - result->GetElements()[0] += first; } if (flagRem == 1) { algo->ModReduceInternalInPlace(result, compositeDegree); - // computes the NTTs for each CRT limb (for the hoisted automorphisms used later on) auto digits = cc->EvalFastRotationPrecompute(result); std::vector> fastRotationRem(p.gRem); #pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(p.gRem)) for (uint32_t j = 0; j < p.gRem; ++j) - fastRotationRem[j] = (rot_in[stop][j] != 0) ? - cc->EvalFastRotationExt(result, rot_in[stop][j], digits, true) : - cc->KeySwitchExt(result, true); + fastRotationRem[j] = + (j != 0) ? cc->EvalFastRotationExt(result, ReduceRotation(j, reduceMod), digits, true) : + cc->KeySwitchExt(result, true); + + // Remainder scale is 1 in CtS; Horner giant stride = gRem + const int32_t tRem = ReduceRotation(p.gRem, reduceMod); Ciphertext outer; - DCRTPoly first; - for (uint32_t i = 0; i < p.bRem; ++i) { - // for the first iteration with j=0: - int32_t GRem = p.gRem * i; - auto inner = EvalMultExt(fastRotationRem[0], A[stop][GRem]); - // continue the loop + for (int32_t i = static_cast(p.bRem) - 1; i >= 0; --i) { + uint32_t GRem = p.gRem * static_cast(i); + auto inner = EvalMultExt(fastRotationRem[0], A[stop][GRem]); for (uint32_t j = 1; j < p.gRem; ++j) { if ((GRem + j) != p.numRotationsRem) EvalAddExtInPlace(inner, EvalMultExt(fastRotationRem[j], A[stop][GRem + j])); } - if (i == 0) { - first = cc->KeySwitchDownFirstElement(inner); + if (i == static_cast(p.bRem) - 1) { outer = std::move(inner); - outer->GetElements()[0].SetValuesToZero(); } else { - if (rot_out[stop][i] != 0) { - inner = cc->KeySwitchDown(inner); - // Find the automorphism index that corresponds to rotation index index. - uint32_t autoIndex = FindAutomorphismIndex2nComplex(rot_out[stop][i], cc->GetCyclotomicOrder()); - PrecomputeAutoMap(N, autoIndex, &map); - first += inner->GetElements()[0].AutomorphismTransform(autoIndex, map); - auto&& innerDigits = cc->EvalFastRotationPrecompute(inner); - EvalAddExtInPlace(outer, cc->EvalFastRotationExt(inner, rot_out[stop][i], innerDigits, false)); - } - else { - first += cc->KeySwitchDownFirstElement(inner); - auto elements = inner->GetElements(); - elements[0].SetValuesToZero(); - inner->SetElements(std::move(elements)); - EvalAddExtInPlace(outer, inner); + if (tRem != 0) { + auto outerStd = cc->KeySwitchDown(outer); + auto outerDigits = cc->EvalFastRotationPrecompute(outerStd); + outer = cc->EvalFastRotationExt(outerStd, tRem, outerDigits, true); } + EvalAddExtInPlace(outer, inner); } } result = cc->KeySwitchDown(outer); - result->GetElements()[0] += first; } + + // undo the per-level accumulated rotation so the output is in standard slot ordering. + const int32_t CtSDelta = ReduceRotation(-static_cast(slots - 1), reduceMod); + if (CtSDelta != 0) + result = cc->EvalAtIndex(result, CtSDelta); + return result; } @@ -2068,155 +1999,104 @@ Ciphertext FHECKKSRNS::EvalSlotsToCoeffs(const std::vector> rot_out(p.lvlb, std::vector(p.b + p.bRem)); - std::vector> rot_in(p.lvlb, std::vector(p.numRotations + 1)); const int32_t flagRem = (p.remCollapse == 0) ? 0 : 1; - if (flagRem == 1) { - // remainder corresponds to index 0 in encoding and to last index in decoding - rot_in[p.lvlb - 1].resize(p.numRotationsRem + 1); - } auto cc = ctxt->GetCryptoContext(); - const uint32_t M4 = cc->GetCyclotomicOrder() / 4; - const int32_t smax = p.lvlb - flagRem; - const int32_t offset = static_cast((p.numRotations + 1) / 2) - 1; - for (int32_t s = 0; s < smax; ++s) { - const int32_t scale = 1 << (s * p.layersCollapse); - for (uint32_t j = 0; j < p.g; ++j) - rot_in[s][j] = ReduceRotation((j - offset) * scale, M4); - for (uint32_t i = 0; i < p.b; ++i) - rot_out[s][i] = ReduceRotation((p.g * i) * scale, M4); - } - - if (flagRem == 1) { - const int32_t scaleRem = 1 << (smax * p.layersCollapse); - const int32_t offsetRem = static_cast((p.numRotationsRem + 1) / 2) - 1; - for (uint32_t j = 0; j < p.gRem; ++j) - rot_in[smax][j] = ReduceRotation((j - offsetRem) * scaleRem, M4); - for (uint32_t i = 0; i < p.bRem; ++i) - rot_out[smax][i] = ReduceRotation((p.gRem * i) * scaleRem, M4); - } + const uint32_t M4 = cc->GetCyclotomicOrder() / 4; + const uint32_t reduceMod = std::min(M4, 2 * slots); + const int32_t smax = p.lvlb - flagRem; - // No need for Encrypted Bit Reverse auto result = ctxt->Clone(); - uint32_t N = cc->GetRingDimension(); - std::vector map(N); - auto algo = cc->GetScheme(); const auto cryptoParams = std::dynamic_pointer_cast(cc->GetCryptoParameters()); uint32_t compositeDegree = cryptoParams->GetCompositeDegree(); - // hoisted automorphisms for (int32_t s = 0; s < smax; ++s) { if (s != 0) algo->ModReduceInternalInPlace(result, compositeDegree); - // computes the NTTs for each CRT limb (for the hoisted automorphisms used later on) + const int32_t scale = 1 << (s * p.layersCollapse); + auto digits = cc->EvalFastRotationPrecompute(result); std::vector> fastRotation(p.g); #pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(p.g)) for (uint32_t j = 0; j < p.g; ++j) - fastRotation[j] = (rot_in[s][j] != 0) ? cc->EvalFastRotationExt(result, rot_in[s][j], digits, true) : - cc->KeySwitchExt(result, true); + fastRotation[j] = (j != 0) ? + cc->EvalFastRotationExt(result, ReduceRotation(scale * j, reduceMod), digits, true) : + cc->KeySwitchExt(result, true); + + // Horner backward accumulation with single giant stride t = g * scale + const int32_t t = ReduceRotation(scale * p.g, reduceMod); Ciphertext outer; - DCRTPoly first; - for (uint32_t i = 0; i < p.b; ++i) { - // for the first iteration with j=0: - uint32_t G = i * p.g; + for (int32_t i = static_cast(p.b) - 1; i >= 0; --i) { + uint32_t G = static_cast(i) * p.g; auto inner = EvalMultExt(fastRotation[0], A[s][G]); - // continue the loop for (uint32_t j = 1; j < p.g; ++j) { if ((G + j) != p.numRotations) EvalAddExtInPlace(inner, EvalMultExt(fastRotation[j], A[s][G + j])); } - if (i == 0) { - first = cc->KeySwitchDownFirstElement(inner); - auto elements = inner->GetElements(); - elements[0].SetValuesToZero(); - inner->SetElements(std::move(elements)); + if (i == static_cast(p.b) - 1) { outer = std::move(inner); } else { - if (rot_out[s][i] != 0) { - inner = cc->KeySwitchDown(inner); - // Find the automorphism index that corresponds to rotation index index. - auto autoIndex = FindAutomorphismIndex2nComplex(rot_out[s][i], cc->GetCyclotomicOrder()); - PrecomputeAutoMap(N, autoIndex, &map); - first += inner->GetElements()[0].AutomorphismTransform(autoIndex, map); - auto&& innerDigits = cc->EvalFastRotationPrecompute(inner); - EvalAddExtInPlace(outer, cc->EvalFastRotationExt(inner, rot_out[s][i], innerDigits, false)); - } - else { - first += cc->KeySwitchDownFirstElement(inner); - auto elements = inner->GetElements(); - elements[0].SetValuesToZero(); - inner->SetElements(std::move(elements)); - EvalAddExtInPlace(outer, inner); + if (t != 0) { + auto outerStd = cc->KeySwitchDown(outer); + auto outerDigits = cc->EvalFastRotationPrecompute(outerStd); + outer = cc->EvalFastRotationExt(outerStd, t, outerDigits, true); } + EvalAddExtInPlace(outer, inner); } } result = cc->KeySwitchDown(outer); - result->GetElements()[0] += first; } if (flagRem == 1) { algo->ModReduceInternalInPlace(result, compositeDegree); - // computes the NTTs for each CRT limb (for the hoisted automorphisms used later on) + const int32_t scaleRem = 1 << (smax * p.layersCollapse); + auto digits = cc->EvalFastRotationPrecompute(result); std::vector> fastRotationRem(p.gRem); #pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(p.gRem)) for (uint32_t j = 0; j < p.gRem; ++j) - fastRotationRem[j] = (rot_in[smax][j] != 0) ? - cc->EvalFastRotationExt(result, rot_in[smax][j], digits, true) : - cc->KeySwitchExt(result, true); + fastRotationRem[j] = (j != 0) ? cc->EvalFastRotationExt(result, ReduceRotation(scaleRem * j, reduceMod), digits, true) : cc->KeySwitchExt(result, true); + + // Horner giant stride for remainder: g_rem * scale_rem + const int32_t tRem = ReduceRotation(scaleRem * p.gRem, reduceMod); Ciphertext outer; - DCRTPoly first; - for (uint32_t i = 0; i < p.bRem; ++i) { - // for the first iteration with j=0: - uint32_t GRem = i * p.gRem; + for (int32_t i = static_cast(p.bRem) - 1; i >= 0; --i) { + uint32_t GRem = static_cast(i) * p.gRem; auto inner = EvalMultExt(fastRotationRem[0], A[smax][GRem]); - // continue the loop for (uint32_t j = 1; j < p.gRem; ++j) { if ((GRem + j) != p.numRotationsRem) EvalAddExtInPlace(inner, EvalMultExt(fastRotationRem[j], A[smax][GRem + j])); } - if (i == 0) { - first = cc->KeySwitchDownFirstElement(inner); - auto elements = inner->GetElements(); - elements[0].SetValuesToZero(); - inner->SetElements(std::move(elements)); + if (i == static_cast(p.bRem) - 1) { outer = std::move(inner); } else { - if (rot_out[smax][i] != 0) { - inner = cc->KeySwitchDown(inner); - // Find the automorphism index that corresponds to rotation index index. - auto autoIndex = FindAutomorphismIndex2nComplex(rot_out[smax][i], cc->GetCyclotomicOrder()); - PrecomputeAutoMap(N, autoIndex, &map); - first += inner->GetElements()[0].AutomorphismTransform(autoIndex, map); - auto innerDigits = cc->EvalFastRotationPrecompute(inner); - EvalAddExtInPlace(outer, cc->EvalFastRotationExt(inner, rot_out[smax][i], innerDigits, false)); - } - else { - first += cc->KeySwitchDownFirstElement(inner); - auto elements = inner->GetElements(); - elements[0].SetValuesToZero(); - inner->SetElements(std::move(elements)); - EvalAddExtInPlace(outer, inner); + if (tRem != 0) { + auto outerStd = cc->KeySwitchDown(outer); + auto outerDigits = cc->EvalFastRotationPrecompute(outerStd); + outer = cc->EvalFastRotationExt(outerStd, tRem, outerDigits, true); } + EvalAddExtInPlace(outer, inner); } } result = cc->KeySwitchDown(outer); - result->GetElements()[0] += first; } + + // apply StC accumulated correction Aut_{-Delta} where Delta = slots-1 + const int32_t Delta = ReduceRotation(-static_cast(slots - 1), reduceMod); + if (Delta != 0) + result = cc->EvalAtIndex(result, Delta); + return result; } @@ -2322,9 +2202,12 @@ void FHECKKSRNS::ExtendCiphertext(std::vector& ctxtDCRTs, const Crypto const auto cryptoParams = std::dynamic_pointer_cast(cc.GetCryptoParameters()); uint32_t compositeDegree = cryptoParams->GetCompositeDegree(); + const auto& raisedParams = elementParamsRaisedPtr->GetParams(); + const uint32_t numRaisedTowers = raisedParams.size(); + std::vector qj(compositeDegree); for (uint32_t j = 0; j < compositeDegree; ++j) { - qj[j] = elementParamsRaisedPtr->GetParams()[j]->GetModulus().ConvertToInt(); + qj[j] = raisedParams[j]->GetModulus().ConvertToInt(); } std::vector qhat_modqj(compositeDegree); @@ -2352,10 +2235,12 @@ void FHECKKSRNS::ExtendCiphertext(std::vector& ctxtDCRTs, const Crypto std::vector tmp(compositeDegree + 1, DCRTPoly(elementParamsRaisedPtr, COEFFICIENT)); std::vector ctxtDCRTs_modq(compositeDegree, DCRTPoly(elementParamsRaisedPtr, COEFFICIENT)); -#pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(dcrt.GetNumOfElements())) - for (size_t j = 0; j < dcrt.GetNumOfElements(); ++j) { + const size_t numEl = dcrt.GetNumOfElements(); +#pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(numEl)) + for (size_t j = 0; j < numEl; ++j) { + const auto& ej = dcrt.GetElementAtIndex(j); for (uint32_t k = 0; k < compositeDegree; ++k) - ctxtDCRTs_modq[k].SetElementAtIndex(j, dcrt.GetElementAtIndex(j) * qhat_inv_modqj[k]); + ctxtDCRTs_modq[k].SetElementAtIndex(j, ej * qhat_inv_modqj[k]); } tmp[0] = ctxtDCRTs_modq[0].GetElementAtIndex(0); @@ -2375,7 +2260,7 @@ void FHECKKSRNS::ExtendCiphertext(std::vector& ctxtDCRTs, const Crypto } } - for (uint32_t j = compositeDegree; j < elementParamsRaisedPtr->GetParams().size(); ++j) + for (uint32_t j = compositeDegree; j < numRaisedTowers; ++j) tmp[d].SetElementAtIndex(j, tmp[init_element_index].GetElementAtIndex(j) * qjProductD); tmp[d].SetElementAtIndex(d, tmp[init_element_index].GetElementAtIndex(d) * qjProductD); @@ -2497,7 +2382,7 @@ Plaintext FHECKKSRNS::MakeAuxPlaintext(const CryptoContextImpl& cc, co im = im64 >> (-pRemaining); } else { - int128_t pPowRemaining = (static_cast(1)) << pRemaining; + int128_t pPowRemaining = static_cast(1) << pRemaining; im = pPowRemaining * im64; } @@ -2921,7 +2806,7 @@ void FHECKKSRNS::EvalFBTSetupInternal(const CryptoContextImpl& cc, con BigInteger q = cryptoParams->GetElementParams()->GetParams()[0]->GetModulus().ConvertToInt(); auto qDouble = q.ConvertToLongDouble(); - double factor = static_cast(1) << static_cast(std::round(std::log2(qDouble))); + double factor = std::ldexp(1.0, static_cast(std::round(std::log2(qDouble)))); double pre = qDouble / factor; double scaleEnc = pre / k; double scaleMod = QPrime.ConvertToLongDouble() / (Bigq.ConvertToLongDouble() * POut.ConvertToDouble()); @@ -3038,8 +2923,8 @@ std::shared_ptr> FHECKKSRNS::EvalMVBPrecomputeInternal( if (cryptoParams->GetKeySwitchTechnique() != HYBRID) OPENFHE_THROW("CKKS Bootstrapping is only supported for the Hybrid key switching method."); - auto paramsQ = cryptoParams->GetElementParams()->GetParams(); - uint32_t sizeQ = paramsQ.size(); + const auto& paramsQ = cryptoParams->GetElementParams()->GetParams(); + uint32_t sizeQ = paramsQ.size(); std::vector moduli(sizeQ); std::vector roots(sizeQ); for (uint32_t i = 0; i < sizeQ; ++i) { From cd91d8b4eb1127c7066a9674d6e57b189b05d697 Mon Sep 17 00:00:00 2001 From: Carlo Pascoe Date: Wed, 24 Jun 2026 16:24:48 -0400 Subject: [PATCH 2/4] Apply single-giant-step Horner BSGS to CKKS/FHEW scheme-switching transforms --- src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp | 93 +++++++++++------- .../ckksrns/ckksrns-schemeswitching.cpp | 96 ++++++++----------- src/pke/lib/scheme/ckksrns/ckksrns-utils.cpp | 20 ++-- 3 files changed, 108 insertions(+), 101 deletions(-) diff --git a/src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp b/src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp index 09d51012b..420848954 100644 --- a/src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp +++ b/src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp @@ -1588,9 +1588,11 @@ std::vector> FHECKKSRNS::EvalCoeffsToSlotsPrecomp c *= scale; } - auto rot = Rotate(coeff[s][ij], ReduceRotation(-rotScale * (ij / p.g) + offset * shiftScale + runningAcc, slots)); + auto rot = Rotate(coeff[s][ij], + ReduceRotation(-rotScale * (ij / p.g) + offset * shiftScale + runningAcc, slots)); - result[s][ij] = MakeAuxPlaintext(cc, paramsVector[s - stop], rot, 1, level0 - compositeDegree * s, rot.size()); + result[s][ij] = + MakeAuxPlaintext(cc, paramsVector[s - stop], rot, 1, level0 - compositeDegree * s, rot.size()); } } runningAcc += offset * shiftScale; @@ -1606,7 +1608,9 @@ std::vector> FHECKKSRNS::EvalCoeffsToSlotsPrecomp for (auto& c : coeff[stop][ij]) c *= scale; - auto rot = Rotate(coeff[stop][ij], ReduceRotation(-static_cast(p.gRem) * (ij / p.gRem) + offsetRem + runningAcc, slots)); + auto rot = Rotate( + coeff[stop][ij], + ReduceRotation(-static_cast(p.gRem) * (ij / p.gRem) + offsetRem + runningAcc, slots)); result[stop][ij] = MakeAuxPlaintext(cc, paramsVector[0], rot, 1, level0, rot.size()); } @@ -1639,9 +1643,11 @@ std::vector> FHECKKSRNS::EvalCoeffsToSlotsPrecomp c *= scale; } - auto rot = Rotate(clearTmp, ReduceRotation(-rotScale * (ij / p.g) + offset * shiftScale + runningAcc, 2 * slots)); + auto rot = Rotate( + clearTmp, ReduceRotation(-rotScale * (ij / p.g) + offset * shiftScale + runningAcc, 2 * slots)); - result[s][ij] = MakeAuxPlaintext(cc, paramsVector[s - stop], rot, 1, level0 - compositeDegree * s, rot.size()); + result[s][ij] = + MakeAuxPlaintext(cc, paramsVector[s - stop], rot, 1, level0 - compositeDegree * s, rot.size()); } } runningAcc += offset * shiftScale; @@ -1661,7 +1667,8 @@ std::vector> FHECKKSRNS::EvalCoeffsToSlotsPrecomp for (auto& c : clearTmp) c *= scale; - auto rot = Rotate(clearTmp, ReduceRotation(-p.gRem * (ij / p.gRem) + offsetRem + runningAcc, 2 * slots)); + auto rot = + Rotate(clearTmp, ReduceRotation(-p.gRem * (ij / p.gRem) + offsetRem + runningAcc, 2 * slots)); result[stop][ij] = MakeAuxPlaintext(cc, paramsVector[0], rot, 1, level0, rot.size()); } @@ -1747,9 +1754,11 @@ std::vector> FHECKKSRNS::EvalSlotsToCoeffsPrecomp c *= scale; } - auto rot = Rotate(coeff[s][ij], ReduceRotation(-rotScale * (ij / p.g) + offset * shiftScale + runningAcc, slots)); + auto rot = Rotate(coeff[s][ij], + ReduceRotation(-rotScale * (ij / p.g) + offset * shiftScale + runningAcc, slots)); - result[s][ij] = MakeAuxPlaintext(cc, paramsVector[s], rot, 1, towersToDrop + compositeDegree * s, rot.size()); + result[s][ij] = + MakeAuxPlaintext(cc, paramsVector[s], rot, 1, towersToDrop + compositeDegree * s, rot.size()); } } runningAcc += offset * shiftScale; @@ -1767,9 +1776,12 @@ std::vector> FHECKKSRNS::EvalSlotsToCoeffsPrecomp for (auto& c : coeff[smax][ij]) c *= scale; - auto rot = Rotate(coeff[smax][ij], ReduceRotation(-rotScale * (ij / p.gRem) + offsetRem * shiftScaleRem + runningAcc, slots)); + auto rot = Rotate( + coeff[smax][ij], + ReduceRotation(-rotScale * (ij / p.gRem) + offsetRem * shiftScaleRem + runningAcc, slots)); - result[smax][ij] = MakeAuxPlaintext(cc, paramsVector[smax], rot, 1, towersToDrop + compositeDegree * smax, rot.size()); + result[smax][ij] = MakeAuxPlaintext(cc, paramsVector[smax], rot, 1, + towersToDrop + compositeDegree * smax, rot.size()); } } } @@ -1802,9 +1814,11 @@ std::vector> FHECKKSRNS::EvalSlotsToCoeffsPrecomp c *= scale; } - auto rot = Rotate(clearTmp, ReduceRotation(-rotScale * (ij / p.g) + offset * shiftScale + runningAcc, 2 * slots)); + auto rot = Rotate( + clearTmp, ReduceRotation(-rotScale * (ij / p.g) + offset * shiftScale + runningAcc, 2 * slots)); - result[s][ij] = MakeAuxPlaintext(cc, paramsVector[s], rot, 1, towersToDrop + compositeDegree * s, rot.size()); + result[s][ij] = + MakeAuxPlaintext(cc, paramsVector[s], rot, 1, towersToDrop + compositeDegree * s, rot.size()); } } runningAcc += offset * shiftScale; @@ -1826,9 +1840,12 @@ std::vector> FHECKKSRNS::EvalSlotsToCoeffsPrecomp for (auto& c : clearTmp) c *= scale; - auto rot = Rotate(clearTmp, ReduceRotation(-rotScale * (ij / p.gRem) + offsetRem * shiftScaleRem + runningAcc, 2 * slots)); + auto rot = Rotate( + clearTmp, + ReduceRotation(-rotScale * (ij / p.gRem) + offsetRem * shiftScaleRem + runningAcc, 2 * slots)); - result[smax][ij] = MakeAuxPlaintext(cc, paramsVector[smax], rot, 1, towersToDrop + compositeDegree * smax, rot.size()); + result[smax][ij] = MakeAuxPlaintext(cc, paramsVector[smax], rot, 1, + towersToDrop + compositeDegree * smax, rot.size()); } } } @@ -1917,22 +1934,23 @@ Ciphertext FHECKKSRNS::EvalCoeffsToSlots(const std::vectorEvalFastRotationExt(result, ReduceRotation(scale * j, reduceMod), digits, true) : - cc->KeySwitchExt(result, true); + cc->EvalFastRotationExt(result, ReduceRotation(scale * j, reduceMod), digits, true) : + cc->KeySwitchExt(result, true); // Horner backward accumulation with single giant stride t = g * scale const int32_t t = ReduceRotation(scale * p.g, reduceMod); + const int32_t bLast = static_cast(p.b) - 1; Ciphertext outer; - for (int32_t i = static_cast(p.b) - 1; i >= 0; --i) { - uint32_t G = p.g * static_cast(i); + for (int32_t i = bLast; i >= 0; --i) { + uint32_t G = p.g * i; auto inner = EvalMultExt(fastRotation[0], A[s][G]); for (uint32_t j = 1; j < p.g; ++j) { if ((G + j) != p.numRotations) EvalAddExtInPlace(inner, EvalMultExt(fastRotation[j], A[s][G + j])); } - if (i == static_cast(p.b) - 1) { + if (i == bLast) { outer = std::move(inner); } else { @@ -1954,23 +1972,24 @@ Ciphertext FHECKKSRNS::EvalCoeffsToSlots(const std::vector> fastRotationRem(p.gRem); #pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(p.gRem)) for (uint32_t j = 0; j < p.gRem; ++j) - fastRotationRem[j] = - (j != 0) ? cc->EvalFastRotationExt(result, ReduceRotation(j, reduceMod), digits, true) : - cc->KeySwitchExt(result, true); + fastRotationRem[j] = (j != 0) ? + cc->EvalFastRotationExt(result, ReduceRotation(j, reduceMod), digits, true) : + cc->KeySwitchExt(result, true); // Remainder scale is 1 in CtS; Horner giant stride = gRem const int32_t tRem = ReduceRotation(p.gRem, reduceMod); + const int32_t bLast = static_cast(p.bRem) - 1; Ciphertext outer; - for (int32_t i = static_cast(p.bRem) - 1; i >= 0; --i) { - uint32_t GRem = p.gRem * static_cast(i); + for (int32_t i = bLast; i >= 0; --i) { + uint32_t GRem = p.gRem * i; auto inner = EvalMultExt(fastRotationRem[0], A[stop][GRem]); for (uint32_t j = 1; j < p.gRem; ++j) { if ((GRem + j) != p.numRotationsRem) EvalAddExtInPlace(inner, EvalMultExt(fastRotationRem[j], A[stop][GRem + j])); } - if (i == static_cast(p.bRem) - 1) { + if (i == bLast) { outer = std::move(inner); } else { @@ -2024,22 +2043,23 @@ Ciphertext FHECKKSRNS::EvalSlotsToCoeffs(const std::vectorEvalFastRotationExt(result, ReduceRotation(scale * j, reduceMod), digits, true) : - cc->KeySwitchExt(result, true); + cc->EvalFastRotationExt(result, ReduceRotation(scale * j, reduceMod), digits, true) : + cc->KeySwitchExt(result, true); // Horner backward accumulation with single giant stride t = g * scale const int32_t t = ReduceRotation(scale * p.g, reduceMod); + const int32_t bLast = static_cast(p.b) - 1; Ciphertext outer; - for (int32_t i = static_cast(p.b) - 1; i >= 0; --i) { - uint32_t G = static_cast(i) * p.g; + for (int32_t i = bLast; i >= 0; --i) { + uint32_t G = p.g * i; auto inner = EvalMultExt(fastRotation[0], A[s][G]); for (uint32_t j = 1; j < p.g; ++j) { if ((G + j) != p.numRotations) EvalAddExtInPlace(inner, EvalMultExt(fastRotation[j], A[s][G + j])); } - if (i == static_cast(p.b) - 1) { + if (i == bLast) { outer = std::move(inner); } else { @@ -2063,21 +2083,24 @@ Ciphertext FHECKKSRNS::EvalSlotsToCoeffs(const std::vector> fastRotationRem(p.gRem); #pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(p.gRem)) for (uint32_t j = 0; j < p.gRem; ++j) - fastRotationRem[j] = (j != 0) ? cc->EvalFastRotationExt(result, ReduceRotation(scaleRem * j, reduceMod), digits, true) : cc->KeySwitchExt(result, true); + fastRotationRem[j] = + (j != 0) ? cc->EvalFastRotationExt(result, ReduceRotation(scaleRem * j, reduceMod), digits, true) : + cc->KeySwitchExt(result, true); // Horner giant stride for remainder: g_rem * scale_rem const int32_t tRem = ReduceRotation(scaleRem * p.gRem, reduceMod); + const int32_t bLast = static_cast(p.bRem) - 1; Ciphertext outer; - for (int32_t i = static_cast(p.bRem) - 1; i >= 0; --i) { - uint32_t GRem = static_cast(i) * p.gRem; + for (int32_t i = bLast; i >= 0; --i) { + uint32_t GRem = p.gRem * i; auto inner = EvalMultExt(fastRotationRem[0], A[smax][GRem]); for (uint32_t j = 1; j < p.gRem; ++j) { if ((GRem + j) != p.numRotationsRem) EvalAddExtInPlace(inner, EvalMultExt(fastRotationRem[j], A[smax][GRem + j])); } - if (i == static_cast(p.bRem) - 1) { + if (i == bLast) { outer = std::move(inner); } else { @@ -2202,7 +2225,7 @@ void FHECKKSRNS::ExtendCiphertext(std::vector& ctxtDCRTs, const Crypto const auto cryptoParams = std::dynamic_pointer_cast(cc.GetCryptoParameters()); uint32_t compositeDegree = cryptoParams->GetCompositeDegree(); - const auto& raisedParams = elementParamsRaisedPtr->GetParams(); + const auto& raisedParams = elementParamsRaisedPtr->GetParams(); const uint32_t numRaisedTowers = raisedParams.size(); std::vector qj(compositeDegree); diff --git a/src/pke/lib/scheme/ckksrns/ckksrns-schemeswitching.cpp b/src/pke/lib/scheme/ckksrns/ckksrns-schemeswitching.cpp index fac5e5b18..a3c99e640 100644 --- a/src/pke/lib/scheme/ckksrns/ckksrns-schemeswitching.cpp +++ b/src/pke/lib/scheme/ckksrns/ckksrns-schemeswitching.cpp @@ -501,7 +501,7 @@ Ciphertext SWITCHCKKSRNS::EvalLTWithPrecomputeSwitch(const CryptoConte // Computing the baby-step bStep and the giant-step gStep uint32_t slots = A.size(); uint32_t bStep = dim1; - uint32_t gStep = std::ceil(static_cast(slots) / bStep); + int32_t gStep = std::ceil(static_cast(slots) / bStep); // Computes the NTTs for each CRT limb (for the hoisted automorphisms used later on) auto digits = cc.EvalFastRotationPrecompute(ctxt); @@ -512,41 +512,31 @@ Ciphertext SWITCHCKKSRNS::EvalLTWithPrecomputeSwitch(const CryptoConte for (uint32_t j = 1; j < bStep; ++j) fastRotation[j - 1] = cc.EvalFastRotationExt(ctxt, j, digits, true); - uint32_t M = cc.GetCyclotomicOrder(); - uint32_t N = cc.GetRingDimension(); - std::vector map(N); - Ciphertext result; - DCRTPoly first; - - for (uint32_t j = 0; j < gStep; ++j) { - auto inner = FHECKKSRNS::EvalMultExt(cc.KeySwitchExt(ctxt, true), A[bStep * j]); - for (uint32_t i = 1; i < bStep; ++i) { - if (bStep * j + i < slots) - FHECKKSRNS::EvalAddExtInPlace(inner, FHECKKSRNS::EvalMultExt(fastRotation[i - 1], A[bStep * j + i])); - } + auto ctxtExt = cc.KeySwitchExt(ctxt, true); - if (j == 0) { - first = cc.KeySwitchDownFirstElement(inner); - auto elements = inner->GetElements(); - elements[0].SetValuesToZero(); - inner->SetElements(std::move(elements)); + // Horner backward accumulation with a single giant stride t = bStep: + // result = block_0 + Aut_bStep(block_1 + Aut_bStep(...)) + // so only one giant-step key (bStep) is needed instead of the distinct keys + // {2*bStep, ..., (gStep-1)*bStep} the forward form required. + Ciphertext result; + for (int32_t j = gStep - 1; j >= 0; --j) { + const int32_t G = bStep * j; + auto inner = FHECKKSRNS::EvalMultExt(ctxtExt, A[G]); + for (uint32_t i = 1; i < bStep; ++i) + if (G + i < slots) + FHECKKSRNS::EvalAddExtInPlace(inner, FHECKKSRNS::EvalMultExt(fastRotation[i - 1], A[G + i])); + + if (j == gStep - 1) { result = std::move(inner); } else { - inner = cc.KeySwitchDown(inner); - // Find the automorphism index that corresponds to the rotation index. - uint32_t autoIndex = FindAutomorphismIndex2nComplex(bStep * j, M); - PrecomputeAutoMap(N, autoIndex, &map); - first += inner->GetElements()[0].AutomorphismTransform(autoIndex, map); - - auto&& innerDigits = cc.EvalFastRotationPrecompute(inner); - FHECKKSRNS::EvalAddExtInPlace(result, cc.EvalFastRotationExt(inner, bStep * j, innerDigits, false)); + auto resStd = cc.KeySwitchDown(result); + auto resDigits = cc.EvalFastRotationPrecompute(resStd); + result = cc.EvalFastRotationExt(resStd, bStep, resDigits, true); + FHECKKSRNS::EvalAddExtInPlace(result, inner); } } - - result = cc.KeySwitchDown(result); - result->GetElements()[0] += first; - return result; + return cc.KeySwitchDown(result); } Ciphertext SWITCHCKKSRNS::EvalLTRectWithPrecomputeSwitch( @@ -556,7 +546,7 @@ Ciphertext SWITCHCKKSRNS::EvalLTRectWithPrecomputeSwitch( // Computing the baby-step bStep and the giant-step gStep uint32_t bStep = (dim1 == 0) ? getRatioBSGSLT(n) : dim1; - uint32_t gStep = std::ceil(static_cast(n) / bStep); + int32_t gStep = std::ceil(static_cast(n) / bStep); uint32_t M = cc.GetCyclotomicOrder(); uint32_t N = cc.GetRingDimension(); @@ -608,43 +598,39 @@ Ciphertext SWITCHCKKSRNS::EvalLTRectWithPrecomputeSwitch( for (uint32_t j = 1; j < bStep; ++j) fastRotation[j - 1] = cc.EvalFastRotationExt(ct, j, digits, true); - std::vector map(N); Ciphertext result; - DCRTPoly first; - - for (uint32_t j = 0; j < gStep; j++) { - int32_t offset = (j == 0) ? 0 : -static_cast(bStep * j); - auto temp = cc.MakeCKKSPackedPlaintext(Rotate(Fill(A[bStep * j], N / 2), offset), 1, towersToDrop, - elementParamsPtr2, N / 2); - auto inner = FHECKKSRNS::EvalMultExt(cc.KeySwitchExt(ct, true), temp); + // ct raised to the extended (P*Q) basis -- the unrotated (i=0) baby term. + auto ctExt = cc.KeySwitchExt(ct, true); + + // Horner backward accumulation with a single giant stride t = bStep (one giant + // key instead of {2*bStep, ..., (gStep-1)*bStep}). Per-block plaintexts stay + // pre-rotated by -bStep*j; that is independent of the outer-rotation Horner form. + for (int32_t j = gStep - 1; j >= 0; --j) { + const int32_t G = bStep * j; + int32_t offset = (j == 0) ? 0 : -G; + auto temp = + cc.MakeCKKSPackedPlaintext(Rotate(Fill(A[G], N / 2), offset), 1, towersToDrop, elementParamsPtr2, N / 2); + auto inner = FHECKKSRNS::EvalMultExt(ctExt, temp); for (uint32_t i = 1; i < bStep; i++) { - if (bStep * j + i < n) { - auto tempi = cc.MakeCKKSPackedPlaintext(Rotate(Fill(A[bStep * j + i], N / 2), offset), 1, towersToDrop, + if (G + i < n) { + auto tempi = cc.MakeCKKSPackedPlaintext(Rotate(Fill(A[G + i], N / 2), offset), 1, towersToDrop, elementParamsPtr2, N / 2); FHECKKSRNS::EvalAddExtInPlace(inner, FHECKKSRNS::EvalMultExt(fastRotation[i - 1], tempi)); } } - if (j == 0) { - first = cc.KeySwitchDownFirstElement(inner); - auto elements = inner->GetElements(); - elements[0].SetValuesToZero(); - inner->SetElements(std::move(elements)); + if (j == gStep - 1) { result = std::move(inner); } else { - inner = cc.KeySwitchDown(inner); - // Find the automorphism index that corresponds to rotation index index. - uint32_t autoIndex = FindAutomorphismIndex2nComplex(bStep * j, M); - PrecomputeAutoMap(N, autoIndex, &map); - first += inner->GetElements()[0].AutomorphismTransform(autoIndex, map); - auto&& innerDigits = cc.EvalFastRotationPrecompute(inner); - FHECKKSRNS::EvalAddExtInPlace(result, cc.EvalFastRotationExt(inner, bStep * j, innerDigits, false)); + auto resStd = cc.KeySwitchDown(result); + auto resDigits = cc.EvalFastRotationPrecompute(resStd); + result = cc.EvalFastRotationExt(resStd, bStep, resDigits, true); + FHECKKSRNS::EvalAddExtInPlace(result, inner); } } result = cc.KeySwitchDown(result); - result->GetElements()[0] += first; // A represents the diagonals, which lose the information whether the initial matrix is tall or wide if (wide) { @@ -920,7 +906,7 @@ std::vector> SWITCHCKKSRNS::EvalCKKStoFHEW(Co #pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(numCtxts)) for (uint32_t i = 0; i < numCtxts; ++i) { auto& original_a = LWEciphertexts[i]->GetA(); - auto original_b = LWEciphertexts[i]->GetB(); + auto original_b = LWEciphertexts[i]->GetB(); // multiply by Q_LWE/Q' and round to Q_LWE NativeVector a_round(n, m_modulus_LWE); for (uint32_t j = 0; j < n; ++j) diff --git a/src/pke/lib/scheme/ckksrns/ckksrns-utils.cpp b/src/pke/lib/scheme/ckksrns/ckksrns-utils.cpp index c9df449a1..1ea85d6c9 100644 --- a/src/pke/lib/scheme/ckksrns/ckksrns-utils.cpp +++ b/src/pke/lib/scheme/ckksrns/ckksrns-utils.cpp @@ -741,17 +741,16 @@ std::vector FindLTRotationIndicesSwitch(uint32_t dim1, uint32_t m, uint else slots = blockDimension; - // Computing the baby-step g and the giant-step h + // Computing the baby-step bStep uint32_t bStep = (dim1 == 0) ? getRatioBSGSLT(slots) : dim1; - uint32_t gStep = static_cast(std::ceil(static_cast(slots) / bStep)); - // Computing all indices for baby-step giant-step procedure + // The Horner form of the LT transforms uses a single giant stride (bStep, + // already contained in {1..bStep}), so the giant-step keys + // {2*bStep, ..., (gStep-1)*bStep} the forward form needed are no longer generated. std::vector indexList; - indexList.reserve(bStep + gStep - 2); + indexList.reserve(bStep); for (uint32_t i = 1; i <= bStep; i++) indexList.emplace_back(i); - for (uint32_t i = 2; i < gStep; i++) - indexList.emplace_back(bStep * i); // Remove possible duplicates sort(indexList.begin(), indexList.end()); @@ -782,11 +781,11 @@ std::vector FindLTRotationIndicesSwitchArgmin(uint32_t m, uint32_t bloc indexList.reserve(bStep + gStep + cols); while (slots >= 1) { - // Computing all indices for baby-step giant-step procedure + // Computing all indices for baby-step giant-step procedure. The Horner form + // uses a single giant stride bStep (already in {1..bStep}), so the giant-step + // keys {2*bStep, ..., (gStep-1)*bStep} are no longer generated. for (uint32_t i = 1; i <= bStep; i++) indexList.emplace_back(i); - for (uint32_t i = 2; i < gStep; i++) - indexList.emplace_back(bStep * i); // If the linear transform is wide instead of tall, we need extra rotations if (slots < cols) { @@ -799,9 +798,8 @@ std::vector FindLTRotationIndicesSwitchArgmin(uint32_t m, uint32_t bloc // Go deeper into the binary tree slots /= 2; - // Computing the baby-step g and the giant-step h + // Recompute the baby-step for the next (halved) tree level bStep = getRatioBSGSLT(slots); - gStep = std::ceil(static_cast(slots) / bStep); } // Remove possible duplicates From 2eb89dc428e09bb5324d9f2dc80767a5b0a5d16c Mon Sep 17 00:00:00 2001 From: Carlo Pascoe Date: Fri, 26 Jun 2026 12:21:42 -0400 Subject: [PATCH 3/4] map-hoist + seed-peel of the Horner giant-step loop --- src/pke/include/scheme/ckksrns/ckksrns-fhe.h | 13 + src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp | 244 +++++++++++------- .../ckksrns/ckksrns-schemeswitching.cpp | 88 ++++--- 3 files changed, 210 insertions(+), 135 deletions(-) diff --git a/src/pke/include/scheme/ckksrns/ckksrns-fhe.h b/src/pke/include/scheme/ckksrns/ckksrns-fhe.h index ff9c1c8eb..acc0af7af 100644 --- a/src/pke/include/scheme/ckksrns/ckksrns-fhe.h +++ b/src/pke/include/scheme/ckksrns/ckksrns-fhe.h @@ -352,6 +352,19 @@ class FHECKKSRNS : public FHERNS { static Ciphertext EvalAddExt(ConstCiphertext ciphertext1, ConstCiphertext ciphertext2); + // Loads the loop-invariant automorphism key, index and O(N) permutation map for a constant + // Horner giant stride, so they are built once instead of re-derived inside the BSGS + // accumulation loop. Shared by the bootstrapping and scheme-switching linear transforms. + static EvalKey GetGiantStepRotation(ConstCiphertext ct, int32_t stride, uint32_t& autoIndex, + std::vector& map); + + // Inlined giant-step rotation for the Horner accumulation: equivalent to + // EvalFastRotationExt(KeySwitchDown(outer), stride, precompute(.), addFirst=true), reusing the + // caller-supplied loop-invariant (autoIndex, map, giantKey) from GetGiantStepRotation. + static Ciphertext EvalHornerGiantRotate(ConstCiphertext outer, uint32_t autoIndex, + const std::vector& map, + const EvalKey& giantKey); + static EvalKey ConjugateKeyGen(const PrivateKey privateKey); static Ciphertext Conjugate(ConstCiphertext ciphertext, diff --git a/src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp b/src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp index 420848954..6d7b5c8ae 100644 --- a/src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp +++ b/src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp @@ -76,6 +76,58 @@ double GetBigModulus(const std::shared_ptr& c namespace lbcrypto { +// Looks up the automorphism key, index and O(N) permutation map for a constant Horner giant +// stride. Hoisted out of the accumulation loop so these loop-invariant quantities are built once +// per stride instead of being re-derived by every EvalFastRotationExt call. +EvalKey FHECKKSRNS::GetGiantStepRotation(ConstCiphertext ct, int32_t stride, uint32_t& autoIndex, + std::vector& map) { + auto cc = ct->GetCryptoContext(); + const uint32_t M = cc->GetCyclotomicOrder(); + const uint32_t N = cc->GetRingDimension(); + autoIndex = FindAutomorphismIndex2nComplex(stride, M); + auto evalKeyMap = cc->GetEvalAutomorphismKeyMap(ct->GetKeyTag()); + auto evalKeyIt = evalKeyMap.find(autoIndex); + if (evalKeyIt == evalKeyMap.end()) + OPENFHE_THROW("EvalKey for index [" + std::to_string(autoIndex) + "] is not found."); + map.resize(N); + PrecomputeAutoMap(N, autoIndex, &map); + return evalKeyIt->second; +} + +// Inlined giant-step rotation for the Horner accumulation: equivalent to +// EvalFastRotationExt(KeySwitchDown(outer), stride, precompute(.), addFirst=true) but reusing the +// caller-supplied loop-invariant (autoIndex, map, giantKey) instead of having EvalFastRotationExt +// re-derive the automorphism index and O(N) map every call. +Ciphertext FHECKKSRNS::EvalHornerGiantRotate(ConstCiphertext outer, uint32_t autoIndex, + const std::vector& map, + const EvalKey& giantKey) { + auto cc = outer->GetCryptoContext(); + auto algo = cc->GetScheme(); + const auto cryptoParams = std::dynamic_pointer_cast(outer->GetCryptoParameters()); + + auto outerStd = cc->KeySwitchDown(outer); + const auto& cvStd = outerStd->GetElements(); + const auto paramsQl = cvStd[0].GetParams(); + + auto digits = algo->EvalKeySwitchPrecomputeCore(cvStd[1], cryptoParams); + auto cTilda = algo->EvalFastKeySwitchCoreExt(digits, giantKey, paramsQl); + + // addFirst: lift element 0 into the extended basis (c0 * P) and fold it in. + DCRTPoly psiC0(cTilda[0].GetParams(), Format::EVALUATION, true); + auto cMult = cvStd[0].TimesNoCheck(cryptoParams->GetPModq()); + const uint32_t sizeQl = paramsQl->GetParams().size(); + for (uint32_t i = 0; i < sizeQl; ++i) + psiC0.SetElementAtIndex(i, std::move(cMult.GetElementAtIndex(i))); + cTilda[0] += psiC0; + + cTilda[0] = cTilda[0].AutomorphismTransform(autoIndex, map); + cTilda[1] = cTilda[1].AutomorphismTransform(autoIndex, map); + + auto rotated = outerStd->CloneEmpty(); + rotated->SetElements(std::move(cTilda)); + return rotated; +} + //------------------------------------------------------------------------------ // Bootstrap Wrapper //------------------------------------------------------------------------------ @@ -165,7 +217,7 @@ void FHECKKSRNS::EvalBootstrapSetup(const CryptoContextImpl& cc, std:: // computes all powers of a primitive root of unity exp(2 * M_PI/m) std::vector> ksiPows(m + 1); double ak = 2 * M_PI / m; -#pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(4)) +#pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(2)) for (uint32_t j = 0; j < m; ++j) { double angle = ak * j; ksiPows[j].real(std::cos(angle)); @@ -355,7 +407,7 @@ void FHECKKSRNS::EvalBootstrapPrecompute(const CryptoContextImpl& cc, // computes all powers of a primitive root of unity exp(2 * M_PI/m) std::vector> ksiPows(m + 1); double ak = 2 * M_PI / m; -#pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(4)) +#pragma omp parallel for num_threads(OpenFHEParallelControls.GetThreadLimit(2)) for (uint32_t j = 0; j < m; ++j) { double angle = ak * j; ksiPows[j].real(std::cos(angle)); @@ -1875,24 +1927,26 @@ Ciphertext FHECKKSRNS::EvalLinearTransform(const std::vectorKeySwitchExt(ct, true); - // Horner backward accumulation with a single giant stride t = bStep: - Ciphertext outer; - for (int32_t j = gStep - 1; j >= 0; --j) { + uint32_t autoIndex = 0; + std::vector map; + EvalKey giantKey; + if (gStep > 1) + giantKey = GetGiantStepRotation(ct, bStep, autoIndex, map); + + // Horner backward accumulation with a single giant stride t = bStep. + const int32_t Gtop = (gStep - 1) * bStep; + Ciphertext outer = EvalMultExt(ctExt, A[Gtop]); + for (int32_t i = 1; i < bStep; ++i) + if (Gtop + i < slots) + EvalAddExtInPlace(outer, EvalMultExt(fastRotation[i - 1], A[Gtop + i])); + + for (int32_t j = gStep - 2; j >= 0; --j) { + outer = EvalHornerGiantRotate(outer, autoIndex, map, giantKey); const int32_t G = j * bStep; auto inner = EvalMultExt(ctExt, A[G]); for (int32_t i = 1; i < bStep; ++i) - if (G + i < slots) - EvalAddExtInPlace(inner, EvalMultExt(fastRotation[i - 1], A[G + i])); - - if (j == gStep - 1) { - outer = std::move(inner); - } - else { - auto outerStd = cc->KeySwitchDown(outer); - auto outerDigits = cc->EvalFastRotationPrecompute(outerStd); - outer = cc->EvalFastRotationExt(outerStd, bStep, outerDigits, true); - EvalAddExtInPlace(outer, inner); - } + EvalAddExtInPlace(inner, EvalMultExt(fastRotation[i - 1], A[G + i])); + EvalAddExtInPlace(outer, inner); } return cc->KeySwitchDown(outer); } @@ -1941,26 +1995,27 @@ Ciphertext FHECKKSRNS::EvalCoeffsToSlots(const std::vector(p.b) - 1; - Ciphertext outer; - for (int32_t i = bLast; i >= 0; --i) { - uint32_t G = p.g * i; - auto inner = EvalMultExt(fastRotation[0], A[s][G]); - for (uint32_t j = 1; j < p.g; ++j) { - if ((G + j) != p.numRotations) - EvalAddExtInPlace(inner, EvalMultExt(fastRotation[j], A[s][G + j])); - } - if (i == bLast) { - outer = std::move(inner); - } - else { - if (t != 0) { - auto outerStd = cc->KeySwitchDown(outer); - auto outerDigits = cc->EvalFastRotationPrecompute(outerStd); - outer = cc->EvalFastRotationExt(outerStd, t, outerDigits, true); - } - EvalAddExtInPlace(outer, inner); - } + uint32_t giantAutoIndex = 0; + std::vector giantMap; + EvalKey giantKey; + if (t != 0 && bLast > 0) + giantKey = GetGiantStepRotation(result, t, giantAutoIndex, giantMap); + + const uint32_t Gtop = p.g * bLast; + Ciphertext outer = EvalMultExt(fastRotation[0], A[s][Gtop]); + for (uint32_t j = 1; j < p.g; ++j) + if ((Gtop + j) != p.numRotations) + EvalAddExtInPlace(outer, EvalMultExt(fastRotation[j], A[s][Gtop + j])); + + for (int32_t i = bLast - 1; i >= 0; --i) { + if (t != 0) + outer = EvalHornerGiantRotate(outer, giantAutoIndex, giantMap, giantKey); + const uint32_t G = p.g * i; + auto inner = EvalMultExt(fastRotation[0], A[s][G]); + for (uint32_t j = 1; j < p.g; ++j) + EvalAddExtInPlace(inner, EvalMultExt(fastRotation[j], A[s][G + j])); + EvalAddExtInPlace(outer, inner); } result = cc->KeySwitchDown(outer); } @@ -1980,26 +2035,27 @@ Ciphertext FHECKKSRNS::EvalCoeffsToSlots(const std::vector(p.bRem) - 1; - Ciphertext outer; - for (int32_t i = bLast; i >= 0; --i) { - uint32_t GRem = p.gRem * i; - auto inner = EvalMultExt(fastRotationRem[0], A[stop][GRem]); - for (uint32_t j = 1; j < p.gRem; ++j) { - if ((GRem + j) != p.numRotationsRem) - EvalAddExtInPlace(inner, EvalMultExt(fastRotationRem[j], A[stop][GRem + j])); - } - if (i == bLast) { - outer = std::move(inner); - } - else { - if (tRem != 0) { - auto outerStd = cc->KeySwitchDown(outer); - auto outerDigits = cc->EvalFastRotationPrecompute(outerStd); - outer = cc->EvalFastRotationExt(outerStd, tRem, outerDigits, true); - } - EvalAddExtInPlace(outer, inner); - } + uint32_t giantAutoIndex = 0; + std::vector giantMap; + EvalKey giantKey; + if (tRem != 0 && bLast > 0) + giantKey = GetGiantStepRotation(result, tRem, giantAutoIndex, giantMap); + + const uint32_t GtopRem = p.gRem * bLast; + Ciphertext outer = EvalMultExt(fastRotationRem[0], A[stop][GtopRem]); + for (uint32_t j = 1; j < p.gRem; ++j) + if ((GtopRem + j) != p.numRotationsRem) + EvalAddExtInPlace(outer, EvalMultExt(fastRotationRem[j], A[stop][GtopRem + j])); + + for (int32_t i = bLast - 1; i >= 0; --i) { + if (tRem != 0) + outer = EvalHornerGiantRotate(outer, giantAutoIndex, giantMap, giantKey); + const uint32_t GRem = p.gRem * i; + auto inner = EvalMultExt(fastRotationRem[0], A[stop][GRem]); + for (uint32_t j = 1; j < p.gRem; ++j) + EvalAddExtInPlace(inner, EvalMultExt(fastRotationRem[j], A[stop][GRem + j])); + EvalAddExtInPlace(outer, inner); } result = cc->KeySwitchDown(outer); } @@ -2050,26 +2106,27 @@ Ciphertext FHECKKSRNS::EvalSlotsToCoeffs(const std::vector(p.b) - 1; - Ciphertext outer; - for (int32_t i = bLast; i >= 0; --i) { - uint32_t G = p.g * i; - auto inner = EvalMultExt(fastRotation[0], A[s][G]); - for (uint32_t j = 1; j < p.g; ++j) { - if ((G + j) != p.numRotations) - EvalAddExtInPlace(inner, EvalMultExt(fastRotation[j], A[s][G + j])); - } - if (i == bLast) { - outer = std::move(inner); - } - else { - if (t != 0) { - auto outerStd = cc->KeySwitchDown(outer); - auto outerDigits = cc->EvalFastRotationPrecompute(outerStd); - outer = cc->EvalFastRotationExt(outerStd, t, outerDigits, true); - } - EvalAddExtInPlace(outer, inner); - } + uint32_t giantAutoIndex = 0; + std::vector giantMap; + EvalKey giantKey; + if (t != 0 && bLast > 0) + giantKey = GetGiantStepRotation(result, t, giantAutoIndex, giantMap); + + const uint32_t Gtop = p.g * bLast; + Ciphertext outer = EvalMultExt(fastRotation[0], A[s][Gtop]); + for (uint32_t j = 1; j < p.g; ++j) + if ((Gtop + j) != p.numRotations) + EvalAddExtInPlace(outer, EvalMultExt(fastRotation[j], A[s][Gtop + j])); + + for (int32_t i = bLast - 1; i >= 0; --i) { + if (t != 0) + outer = EvalHornerGiantRotate(outer, giantAutoIndex, giantMap, giantKey); + const uint32_t G = p.g * i; + auto inner = EvalMultExt(fastRotation[0], A[s][G]); + for (uint32_t j = 1; j < p.g; ++j) + EvalAddExtInPlace(inner, EvalMultExt(fastRotation[j], A[s][G + j])); + EvalAddExtInPlace(outer, inner); } result = cc->KeySwitchDown(outer); } @@ -2091,26 +2148,27 @@ Ciphertext FHECKKSRNS::EvalSlotsToCoeffs(const std::vector(p.bRem) - 1; - Ciphertext outer; - for (int32_t i = bLast; i >= 0; --i) { - uint32_t GRem = p.gRem * i; - auto inner = EvalMultExt(fastRotationRem[0], A[smax][GRem]); - for (uint32_t j = 1; j < p.gRem; ++j) { - if ((GRem + j) != p.numRotationsRem) - EvalAddExtInPlace(inner, EvalMultExt(fastRotationRem[j], A[smax][GRem + j])); - } - if (i == bLast) { - outer = std::move(inner); - } - else { - if (tRem != 0) { - auto outerStd = cc->KeySwitchDown(outer); - auto outerDigits = cc->EvalFastRotationPrecompute(outerStd); - outer = cc->EvalFastRotationExt(outerStd, tRem, outerDigits, true); - } - EvalAddExtInPlace(outer, inner); - } + uint32_t giantAutoIndex = 0; + std::vector giantMap; + EvalKey giantKey; + if (tRem != 0 && bLast > 0) + giantKey = GetGiantStepRotation(result, tRem, giantAutoIndex, giantMap); + + const uint32_t GtopRem = p.gRem * bLast; + Ciphertext outer = EvalMultExt(fastRotationRem[0], A[smax][GtopRem]); + for (uint32_t j = 1; j < p.gRem; ++j) + if ((GtopRem + j) != p.numRotationsRem) + EvalAddExtInPlace(outer, EvalMultExt(fastRotationRem[j], A[smax][GtopRem + j])); + + for (int32_t i = bLast - 1; i >= 0; --i) { + if (tRem != 0) + outer = EvalHornerGiantRotate(outer, giantAutoIndex, giantMap, giantKey); + const uint32_t GRem = p.gRem * i; + auto inner = EvalMultExt(fastRotationRem[0], A[smax][GRem]); + for (uint32_t j = 1; j < p.gRem; ++j) + EvalAddExtInPlace(inner, EvalMultExt(fastRotationRem[j], A[smax][GRem + j])); + EvalAddExtInPlace(outer, inner); } result = cc->KeySwitchDown(outer); } diff --git a/src/pke/lib/scheme/ckksrns/ckksrns-schemeswitching.cpp b/src/pke/lib/scheme/ckksrns/ckksrns-schemeswitching.cpp index a3c99e640..799fe2512 100644 --- a/src/pke/lib/scheme/ckksrns/ckksrns-schemeswitching.cpp +++ b/src/pke/lib/scheme/ckksrns/ckksrns-schemeswitching.cpp @@ -514,27 +514,26 @@ Ciphertext SWITCHCKKSRNS::EvalLTWithPrecomputeSwitch(const CryptoConte auto ctxtExt = cc.KeySwitchExt(ctxt, true); - // Horner backward accumulation with a single giant stride t = bStep: - // result = block_0 + Aut_bStep(block_1 + Aut_bStep(...)) - // so only one giant-step key (bStep) is needed instead of the distinct keys - // {2*bStep, ..., (gStep-1)*bStep} the forward form required. - Ciphertext result; - for (int32_t j = gStep - 1; j >= 0; --j) { + // Horner backward accumulation with a single giant stride t = bStep + uint32_t autoIndex = 0; + std::vector map; + EvalKey giantKey; + if (gStep > 1) + giantKey = FHECKKSRNS::GetGiantStepRotation(ctxt, bStep, autoIndex, map); + + const int32_t Gtop = bStep * (gStep - 1); + Ciphertext result = FHECKKSRNS::EvalMultExt(ctxtExt, A[Gtop]); + for (uint32_t i = 1; i < bStep; ++i) + if (Gtop + i < slots) + FHECKKSRNS::EvalAddExtInPlace(result, FHECKKSRNS::EvalMultExt(fastRotation[i - 1], A[Gtop + i])); + + for (int32_t j = gStep - 2; j >= 0; --j) { + result = FHECKKSRNS::EvalHornerGiantRotate(result, autoIndex, map, giantKey); const int32_t G = bStep * j; auto inner = FHECKKSRNS::EvalMultExt(ctxtExt, A[G]); for (uint32_t i = 1; i < bStep; ++i) - if (G + i < slots) - FHECKKSRNS::EvalAddExtInPlace(inner, FHECKKSRNS::EvalMultExt(fastRotation[i - 1], A[G + i])); - - if (j == gStep - 1) { - result = std::move(inner); - } - else { - auto resStd = cc.KeySwitchDown(result); - auto resDigits = cc.EvalFastRotationPrecompute(resStd); - result = cc.EvalFastRotationExt(resStd, bStep, resDigits, true); - FHECKKSRNS::EvalAddExtInPlace(result, inner); - } + FHECKKSRNS::EvalAddExtInPlace(inner, FHECKKSRNS::EvalMultExt(fastRotation[i - 1], A[G + i])); + FHECKKSRNS::EvalAddExtInPlace(result, inner); } return cc.KeySwitchDown(result); } @@ -598,37 +597,42 @@ Ciphertext SWITCHCKKSRNS::EvalLTRectWithPrecomputeSwitch( for (uint32_t j = 1; j < bStep; ++j) fastRotation[j - 1] = cc.EvalFastRotationExt(ct, j, digits, true); - Ciphertext result; - // ct raised to the extended (P*Q) basis -- the unrotated (i=0) baby term. auto ctExt = cc.KeySwitchExt(ct, true); - // Horner backward accumulation with a single giant stride t = bStep (one giant - // key instead of {2*bStep, ..., (gStep-1)*bStep}). Per-block plaintexts stay - // pre-rotated by -bStep*j; that is independent of the outer-rotation Horner form. - for (int32_t j = gStep - 1; j >= 0; --j) { - const int32_t G = bStep * j; - int32_t offset = (j == 0) ? 0 : -G; + // Horner backward accumulation with a single giant stride t = bStep + uint32_t autoIndex = 0; + std::vector map; + EvalKey giantKey; + if (gStep > 1) + giantKey = FHECKKSRNS::GetGiantStepRotation(ct, bStep, autoIndex, map); + + const int32_t Gtop = bStep * (gStep - 1); + const int32_t offTop = (gStep == 1) ? 0 : -Gtop; + auto tempTop = + cc.MakeCKKSPackedPlaintext(Rotate(Fill(A[Gtop], N / 2), offTop), 1, towersToDrop, elementParamsPtr2, N / 2); + Ciphertext result = FHECKKSRNS::EvalMultExt(ctExt, tempTop); + for (uint32_t i = 1; i < bStep; ++i) { + if (Gtop + i < n) { + auto tempi = cc.MakeCKKSPackedPlaintext(Rotate(Fill(A[Gtop + i], N / 2), offTop), 1, towersToDrop, + elementParamsPtr2, N / 2); + FHECKKSRNS::EvalAddExtInPlace(result, FHECKKSRNS::EvalMultExt(fastRotation[i - 1], tempi)); + } + } + + for (int32_t j = gStep - 2; j >= 0; --j) { + result = FHECKKSRNS::EvalHornerGiantRotate(result, autoIndex, map, giantKey); + const int32_t G = bStep * j; + const int32_t offset = (j == 0) ? 0 : -G; auto temp = cc.MakeCKKSPackedPlaintext(Rotate(Fill(A[G], N / 2), offset), 1, towersToDrop, elementParamsPtr2, N / 2); auto inner = FHECKKSRNS::EvalMultExt(ctExt, temp); - for (uint32_t i = 1; i < bStep; i++) { - if (G + i < n) { - auto tempi = cc.MakeCKKSPackedPlaintext(Rotate(Fill(A[G + i], N / 2), offset), 1, towersToDrop, - elementParamsPtr2, N / 2); - FHECKKSRNS::EvalAddExtInPlace(inner, FHECKKSRNS::EvalMultExt(fastRotation[i - 1], tempi)); - } - } - - if (j == gStep - 1) { - result = std::move(inner); - } - else { - auto resStd = cc.KeySwitchDown(result); - auto resDigits = cc.EvalFastRotationPrecompute(resStd); - result = cc.EvalFastRotationExt(resStd, bStep, resDigits, true); - FHECKKSRNS::EvalAddExtInPlace(result, inner); + for (uint32_t i = 1; i < bStep; ++i) { + auto tempi = cc.MakeCKKSPackedPlaintext(Rotate(Fill(A[G + i], N / 2), offset), 1, towersToDrop, + elementParamsPtr2, N / 2); + FHECKKSRNS::EvalAddExtInPlace(inner, FHECKKSRNS::EvalMultExt(fastRotation[i - 1], tempi)); } + FHECKKSRNS::EvalAddExtInPlace(result, inner); } result = cc.KeySwitchDown(result); From 25e7d935e664e3b973ead4cdc66bf5dfef0ad717 Mon Sep 17 00:00:00 2001 From: Carlo Pascoe Date: Thu, 2 Jul 2026 11:46:58 -0400 Subject: [PATCH 4/4] updated comments and simplified scheme-switching LT rotation-index generation to closed form --- src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp | 5 ++ src/pke/lib/scheme/ckksrns/ckksrns-utils.cpp | 84 +++++--------------- 2 files changed, 26 insertions(+), 63 deletions(-) diff --git a/src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp b/src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp index 6d7b5c8ae..e2223e8a7 100644 --- a/src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp +++ b/src/pke/lib/scheme/ckksrns/ckksrns-fhe.cpp @@ -1333,6 +1333,7 @@ std::vector FHECKKSRNS::FindLinearTransformRotationIndices(uint32_t sl const auto& p = GetBootPrecom(slots); const uint32_t g = (p.m_paramsEnc.g == 0) ? std::ceil(std::sqrt(slots)) : p.m_paramsEnc.g; + // +32 is fixed upper bound for the log2((M/4)/slots) sparse-packing indices std::vector indexList; indexList.reserve(g + 32); @@ -1354,6 +1355,8 @@ std::vector FHECKKSRNS::FindLinearTransformRotationIndices(uint32_t sl std::vector FHECKKSRNS::FindCoeffsToSlotsRotationIndices(uint32_t slots, uint32_t M) { const auto& p = GetBootPrecom(slots).m_paramsEnc; + // lvlb*(g+1)+gRem sizes the per-level baby/giant rotations + // +32 is fixed upper bound for the log2((M/4)/slots) sparse-packing indices + correction key std::vector indexList; indexList.reserve(p.lvlb * (p.g + 1) + p.gRem + 32); @@ -1391,6 +1394,8 @@ std::vector FHECKKSRNS::FindCoeffsToSlotsRotationIndices(uint32_t slot std::vector FHECKKSRNS::FindSlotsToCoeffsRotationIndices(uint32_t slots, uint32_t M) { const auto& p = GetBootPrecom(slots).m_paramsDec; + // lvlb*(g+1)+gRem sizes the per-level baby/giant rotations + // +32 is fixed upper bound for the log2((M/4)/slots) sparse-packing indices std::vector indexList; indexList.reserve(p.lvlb * (p.g + 1) + p.gRem + 32); diff --git a/src/pke/lib/scheme/ckksrns/ckksrns-utils.cpp b/src/pke/lib/scheme/ckksrns/ckksrns-utils.cpp index 1ea85d6c9..6e12d80e5 100644 --- a/src/pke/lib/scheme/ckksrns/ckksrns-utils.cpp +++ b/src/pke/lib/scheme/ckksrns/ckksrns-utils.cpp @@ -734,82 +734,40 @@ uint32_t getRatioBSGSLT(uint32_t slots) { // returns powers of two } std::vector FindLTRotationIndicesSwitch(uint32_t dim1, uint32_t m, uint32_t blockDimension) { - uint32_t slots; // Set slots depending on packing mode (fully-packed or sparsely-packed) - if ((blockDimension == 0) || (blockDimension == m / 4)) - slots = m / 4; - else - slots = blockDimension; + const uint32_t slots = ((blockDimension == 0) || (blockDimension == m / 4)) ? m / 4 : blockDimension; // Computing the baby-step bStep - uint32_t bStep = (dim1 == 0) ? getRatioBSGSLT(slots) : dim1; + const uint32_t bStep = (dim1 == 0) ? getRatioBSGSLT(slots) : dim1; - // The Horner form of the LT transforms uses a single giant stride (bStep, - // already contained in {1..bStep}), so the giant-step keys - // {2*bStep, ..., (gStep-1)*bStep} the forward form needed are no longer generated. + // The Horner form of the LT transforms uses a single giant stride bStep (already contained in {1..bStep}) std::vector indexList; indexList.reserve(bStep); - for (uint32_t i = 1; i <= bStep; i++) + for (uint32_t i = 1; i <= bStep; ++i) indexList.emplace_back(i); - - // Remove possible duplicates - sort(indexList.begin(), indexList.end()); - indexList.erase(unique(indexList.begin(), indexList.end()), indexList.end()); - - // Remove automorphisms corresponding to 0 - auto it = std::find(indexList.begin(), indexList.end(), 0); - if (it != indexList.end()) - indexList.erase(it); return indexList; } std::vector FindLTRotationIndicesSwitchArgmin(uint32_t m, uint32_t blockDimension, uint32_t cols) { - uint32_t slots; // Set slots depending on packing mode (fully-packed or sparsely-packed) - if ((blockDimension == 0) || (blockDimension == m / 4)) - slots = m / 4; - else - slots = blockDimension; - - // Computing the baby-step g and the giant-step h - uint32_t bStep = getRatioBSGSLT(slots); - uint32_t gStep = std::ceil(static_cast(slots) / bStep); - uint32_t logl = std::log2(cols / slots); // These are powers of two, so log(l) is integer - - // There will be a lot of intersection between the rotations, provide an upper bound + const uint32_t slots = ((blockDimension == 0) || (blockDimension == m / 4)) ? m / 4 : blockDimension; + + // The argmin routine descends a binary tree, halving the slot count at each level. Each level needs + // the baby steps {1..getRatioBSGSLT(level slots)}; getRatioBSGSLT is monotonic non-decreasing in + // slots, so the union of the baby steps over every level collapses to {1, 2, ..., bStep} for bStep + // at the (largest) top-level slot count. The Horner form emits no giant-step keys. + const uint32_t bStep = getRatioBSGSLT(slots); + + // Wide (cols > slots) transforms additionally need power-of-two rotations up to cols/2; summed over + // the halving levels these collapse to exactly {1, 2, 4, ..., cols/2}. The powers <= bStep are + // already in the baby-step range, so only those in (bStep, cols/2] are appended -- and since bStep + // is itself a power of two, appending 2*bStep, 4*bStep, ... keeps the list sorted with no dedup pass. std::vector indexList; - indexList.reserve(bStep + gStep + cols); - - while (slots >= 1) { - // Computing all indices for baby-step giant-step procedure. The Horner form - // uses a single giant stride bStep (already in {1..bStep}), so the giant-step - // keys {2*bStep, ..., (gStep-1)*bStep} are no longer generated. - for (uint32_t i = 1; i <= bStep; i++) - indexList.emplace_back(i); - - // If the linear transform is wide instead of tall, we need extra rotations - if (slots < cols) { - logl = std::log2(cols / slots); // These are powers of two, so log(l) is integer - for (uint32_t j = 1; j <= logl; ++j) { - indexList.emplace_back(slots * (1U << (j - 1))); - } - } - - // Go deeper into the binary tree - slots /= 2; - - // Recompute the baby-step for the next (halved) tree level - bStep = getRatioBSGSLT(slots); - } - - // Remove possible duplicates - sort(indexList.begin(), indexList.end()); - indexList.erase(unique(indexList.begin(), indexList.end()), indexList.end()); - - // Remove automorphisms corresponding to 0 - auto it = std::find(indexList.begin(), indexList.end(), 0); - if (it != indexList.end()) - indexList.erase(it); + indexList.reserve(bStep + 32); // +32 bounds the (<= log2(cols)) wide power-of-two rotations + for (uint32_t i = 1; i <= bStep; ++i) + indexList.emplace_back(i); + for (uint32_t p = bStep << 1; p <= cols / 2; p <<= 1) + indexList.emplace_back(static_cast(p)); return indexList; }