g_param.c

这是一个C的源代码

    // cm1 = 0
    // C = (2y)^-1
    element_set0(cm1);
    element_invert(C, c1);

    element_set1(dm1);
    element_set1(d0);

    element_t sm2, sm1;
    element_t s0, s1, s2, s3;
    element_t tm2, tm1;
    element_t t0, t1, t2, t3;
    element_t e0, e1;
    element_t u, v;

    element_init_same_as(sm2, x);
    element_init_same_as(sm1, x);
    element_init_same_as(s0, x);
    element_init_same_as(s1, x);
    element_init_same_as(s2, x);
    element_init_same_as(s3, x);

    element_init_same_as(tm2, x);
    element_init_same_as(tm1, x);
    element_init_same_as(t0, x);
    element_init_same_as(t1, x);
    element_init_same_as(t2, x);
    element_init_same_as(t3, x);

    element_init_same_as(e0, x);
    element_init_same_as(e1, x);

    element_init_same_as(u, d0);
    element_init_same_as(v, d0);

    int m = mpz_sizeinbase(pairing->r, 2) - 2;
    for (;;) {
	element_square(sm2, cm2);
	element_square(sm1, cm1);
	element_square(s0, c0);
	element_square(s1, c1);
	element_square(s2, c2);
	element_square(s3, c3);

	element_mul(tm2, cm3, cm1);
	element_mul(tm1, cm2, c0);
	element_mul(t0, cm1, c1);
	element_mul(t1, c0, c2);
	element_mul(t2, c1, c3);
	element_mul(t3, c2, c4);

	element_square(u, d0);
	element_mul(v, dm1, d1);

	if (mpz_tstbit(pairing->r, m)) {
	    //double-and-add
	    element_mul(e0, t0, sm2);
	    element_mul(e1, tm2, s0);
	    element_sub(cm3, e0, e1);
	    element_mul(cm3, cm3, C);

	    element_mul(e0, t0, sm1);
	    element_mul(e1, tm1, s0);
	    element_sub(cm2, e0, e1);

	    element_mul(e0, t1, sm1);
	    element_mul(e1, tm1, s1);
	    element_sub(cm1, e0, e1);
	    element_mul(cm1, cm1, C);

	    element_mul(e0, t1, s0);
	    element_mul(e1, t0, s1);
	    element_sub(c0, e0, e1);

	    element_mul(e0, t2, s0);
	    element_mul(e1, t0, s2);
	    element_sub(c1, e0, e1);
	    element_mul(c1, c1, C);

	    element_mul(e0, t2, s1);
	    element_mul(e1, t1, s2);
	    element_sub(c2, e0, e1);

	    element_mul(e0, t3, s1);
	    element_mul(e1, t1, s3);
	    element_sub(c3, e0, e1);
	    element_mul(c3, c3, C);

	    element_mul(e0, t3, s2);
	    element_mul(e1, t2, s3);
	    element_sub(c4, e0, e1);

	    polymod_const_mul(fi_re(out), t0,  fi_re(u));
	    polymod_const_mul(fi_im(out), t0,  fi_im(u));
	    polymod_const_mul(fi_re(dm1), s0,  fi_re(v));
	    polymod_const_mul(fi_im(dm1), s0,  fi_im(v));
	    element_sub(dm1, dm1, out);

	    polymod_const_mul(fi_re(out), t1, fi_re(u));
	    polymod_const_mul(fi_im(out), t1, fi_im(u));
	    polymod_const_mul(fi_re(d0), s1, fi_re(v));
	    polymod_const_mul(fi_im(d0), s1, fi_im(v));
	    element_sub(d0, d0, out);
	    element_mul(d0, d0, A);

	    polymod_const_mul(fi_re(out), t2, fi_re(u));
	    polymod_const_mul(fi_im(out), t2, fi_im(u));
	    polymod_const_mul(fi_re(d1), s2, fi_re(v));
	    polymod_const_mul(fi_im(d1), s2, fi_im(v));
	    element_sub(d1, d1, out);
	    element_mul(d1, d1, B);
	} else {
	    //double
	    element_mul(e0, tm1, sm2);
	    element_mul(e1, tm2, sm1);
	    element_sub(cm3, e0, e1);

	    element_mul(e0, t0, sm2);
	    element_mul(e1, tm2, s0);
	    element_sub(cm2, e0, e1);
	    element_mul(cm2, cm2, C);

	    element_mul(e0, t0, sm1);
	    element_mul(e1, tm1, s0);
	    element_sub(cm1, e0, e1);

	    element_mul(e0, t1, sm1);
	    element_mul(e1, tm1, s1);
	    element_sub(c0, e0, e1);
	    element_mul(c0, c0, C);

	    element_mul(e0, t1, s0);
	    element_mul(e1, t0, s1);
	    element_sub(c1, e0, e1);

	    element_mul(e0, t2, s0);
	    element_mul(e1, t0, s2);
	    element_sub(c2, e0, e1);
	    element_mul(c2, c2, C);

	    element_mul(e0, t2, s1);
	    element_mul(e1, t1, s2);
	    element_sub(c3, e0, e1);

	    element_mul(e0, t3, s1);
	    element_mul(e1, t1, s3);
	    element_sub(c4, e0, e1);
	    element_mul(c4, c4, C);

	    polymod_const_mul(fi_re(out), tm1, fi_re(u));
	    polymod_const_mul(fi_im(out), tm1, fi_im(u));
	    polymod_const_mul(fi_re(dm1), sm1, fi_re(v));
	    polymod_const_mul(fi_im(dm1), sm1, fi_im(v));
	    element_sub(dm1, dm1, out);

	    polymod_const_mul(fi_re(out), t0, fi_re(u));
	    polymod_const_mul(fi_im(out), t0, fi_im(u));
	    polymod_const_mul(fi_re(d0), s0, fi_re(v));
	    polymod_const_mul(fi_im(d0), s0, fi_im(v));
	    element_sub(d0, d0, out);

	    polymod_const_mul(fi_re(out), t1, fi_re(u));
	    polymod_const_mul(fi_im(out), t1, fi_im(u));
	    polymod_const_mul(fi_re(d1), s1, fi_re(v));
	    polymod_const_mul(fi_im(d1), s1, fi_im(v));
	    element_sub(d1, d1, out);
	    element_mul(d1, d1, A);
	}
	if (!m) break;
	m--;
    }
    // since c_k lies base field
    // it gets killed by the final powering
    //element_invert(c1, c1);
    //element_mul(fi_re(d1), fi_re(d1), c1);
    //element_mul(fi_im(d1), fi_im(d1), c1);

    tatepower10(out, d1, pairing);

    element_clear(dm1);
    element_clear(d0);
    element_clear(d1);

    element_clear(cm3);
    element_clear(cm2);
    element_clear(cm1);
    element_clear(c0);
    element_clear(c1);
    element_clear(c2);
    element_clear(c3);
    element_clear(c4);

    element_clear(sm2);
    element_clear(sm1);
    element_clear(s0);
    element_clear(s1);
    element_clear(s2);
    element_clear(s3);

    element_clear(tm2);
    element_clear(tm1);
    element_clear(t0);
    element_clear(t1);
    element_clear(t2);
    element_clear(t3);

    element_clear(e0);
    element_clear(e1);
    element_clear(A);
    element_clear(B);
    element_clear(C);
    element_clear(u);
    element_clear(v);
}

