From 2901247ecd5254879ea295c58eadd0b1f6da92fe Mon Sep 17 00:00:00 2001 From: armfazh Date: Fri, 29 Jul 2022 18:11:00 -0700 Subject: [PATCH] Adds polynomials and Lagrange polynomials. --- math/polynomial/polynomial.go | 123 +++++++++++++++++++++++++++++ math/polynomial/polynomial_test.go | 101 +++++++++++++++++++++++ 2 files changed, 224 insertions(+) create mode 100644 math/polynomial/polynomial.go create mode 100644 math/polynomial/polynomial_test.go diff --git a/math/polynomial/polynomial.go b/math/polynomial/polynomial.go new file mode 100644 index 000000000..db2a1d579 --- /dev/null +++ b/math/polynomial/polynomial.go @@ -0,0 +1,123 @@ +// Package polynomial provides representations of polynomials over the scalars +// of a group. +package polynomial + +import ( + "fmt" + + "github.com/cloudflare/circl/group" +) + +// Polynomial stores a polynomial over the set of scalars of a group. +type Polynomial struct { + // Internal representation is in polynomial basis: + // Thus, + // p(x) = \sum_i^k c[i] x^i, + // where k = len(c)-1 is the degree of the polynomial. + c []group.Scalar +} + +// New creates a new polynomial given its coefficients in ascending order. +// Thus, +// p(x) = \sum_i^k c[i] x^i, +// where k = len(c)-1 is the degree of the polynomial. +func New(coeffs []group.Scalar) Polynomial { + if len(coeffs) == 0 { + panic("polynomial: invalid length") + } + return Polynomial{coeffs} +} + +func (p Polynomial) Degree() uint { return uint(len(p.c) - 1) } + +func (p Polynomial) Evaluate(x group.Scalar) group.Scalar { + l := len(p.c) + px := p.c[l-1].Copy() + for i := l - 2; i >= 0; i-- { + px.Mul(px, x) + px.Add(px, p.c[i]) + } + return px +} + +// LagrangePolynomial stores a Lagrange polynomial over the set of scalars of a group. +type LagrangePolynomial struct { + // Internal representation is in Lagrange basis: + // Thus, + // p(x) = \sum_i^k y[i] L_j(x), where k is the degree of the polynomial, + // L_j(x) = \prod_i^k (x-x[i])/(x[j]-x[i]), + // y[i] = p(x[i]), and + // all x[i] are different. + g group.Group + x, y []group.Scalar +} + +// NewLagrangePolynomial creates a polynomial in Lagrange basis given a list +// of nodes (x) and values (y), such that: +// p(x) = \sum_i^k y[i] L_j(x), where k is the degree of the polynomial, +// L_j(x) = \prod_i^k (x-x[i])/(x[j]-x[i]), +// y[i] = p(x[i]), and +// all x[i] are different. +// It panics if one of these conditions does not hold. +func NewLagrangePolynomial(g group.Group, x, y []group.Scalar) LagrangePolynomial { + if len(x) != len(y) { + panic("lagrange: invalid length") + } + if !isAllDifferent(x) { + panic("lagrange: nodes must be different") + } + + return LagrangePolynomial{g, x, y} +} + +func (l LagrangePolynomial) Degree() uint { return uint(len(l.x) - 1) } + +func (l LagrangePolynomial) Evaluate(x group.Scalar) group.Scalar { + px := l.g.NewScalar() + tmp := l.g.NewScalar() + for i := range l.x { + LjX := baseRatio(uint(i), l.x, x) + tmp.Mul(l.y[i], LjX) + px.Add(px, tmp) + } + + return px +} + +// LagrangeBase returns the j-th Lagrange polynomial base evaluated at x. +// Thus, L_j(x) = \prod (x - x[i]) / (x[j] - x[i]) for 0 <= i < k, and i != j. +func LagrangeBase(jth uint, xi []group.Scalar, x group.Scalar) group.Scalar { + if jth >= uint(len(xi)) { + panic("lagrange: invalid index") + } + return baseRatio(jth, xi, x) +} + +func baseRatio(jth uint, xi []group.Scalar, x group.Scalar) group.Scalar { + num := x.Copy() + num.SetUint64(1) + den := x.Copy() + den.SetUint64(1) + + tmp := x.Copy() + for i := range xi { + if uint(i) != jth { + num.Mul(num, tmp.Sub(x, xi[i])) + den.Mul(den, tmp.Sub(xi[jth], xi[i])) + } + } + + return num.Mul(num, den.Inv(den)) +} + +func isAllDifferent(x []group.Scalar) bool { + m := make(map[string]struct{}) + for i := range x { + k := x[i].(fmt.Stringer).String() + if _, exists := m[k]; exists { + return false + } + m[k] = struct{}{} + } + return true +} diff --git a/math/polynomial/polynomial_test.go b/math/polynomial/polynomial_test.go new file mode 100644 index 000000000..8784e5299 --- /dev/null +++ b/math/polynomial/polynomial_test.go @@ -0,0 +1,101 @@ +package polynomial_test + +import ( + "testing" + + "github.com/cloudflare/circl/group" + "github.com/cloudflare/circl/internal/test" + "github.com/cloudflare/circl/math/polynomial" +) + +func TestPolyEval(t *testing.T) { + g := group.P256 + c := []group.Scalar{ + g.NewScalar(), + g.NewScalar(), + g.NewScalar(), + } + c[0].SetUint64(5) + c[1].SetUint64(5) + c[2].SetUint64(2) + p := polynomial.New(c) + + x := g.NewScalar() + x.SetUint64(10) + + got := p.Evaluate(x) + want := g.NewScalar() + want.SetUint64(255) + if !got.IsEqual(want) { + test.ReportError(t, got, want) + } +} + +func TestLagrange(t *testing.T) { + g := group.P256 + c := []group.Scalar{ + g.NewScalar(), + g.NewScalar(), + g.NewScalar(), + } + c[0].SetUint64(1234) + c[1].SetUint64(166) + c[2].SetUint64(94) + p := polynomial.New(c) + + x := []group.Scalar{g.NewScalar(), g.NewScalar(), g.NewScalar()} + x[0].SetUint64(2) + x[1].SetUint64(4) + x[2].SetUint64(5) + + y := []group.Scalar{} + for i := range x { + y = append(y, p.Evaluate(x[i])) + } + + zero := g.NewScalar() + l := polynomial.NewLagrangePolynomial(g, x, y) + test.CheckOk(l.Degree() == p.Degree(), "bad degree", t) + + got := l.Evaluate(zero) + want := p.Evaluate(zero) + + if !got.IsEqual(want) { + test.ReportError(t, got, want) + } + + // Test Kronecker's delta of LagrangeBase. + // Thus: + // L_j(x[i]) = { 1, if i == j; + // { 0, otherwise. + one := g.NewScalar() + one.SetUint64(1) + for j := range x { + for i := range x { + got := polynomial.LagrangeBase(uint(j), x, x[i]) + + if i == j { + want = one + } else { + want = zero + } + + if !got.IsEqual(want) { + test.ReportError(t, got, want) + } + } + } + + // Test that inputs are different length + err := test.CheckPanic(func() { polynomial.NewLagrangePolynomial(g, x, y[:1]) }) + test.CheckNoErr(t, err, "should panic") + + // Test that nodes must be different. + x[0].Set(x[1]) + err = test.CheckPanic(func() { polynomial.NewLagrangePolynomial(g, x, y) }) + test.CheckNoErr(t, err, "should panic") + + // Test LagrangeBase wrong index + err = test.CheckPanic(func() { polynomial.LagrangeBase(10, x, zero) }) + test.CheckNoErr(t, err, "should panic") +}