Reproducing "Optimization of Homomorphic Comparison Algorithm on RNS-CKKS Scheme" using GenMinimaxCompositePolynomial

[Time0o-Lyceum]

The paper "Optimization of Homomorphic Comparison Algorithm on RNS-CKKS Scheme" (https://eprint.iacr.org/2021/1215) by Lee et al. uses the multi-interval Remez algorithm in its implementation of Algorithm 5: MinimaxComp (page 5). Lattigo provides the only open source implementation of this algorithm I have been able to find, based on a paper by the same authors.

The authors provide code for calculating optimal polynomial degrees for the comparison operator here. Running their default example (with alpha == 20, epsilon == 0.2002716064453125e-4 and is_comp == true) yields these degrees: 5 3 7 15 15 59. In their paper they then successively run the multi-interval Remez algorithm for the sign function and the intervals [-1, -epsilon], [epsilon, 1].

It seems like GenMinimaxCompositePolynomial is not able to reproduce this. Even just running it for degree 5 alone with logalpha == 16 (corresponding to the epsilon above) and logerr set to some value >= logalpha always results in panic: slope 0 occurred: consider increasing the precision, no matter how large a precision value is chosen (or at least for those that finish in a matter of minutes), e.g.:

package main                                                                        
                                                                                    
import (                                                                            
    "fmt"                                                                           
    "math/big"                                                                      
    "github.com/tuneinsight/lattigo/v6/circuits/ckks/minimax"                       
)                                                                                   
                                                                                    
func main() {                                                                       
    prec := uint(10000)                                                             
    logerr := 16                                                                    
    logalpha := 16                                                                  
                                                                                    
    deg := []int{5}                                                                 
                                                                                    
    f := func(x *big.Float) *big.Float {                                            
        if x.Cmp(big.NewFloat(0)) > 0 {                                             
            return big.NewFloat(1)                                                  
        }                                                                           
        return big.NewFloat(-1)                                                     
    }                                                                               
                                                                                    
    coeffs:= minimax.GenMinimaxCompositePolynomial(prec, logalpha, logerr, deg, f)
                                                                                    
    fmt.Println("Generated Coefficients:", coeffs)                                  
}  

Am I misusing the function?

[Pro7ech]

My original implementation of the multi-interval Remez algorithm, currently present in Lattigo, is unfortunately unstable and suffers from many issues, notably during the steps that tries to find zeroes of the error function.

I re-implemented the Remez algorithm with a different approach, which fixes all these issue. It can be found here. And you can also find the API to compute the minimax polynomials here.

[Pro7ech]

So far, I haven't been able to break it, and it performs well, even on functions with dozens of intervals. I used it to compute 128 different lookup tables, each of 25 entries for the 2024 iDash challenge

[Time0o-Lyceum]

@Pro7ech Ah yeah, that is the paper I was referring to, particularly the "Obtaining Extreme Points" section. I did just get around to trying this and unfortunately it still seems to have problems with the sign function, i.e. the slop panic for prec up to 10.000, with smaller logalpha as well:

prec := uint(10000)
logerr := 8                                
logalpha := 8
 
deg := []int{5} 

I ran it this time by checking out commit cf329e68cfab1d4c3dac5a423598a42aa63b5105 in your repo and then building against that.

[Pro7ech]

It works fine on my side:

package main

import(
	"github.com/Pro7ech/lattigo/he/hefloat"
)

func main(){
	hefloat.GenMinimaxCompositePolynomialForSign(64, 8, 8, []int{5});
}

output

P[0]
Interval: [-1.003906250000, -0.003906250000] U [0.003906250000, 1.003906250000]
Iteration:  0 - 0.718227150335 1.710085200282
Iteration:  1 - 0.967369969686 1.000880600535
Iteration:  2 - 0.967527672271 0.967565834124
MaxErr: 0.967565834124
MinErr: 0.967527672271

Output Precision: -0.047568

COEFFICIENTS:
{
{"0", "1.282361759434", "0", "-0.454696199003", "0", "1.139862111840", },
},

Also works fine with 10k bit-precision (although way overkill).

[Time0o-Lyceum]

Ah I see, I did not realize there was a special purpose function for sign. That works for me again, thanks again and I think that answers my question.

[Pro7ech]

I suspect your issue is that you weren't giving an arbitrary precision f(x).