more ergonomic options management

cbouilla · Sep 27, 2023 · 0bd5b41 · 0bd5b41
1 parent 54bcf2c
commit 0bd5b41
Show file tree

Hide file tree

Showing 7 changed files with 297 additions and 89 deletions.
diff --git a/README.md b/README.md
@@ -10,6 +10,7 @@ Features
 --------
 
 The core of the library is a multithreaded function that computes the row echelon form of a sparse matrix modulo a word-sized prime. 
+
 This enables several of useful computations on sparse matrices modulo _p_:
   * Basis of the row space
   * Basis of the kernel
@@ -18,7 +19,8 @@ This enables several of useful computations on sparse matrices modulo _p_:
 
 In addition, SpaSM is capable of computing a full PLUQ factorization, but this is slower than simple echelonization. This enables:
   * Solution of linear systems
-  * Rank certificates
+  * Extraction of a square submatrix of maximal rank
+  * Production of rank certificates
 
 SpaSM works with all odd 32-bit prime moduli.  While it would be possible to work with _p = 2_, this would require non-trivial tweaks to the code. 
 

diff --git a/src/spasm.h b/src/spasm.h
@@ -15,6 +15,10 @@ typedef int32_t i32;
 #include <omp.h>
 #endif
 
+#define SPASM_VERSION "1.3"
+#define SPASM_BUG_ADDRESS "<charles.bouillaguet@lip6.fr>"
+
+
 /* --- primary SpaSM routines and data structures --- */
 
 // unfortunately we use "n" for #rows and "m" for #columns whereas the rest of the world (BLAS...)

diff --git a/tools/CMakeLists.txt b/tools/CMakeLists.txt
@@ -1,5 +1,7 @@
 include_directories (../src)
 
+add_library(spasmtools common.c)
+
 ############### tools
 
 add_executable(bitmap bitmap.c)
@@ -24,10 +26,10 @@ add_executable(dm dm.c)
 target_link_libraries(dm PUBLIC spasm)
 
 add_executable(rank rank.c)
-target_link_libraries(rank PUBLIC spasm)
+target_link_libraries(rank PUBLIC spasm spasmtools)
 
 add_executable(echelonize echelonize.c)
-target_link_libraries(echelonize PUBLIC spasm)
+target_link_libraries(echelonize PUBLIC spasm spasmtools)
 
 add_executable(kernel kernel.c)
