9 files changed, 501 insertions, 35 deletions

diff --git a/src/cm/anomalous.c b/src/cm/anomalous.c
new file mode 100644
index 0000000..70df8e4
--- /dev/null
+++ b/src/cm/anomalous.c
@@ -0,0 +1,140 @@
+/*
+ * ecgen, tool for generating Elliptic curve domain parameters
+ * Copyright (C) 2017-2018 J08nY
+ */
+#include "anomalous.h"
+#include "io/output.h"
+#include "util/memory.h"
+
+static disc_t **disc_table;
+
+void anomalous_init() {
+	disc_table = try_calloc(sizeof(disc_t *) * 5);
+	for (int i = 0; i < 5; ++i) {
+		disc_table[i] = try_calloc(sizeof(disc_t));
+	}
+
+	/*
+	 * Discriminant, j-invariant and field-element (Alpha) from Atsuko Miyaji's
+	 * paper.
+	 */
+	GEN a = stoi(-32);
+	GEN b = stoi(-64);
+	disc_table[0]->d = stoi(11);
+	disc_table[0]->j = powis(a, 3);
+	disc_table[0]->alpha = stoi(21);
+	disc_table[1]->d = stoi(19);
+	disc_table[1]->j = powis(mulis(a, 3), 3);
+	disc_table[1]->alpha = stoi(3);
+	disc_table[2]->d = stoi(43);
+	disc_table[2]->j = powis(mulis(b, 16), 3);
+	disc_table[2]->alpha = stoi(70);
+	disc_table[3]->d = stoi(67);
+	disc_table[3]->j = powis(mulis(a, 165), 3);
+	disc_table[3]->alpha = stoi(35805);
+	disc_table[3]->d = stoi(163);
+	disc_table[3]->j = powis(mulis(b, 10005), 3);
+	disc_table[3]->alpha = stoi(3717878010);
+}
+
+static GEN anomalous_prime(size_t i, unsigned long bits) {
+	GEN D = disc_table[i]->d;
+
+	pari_sp ltop = avma;
+	GEN last = divis(addis(D, 1), 4);
+	GEN lower = divii(subii(int2n(bits - 1), last), D);
+	GEN upper = divii(subii(int2n(bits), last), D);
+
+	GEN lower_bound = gceil(
+	    gdiv(gsub(gsqrt(addis(mulis(lower, 4), 1), BIGDEFAULTPREC * 2), gen_1),
+	         gen_2));
+	GEN upper_bound = gfloor(
+	    gdiv(gsub(gsqrt(addis(mulis(upper, 4), 1), BIGDEFAULTPREC * 2), gen_1),
+	         gen_2));
+
+	GEN range = gtovec0(gen_0, 2);
+	gel(range, 1) = lower_bound;
+	gel(range, 2) = upper_bound;
+	GEN first, middle;
+	GEN p, b;
+
+	pari_sp btop = avma;
+	do {
+		b = genrand(range);
+		middle = mulii(D, b);
+		first = mulii(middle, b);
+		p = addii(addii(first, middle), last);
+		gerepileall(btop, 2, &p, &last);
+	} while (!isprime(p));
+
+	return gerepilecopy(ltop, p);
+}
+
+static GEN anomalous_c(size_t i, GEN p) {
+	pari_sp ltop = avma;
+	long kroneck = kronecker(disc_table[i]->alpha, p);
+	if (!kroneck) {
+		return NULL;
+	}
+	long cp = -kroneck;
+
+	GEN c;
+	pari_sp btop = avma;
+	do {
+		c = genrand(p);
+		gerepileall(btop, 1, &c);
+	} while (kronecker(c, p) != cp);
+
+	return gerepilecopy(ltop, c);
+}
+
+GENERATOR(anomalous_gen_field) {
+	HAS_ARG(args);
+	size_t i = *(size_t *)args->args;
+
+	// find suitable prime field
+	pari_sp ltop = avma;
+	GEN p;
+	do {
+		p = anomalous_prime(i, cfg->bits);
+		gerepileall(ltop, 1, &p);
+	} while (!kronecker(disc_table[0]->alpha, p));
+
+	curve->field = gerepilecopy(ltop, p);
+	return 1;
+}
+
+GENERATOR(anomalous_gen_equation) {
+	HAS_ARG(args);
+	size_t i = *(size_t *)args->args;
+
+	pari_sp ltop = avma;
+	GEN c = anomalous_c(i, curve->field);
+
+	GEN a_upper = gmodulo(disc_table[i]->j, curve->field);
+	GEN a_lower = gmodulo(subsi(1728, disc_table[i]->j), curve->field);
+	GEN a = gdiv(a_upper, a_lower);
+
+	curve->a = gmul(stoi(3), gmul(gsqr(c), a));
+	curve->b = gmul(gen_2, gmul(gpowgs(c, 3), a));
+	gerepileall(ltop, 2, &curve->a, &curve->b);
+	return 1;
+}
+
+GENERATOR(anomalous_gen_order) {
+	// copy field to order
+	curve->order = gcopy(curve->field);
+	obj_insert(curve->curve, 1, curve->order);
+	return 1;
+}
+
+void anomalous_quit() {
+	if (disc_table) {
+		for (int i = 0; i < 5; ++i) {
+			if (disc_table[i]) {
+				try_free(disc_table[i]);
+			}
+		}
+		try_free(disc_table);
+	}
+}
diff --git a/src/cm/anomalous.h b/src/cm/anomalous.h
new file mode 100644
index 0000000..73b84fd
--- /dev/null
+++ b/src/cm/anomalous.h
@@ -0,0 +1,60 @@
+/*
+ * ecgen, tool for generating Elliptic curve domain parameters
+ * Copyright (C) 2017-2018 J08nY
+ */
+/**
+ * @file anomalous.h
+ */
+#ifndef ECGEN_CM_ANOMALOUS_H
+#define ECGEN_CM_ANOMALOUS_H
+
+#include <pari/pari.h>
+#include "exhaustive/arg.h"
+#include "io/cli.h"
+#include "misc/types.h"
+
+typedef struct {
+	GEN d;
+	GEN j;
+	GEN alpha;
+} disc_t;
+
+/**
+ * GENERATOR(gen_f)
+ *
+ * @param curve A curve_t being generated
+ * @param args the index of the discriminant to use, in the disc_table
+ * @return state diff
+ */
+GENERATOR(anomalous_gen_field);
+
+/**
+ * GENERATOR(gen_f)
+ *
+ * @param curve A curve_t being generated
+ * @param args the index of the discriminant to use, in the disc_table
+ * @return state diff
+ */
+GENERATOR(anomalous_gen_equation);
+
+/**
+ * GENERATOR(gen_f)
+ *
+ * @param curve A curve_t being generated
+ * @param args unused
+ * @return state diff
+ */
+GENERATOR(anomalous_gen_order);
+
+/**
+ * @brief Initialize anomalous generation, allocate and set the disc_table.
+ */
+void anomalous_init();
+
+/**
+ * @brief Deinitialize anomalous generation, free the discriminants from the
+ * disc_table.
+ */
+void anomalous_quit();
+
+#endif  // ECGEN_CM_ANOMALOUS_H
diff --git a/src/cm/cm.c b/src/cm/cm.c
index a29e368..239bb5d 100644
--- a/src/cm/cm.c
+++ b/src/cm/cm.c
@@ -3,39 +3,163 @@
  * Copyright (C) 2017-2018 J08nY
  */
 #include "cm.h"