void g_pairing_clear(pairing_t pairing)
{
    mnt_pairing_data_ptr p = pairing->data;

    element_clear(p->xpowq);
    element_clear(p->xpowq2);
    element_clear(p->xpowq3);
    element_clear(p->xpowq4);
    mpz_clear(p->tateexp);
    mpz_clear(pairing->phikonr);

    field_clear(p->Etwist);
    field_clear(p->Eq);
    element_clear(p->nqrinv);
    element_clear(p->nqrinv2);
    field_clear(p->Fqk);
    field_clear(p->Fqd);
    field_clear(p->Fqx);
    field_clear(p->Fq);
    field_clear(pairing->Zr);
    mpz_clear(pairing->r);
    pbc_free(p);
}

static void g_pairing_option_set(pairing_t pairing, char *key, char *value)
{
    UNUSED_VAR(pairing);
    if (!strcmp(key, "method")) {
	if (!strcmp(value, "miller")) {
	    cc_miller_no_denom_fn = cc_miller_no_denom_proj;
	} else if (!strcmp(value, "miller-affine")) {
	    cc_miller_no_denom_fn = cc_miller_no_denom_affine;
	} else if (!strcmp(value, "shipsey-stange")) {
	    pairing->map = g_pairing_ellnet;
	}
    }
}

