(Parallel) timing code
This commit is contained in:
Miguel M
parent
c7fe98e620
commit
3ef2e95c69
|
|
@ -0,0 +1,171 @@
|
|||
// Computes an integer power by modular exponentiation.
|
||||
//
|
||||
// From https://gist.github.com/orlp/3551590, in turn optimized from
|
||||
// https://stackoverflow.com/a/101613.
|
||||
//
|
||||
// Adapted.
|
||||
#ifndef ACED_CONTRIB_ORLP_IPOW_H_
|
||||
#define ACED_CONTRIB_ORLP_IPOW_H_
|
||||
|
||||
#include <stdint.h>
|
||||
|
||||
int64_t ipow(int64_t base, uint8_t exp) {
|
||||
static const uint8_t highest_bit_set[] = {
|
||||
0, 1, 2, 2, 3, 3, 3, 3,
|
||||
4, 4, 4, 4, 4, 4, 4, 4,
|
||||
5, 5, 5, 5, 5, 5, 5, 5,
|
||||
5, 5, 5, 5, 5, 5, 5, 5,
|
||||
6, 6, 6, 6, 6, 6, 6, 6,
|
||||
6, 6, 6, 6, 6, 6, 6, 6,
|
||||
6, 6, 6, 6, 6, 6, 6, 6,
|
||||
6, 6, 6, 6, 6, 6, 6, 255, // anything past 63 is a guaranteed overflow with base > 1
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
};
|
||||
|
||||
int64_t result = 1;
|
||||
|
||||
switch (highest_bit_set[exp]) {
|
||||
case 255: // we use 255 as an overflow marker and return 0 on overflow/underflow
|
||||
if (base == 1) {
|
||||
return 1;
|
||||
}
|
||||
|
||||
if (base == -1) {
|
||||
return 1 - 2 * (exp & 1);
|
||||
}
|
||||
|
||||
return 0;
|
||||
case 6:
|
||||
if (exp & 1) result *= base;
|
||||
exp >>= 1;
|
||||
base *= base;
|
||||
// fall through
|
||||
case 5:
|
||||
if (exp & 1) result *= base;
|
||||
exp >>= 1;
|
||||
base *= base;
|
||||
// fall through
|
||||
case 4:
|
||||
if (exp & 1) result *= base;
|
||||
exp >>= 1;
|
||||
base *= base;
|
||||
// fall through
|
||||
case 3:
|
||||
if (exp & 1) result *= base;
|
||||
exp >>= 1;
|
||||
base *= base;
|
||||
// fall through
|
||||
case 2:
|
||||
if (exp & 1) result *= base;
|
||||
exp >>= 1;
|
||||
base *= base;
|
||||
// fall through
|
||||
case 1:
|
||||
if (exp & 1) result *= base;
|
||||
// fall through
|
||||
default:
|
||||
return result;
|
||||
}
|
||||
}
|
||||
|
||||
uint64_t upow(uint64_t base, uint8_t exp) {
|
||||
static const uint8_t highest_bit_set[] = {
|
||||
0, 1, 2, 2, 3, 3, 3, 3,
|
||||
4, 4, 4, 4, 4, 4, 4, 4,
|
||||
5, 5, 5, 5, 5, 5, 5, 5,
|
||||
5, 5, 5, 5, 5, 5, 5, 5,
|
||||
6, 6, 6, 6, 6, 6, 6, 6,
|
||||
6, 6, 6, 6, 6, 6, 6, 6,
|
||||
6, 6, 6, 6, 6, 6, 6, 6,
|
||||
6, 6, 6, 6, 6, 6, 6, 255, // anything past 63 is a guaranteed overflow with base > 1
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
255, 255, 255, 255, 255, 255, 255, 255,
|
||||
};
|
||||
|
||||
uint64_t result = 1;
|
||||
|
||||
switch (highest_bit_set[exp]) {
|
||||
case 255: // we use 255 as an overflow marker and return 0 on overflow/underflow
|
||||
if (base == 1) {
|
||||
return 1;
|
||||
}
|
||||
return 0;
|
||||
case 6:
|
||||
if (exp & 1) result *= base;
|
||||
exp >>= 1;
|
||||
base *= base;
|
||||
// fall through
|
||||
case 5:
|
||||
if (exp & 1) result *= base;
|
||||
exp >>= 1;
|
||||
base *= base;
|
||||
// fall through
|
||||
case 4:
|
||||
if (exp & 1) result *= base;
|
||||
exp >>= 1;
|
||||
base *= base;
|
||||
// fall through
|
||||
case 3:
|
||||
if (exp & 1) result *= base;
|
||||
exp >>= 1;
|
||||
base *= base;
|
||||
// fall through
|
||||
case 2:
|
||||
if (exp & 1) result *= base;
|
||||
exp >>= 1;
|
||||
base *= base;
|
||||
// fall through
|
||||
case 1:
|
||||
if (exp & 1) result *= base;
|
||||
// fall through
|
||||
default:
|
||||
return result;
|
||||
}
|
||||
}
|
||||
|
||||
#endif
|
||||
File diff suppressed because it is too large
Load Diff
532
src/main.c
532
src/main.c
|
|
@ -1,17 +1,80 @@
|
|||
#include <alloca.h>
|
||||
#include <stdbool.h>
|
||||
#include <stdio.h>
|
||||
#include <stdlib.h>
|
||||
#include <string.h>
|
||||
#include <time.h>
|
||||
#include <unistd.h>
|
||||
|
||||
#define STB_DS_IMPLEMENTATION
|
||||
#define STBDS_NO_SHORT_NAMES
|
||||
#include "../contrib/orlp/ipow.h"
|
||||
#include "../contrib/stb/stb_ds.h"
|
||||
#include "../includes/cg.h"
|
||||
#include "../includes/matrix.h"
|
||||
#include "../includes/permutation.h"
|
||||
|
||||
static void ReadArguments(char *from, size_t *into, const char *varname) {
|
||||
if (!sscanf(from, "%zu", into)) {
|
||||
fprintf(stderr, "Failed to read %s.", varname);
|
||||
exit(EXIT_FAILURE);
|
||||
static uintmax_t _timing_last_count = 0;
|
||||
static clock_t _timing_last_time = 0.;
|
||||
static float _timing_rate = 0.;
|
||||
static float _timing_exp_backoff = 1.;
|
||||
|
||||
static void ReportProgress(uintmax_t current, uintmax_t total, clock_t now) {
|
||||
_Bool updated = 0;
|
||||
#pragma omp critical
|
||||
{
|
||||
if (current >= _timing_last_count) {
|
||||
float inferred_rate = (float)(current - _timing_last_count) /
|
||||
((float)(now - _timing_last_time) / CLOCKS_PER_SEC);
|
||||
_timing_rate = _timing_rate * (1. - _timing_exp_backoff) +
|
||||
inferred_rate * _timing_exp_backoff;
|
||||
if (_timing_exp_backoff > 0.0001) {
|
||||
_timing_exp_backoff *= 0.99;
|
||||
}
|
||||
_timing_last_count = current;
|
||||
_timing_last_time = now;
|
||||
updated = 1;
|
||||
}
|
||||
}
|
||||
|
||||
if (updated) {
|
||||
uintmax_t time_left = (total - current) / _timing_rate;
|
||||
|
||||
uintmax_t days = time_left / 60 / 60 / 24;
|
||||
if (days > 0) {
|
||||
fprintf(stderr,
|
||||
"\33[2K\r%zu/%zu (%" PRIuMAX " iter/sec | About %" PRIuMAX
|
||||
" days left.)",
|
||||
current, total, (uintmax_t)_timing_rate, days);
|
||||
} else {
|
||||
uint8_t hours = (time_left % (60 * 60 * 24)) / 60;
|
||||
uint8_t minutes = (time_left % (60 * 60)) / 60;
|
||||
uint8_t seconds = time_left % 60;
|
||||
|
||||
fprintf(stderr,
|
||||
"\33[2K\r%zu/%zu (%" PRIuMAX " iter/sec | About %" PRIu8
|
||||
"h%" PRIu8 "m%" PRIu8 "s left.)",
|
||||
current, total, (uintmax_t)_timing_rate, hours, minutes, seconds);
|
||||
}
|
||||
|
||||
fflush(stderr);
|
||||
}
|
||||
}
|
||||
|
||||
static uintmax_t Factorial(uint8_t value) {
|
||||
uintmax_t acc = 1;
|
||||
for (uint8_t i = 2; i <= value; i++) {
|
||||
acc *= (uintmax_t)value;
|
||||
}
|
||||
return acc;
|
||||
}
|
||||
|
||||
static _Bool ReadArguments(char *from, size_t *into, const char *varname) {
|
||||
if (!sscanf(from, "%zu", into)) {
|
||||
fprintf(stderr, "Failed to parse %s as a positive integer.\n", varname);
|
||||
return 0;
|
||||
}
|
||||
return 1;
|
||||
}
|
||||
|
||||
static inline void AdvanceConditionalPermutations(
|
||||
|
|
@ -40,220 +103,325 @@ typedef enum {
|
|||
kPRawFmt,
|
||||
} arg_format_t;
|
||||
|
||||
int main(int argc, char *argv[]) {
|
||||
if (argc < 5) {
|
||||
fprintf(stderr,
|
||||
"Usage:\n ./main.exe <A outputs> <B outputs> <A inputs> <B "
|
||||
"inputs> [--cg|--p|--p-raw]");
|
||||
exit(EXIT_FAILURE);
|
||||
typedef struct {
|
||||
size_t a_out;
|
||||
size_t b_out;
|
||||
size_t a_in;
|
||||
size_t b_in;
|
||||
arg_format_t formatting;
|
||||
} args_t;
|
||||
|
||||
static args_t ParseArgs(size_t argc, char *argv[]) {
|
||||
char **positionals = NULL;
|
||||
char **flags = NULL;
|
||||
|
||||
for (size_t i = 1; i < argc; i++) {
|
||||
char *arg = argv[i];
|
||||
|
||||
if (*arg == '-' && *(arg + 1) == '-') {
|
||||
stbds_arrput(flags, arg + 2);
|
||||
} else {
|
||||
stbds_arrput(positionals, arg);
|
||||
}
|
||||
}
|
||||
|
||||
size_t a_out, b_out, a_in, b_in;
|
||||
ReadArguments(argv[1], &a_out, "A outputs");
|
||||
ReadArguments(argv[2], &b_out, "B outputs");
|
||||
ReadArguments(argv[3], &a_in, "A inputs");
|
||||
ReadArguments(argv[4], &b_in, "B inputs");
|
||||
arg_format_t formatting;
|
||||
|
||||
int cg_format = kRowFmt;
|
||||
if (argc > 5) {
|
||||
if (strcmp(argv[5], "--cg") == 0) {
|
||||
cg_format = kCgFmt;
|
||||
} else if (strcmp(argv[5], "--p") == 0) {
|
||||
cg_format = kPFmt;
|
||||
} else if (strcmp(argv[5], "--p-raw") == 0) {
|
||||
cg_format = kPRawFmt;
|
||||
if (stbds_arrlen(positionals) != 4) {
|
||||
goto arg_error;
|
||||
}
|
||||
|
||||
if (!ReadArguments(positionals[0], &a_out, "A outputs") ||
|
||||
!ReadArguments(positionals[1], &b_out, "B outputs") ||
|
||||
!ReadArguments(positionals[2], &a_in, "A inputs") ||
|
||||
!ReadArguments(positionals[3], &b_in, "B inputs")) {
|
||||
goto arg_error;
|
||||
}
|
||||
|
||||
formatting = kRowFmt;
|
||||
for (size_t i = 0; i < stbds_arrlenu(flags); i++) {
|
||||
char *flag = flags[i];
|
||||
if (strcmp(flag, "cg") == 0) {
|
||||
formatting = kCgFmt;
|
||||
} else if (strcmp(flag, "p") == 0) {
|
||||
formatting = kPFmt;
|
||||
} else if (strcmp(flag, "p-raw") == 0) {
|
||||
formatting = kPRawFmt;
|
||||
}
|
||||
}
|
||||
|
||||
args_t result = {
|
||||
.a_in = a_in,
|
||||
.a_out = a_out,
|
||||
.b_in = b_in,
|
||||
.b_out = b_out,
|
||||
.formatting = formatting,
|
||||
};
|
||||
|
||||
stbds_arrfree(positionals);
|
||||
stbds_arrfree(flags);
|
||||
|
||||
return result;
|
||||
|
||||
arg_error:
|
||||
stbds_arrfree(positionals);
|
||||
stbds_arrfree(flags);
|
||||
|
||||
fprintf(stderr, "Usage:\n ");
|
||||
if (argc > 0) {
|
||||
fprintf(stderr, "%s", argv[0]);
|
||||
} else {
|
||||
fprintf(stderr, "./exe");
|
||||
}
|
||||
fprintf(stderr,
|
||||
" [--cg|--p|--p-raw] <A outputs> <B outputs> <A inputs> <B inputs>");
|
||||
exit(EXIT_FAILURE);
|
||||
}
|
||||
|
||||
int main(int argc, char *argv[]) {
|
||||
args_t args = ParseArgs(argc, argv);
|
||||
size_t a_out, a_in, b_out, b_in;
|
||||
a_out = args.a_out;
|
||||
b_out = args.b_out;
|
||||
a_in = args.a_in;
|
||||
b_in = args.b_in;
|
||||
|
||||
if (a_out < 2 || b_out < 2 || a_in < 2 || b_in < 2) {
|
||||
fprintf(stderr,
|
||||
"This program is not prepared to deal with input or output label "
|
||||
"counts smaller than 2.");
|
||||
exit(EXIT_FAILURE);
|
||||
}
|
||||
|
||||
size_t row_len = (a_out - 1) * (b_out - 1) * a_in * b_in +
|
||||
(a_out - 1) * a_in + (b_out - 1) * b_in + 1;
|
||||
matrix_t matrix = ParseMatrix(row_len);
|
||||
size_t row_count = matrix.len / matrix.row_len;
|
||||
|
||||
data_t *p_buf = malloc(a_out * b_out * a_in * b_in * sizeof(data_t));
|
||||
data_t *cg_buf = malloc((((a_out - 1) * (b_out - 1) * a_in * b_in +
|
||||
(a_out - 1) * a_in + (b_out - 1) * b_in) +
|
||||
1) *
|
||||
sizeof(data_t));
|
||||
_Bool *seen = calloc(row_count, sizeof(_Bool));
|
||||
|
||||
if (p_buf == NULL) {
|
||||
fprintf(stderr, "Failed to allocate P notation buffer. Aborting.");
|
||||
exit(EXIT_FAILURE);
|
||||
}
|
||||
if (cg_buf == NULL) {
|
||||
fprintf(stderr, "Failed to allocate CG notation buffer. Aborting.");
|
||||
exit(EXIT_FAILURE);
|
||||
}
|
||||
if (seen == NULL) {
|
||||
fprintf(stderr, "Failed to allocate equivalence flag string. Aborting.");
|
||||
exit(EXIT_FAILURE);
|
||||
}
|
||||
|
||||
permutation_generator_t a_in_perm = PermutationNewGenerator(a_in);
|
||||
permutation_generator_t b_in_perm = PermutationNewGenerator(b_in);
|
||||
permutation_generator_t *a_out_perms =
|
||||
malloc(a_in * sizeof(permutation_generator_t));
|
||||
permutation_generator_t *b_out_perms =
|
||||
malloc(b_in * sizeof(permutation_generator_t));
|
||||
|
||||
if (a_out_perms == NULL || b_out_perms == NULL) {
|
||||
fprintf(stderr,
|
||||
"Failed to allocate space for conditioned permutation generators. "
|
||||
"Aborting.");
|
||||
exit(EXIT_FAILURE);
|
||||
}
|
||||
|
||||
for (size_t i = 0; i < a_in; i++) {
|
||||
a_out_perms[i] = PermutationNewGenerator(a_out);
|
||||
PermutationReset(a_out_perms + i);
|
||||
}
|
||||
for (size_t i = 0; i < b_in; i++) {
|
||||
b_out_perms[i] = PermutationNewGenerator(b_out);
|
||||
PermutationReset(b_out_perms + i);
|
||||
}
|
||||
|
||||
// Start testing for equivalences.
|
||||
const uintmax_t perm_iterations =
|
||||
Factorial(a_in) * Factorial(b_in) *
|
||||
(uintmax_t)upow((uint64_t)(Factorial(a_in)), (uint8_t)a_out) *
|
||||
(uintmax_t)upow((uint64_t)(Factorial(b_in)), (uint8_t)b_out);
|
||||
const uintmax_t total =
|
||||
(uintmax_t)row_count * (uintmax_t)row_count * perm_iterations;
|
||||
uintmax_t tested = 0;
|
||||
#pragma omp parallel default(none) \
|
||||
shared(stderr, args, a_out, b_out, a_in, b_in, row_len, matrix, row_count, \
|
||||
seen, perm_iterations, total, tested)
|
||||
{
|
||||
data_t *p_buf = malloc(a_out * b_out * a_in * b_in * sizeof(data_t));
|
||||
data_t *cg_buf = malloc((((a_out - 1) * (b_out - 1) * a_in * b_in +
|
||||
(a_out - 1) * a_in + (b_out - 1) * b_in) +
|
||||
1) *
|
||||
sizeof(data_t));
|
||||
|
||||
for (size_t lhs_i = 0; lhs_i < row_count - 1; lhs_i++) {
|
||||
if (seen[lhs_i]) {
|
||||
continue;
|
||||
if (p_buf == NULL) {
|
||||
fprintf(stderr, "Failed to allocate P notation buffer. Aborting.");
|
||||
exit(EXIT_FAILURE);
|
||||
}
|
||||
if (cg_buf == NULL) {
|
||||
fprintf(stderr, "Failed to allocate CG notation buffer. Aborting.");
|
||||
exit(EXIT_FAILURE);
|
||||
}
|
||||
|
||||
fprintf(stderr, "%zu/%zu \033[G", lhs_i, row_count - 1);
|
||||
fflush(stderr);
|
||||
permutation_generator_t a_in_perm = PermutationNewGenerator(a_in);
|
||||
permutation_generator_t b_in_perm = PermutationNewGenerator(b_in);
|
||||
permutation_generator_t *a_out_perms =
|
||||
malloc(a_in * sizeof(permutation_generator_t));
|
||||
permutation_generator_t *b_out_perms =
|
||||
malloc(b_in * sizeof(permutation_generator_t));
|
||||
|
||||
data_t *lhs = matrix.head + lhs_i * matrix.row_len;
|
||||
if (a_out_perms == NULL || b_out_perms == NULL) {
|
||||
fprintf(
|
||||
stderr,
|
||||
"Failed to allocate space for conditioned permutation generators. "
|
||||
"Aborting.");
|
||||
exit(EXIT_FAILURE);
|
||||
}
|
||||
|
||||
for (size_t rhs_i = lhs_i + 1; rhs_i < row_count; rhs_i++) {
|
||||
if (seen[rhs_i]) {
|
||||
continue;
|
||||
}
|
||||
|
||||
data_t *rhs = matrix.head + rhs_i * matrix.row_len;
|
||||
|
||||
PermutationReset(&a_in_perm);
|
||||
while (!a_in_perm.exhausted) {
|
||||
PermutationReset(&b_in_perm);
|
||||
while (!b_in_perm.exhausted) {
|
||||
ResetConditionalPermutations(a_out_perms, a_out);
|
||||
while (!a_out_perms[a_out - 1].exhausted) {
|
||||
ResetConditionalPermutations(b_out_perms, b_out);
|
||||
while (!b_out_perms[b_out - 1].exhausted) {
|
||||
// Compare the two rows
|
||||
FromCgToP(a_out, b_out, a_in, b_in, rhs, p_buf,
|
||||
a_in_perm.permutation, b_in_perm.permutation,
|
||||
a_out_perms, b_out_perms);
|
||||
|
||||
{
|
||||
FromPToCg(a_out, b_out, a_in, b_in, cg_buf, p_buf);
|
||||
_Bool equivalent = 1;
|
||||
for (size_t i = row_len; i > 0; i--) {
|
||||
if (lhs[i - 1] != cg_buf[i - 1]) {
|
||||
equivalent = 0;
|
||||
break;
|
||||
}
|
||||
}
|
||||
|
||||
if (equivalent) {
|
||||
seen[rhs_i] = 1;
|
||||
goto skip_permutations;
|
||||
}
|
||||
}
|
||||
|
||||
// If the number of output labels are the same, and the number of
|
||||
// input labels is also the same, we can also check for equality
|
||||
// under party swapping.
|
||||
// I don't expect this conditional to be very penalizing because
|
||||
// it's very predictable.
|
||||
if (a_in == b_in && a_out == b_out) {
|
||||
// Use the results in p_buf
|
||||
PSwapParties(a_out, a_in, p_buf);
|
||||
FromPToCg(a_out, b_out, a_in, b_in, cg_buf, p_buf);
|
||||
_Bool equivalent = 1;
|
||||
for (size_t i = row_len; i > 0; i--) {
|
||||
if (lhs[i - 1] != cg_buf[i - 1]) {
|
||||
equivalent = 0;
|
||||
break;
|
||||
}
|
||||
}
|
||||
|
||||
if (equivalent) {
|
||||
seen[rhs_i] = 1;
|
||||
goto skip_permutations;
|
||||
}
|
||||
}
|
||||
|
||||
AdvanceConditionalPermutations(b_out_perms, b_out);
|
||||
}
|
||||
|
||||
AdvanceConditionalPermutations(a_out_perms, a_out);
|
||||
}
|
||||
|
||||
PermutationNext(&b_in_perm);
|
||||
}
|
||||
|
||||
PermutationNext(&a_in_perm);
|
||||
}
|
||||
|
||||
skip_permutations:;
|
||||
} // For loop over rhs_i
|
||||
} // For loop over lhs_i
|
||||
|
||||
if (cg_format == kPFmt) {
|
||||
PermutationReset(&a_in_perm);
|
||||
PermutationReset(&b_in_perm);
|
||||
for (size_t i = 0; i < a_in; i++) {
|
||||
a_out_perms[i] = PermutationNewGenerator(a_out);
|
||||
PermutationReset(a_out_perms + i);
|
||||
}
|
||||
for (size_t i = 0; i < b_in; i++) {
|
||||
b_out_perms[i] = PermutationNewGenerator(b_out);
|
||||
PermutationReset(b_out_perms + i);
|
||||
}
|
||||
}
|
||||
|
||||
// Print every unique row
|
||||
for (size_t i = 0; i < row_count; i++) {
|
||||
if (!seen[i]) {
|
||||
switch (cg_format) {
|
||||
default:
|
||||
case kRowFmt: {
|
||||
PrintMatrixRow(&matrix, i);
|
||||
} break;
|
||||
case kCgFmt: {
|
||||
printf("%zu: ", i);
|
||||
PrintCg(a_out, b_out, a_in, b_in, matrix.head + i * matrix.row_len);
|
||||
} break;
|
||||
case kPFmt: {
|
||||
printf("%zu: ", i);
|
||||
FromCgToP(a_out, b_out, a_in, b_in, matrix.head + i * matrix.row_len,
|
||||
p_buf, a_in_perm.permutation, b_in_perm.permutation,
|
||||
a_out_perms, b_out_perms);
|
||||
PrintP(a_out, b_out, a_in, b_in, p_buf);
|
||||
} break;
|
||||
case kPRawFmt: {
|
||||
FromCgToP(a_out, b_out, a_in, b_in, matrix.head + i * matrix.row_len,
|
||||
p_buf, a_in_perm.permutation, b_in_perm.permutation,
|
||||
a_out_perms, b_out_perms);
|
||||
PrintPRaw(a_out, b_out, a_in, b_in, p_buf);
|
||||
} break;
|
||||
#pragma omp for schedule(guided)
|
||||
for (size_t lhs_i = 0; lhs_i < row_count - 1; lhs_i++) {
|
||||
if (seen[lhs_i]) {
|
||||
const uintmax_t remaining = (uintmax_t)row_count * perm_iterations;
|
||||
#pragma omp atomic update
|
||||
tested += remaining;
|
||||
continue;
|
||||
}
|
||||
|
||||
data_t *lhs = matrix.head + lhs_i * matrix.row_len;
|
||||
|
||||
for (size_t rhs_i = lhs_i + 1; rhs_i < row_count; rhs_i++) {
|
||||
if (seen[lhs_i]) {
|
||||
const uintmax_t remaining =
|
||||
(uintmax_t)(row_count - lhs_i) * perm_iterations;
|
||||
#pragma omp atomic update
|
||||
tested += remaining;
|
||||
break;
|
||||
}
|
||||
|
||||
if (seen[rhs_i]) {
|
||||
#pragma omp atomic update
|
||||
tested += perm_iterations;
|
||||
continue;
|
||||
}
|
||||
|
||||
data_t *rhs = matrix.head + rhs_i * matrix.row_len;
|
||||
|
||||
PermutationReset(&a_in_perm);
|
||||
while (!a_in_perm.exhausted) {
|
||||
PermutationReset(&b_in_perm);
|
||||
while (!b_in_perm.exhausted) {
|
||||
ResetConditionalPermutations(a_out_perms, a_out);
|
||||
while (!a_out_perms[a_out - 1].exhausted) {
|
||||
ResetConditionalPermutations(b_out_perms, b_out);
|
||||
while (!b_out_perms[b_out - 1].exhausted) {
|
||||
// Compare the two rows
|
||||
FromCgToP(a_out, b_out, a_in, b_in, rhs, p_buf,
|
||||
a_in_perm.permutation, b_in_perm.permutation,
|
||||
a_out_perms, b_out_perms);
|
||||
|
||||
{
|
||||
FromPToCg(a_out, b_out, a_in, b_in, cg_buf, p_buf);
|
||||
_Bool equivalent = 1;
|
||||
for (size_t i = row_len; i > 0; i--) {
|
||||
if (lhs[i - 1] != cg_buf[i - 1]) {
|
||||
equivalent = 0;
|
||||
break;
|
||||
}
|
||||
}
|
||||
|
||||
if (equivalent) {
|
||||
seen[rhs_i] = 1;
|
||||
goto skip_permutations;
|
||||
}
|
||||
}
|
||||
|
||||
// If the number of output labels are the same, and the number
|
||||
// of input labels is also the same, we can also check for
|
||||
// equality under party swapping. I don't expect this
|
||||
// conditional to be very penalizing because it's very
|
||||
// predictable.
|
||||
if (a_in == b_in && a_out == b_out) {
|
||||
// Use the results in p_buf
|
||||
PSwapParties(a_out, a_in, p_buf);
|
||||
FromPToCg(a_out, b_out, a_in, b_in, cg_buf, p_buf);
|
||||
_Bool equivalent = 1;
|
||||
for (size_t i = row_len; i > 0; i--) {
|
||||
if (lhs[i - 1] != cg_buf[i - 1]) {
|
||||
equivalent = 0;
|
||||
break;
|
||||
}
|
||||
}
|
||||
|
||||
if (equivalent) {
|
||||
seen[rhs_i] = 1;
|
||||
goto skip_permutations;
|
||||
}
|
||||
}
|
||||
|
||||
AdvanceConditionalPermutations(b_out_perms, b_out);
|
||||
}
|
||||
|
||||
AdvanceConditionalPermutations(a_out_perms, a_out);
|
||||
}
|
||||
|
||||
PermutationNext(&b_in_perm);
|
||||
}
|
||||
|
||||
PermutationNext(&a_in_perm);
|
||||
}
|
||||
|
||||
skip_permutations:;
|
||||
#pragma omp atomic update
|
||||
tested += perm_iterations;
|
||||
{
|
||||
clock_t now = clock();
|
||||
#pragma omp task default(none) shared(tested, total) firstprivate(now)
|
||||
ReportProgress(tested, total, now);
|
||||
}
|
||||
|
||||
} // For loop over rhs_i
|
||||
} // For loop over lhs_i
|
||||
|
||||
#pragma omp single
|
||||
{
|
||||
if (args.formatting == kPFmt) {
|
||||
PermutationReset(&a_in_perm);
|
||||
PermutationReset(&b_in_perm);
|
||||
for (size_t i = 0; i < a_in; i++) {
|
||||
PermutationReset(a_out_perms + i);
|
||||
}
|
||||
for (size_t i = 0; i < b_in; i++) {
|
||||
PermutationReset(b_out_perms + i);
|
||||
}
|
||||
}
|
||||
|
||||
// Print every unique row
|
||||
for (size_t i = 0; i < row_count; i++) {
|
||||
if (!seen[i]) {
|
||||
switch (args.formatting) {
|
||||
default:
|
||||
case kRowFmt: {
|
||||
PrintMatrixRow(&matrix, i);
|
||||
} break;
|
||||
case kCgFmt: {
|
||||
printf("%zu: ", i);
|
||||
PrintCg(a_out, b_out, a_in, b_in,
|
||||
matrix.head + i * matrix.row_len);
|
||||
} break;
|
||||
case kPFmt: {
|
||||
printf("%zu: ", i);
|
||||
FromCgToP(a_out, b_out, a_in, b_in,
|
||||
matrix.head + i * matrix.row_len, p_buf,
|
||||
a_in_perm.permutation, b_in_perm.permutation,
|
||||
a_out_perms, b_out_perms);
|
||||
PrintP(a_out, b_out, a_in, b_in, p_buf);
|
||||
} break;
|
||||
case kPRawFmt: {
|
||||
FromCgToP(a_out, b_out, a_in, b_in,
|
||||
matrix.head + i * matrix.row_len, p_buf,
|
||||
a_in_perm.permutation, b_in_perm.permutation,
|
||||
a_out_perms, b_out_perms);
|
||||
PrintPRaw(a_out, b_out, a_in, b_in, p_buf);
|
||||
} break;
|
||||
}
|
||||
}
|
||||
}
|
||||
} // End of omp single region
|
||||
|
||||
PermutationFree(&a_in_perm);
|
||||
PermutationFree(&b_in_perm);
|
||||
for (size_t i = 0; i < a_in; i++) {
|
||||
PermutationFree(a_out_perms + i);
|
||||
}
|
||||
}
|
||||
free(a_out_perms);
|
||||
for (size_t i = 0; i < b_in; i++) {
|
||||
PermutationFree(b_out_perms + i);
|
||||
}
|
||||
free(b_out_perms);
|
||||
free(cg_buf);
|
||||
free(p_buf);
|
||||
} // End of omp parallel region
|
||||
|
||||
// Free all the memory so the sanitizer is happy.
|
||||
free(seen);
|
||||
free(cg_buf);
|
||||
free(p_buf);
|
||||
PermutationFree(&a_in_perm);
|
||||
PermutationFree(&b_in_perm);
|
||||
for (size_t i = 0; i < a_in; i++) {
|
||||
PermutationFree(a_out_perms + i);
|
||||
}
|
||||
free(a_out_perms);
|
||||
for (size_t i = 0; i < b_in; i++) {
|
||||
PermutationFree(b_out_perms + i);
|
||||
}
|
||||
free(b_out_perms);
|
||||
MatrixFree(&matrix);
|
||||
return EXIT_SUCCESS;
|
||||
}
|
||||
|
|
|
|||
|
|
@ -84,14 +84,14 @@ void PermutationFree(permutation_generator_t *permutation_generator) {
|
|||
}
|
||||
|
||||
void PrintPermutation(permutation_generator_t *permutation_generator) {
|
||||
printf("Permutation{");
|
||||
printf("(");
|
||||
for (size_t i = 0; i < permutation_generator->len - 1; i++) {
|
||||
printf("%zu, ", permutation_generator->permutation[i]);
|
||||
}
|
||||
if (permutation_generator->len > 0) {
|
||||
printf("%zu}\n",
|
||||
printf("%zu)",
|
||||
permutation_generator->permutation[permutation_generator->len - 1]);
|
||||
} else {
|
||||
printf("}\n");
|
||||
printf(")");
|
||||
}
|
||||
}
|
||||
Loading…
Reference in New Issue