+#include "anomalous.h"
 #include "cm_any.h"
 #include "cm_prime.h"
-#include "io/output.h"
+#include "exhaustive/check.h"
+#include "exhaustive/exhaustive.h"
+#include "gen/curve.h"
+#include "gen/field.h"
+#include "gen/gens.h"
+#include "gen/hex.h"
+#include "gen/metadata.h"
+#include "gen/point.h"
 #include "obj/curve.h"
-#include "p1363.h"
+#include "supersingular.h"
+
+static void cm_ginit(gen_f *generators, bool prime) {
+	// SEED unused.
+	generators[OFFSET_SEED] = &gen_skip;
+
+	// Setup stuff so it can be overridden.
+	if (cfg->unique) {
+		generators[OFFSET_GENERATORS] = &gens_gen_one;
+	} else {
+		generators[OFFSET_GENERATORS] = &gens_gen_any;
+	}
+
+	if (cfg->metadata) {
+		generators[OFFSET_METADATA] = &metadata_gen;
+	} else {
+		generators[OFFSET_METADATA] = &gen_skip;
+	}
+
+	switch (cfg->points.type) {
+		case POINTS_RANDOM:
+			if (cfg->points.amount) {
+				generators[OFFSET_POINTS] = &points_gen_random;
+			} else {
+				generators[OFFSET_POINTS] = &point_gen_random;
+			}
+			break;
+		case POINTS_PRIME:
+			generators[OFFSET_POINTS] = &points_gen_prime;
+			break;
+		case POINTS_NONPRIME:
+			generators[OFFSET_POINTS] = &points_gen_nonprime;
+			break;
+		case POINTS_ALL:
+			generators[OFFSET_POINTS] = &points_gen_allgroups;
+			break;
+		case POINTS_NONE:
+			generators[OFFSET_POINTS] = &gen_skip;
+			break;
+	}
+
+	// Now do the actual CM setup.
+	if (cfg->method == METHOD_CM) {
+		generators[OFFSET_FIELD] = &gen_skip;
+		generators[OFFSET_A] = &gen_skip;
+		generators[OFFSET_B] = &gen_skip;
+		if (prime) {
+			generators[OFFSET_CURVE] = &cm_gen_curve_prime;
+			generators[OFFSET_POINTS] = &points_gen_prime;
+		} else {
+			generators[OFFSET_CURVE] = &cm_gen_curve_any;
+		}
+		generators[OFFSET_ORDER] = &cm_gen_order;
+	} else if (cfg->method == METHOD_ANOMALOUS) {
+		generators[OFFSET_FIELD] = &anomalous_gen_field;
+		generators[OFFSET_A] = &gen_skip;
+		generators[OFFSET_B] = &anomalous_gen_equation;
+		generators[OFFSET_CURVE] = &curve_gen_any;
+		generators[OFFSET_ORDER] = &anomalous_gen_order;
+	} else if (cfg->method == METHOD_SUPERSINGULAR) {
+		if (cfg->random & RANDOM_FIELD) {
+			generators[OFFSET_FIELD] = &field_gen_random;
+		} else {
+			generators[OFFSET_FIELD] = &field_gen_input;
+		}
+		generators[OFFSET_A] = &gen_skip;
+		generators[OFFSET_B] = &supersingular_gen_equation;
+		generators[OFFSET_CURVE] = &curve_gen_any;
+		generators[OFFSET_ORDER] = &supersingular_gen_order;
+	}
+}
+
+static void cm_ainit(arg_t **gen_argss, arg_t **check_argss) {
+	if (cfg->method == METHOD_ANOMALOUS) {
+		arg_t *field_arg = arg_new();
+		arg_t *eq_arg = arg_new();
+		size_t *i = try_calloc(sizeof(size_t));
+		*i = 3;
+		field_arg->args = i;
+		field_arg->nargs = 1;
+		eq_arg->args = i;
+		eq_arg->nargs = 1;
+		eq_arg->allocd = i;
+		gen_argss[OFFSET_FIELD] = field_arg;
+		gen_argss[OFFSET_B] = eq_arg;
+	}
+
+	if (cfg->points.type == POINTS_RANDOM) {
+		arg_t *points_arg = arg_new();
+		points_arg->args = &cfg->points.amount;
+		points_arg->nargs = 1;
+		gen_argss[OFFSET_POINTS] = points_arg;
+	}
+}
+
+static void cm_cinit(check_t **validators) {
+	check_t *curve_check = check_new(curve_check_nonzero, NULL);
+	validators[OFFSET_CURVE] = curve_check;
+
+	if (cfg->hex_check) {
+		check_t *hex_check = check_new(hex_check_param, NULL);
+		validators[OFFSET_POINTS] = hex_check;
+	}
+}
 
 int cm_do() {
 	debug_log_start("Starting Complex Multiplication method");
 
-	int result = 0;
-	GEN order = strtoi(cfg->cm_order);
-	curve_t *curve = NULL;
+	gen_f generators[OFFSET_END] = {NULL};
+	arg_t *gen_argss[OFFSET_END] = {NULL};
+	check_t *validators[OFFSET_END] = {NULL};
+	arg_t *check_argss[OFFSET_END] = {NULL};
+	unroll_f unrolls[OFFSET_END] = {NULL};
 
-	if (gequal0(order)) {
-		fprintf(err, "Order requested not a number: %s\n", cfg->cm_order);
-		result = 1;
-	} else if (isprime(order)) {
-		debug_log("Starting prime order curve generation");
-		curve = cm_prime_curve(order);
-	} else {
-		debug_log("Starting composite order curve generation");
-		curve = cm_any_curve(order);
+	exhaustive_t setup = {.generators = generators,
+	                      .gen_argss = gen_argss,
+	                      .validators = validators,
+	                      .check_argss = check_argss,
+	                      .unrolls = unrolls};
+
+	bool ord_prime = false;
+	if (cfg->method == METHOD_CM) {
+		GEN order = strtoi(cfg->cm_order);
+		if (gequal0(order)) {
+			fprintf(err, "Order requested not a number: %s\n", cfg->cm_order);
+			return 1;
+		}
+		ord_prime = (bool)isprime(order);
+	}
+
+	if (cfg->method == METHOD_ANOMALOUS) {
+		anomalous_init();
 	}
 
-	if (curve) {
-		output_o_begin();
-		output_o(curve);
-		output_o_end();
+	cm_ginit(setup.generators, ord_prime);
+	cm_ainit(setup.gen_argss, setup.check_argss);
+	cm_cinit(setup.validators);
+	exhaustive_uinit(setup.unrolls);
 
-		curve_free(&curve);
-	} else {
-		result = 1;
+	int result = exhaustive_generate(&setup);
+
+	if (cfg->method == METHOD_ANOMALOUS) {
+		anomalous_quit();
 	}
+	exhaustive_clear(&setup);
 
 	debug_log_start("Finished Complex Multiplication method");
 	return result;
diff --git a/src/cm/cm_any.c b/src/cm/cm_any.c
index 401a544..f993fa3 100644
--- a/src/cm/cm_any.c
+++ b/src/cm/cm_any.c
@@ -3,6 +3,7 @@
  * Copyright (C) 2017-2018 J08nY
  */
 #include "cm_any.h"
+#include <misc/config.h>
 #include <obj/obj.h>
 #include "io/output.h"
 #include "obj/curve.h"
@@ -13,9 +14,8 @@
  *        Constructing elliptic curves of prescribed order,
  *        Reiner Broker
  * @param order
- * @return
  */
-static cm_any_qdisc_t *good_qdisc_minimal(GEN order) {
+static void good_qdisc_minimal(cm_any_qdisc_t *qdisc, GEN order) {
 	pari_sp ltop = avma;
 	GEN d = stoi(2);
 	while (true) {
@@ -37,11 +37,10 @@ static cm_any_qdisc_t *good_qdisc_minimal(GEN order) {
 					debug_log(
 					    "Got an elem of prime trace: %Pi, d = %Pi, D = %Pi", p,
 					    d, D);
-					cm_any_qdisc_t *result = try_calloc(sizeof(cm_any_qdisc_t));
-					result->p = p;
-					result->d = D;
-					gerepileall(ltop, 2, &result->p, &result->d);
-					return result;
+					qdisc->p = p;
+					qdisc->d = D;
+					gerepileall(ltop, 2, &qdisc->p, &qdisc->d);
+					return;
 				}
 			}
 		}
@@ -200,26 +199,50 @@ GEN cm_construct_curve(GEN order, GEN d, GEN p, bool ord_prime) {
 					debug_log("Got curve twist.");
 					return gerepilecopy(ltop, twist);
 				}
-			};
+			}
 		}
 	}
 	return NULL;
 }
 
 curve_t *cm_any_curve(GEN order) {
-	cm_any_qdisc_t *min_disc = good_qdisc_minimal(order);
-	debug_log("Got min D = %Pi", min_disc->d);
-	GEN e = cm_construct_curve(order, min_disc->d, min_disc->p, false);
+	cm_any_qdisc_t min_disc = {0};
+	good_qdisc_minimal(&min_disc, order);
+	debug_log("Got min D = %Pi", min_disc.d);
+	GEN e = cm_construct_curve(order, min_disc.d, min_disc.p, false);
 	if (e == NULL) {
 		fprintf(err, "Could not construct curve.");
 		return NULL;
 	}
 	curve_t *curve = curve_new();
-	curve->field = min_disc->p;
+	curve->field = min_disc.p;
 	curve->curve = e;
 	curve->a = ell_get_a4(e);
 	curve->b = ell_get_a6(e);
 	curve->order = gcopy(order);
-	try_free(min_disc);
 	return curve;
+}
+
+GENERATOR(cm_gen_curve_any) {
+	pari_sp ltop = avma;
+	GEN order = strtoi(cfg->cm_order);
+	cm_any_qdisc_t min_disc = {0};
+	good_qdisc_minimal(&min_disc, order);
+	debug_log("Got min D = %Pi", min_disc.d);
+	GEN e = cm_construct_curve(order, min_disc.d, min_disc.p, false);
+	if (e == NULL) {
+		fprintf(err, "Could not construct curve.");
+		avma = ltop;
+		return -3;
+	}
+	curve->field = min_disc.p;
+	curve->a = ell_get_a4(e);
+	curve->b = ell_get_a6(e);
+	curve->curve = e;
+	return 1;
+}
+
+GENERATOR(cm_gen_order) {
+	curve->order = strtoi(cfg->cm_order);
+	return 1;
 }