void pairing_init_g_param(pairing_t pairing, g_param_t param)
{
    mnt_pairing_data_ptr p;
    element_t a, b;
    element_t irred;
    int i;

    mpz_init(pairing->r);
    mpz_set(pairing->r, param->r);
    field_init_fp(pairing->Zr, pairing->r);
    pairing->map = cc_pairing;
    pairing->is_almost_coddh = cc_is_almost_coddh;

    p =	pairing->data = pbc_malloc(sizeof(mnt_pairing_data_t));
    field_init_fp(p->Fq, param->q);
    element_init(a, p->Fq);
    element_init(b, p->Fq);
    element_set_mpz(a, param->a);
    element_set_mpz(b, param->b);
    field_init_curve_ab(p->Eq, a, b, pairing->r, param->h);

    field_init_poly(p->Fqx, p->Fq);
    element_init(irred, p->Fqx);
    poly_alloc(irred, 5 + 1);
    for (i=0; i<5; i++) {
	element_set_mpz(poly_coeff(irred, i), param->coeff[i]);
    }
    element_set1(poly_coeff(irred, 5));

    field_init_polymod(p->Fqd, irred);
    element_clear(irred);

    p->Fqd->nqr = pbc_malloc(sizeof(element_t));
    element_init(p->Fqd->nqr, p->Fqd);
    element_set_mpz(((element_t *) p->Fqd->nqr->data)[0], param->nqr);

    field_init_quadratic(p->Fqk, p->Fqd);

    { //compute phi(k)/r = (q^4 - q^3 + ... + 1)/r 
	element_ptr e = p->xpowq;
	mpz_t z0;
	mpz_ptr q = param->q;
	mpz_ptr z = pairing->phikonr;
	mpz_init(z);
	mpz_init(z0);
	mpz_set_ui(z, 1);
	mpz_sub(z, z, q);
	mpz_mul(z0, q, q);
	mpz_add(z, z, z0);
	mpz_mul(z0, z0, q);
	mpz_sub(z, z, z0);
	mpz_mul(z0, z0, q);
	mpz_add(z, z, z0);
	mpz_clear(z0);
	mpz_divexact(z, z, pairing->r);

	element_init(e, p->Fqd);
	element_init(p->xpowq2, p->Fqd);
	element_init(p->xpowq3, p->Fqd);
	element_init(p->xpowq4, p->Fqd);
	element_set1(((element_t *) e->data)[1]);
	element_pow_mpz(e, e, q);

	element_square(p->xpowq2, p->xpowq);
	element_square(p->xpowq4, p->xpowq2);
	element_mul(p->xpowq3, p->xpowq2, p->xpowq);
    }
    mpz_init(p->tateexp);
    mpz_sub_ui(p->tateexp, p->Fqk->order, 1);
    mpz_divexact(p->tateexp, p->tateexp, pairing->r);

    field_init_curve_ab_map(p->Etwist, p->Eq, element_field_to_polymod, p->Fqd, pairing->r, NULL);

    twist_curve(p->Etwist);

    element_init(p->nqrinv, p->Fqd);
    element_invert(p->nqrinv, field_get_nqr(p->Fqd));
    element_init(p->nqrinv2, p->Fqd);
    element_square(p->nqrinv2, p->nqrinv);

    pairing->G1 = p->Eq;
    pairing->G2 = p->Etwist;

    pairing->GT = p->Fqk;

    cc_miller_no_denom_fn = cc_miller_no_denom_affine;
    pairing->option_set = g_pairing_option_set;
    pairing->pp_init = g_pairing_pp_init;
    pairing->pp_clear = g_pairing_pp_clear;
    pairing->pp_apply = g_pairing_pp_apply;

    pairing->clear_func = g_pairing_clear;

    element_clear(a);
    element_clear(b);
}