-target_link_libraries(kernel PUBLIC spasm)
+target_link_libraries(kernel PUBLIC spasm spasmtools)
diff --git a/tools/common.c b/tools/common.c
@@ -0,0 +1,134 @@
+#include <stdio.h>
+#include <stdlib.h>
+#include <argp.h>
+#include <err.h>
+
+#include "spasm.h"
+#include "common.h"
+
+const char *argp_program_version = SPASM_VERSION;
+const char *argp_program_bug_address = SPASM_BUG_ADDRESS;
+
+struct argp_option input_options[] = {
+	{0,               0,   0,     0, "Input problem", 1 },
+	{"matrix",       'm', "FILE", 0, "Read the input matrix from FILE", 1 },
+	{"modulus",      'p', "P",    0, "Perform arithmetic modulo P", 1 },
+	{0}
+};
+
+/* Parse a single option --- rank program. */
+error_t parse_input_opt(int key, char *arg, struct argp_state *state)
+{
+	struct input_matrix *arguments = state->input;
+	switch (key) {
+	case ARGP_KEY_INIT:  /* set defaults */
+		arguments->filename = NULL;
+		arguments->prime = 42013;
+		break;
+	case 'm':
+		arguments->filename = arg;
+		break;
+	case 'p':
+		arguments->prime = atoll(arg);
+		break;
+	default:
+		return ARGP_ERR_UNKNOWN;
+	}
+	return 0;
+}
+
+struct argp input_argp = { input_options, parse_input_opt, NULL, NULL, NULL, NULL, NULL};
+
+
+/* The options of the echelonization code */
+enum ech_opt_key {
+	NO_LOW_RANK, NO_DENSE, NO_GPLU,
+	MAX_ITER, DENSE_THR, MIN_PIV_RATIO,
+	DENSE_BLKSZ, MIN_RANK_RATIO, MAX_ASPECT_RATIO
+};
+
+struct argp_option echelonize_options[] = {
+	{0,                     0,                 0, 0, "Echelonization sub-algorithms", -2},
+	{"no-low-rank-mode",    NO_LOW_RANK,       0, 0, "Disable the (dense) low-rank mode", -2 },
+	{"no-dense-mode",       NO_DENSE,          0, 0, "Don't use FFPACK", -2 },
+	{"no-GPLU",             NO_GPLU,           0, 0, "Don't use GPLU", -2 },
+
+	{0,                     0,                 0,  0, "Main echelonization options", -3 },
+	{"max-iterations",      MAX_ITER,         "N", 0, "Compute at most N sparse Schur complements ", -3},
+	{"dense-threshold",     DENSE_THR,        "D", 0, "Use FFPACK when the density is greater than D", -3},
+	{"min-pivot-proportion", MIN_PIV_RATIO,   "P", 0, "Stop if the proportion of pivots found is less than P", -3},
+
+	{0,                     0,                 0,  0, "Dense code options", -4},
+	{"dense-block-size",    DENSE_BLKSZ,      "N", 0, "Use dense matrices of at most N rows", -4},
+	{"min-rank-ratio",      MIN_RANK_RATIO,   "R", 0, "Use low-rank mode if k rows have rank <= k * X", -4},
+	{"max-aspect-ratio",    MAX_ASPECT_RATIO, "R", 0, "Use low-rank mode if #rows / #columns >= R", -4},
+
+	{ 0 }
+};
+
+/* Parse a single option --- rank program. */
+error_t parse_ech_opt(int key, char *arg, struct argp_state *state)
+{
+	struct echelonize_opts *opts = state->input;
+	switch (key) {
+	case ARGP_KEY_INIT:  /* set defaults */
+		spasm_echelonize_init_opts(opts);
+		break;
+	case NO_LOW_RANK:
+		opts->enable_tall_and_skinny = 0;
+		break;
+	case NO_DENSE:
+		opts->enable_dense = 0;
+		break;
+	case NO_GPLU:
+		opts->enable_GPLU = 0;
+		break;
+	case MAX_ITER:
+		opts->max_round = atoi(arg);
+		break;
+	case DENSE_THR:
+		opts->sparsity_threshold = atof(arg);
+		break;
+	case MIN_PIV_RATIO:
+		opts->min_pivot_proportion = atof(arg);
+		break;
+	case DENSE_BLKSZ:
+		opts->dense_block_size = atoi(arg);
+		break;
+	case MIN_RANK_RATIO:
+		opts->low_rank_ratio = atof(arg);
+		break;
+	case MAX_ASPECT_RATIO:
+		opts->tall_and_skinny_ratio = atof(arg);
+		break;
+	default:
+		return ARGP_ERR_UNKNOWN;
+	}
+	return 0;
+}
+
+struct argp echelonize_argp = { echelonize_options, parse_ech_opt, NULL, NULL, NULL, NULL, NULL};
+
+
+spasm_triplet * load_input_matrix(struct input_matrix *in)
+{
+	if (in->filename == NULL)
+		return spasm_load_sms(stdin, in->prime);    /* load from stdin */
+	FILE *f = fopen(in->filename, "r");
+	if (f == NULL)
+		err(1, "Cannot open %s", in->filename);
+	spasm_triplet *T = spasm_load_sms(f, in->prime);
+	fclose(f);
+	return T;
+}
+
+
+FILE * open_output(const char *filename)
+{
+	if (filename == NULL)
+		return stdout;
+	FILE *f = fopen(filename, "w");
+	if (f == NULL)
+		err(1, "Cannot open %s", filename);
+	return f;
+}
diff --git a/tools/common.h b/tools/common.h
@@ -0,0 +1,20 @@
+#ifndef _SPASM_COMMON_H
+#define _SPASM_COMMON_H
+
+#include <stdio.h>
+#include <argp.h>
+
+struct input_matrix {
+	char *filename;
+	i64 prime;
+};
+
+extern struct argp echelonize_argp;
+extern struct argp input_argp;
+extern const char *argp_program_version;
+extern const char *argp_program_bug_address;
+
+spasm_triplet * load_input_matrix(struct input_matrix *in);
+FILE * open_output(const char *filename);
+
+#endif
diff --git a/tools/kernel.c b/tools/kernel.c
@@ -5,72 +5,87 @@
 #include <err.h>
 
 #include "spasm.h"
+#include "common.h"
 
-i64 prime = 42013;
-struct echelonize_opts opts;
-bool left = 0;
+/* Program documentation. */
+char doc[] = "Compute a kernel basis of a sparse matrix";
 