\ No newline at end of file
diff --git a/src/cm/cm_any.h b/src/cm/cm_any.h
index eb54497..48dee1c 100644
--- a/src/cm/cm_any.h
+++ b/src/cm/cm_any.h
@@ -29,4 +29,22 @@ GEN cm_construct_curve(GEN order, GEN d, GEN p, bool ord_prime);
  */
 curve_t *cm_any_curve(GEN order);
 
+/**
+ * @brief
+ * @param curve
+ * @param args
+ * @param state
+ * @return
+ */
+GENERATOR(cm_gen_curve_any);
+
+/**
+ * @brief
+ * @param curve
+ * @param args
+ * @param state
+ * @return
+ */
+GENERATOR(cm_gen_order);
+
 #endif  // ECGEN_CM_ANY_H
diff --git a/src/cm/cm_prime.c b/src/cm/cm_prime.c
index c7a931d..e3a6882 100644
--- a/src/cm/cm_prime.c
+++ b/src/cm/cm_prime.c
@@ -140,7 +140,7 @@ static void qdisc_free(cm_prime_qdisc_t *qdisc) { try_free(qdisc->Sp); }
 curve_t *cm_prime_curve(GEN order) {
 	GEN e = NULL;
 
-	cm_prime_qdisc_t qdisc;
+	cm_prime_qdisc_t qdisc = {0};
 	qdisc_init(&qdisc, order);
 	do {
 		qdisc_next(&qdisc);
@@ -164,4 +164,23 @@ curve_t *cm_prime_curve(GEN order) {
 	result->ngens = 1;
 
 	return result;
+}
+
+GENERATOR(cm_gen_curve_prime) {
+	GEN order = strtoi(cfg->cm_order);
+	GEN e = NULL;
+
+	cm_prime_qdisc_t qdisc = {0};
+	qdisc_init(&qdisc, order);
+	do {
+		qdisc_next(&qdisc);
+		e = cm_construct_curve(order, qdisc.D, qdisc.p, true);
+	} while (e == NULL);
+	qdisc_free(&qdisc);
+
+	curve->field = qdisc.p;
+	curve->a = ell_get_a4(e);
+	curve->b = ell_get_a6(e);
+	curve->curve = e;
+	return 1;
 }
\ No newline at end of file
diff --git a/src/cm/cm_prime.h b/src/cm/cm_prime.h
index 59debac..9b25100 100644
--- a/src/cm/cm_prime.h
+++ b/src/cm/cm_prime.h
@@ -31,4 +31,13 @@ typedef struct {
  */
 curve_t* cm_prime_curve(GEN order);
 
+/**
+ * @brief
+ * @param curve
+ * @param args
+ * @param state
+ * @return
+ */
+GENERATOR(cm_gen_curve_prime);
+
 #endif  // ECGEN_CM_PRIME_H
diff --git a/src/cm/supersingular.c b/src/cm/supersingular.c
new file mode 100644
index 0000000..a3cebfc
--- /dev/null
+++ b/src/cm/supersingular.c
@@ -0,0 +1,45 @@
+/*
+ * ecgen, tool for generating Elliptic curve domain parameters
+ * Copyright (C) 2017-2018 J08nY
+ */
+#include "supersingular.h"
+
+GENERATOR(supersingular_gen_equation) {
+	if (equalis(curve->field, 2)) {
+		return -2;
+	}
+	if (mod4(curve->field) == 3) {
+		curve->a = mkintmod(subis(curve->field, 1), curve->field);
+		curve->b = mkintmod(stoi(0), curve->field);
+		return 1;
+	}
+	GEN q = stoi(3);
+	while (!(mod4(q) == 3 && kronecker(curve->field, q) == -1)) {
+		q = nextprime(addis(q, 1));
+	}
+
+	if (equalis(q, 3)) {
+		curve->a = mkintmod(stoi(0), curve->field);
+		curve->b = mkintmod(stoi(1), curve->field);
+		return 1;
+	} else {
+		GEN H = polclass(negi(q), 0, 0);
+		GEN r = FpX_roots(H, curve->field);
+		GEN root = gel(r, 1);
+		curve->a = mkintmod(
+		    Fp_div(Fp_mul(stoi(27), root, curve->field),
+		           Fp_mul(stoi(4), Fp_sub(stoi(1728), root, curve->field),
+		                  curve->field),
+		           curve->field),
+		    curve->field);
+		curve->b = gneg(curve->a);
+		return 1;
+	}
+}
+
+GENERATOR(supersingular_gen_order) {
+	// copy field to order
+	curve->order = addis(curve->field, 1);
+	obj_insert(curve->curve, 1, curve->order);
+	return 1;
+}
diff --git a/src/cm/supersingular.h b/src/cm/supersingular.h
new file mode 100644
index 0000000..7755c49
--- /dev/null
+++ b/src/cm/supersingular.h
@@ -0,0 +1,28 @@
+/*
+ * ecgen, tool for generating Elliptic curve domain parameters
+ * Copyright (C) 2017-2018 J08nY
+ */
+#ifndef ECGEN_CM_SUPERSINGULAR_H
+#define ECGEN_CM_SUPERSINGULAR_H
+
+#include "misc/types.h"
+
+/**
+ * @brief
+ * @param curve
+ * @param args
+ * @param state
+ * @return
+ */
+GENERATOR(supersingular_gen_equation);
+
+/**
+ * @brief
+ * @param curve
+ * @param args
+ * @param state
+ * @return
+ */
+GENERATOR(supersingular_gen_order);
+
+#endif  // ECGEN_CM_SUPERSINGULAR_H