static void compute_cm_curve(g_param_ptr param, cm_info_ptr cm)
    //computes a curve and sets fp to the field it is defined over
    //using the complex multiplication method, where cm holds
    //the appropriate information (e.g. discriminant, field order)
{
    darray_t coefflist;
    element_t hp, root;
    field_t fp, fpx;
    int i, n;
    field_t cc;

    field_init_fp(fp, cm->q);
    field_init_poly(fpx, fp);
    element_init(hp, fpx);

    darray_init(coefflist);

    hilbert_poly(coefflist, cm->D);

    n = coefflist->count;
    poly_alloc(hp, n);
    for (i=0; i<n; i++) {
	element_set_mpz(poly_coeff(hp, i), coefflist->item[i]);
    }

    hilbert_poly_clear(coefflist);

    darray_clear(coefflist);
    //TODO: remove x = 0, 1728 roots
    //TODO: what if there's no roots?
    //printf("hp ");
    //element_out_str(stdout, 0, hp);
    //printf("\n");

    element_init(root, fp);
    findroot(root, hp);
    //printf("root = ");
    //element_out_str(stdout, 0, root);
    //printf("\n");
    element_clear(hp);
    field_clear(fpx);

    //the root is the j-invariant of our desired curve
    field_init_curve_j(cc, root, cm->n, NULL);
    element_clear(root);

    //we may need to twist it however
    {
	element_t P;

	//pick a random point P and see if it has the right order
	element_init(P, cc);
	element_random(P);
	element_mul_mpz(P, P, cm->n);
	//element_printf("P = %B", P);
	//if not, we twist the curve
	if (!element_is0(P)) {
	    twist_curve(cc);
	}
	element_clear(P);
    }

    mpz_set(param->q, cm->q);
    mpz_set(param->n, cm->n);
    mpz_set(param->h, cm->h);
    mpz_set(param->r, cm->r);
    element_to_mpz(param->a, curve_field_a_coeff(cc));
    element_to_mpz(param->b, curve_field_b_coeff(cc));
    {
	mpz_t z;
	mpz_init(z);
	//compute order of curve in F_q^k
	//n = q - t + 1 hence t = q - n + 1
	mpz_sub(z, param->q, param->n);
	mpz_add_ui(z, z, 1);
	compute_trace_n(z, param->q, z, 10);
	mpz_pow_ui(param->nk, param->q, 10);
	mpz_sub_ui(z, z, 1);
	mpz_sub(param->nk, param->nk, z);
	mpz_mul(z, param->r, param->r);
	mpz_divexact(param->hk, param->nk, z);
	mpz_clear(z);
    }
    field_clear(cc);
    field_clear(fp);
}

void g_param_from_cm(g_param_t param, cm_info_ptr cm)
{
    field_t Fq, Fqx, Fqd;
    element_t irred, nqr;
    int i;

    compute_cm_curve(param, cm);

    field_init_fp(Fq, param->q);
    field_init_poly(Fqx, Fq);
    element_init(irred, Fqx);
    do {
	poly_random_monic(irred, 5);
    } while (!poly_is_irred(irred));
    field_init_polymod(Fqd, irred);

    //find a quadratic nonresidue of Fqd lying in Fq
    element_init(nqr, Fqd);
    do {
	element_random(((element_t *) nqr->data)[0]);
    } while (element_is_sqr(nqr));

    param->coeff = pbc_realloc(param->coeff, sizeof(mpz_t) * 5);

    for (i=0; i<5; i++) {
	mpz_init(param->coeff[i]);
	element_to_mpz(param->coeff[i], poly_coeff(irred, i));
    }
    element_to_mpz(param->nqr, ((element_t *) nqr->data)[0]);

    element_clear(nqr);
    element_clear(irred);

    field_clear(Fqx);
    field_clear(Fqd);
    field_clear(Fq);
}