1// Copyright 2012 The Go Authors. All rights reserved. 2// Use of this source code is governed by a BSD-style 3// license that can be found in the LICENSE file. 4 5package bn256 6 7func lineFunctionAdd(r, p *twistPoint, q *curvePoint, r2 *gfP2, pool *bnPool) (a, b, c *gfP2, rOut *twistPoint) { 8 // See the mixed addition algorithm from "Faster Computation of the 9 // Tate Pairing", http://arxiv.org/pdf/0904.0854v3.pdf 10 11 B := newGFp2(pool).Mul(p.x, r.t, pool) 12 13 D := newGFp2(pool).Add(p.y, r.z) 14 D.Square(D, pool) 15 D.Sub(D, r2) 16 D.Sub(D, r.t) 17 D.Mul(D, r.t, pool) 18 19 H := newGFp2(pool).Sub(B, r.x) 20 I := newGFp2(pool).Square(H, pool) 21 22 E := newGFp2(pool).Add(I, I) 23 E.Add(E, E) 24 25 J := newGFp2(pool).Mul(H, E, pool) 26 27 L1 := newGFp2(pool).Sub(D, r.y) 28 L1.Sub(L1, r.y) 29 30 V := newGFp2(pool).Mul(r.x, E, pool) 31 32 rOut = newTwistPoint(pool) 33 rOut.x.Square(L1, pool) 34 rOut.x.Sub(rOut.x, J) 35 rOut.x.Sub(rOut.x, V) 36 rOut.x.Sub(rOut.x, V) 37 38 rOut.z.Add(r.z, H) 39 rOut.z.Square(rOut.z, pool) 40 rOut.z.Sub(rOut.z, r.t) 41 rOut.z.Sub(rOut.z, I) 42 43 t := newGFp2(pool).Sub(V, rOut.x) 44 t.Mul(t, L1, pool) 45 t2 := newGFp2(pool).Mul(r.y, J, pool) 46 t2.Add(t2, t2) 47 rOut.y.Sub(t, t2) 48 49 rOut.t.Square(rOut.z, pool) 50 51 t.Add(p.y, rOut.z) 52 t.Square(t, pool) 53 t.Sub(t, r2) 54 t.Sub(t, rOut.t) 55 56 t2.Mul(L1, p.x, pool) 57 t2.Add(t2, t2) 58 a = newGFp2(pool) 59 a.Sub(t2, t) 60 61 c = newGFp2(pool) 62 c.MulScalar(rOut.z, q.y) 63 c.Add(c, c) 64 65 b = newGFp2(pool) 66 b.SetZero() 67 b.Sub(b, L1) 68 b.MulScalar(b, q.x) 69 b.Add(b, b) 70 71 B.Put(pool) 72 D.Put(pool) 73 H.Put(pool) 74 I.Put(pool) 75 E.Put(pool) 76 J.Put(pool) 77 L1.Put(pool) 78 V.Put(pool) 79 t.Put(pool) 80 t2.Put(pool) 81 82 return 83} 84 85func lineFunctionDouble(r *twistPoint, q *curvePoint, pool *bnPool) (a, b, c *gfP2, rOut *twistPoint) { 86 // See the doubling algorithm for a=0 from "Faster Computation of the 87 // Tate Pairing", http://arxiv.org/pdf/0904.0854v3.pdf 88 89 A := newGFp2(pool).Square(r.x, pool) 90 B := newGFp2(pool).Square(r.y, pool) 91 C := newGFp2(pool).Square(B, pool) 92 93 D := newGFp2(pool).Add(r.x, B) 94 D.Square(D, pool) 95 D.Sub(D, A) 96 D.Sub(D, C) 97 D.Add(D, D) 98 99 E := newGFp2(pool).Add(A, A) 100 E.Add(E, A) 101 102 G := newGFp2(pool).Square(E, pool) 103 104 rOut = newTwistPoint(pool) 105 rOut.x.Sub(G, D) 106 rOut.x.Sub(rOut.x, D) 107 108 rOut.z.Add(r.y, r.z) 109 rOut.z.Square(rOut.z, pool) 110 rOut.z.Sub(rOut.z, B) 111 rOut.z.Sub(rOut.z, r.t) 112 113 rOut.y.Sub(D, rOut.x) 114 rOut.y.Mul(rOut.y, E, pool) 115 t := newGFp2(pool).Add(C, C) 116 t.Add(t, t) 117 t.Add(t, t) 118 rOut.y.Sub(rOut.y, t) 119 120 rOut.t.Square(rOut.z, pool) 121 122 t.Mul(E, r.t, pool) 123 t.Add(t, t) 124 b = newGFp2(pool) 125 b.SetZero() 126 b.Sub(b, t) 127 b.MulScalar(b, q.x) 128 129 a = newGFp2(pool) 130 a.Add(r.x, E) 131 a.Square(a, pool) 132 a.Sub(a, A) 133 a.Sub(a, G) 134 t.Add(B, B) 135 t.Add(t, t) 136 a.Sub(a, t) 137 138 c = newGFp2(pool) 139 c.Mul(rOut.z, r.t, pool) 140 c.Add(c, c) 141 c.MulScalar(c, q.y) 142 143 A.Put(pool) 144 B.Put(pool) 145 C.Put(pool) 146 D.Put(pool) 147 E.Put(pool) 148 G.Put(pool) 149 t.Put(pool) 150 151 return 152} 153 154func mulLine(ret *gfP12, a, b, c *gfP2, pool *bnPool) { 155 a2 := newGFp6(pool) 156 a2.x.SetZero() 157 a2.y.Set(a) 158 a2.z.Set(b) 159 a2.Mul(a2, ret.x, pool) 160 t3 := newGFp6(pool).MulScalar(ret.y, c, pool) 161 162 t := newGFp2(pool) 163 t.Add(b, c) 164 t2 := newGFp6(pool) 165 t2.x.SetZero() 166 t2.y.Set(a) 167 t2.z.Set(t) 168 ret.x.Add(ret.x, ret.y) 169 170 ret.y.Set(t3) 171 172 ret.x.Mul(ret.x, t2, pool) 173 ret.x.Sub(ret.x, a2) 174 ret.x.Sub(ret.x, ret.y) 175 a2.MulTau(a2, pool) 176 ret.y.Add(ret.y, a2) 177 178 a2.Put(pool) 179 t3.Put(pool) 180 t2.Put(pool) 181 t.Put(pool) 182} 183 184// sixuPlus2NAF is 6u+2 in non-adjacent form. 185var sixuPlus2NAF = []int8{0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 1, 0, 0, 0, -1, 0, 1, 0, 1, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, -1, 0, 1, 0, 0, 0, 1, 0, -1, 0, 0, 0, -1, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, -1, 0, -1, 0, 0, 0, 0, 1, 0, 0, 0, 1} 186 187// miller implements the Miller loop for calculating the Optimal Ate pairing. 188// See algorithm 1 from http://cryptojedi.org/papers/dclxvi-20100714.pdf 189func miller(q *twistPoint, p *curvePoint, pool *bnPool) *gfP12 { 190 ret := newGFp12(pool) 191 ret.SetOne() 192 193 aAffine := newTwistPoint(pool) 194 aAffine.Set(q) 195 aAffine.MakeAffine(pool) 196 197 bAffine := newCurvePoint(pool) 198 bAffine.Set(p) 199 bAffine.MakeAffine(pool) 200 201 minusA := newTwistPoint(pool) 202 minusA.Negative(aAffine, pool) 203 204 r := newTwistPoint(pool) 205 r.Set(aAffine) 206 207 r2 := newGFp2(pool) 208 r2.Square(aAffine.y, pool) 209 210 for i := len(sixuPlus2NAF) - 1; i > 0; i-- { 211 a, b, c, newR := lineFunctionDouble(r, bAffine, pool) 212 if i != len(sixuPlus2NAF)-1 { 213 ret.Square(ret, pool) 214 } 215 216 mulLine(ret, a, b, c, pool) 217 a.Put(pool) 218 b.Put(pool) 219 c.Put(pool) 220 r.Put(pool) 221 r = newR 222 223 switch sixuPlus2NAF[i-1] { 224 case 1: 225 a, b, c, newR = lineFunctionAdd(r, aAffine, bAffine, r2, pool) 226 case -1: 227 a, b, c, newR = lineFunctionAdd(r, minusA, bAffine, r2, pool) 228 default: 229 continue 230 } 231 232 mulLine(ret, a, b, c, pool) 233 a.Put(pool) 234 b.Put(pool) 235 c.Put(pool) 236 r.Put(pool) 237 r = newR 238 } 239 240 // In order to calculate Q1 we have to convert q from the sextic twist 241 // to the full GF(p^12) group, apply the Frobenius there, and convert 242 // back. 243 // 244 // The twist isomorphism is (x', y') -> (xω², yω³). If we consider just 245 // x for a moment, then after applying the Frobenius, we have x̄ω^(2p) 246 // where x̄ is the conjugate of x. If we are going to apply the inverse 247 // isomorphism we need a value with a single coefficient of ω² so we 248 // rewrite this as x̄ω^(2p-2)ω². ξ⁶ = ω and, due to the construction of 249 // p, 2p-2 is a multiple of six. Therefore we can rewrite as 250 // x̄ξ^((p-1)/3)ω² and applying the inverse isomorphism eliminates the 251 // ω². 252 // 253 // A similar argument can be made for the y value. 254 255 q1 := newTwistPoint(pool) 256 q1.x.Conjugate(aAffine.x) 257 q1.x.Mul(q1.x, xiToPMinus1Over3, pool) 258 q1.y.Conjugate(aAffine.y) 259 q1.y.Mul(q1.y, xiToPMinus1Over2, pool) 260 q1.z.SetOne() 261 q1.t.SetOne() 262 263 // For Q2 we are applying the p² Frobenius. The two conjugations cancel 264 // out and we are left only with the factors from the isomorphism. In 265 // the case of x, we end up with a pure number which is why 266 // xiToPSquaredMinus1Over3 is ∈ GF(p). With y we get a factor of -1. We 267 // ignore this to end up with -Q2. 268 269 minusQ2 := newTwistPoint(pool) 270 minusQ2.x.MulScalar(aAffine.x, xiToPSquaredMinus1Over3) 271 minusQ2.y.Set(aAffine.y) 272 minusQ2.z.SetOne() 273 minusQ2.t.SetOne() 274 275 r2.Square(q1.y, pool) 276 a, b, c, newR := lineFunctionAdd(r, q1, bAffine, r2, pool) 277 mulLine(ret, a, b, c, pool) 278 a.Put(pool) 279 b.Put(pool) 280 c.Put(pool) 281 r.Put(pool) 282 r = newR 283 284 r2.Square(minusQ2.y, pool) 285 a, b, c, newR = lineFunctionAdd(r, minusQ2, bAffine, r2, pool) 286 mulLine(ret, a, b, c, pool) 287 a.Put(pool) 288 b.Put(pool) 289 c.Put(pool) 290 r.Put(pool) 291 r = newR 292 293 aAffine.Put(pool) 294 bAffine.Put(pool) 295 minusA.Put(pool) 296 r.Put(pool) 297 r2.Put(pool) 298 299 return ret 300} 301 302// finalExponentiation computes the (p¹²-1)/Order-th power of an element of 303// GF(p¹²) to obtain an element of GT (steps 13-15 of algorithm 1 from 304// http://cryptojedi.org/papers/dclxvi-20100714.pdf) 305func finalExponentiation(in *gfP12, pool *bnPool) *gfP12 { 306 t1 := newGFp12(pool) 307 308 // This is the p^6-Frobenius 309 t1.x.Negative(in.x) 310 t1.y.Set(in.y) 311 312 inv := newGFp12(pool) 313 inv.Invert(in, pool) 314 t1.Mul(t1, inv, pool) 315 316 t2 := newGFp12(pool).FrobeniusP2(t1, pool) 317 t1.Mul(t1, t2, pool) 318 319 fp := newGFp12(pool).Frobenius(t1, pool) 320 fp2 := newGFp12(pool).FrobeniusP2(t1, pool) 321 fp3 := newGFp12(pool).Frobenius(fp2, pool) 322 323 fu, fu2, fu3 := newGFp12(pool), newGFp12(pool), newGFp12(pool) 324 fu.Exp(t1, u, pool) 325 fu2.Exp(fu, u, pool) 326 fu3.Exp(fu2, u, pool) 327 328 y3 := newGFp12(pool).Frobenius(fu, pool) 329 fu2p := newGFp12(pool).Frobenius(fu2, pool) 330 fu3p := newGFp12(pool).Frobenius(fu3, pool) 331 y2 := newGFp12(pool).FrobeniusP2(fu2, pool) 332 333 y0 := newGFp12(pool) 334 y0.Mul(fp, fp2, pool) 335 y0.Mul(y0, fp3, pool) 336 337 y1, y4, y5 := newGFp12(pool), newGFp12(pool), newGFp12(pool) 338 y1.Conjugate(t1) 339 y5.Conjugate(fu2) 340 y3.Conjugate(y3) 341 y4.Mul(fu, fu2p, pool) 342 y4.Conjugate(y4) 343 344 y6 := newGFp12(pool) 345 y6.Mul(fu3, fu3p, pool) 346 y6.Conjugate(y6) 347 348 t0 := newGFp12(pool) 349 t0.Square(y6, pool) 350 t0.Mul(t0, y4, pool) 351 t0.Mul(t0, y5, pool) 352 t1.Mul(y3, y5, pool) 353 t1.Mul(t1, t0, pool) 354 t0.Mul(t0, y2, pool) 355 t1.Square(t1, pool) 356 t1.Mul(t1, t0, pool) 357 t1.Square(t1, pool) 358 t0.Mul(t1, y1, pool) 359 t1.Mul(t1, y0, pool) 360 t0.Square(t0, pool) 361 t0.Mul(t0, t1, pool) 362 363 inv.Put(pool) 364 t1.Put(pool) 365 t2.Put(pool) 366 fp.Put(pool) 367 fp2.Put(pool) 368 fp3.Put(pool) 369 fu.Put(pool) 370 fu2.Put(pool) 371 fu3.Put(pool) 372 fu2p.Put(pool) 373 fu3p.Put(pool) 374 y0.Put(pool) 375 y1.Put(pool) 376 y2.Put(pool) 377 y3.Put(pool) 378 y4.Put(pool) 379 y5.Put(pool) 380 y6.Put(pool) 381 382 return t0 383} 384 385func optimalAte(a *twistPoint, b *curvePoint, pool *bnPool) *gfP12 { 386 e := miller(a, b, pool) 387 ret := finalExponentiation(e, pool) 388 e.Put(pool) 389 390 if a.IsInfinity() || b.IsInfinity() { 391 ret.SetOne() 392 } 393 394 return ret 395} 396