1package r1csqap
  2
  3import (
  4	"math/big"
  5
  6	"github.com/arnaucube/go-snark/fields"
  7)
  8
  9func Transpose(matrix [][]*big.Int) [][]*big.Int {
 10	var r [][]*big.Int
 11	for i := 0; i < len(matrix[0]); i++ {
 12		var row []*big.Int
 13		for j := 0; j < len(matrix); j++ {
 14			row = append(row, matrix[j][i])
 15		}
 16		r = append(r, row)
 17	}
 18	return r
 19}
 20
 21func ArrayOfBigZeros(num int) []*big.Int {
 22	bigZero := big.NewInt(int64(0))
 23	var r []*big.Int
 24	for i := 0; i < num; i++ {
 25		r = append(r, bigZero)
 26	}
 27	return r
 28}
 29
 30type PolynomialField struct {
 31	F fields.Fq
 32}
 33
 34func NewPolynomialField(f fields.Fq) PolynomialField {
 35	return PolynomialField{
 36		f,
 37	}
 38}
 39func (pf PolynomialField) Mul(a, b []*big.Int) []*big.Int {
 40	r := ArrayOfBigZeros(len(a) + len(b) - 1)
 41	for i := 0; i < len(a); i++ {
 42		for j := 0; j < len(b); j++ {
 43			r[i+j] = pf.F.Add(
 44				r[i+j],
 45				pf.F.Mul(a[i], b[j]))
 46		}
 47	}
 48	return r
 49}
 50func (pf PolynomialField) Div(a, b []*big.Int) ([]*big.Int, []*big.Int) {
 51	// https://en.wikipedia.org/wiki/Division_algorithm
 52	r := ArrayOfBigZeros(len(a) - len(b) + 1)
 53	rem := a
 54	for len(rem) >= len(b) {
 55		l := pf.F.Div(rem[len(rem)-1], b[len(b)-1])
 56		pos := len(rem) - len(b)
 57		r[pos] = l
 58		aux := ArrayOfBigZeros(pos)
 59		aux1 := append(aux, l)
 60		aux2 := pf.Sub(rem, pf.Mul(b, aux1))
 61		rem = aux2[:len(aux2)-1]
 62	}
 63	return r, rem
 64}
 65
 66func max(a, b int) int {
 67	if a > b {
 68		return a
 69	}
 70	return b
 71}
 72
 73func (pf PolynomialField) Add(a, b []*big.Int) []*big.Int {
 74	r := ArrayOfBigZeros(max(len(a), len(b)))
 75	for i := 0; i < len(a); i++ {
 76		r[i] = pf.F.Add(r[i], a[i])
 77	}
 78	for i := 0; i < len(b); i++ {
 79		r[i] = pf.F.Add(r[i], b[i])
 80	}
 81	return r
 82}
 83
 84func (pf PolynomialField) Sub(a, b []*big.Int) []*big.Int {
 85	r := ArrayOfBigZeros(max(len(a), len(b)))
 86	for i := 0; i < len(a); i++ {
 87		r[i] = pf.F.Add(r[i], a[i])
 88	}
 89	for i := 0; i < len(b); i++ {
 90		r[i] = pf.F.Sub(r[i], b[i])
 91	}
 92	return r
 93}
 94
 95func (pf PolynomialField) Eval(v []*big.Int, x *big.Int) *big.Int {
 96	r := big.NewInt(int64(0))
 97	for i := 0; i < len(v); i++ {
 98		// xi := FloatPow(x, i)
 99		xi := pf.F.Exp(x, big.NewInt(int64(i)))
100		elem := pf.F.Mul(v[i], xi)
101		r = pf.F.Add(r, elem)
102	}
103	return r
104}
105
106func (pf PolynomialField) NewPolZeroAt(pointPos, totalPoints int, height *big.Int) []*big.Int {
107	fac := 1
108	for i := 1; i < totalPoints+1; i++ {
109		if i != pointPos {
110			fac = fac * (pointPos - i)
111		}
112	}
113	facBig := big.NewInt(int64(fac))
114	hf := pf.F.Div(height, facBig)
115	r := []*big.Int{hf}
116	for i := 1; i < totalPoints+1; i++ {
117		if i != pointPos {
118			ineg := big.NewInt(int64(-i))
119			b1 := big.NewInt(int64(1))
120			r = pf.Mul(r, []*big.Int{ineg, b1})
121		}
122	}
123	return r
124}
125
126func (pf PolynomialField) LagrangeInterpolation(v []*big.Int) []*big.Int {
127	// https://en.wikipedia.org/wiki/Lagrange_polynomial
128	var r []*big.Int
129	for i := 0; i < len(v); i++ {
130		r = pf.Add(r, pf.NewPolZeroAt(i+1, len(v), v[i]))
131	}
132	//
133	return r
134}
135
136func (pf PolynomialField) R1CSToQAP(a, b, c [][]*big.Int) ([][]*big.Int, [][]*big.Int, [][]*big.Int, []*big.Int) {
137	aT := Transpose(a)
138	bT := Transpose(b)
139	cT := Transpose(c)
140	var alphas [][]*big.Int
141	for i := 0; i < len(aT); i++ {
142		alphas = append(alphas, pf.LagrangeInterpolation(aT[i]))
143	}
144	var betas [][]*big.Int
145	for i := 0; i < len(bT); i++ {
146		betas = append(betas, pf.LagrangeInterpolation(bT[i]))
147	}
148	var gammas [][]*big.Int
149	for i := 0; i < len(cT); i++ {
150		gammas = append(gammas, pf.LagrangeInterpolation(cT[i]))
151	}
152	z := []*big.Int{big.NewInt(int64(1))}
153	for i := 1; i < len(aT[0])+1; i++ {
154		ineg := big.NewInt(int64(-i))
155		b1 := big.NewInt(int64(1))
156		z = pf.Mul(z, []*big.Int{ineg, b1})
157	}
158	return alphas, betas, gammas, z
159}
160
161func (pf PolynomialField) CombinePolynomials(r []*big.Int, ap, bp, cp [][]*big.Int) ([]*big.Int, []*big.Int, []*big.Int, []*big.Int) {
162	var alpha []*big.Int
163	for i := 0; i < len(r); i++ {
164		m := pf.Mul([]*big.Int{r[i]}, ap[i])
165		alpha = pf.Add(alpha, m)
166	}
167	var beta []*big.Int
168	for i := 0; i < len(r); i++ {
169		m := pf.Mul([]*big.Int{r[i]}, bp[i])
170		beta = pf.Add(beta, m)
171	}
172	var gamma []*big.Int
173	for i := 0; i < len(r); i++ {
174		m := pf.Mul([]*big.Int{r[i]}, cp[i])
175		gamma = pf.Add(gamma, m)
176	}
177
178	px := pf.Sub(pf.Mul(alpha, beta), gamma)
179	return alpha, beta, gamma, px
180}
181
182func (pf PolynomialField) DivisorPolinomial(px, z []*big.Int) []*big.Int {
183	quo, _ := pf.Div(px, z)
184	return quo
185}