-void parse_command_line_options(int argc, char **argv)
+struct cmdline_args {
+	/* input problem */
+	struct input_matrix input;
+
+	/* common options to all programs that call spasm_echelonize() */
+	struct echelonize_opts opts;
+
+	/* options specific to the kernel program */
+	bool left;
+	char *output_filename;
+};
+
+/* The options we understand. */
+struct argp_option options[] = {
+	{0,               0,  0,      0, "Kernel options", 2 },
+	{"left",         'l', 0,      0, "Compute the left-kernel", 2},
+	{"output",       'o', "FILE", 0, "Write the kernel basis in FILE", 2 },
+	{ 0 }
+};
+
+/* Parse a single option --- rank program. */
+error_t parse_ker_opt(int key, char *arg, struct argp_state *state)
 {
-	struct option longopts[] = {
-		{"modulus", required_argument, NULL, 'p'},
-		{"left", no_argument, NULL, 'L'},
-		{"no-greedy-pivot-search", no_argument, NULL, 'g'},
-		{"no-low-rank-mode", no_argument, NULL, 'l'},
-		{"dense-block-size", required_argument, NULL, 'd'},
-		{"low-rank-start-weight", required_argument, NULL, 'w'},
-		{NULL, 0, NULL, 0}
-	};
-	char ch;
-	while ((ch = getopt_long(argc, argv, "", longopts, NULL)) != -1) {
-		switch (ch) {
-		case 'p':
-			prime = atoll(optarg);
-			break;
-		case 'L':
-			left = 1;
-			break;
-		case 'g':
-			opts.enable_greedy_pivot_search = 0;
-			break;
-		case 'd':
-			opts.dense_block_size = atoi(optarg);
-			fprintf(stderr, "Using dense block size %d\n", opts.dense_block_size);
-			break;
-		case 'l':
-			opts.enable_tall_and_skinny = 0;
-			break;
-		case 'w':
-			opts.low_rank_start_weight = atoi(optarg);
-			break;
-		default:
-			errx(1, "Unknown option\n");
-		}
+	struct cmdline_args *arguments = state->input;
+	switch (key) {
+	case 'l':
+		arguments->left = 1;
+		break;
+	case 'o':
+		arguments->output_filename = arg;
+		break;
+	case ARGP_KEY_INIT:
+		arguments->left = 0;
+		arguments->output_filename = NULL;
+		state->child_inputs[0] = &arguments->input;
+		state->child_inputs[1] = &arguments->opts;
+		break;
+	default:
+		return ARGP_ERR_UNKNOWN;
 	}
+	return 0;
 }
 
+struct argp_child children_parsers[] = {
+	{ &input_argp,      0, 0, 0 },
+	{ &echelonize_argp, 0, 0, 0 },
+	{ 0 }
+};
+
+/* argp parser for the rank program */
+struct argp argp = { options, parse_ker_opt, NULL, doc, children_parsers, NULL, NULL };
+
+
 int main(int argc, char **argv)
 {
-	spasm_echelonize_init_opts(&opts);
-	parse_command_line_options(argc, argv);
+	/* process command-line options */
+	struct cmdline_args args;
+	argp_parse(&argp, argc, argv, 0, 0, &args);
 
-	spasm_triplet *T = spasm_load_sms(stdin, prime);
-	if (left) {
+	/* load input matrix */
+	spasm_triplet *T = load_input_matrix(&args.input);
+	if (args.left) {
 		fprintf(stderr, "Left-kernel, transposing\n");
 		spasm_triplet_transpose(T);
 	}
 	spasm *A = spasm_compress(T);
 	spasm_triplet_free(T);
 
 	/* echelonize A */
-	spasm_lu *fact = spasm_echelonize(A, &opts);
+	spasm_lu *fact = spasm_echelonize(A, &args.opts);
 	spasm_csr_free(A);
 
 	/* kernel basis */
 	spasm *K = spasm_kernel(fact);
-	spasm_save_csr(stdout, K);
 	fprintf(stderr, "Kernel basis matrix is %d x %d with %" PRId64 " nz\n", K->n, K->m, spasm_nnz(K));
-	spasm_lu_free(fact);
-	spasm_csr_free(K);
-	exit(EXIT_SUCCESS);
-}
+
+	FILE *f = open_output(args.output_filename);
+	spasm_save_csr(f, K);
+}