1 /*
2  * Example of a singular curve, similar to 19.c
3  * but the Tate pairing degenerates
4  *
5  * Consider the curve E: y^2 = x^3 + x^2 over F_19:
6  * E_ns(F_19) is a cyclic group of order 18.
7  */
8 
9 #include "pbc.h"
10 #include "pbc_singular.h"
11 #include "pbc_fp.h"
12 
miller(element_t res,element_t P,element_t Q,element_t R,int n)13 static void miller(element_t res, element_t P, element_t Q, element_t R, int n)
14 {
15     //collate divisions
16     int m;
17     element_t v, vd;
18     element_t Z;
19     element_t a, b, c;
20     element_t e0, e1;
21     mpz_t q;
22     element_ptr Zx, Zy;
23     const element_ptr Px = curve_x_coord(P);
24     const element_ptr Py = curve_y_coord(P);
25     const element_ptr numx = curve_x_coord(Q);
26     const element_ptr numy = curve_y_coord(Q);
27     const element_ptr denomx = curve_x_coord(R);
28     const element_ptr denomy = curve_y_coord(R);
29 
30     void do_vertical(element_t e, element_t edenom)
31     {
32         element_sub(e0, numx, Zx);
33         element_mul(e, e, e0);
34 
35         element_sub(e0, denomx, Zx);
36         element_mul(edenom, edenom, e0);
37     }
38 
39     void do_tangent(element_t e, element_t edenom)
40     {
41         //a = -slope_tangent(A.x, A.y);
42         //b = 1;
43         //c = -(A.y + a * A.x);
44         //but we multiply by 2*A.y to avoid division
45 
46         //a = -Ax * (Ax + Ax + Ax + twicea_2) - a_4;
47         //This curve is special:
48         //a = -(3 Ax^2 + 2Ax)
49         //b = 2 * Ay
50         //c = -(2 Ay^2 + a Ax);
51 
52         if (element_is0(Zy)) {
53             do_vertical(e, edenom);
54             return;
55         }
56         element_square(a, Zx);
57         element_mul_si(a, a, 3);
58         element_add(a, a, Zx);
59         element_add(a, a, Zx);
60         element_neg(a, a);
61 
62         element_add(b, Zy, Zy);
63 
64         element_mul(e0, b, Zy);
65         element_mul(c, a, Zx);
66         element_add(c, c, e0);
67         element_neg(c, c);
68 
69         element_mul(e0, a, numx);
70         element_mul(e1, b, numy);
71         element_add(e0, e0, e1);
72         element_add(e0, e0, c);
73         element_mul(e, e, e0);
74 
75         element_mul(e0, a, denomx);
76         element_mul(e1, b, denomy);
77         element_add(e0, e0, e1);
78         element_add(e0, e0, c);
79         element_mul(edenom, edenom, e0);
80     }
81 
82     void do_line(element_ptr e, element_ptr edenom)
83     {
84         if (!element_cmp(Zx, Px)) {
85             if (!element_cmp(Zy, Py)) {
86                 do_tangent(e, edenom);
87             } else {
88                 do_vertical(e, edenom);
89             }
90             return;
91         }
92 
93         element_sub(b, Px, Zx);
94         element_sub(a, Zy, Py);
95         element_mul(c, Zx, Py);
96         element_mul(e0, Zy, Px);
97         element_sub(c, c, e0);
98 
99         element_mul(e0, a, numx);
100         element_mul(e1, b, numy);
101         element_add(e0, e0, e1);
102         element_add(e0, e0, c);
103         element_mul(e, e, e0);
104 
105         element_mul(e0, a, denomx);
106         element_mul(e1, b, denomy);
107         element_add(e0, e0, e1);
108         element_add(e0, e0, c);
109         element_mul(edenom, edenom, e0);
110     }
111 
112     element_init(a, res->field);
113     element_init(b, res->field);
114     element_init(c, res->field);
115     element_init(e0, res->field);
116     element_init(e1, res->field);
117 
118     element_init(v, res->field);
119     element_init(vd, res->field);
120     element_init(Z, P->field);
121 
122     element_set(Z, P);
123     Zx = curve_x_coord(Z);
124     Zy = curve_y_coord(Z);
125 
126     element_set1(v);
127     element_set1(vd);
128 
129     mpz_init(q);
130     mpz_set_ui(q, n);
131     m = mpz_sizeinbase(q, 2) - 2;
132 
133     while(m >= 0) {
134         element_square(v, v);
135         element_square(vd, vd);
136         do_tangent(v, vd);
137         element_double(Z, Z);
138         do_vertical(vd, v);
139 
140         if (mpz_tstbit(q, m)) {
141             do_line(v, vd);
142             element_add(Z, Z, P);
143             if (m) {
144                 do_vertical(vd, v);
145             }
146         }
147         m--;
148     }
149 
150     mpz_clear(q);
151 
152     element_invert(vd, vd);
153     element_mul(res, v, vd);
154 
155     element_clear(v);
156     element_clear(vd);
157     element_clear(Z);
158     element_clear(a);
159     element_clear(b);
160     element_clear(c);
161     element_clear(e0);
162     element_clear(e1);
163 }
164 
tate_3(element_ptr out,element_ptr P,element_ptr Q,element_ptr R)165 static void tate_3(element_ptr out, element_ptr P, element_ptr Q, element_ptr R)
166 {
167     mpz_t six;
168 
169     mpz_init(six);
170     mpz_set_ui(six, 6);
171     element_t QR;
172     element_t e0;
173 
174     element_init(QR, P->field);
175     element_init(e0, out->field);
176 
177     element_add(QR, Q, R);
178 
179     //for subgroup size 3, -2P = P, hence
180     //the tangent line at P has divisor 3(P) - 3(O)
181 
182     miller(out, P, QR, R, 3);
183 
184     element_pow_mpz(out, out, six);
185     element_clear(QR);
186     element_clear(e0);
187     mpz_clear(six);
188 }
189 
tate_9(element_ptr out,element_ptr P,element_ptr Q,element_ptr R)190 static void tate_9(element_ptr out, element_ptr P, element_ptr Q, element_ptr R)
191 {
192     element_t QR;
193     element_init(QR, P->field);
194 
195     element_add(QR, Q, R);
196 
197     miller(out, P, QR, R, 9);
198 
199     element_square(out, out);
200 
201     element_clear(QR);
202 }
203 
main(void)204 int main(void)
205 {
206     field_t c;
207     field_t Z19;
208     element_t P, Q, R;
209     mpz_t q, z;
210     element_t a;
211     int i;
212 
213     mpz_init(q);
214     mpz_init(z);
215 
216     mpz_set_ui(q, 19);
217 
218     field_init_fp(Z19, q);
219     element_init(a, Z19);
220 
221     field_init_curve_singular_with_node(c, Z19);
222 
223     element_init(P, c);
224     element_init(Q, c);
225     element_init(R, c);
226 
227     //(3,+/-6) is a generator
228     //we have an isomorphism from E_ns to F_19^*
229     // (3,6) --> 3
230     //(generally (x,y) --> (y+x)/(y-x)
231 
232     curve_set_si(R, 3, 6);
233 
234     for (i=1; i<=18; i++) {
235         mpz_set_si(z, i);
236         element_mul_mpz(Q, R, z);
237         element_printf("%dR = %B\n", i, Q);
238     }
239 
240     mpz_set_ui(z, 6);
241     element_mul_mpz(P, R, z);
242     //P has order 3
243     element_printf("P = %B\n", P);
244 
245     for (i=1; i<=3; i++) {
246         mpz_set_si(z, i);
247         element_mul_mpz(Q, R, z);
248         tate_3(a, P, Q, R);
249         element_printf("e_3(P,%dP) = %B\n", i, a);
250     }
251 
252     element_double(P, R);
253     //P has order 9
254     element_printf("P = %B\n", P);
255     for (i=1; i<=9; i++) {
256         mpz_set_si(z, i);
257         element_mul_mpz(Q, P, z);
258         tate_9(a, P, Q, R);
259         element_printf("e_9(P,%dP) = %B\n", i, a);
260     }
261 
262     return 0;
263 }
264