diff --git a/DESCRIPTION b/DESCRIPTION index 915afe7..65f7de9 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -10,7 +10,7 @@ Authors@R: c( URL: https://github.com/dusadrian/IEEE BugReports: https://github.com/dusadrian/IEEE/issues Depends: R (>= 3.6.0) -Imports: LogicOpt, Matrix, gurobi +Imports: admisc, LogicOpt, lpSolve, Matrix, gurobi Description: An extensive set of functions to perform Qualitative Comparative Analysis: crisp sets ('csQCA'), temporal ('tQCA'), multi-value ('mvQCA') and fuzzy sets ('fsQCA'), using a GUI - graphical user interface. diff --git a/NAMESPACE b/NAMESPACE index 3faae86..1cc70d9 100755 --- a/NAMESPACE +++ b/NAMESPACE @@ -1,7 +1,10 @@ importFrom("gurobi", "gurobi") +importFrom("lpSolve", "lp") importFrom("Matrix", "Matrix") importFrom("utils", "tail") importFrom("LogicOpt", "logicopt") +importFrom("admisc", "recode", "export", "tryCatchWEM") +export(compare) export(gendat) export(make_infile) export(solvechart) diff --git a/R/compare.R b/R/compare.R new file mode 100644 index 0000000..afd97fd --- /dev/null +++ b/R/compare.R @@ -0,0 +1,90 @@ +compare <- function( + iter = 1:10, + ncols = 21:27, + nrows = 25:35, + file = "results.csv", + espresso = TRUE, + avg = 1 +) { + on.exit({ + suppressWarnings(sink()) + }) + + for (nc in ncols) { + for (nr in nrows) { + for (i in iter) { + dat <- gendat(i, nc, nr) + make_infile(dat) + ## $ sudo port install espresso + # lo <- system("espresso -Dexact -o eqntott infile.esp", intern = TRUE)[1] + + if (espresso) { + etc <- admisc::tryCatchWEM( + es_time <- system.time( + elo <- LogicOpt::logicopt(esp_file = "infile.esp") + ) + ) + } + + ctc <- admisc::tryCatchWEM({ + cc_time <- NULL + for (ii in seq(avg)) { + cc_time <- c( + cc_time, + system.time( + cc <- .Call("CCubes", dat, PACKAGE = "IEEE") + )[3] + ) + } + }) + + if (is.null(ctc)) { + colnames(cc) <- colnames(dat)[-nc] + if (espresso) { + soles <- apply( + elo[[1]][, -nc], 1, + function(x) { + x <- x[x != "-"] + names(x)[x == 0] <- paste("!", names(x)[x == 0], sep = "") + return(c(paste(names(x), collapse = "&"), length(x))) + } + ) + } + + solcc <- apply(cc, 1, function(x) { + x <- x[x > 0] - 1 + names(x)[x == 0] <- paste("!", names(x)[x == 0], sep = "") + return(c(paste(names(x), collapse = "&"), length(x))) + }) + + sink(file, append = TRUE) + + cat( + sprintf( + "%s, %s, %s, %s, %s, %s%s, %s%s\n", + i, nc, nr, + round(es_time[3], 3), + round(min(cc_time), 3), + ifelse( + espresso, + paste0("\"(", paste(soles[1, ], collapse = ") + ("), ")\", "), + "" + ), + paste0("\"(", paste(solcc[1, ], collapse = ") + ("), ")\""), + ifelse( + espresso, + paste(paste(soles[2, ], collapse = ""), ", ", sep = ""), + "" + ), + paste(solcc[2, ], collapse = "") + ) + ) + + sink() + } else { + cat(sprintf("i: %s, nc: %s, nr: %s\n", i, nc, nr)) + } + } + } + } +} diff --git a/R/gendat.R b/R/gendat.R index 033d8b0..c9fa429 100644 --- a/R/gendat.R +++ b/R/gendat.R @@ -1,6 +1,6 @@ -gendat <- function(i, ncols, nrows) { - set.seed(i * nrows * ncols) +gendat <- function(iter, ncols, nrows) { + set.seed(iter * nrows * ncols) # all configurations BUT the outcome have to be unique dat <- unique( diff --git a/R/solvechart.R b/R/solvechart.R index 725df39..70c7631 100644 --- a/R/solvechart.R +++ b/R/solvechart.R @@ -48,7 +48,6 @@ value <- value[seq(PIs)] } - # cat("Saving the PI chart to a file\n") # admisc::export( # as.data.frame(t(x)), # file = paste("pic", sample(5000, 1), ".csv", sep = ""), diff --git a/man/compare.Rd b/man/compare.Rd new file mode 100644 index 0000000..88bd678 --- /dev/null +++ b/man/compare.Rd @@ -0,0 +1,41 @@ +\name{compare} + +\alias{compare} + +\title{Compare CCubes results against Espresso} + +\description{ +Compare CCubes results against Espresso +} + +\usage{ +compare( + iter = 1:10, + ncols = 21:27, + nrows = 25:35, + file = "results.csv", + espresso = TRUE, + avg = FALSE +) +} + +\arguments{ + \item{iter}{Iterator} + \item{ncols}{Number of columns} + \item{nrows}{Number of rows} + \item{file}{Character, name / path to the file where results are written} + \item{espresso}{Boolean, also write Espresso results} + \item{avg}{Numeric, number of runs to average results over} +} + +\author{ +Adrian Dusa +} + +\examples{ +\dontrun{ + compare(1, 7, 15) +} +} + +\keyword{functions} diff --git a/man/gendat.Rd b/man/gendat.Rd new file mode 100644 index 0000000..bd446d9 --- /dev/null +++ b/man/gendat.Rd @@ -0,0 +1,33 @@ +\name{gendat} + +\alias{gendat} + +\title{Generate a binary matrix} + +\description{ +Generate a random binary matrix for testing purposes. +} + +\usage{ +gendat(iter, ncols, nrows) +} + +\arguments{ + \item{iter}{Iterator} + \item{ncols}{Number of columns} + \item{nrows}{Number of rows} +} + +\value{ +A matrix containing binary values. +} + +\author{ +Adrian Dusa +} + +\examples{ +gendat(1, 7, 15) +} + +\keyword{functions} diff --git a/man/make_infile.Rd b/man/make_infile.Rd new file mode 100644 index 0000000..dfbc39f --- /dev/null +++ b/man/make_infile.Rd @@ -0,0 +1,39 @@ +\name{make_infile} + +\alias{make_infile} + +\title{Create an infile for espresso} + +\description{ +Create a text file containing the data for the expresso minimizer. +} + +\usage{ +make_infile(dat, espname = "infile.esp", copy = FALSE, ...) +} + +\arguments{ + \item{dat}{A binary matrix} + \item{espname}{Name of the text file} + \item{copy}{Logical, copy the file with a different name.} + \item{...}{Other arguments.} +} + +\details{ +Different iterations generate different versions of the data. +They are all temporarily saved into \dQuote{infile.esp} unless +the argument \code{copy} is activated, which copies the infile +into a file with a name containing the iterator, the number of +columns and the number of rows. +} + +\author{ +Adrian Dusa +} + +\examples{ +dat <- gendat(1, 7, 15) +make_infile(dat) +} + +\keyword{functions} diff --git a/man/solvechart.Rd b/man/solvechart.Rd new file mode 100644 index 0000000..fae9533 --- /dev/null +++ b/man/solvechart.Rd @@ -0,0 +1,43 @@ +\name{solvechart} + +\alias{solvechart} + +\title{Find the minimum for a prime implicants chart} + +\description{ +Solve a PI chart using a minimum number of prime implicants. +} + +\usage{ +solvechart(x, ...) +} + +\arguments{ + \item{x}{A binary matrix.} + \item{...}{Other arguments.} +} + +\details{ +A PI chart, in this package, is a binary matrix containing the prime implicants on the columns +and the observed positive minterms on the rows. +} + +\value{ +A vector of binary values +} + +\author{ +Adrian Dusa +} + +\examples{ +set.seed(12345) +x <- matrix( + sample(0:1, prob = c(0.65, 0.35), size = 35, replace = TRUE), + nrow = 5 +) + +solvechart(x) +} + +\keyword{functions} diff --git a/src/CCubes.c b/src/CCubes.c index d2996ec..9fd17ed 100644 --- a/src/CCubes.c +++ b/src/CCubes.c @@ -1,6 +1,6 @@ +#include #include #include -#include #include #include #include @@ -24,21 +24,26 @@ #include "config.h" #include "logging.h" +// $ export BITS_PER_WORD=32 in the Terminal, before running R SEXP CCubes(SEXP tt) { - // $ export BITS_PER_WORD=32 in the Terminal, before running R - // 32 bits per word, in bit shifting representation - char *bits_per_word = getenv("BITS_PER_WORD"); // Read from the PATH + char *bits_per_word = getenv("BITS_PER_WORD"); // Read from the PATH int BITS_PER_WORD = bits_per_word ? atoi(bits_per_word) : 32; if (BITS_PER_WORD != 8 && BITS_PER_WORD != 16 && BITS_PER_WORD != 32 && BITS_PER_WORD != 64) { - BITS_PER_WORD = 32; // Default to 32 + BITS_PER_WORD = 32; // Default to 32 } - // $ export PRINT_INFO=1 in the Terminal, before running R - char *print_info = getenv("PRINT_INFO"); // Read from the PATH - Rboolean PRINT_INFO = print_info && print_info[0] == '1'; + // $ export PRINT_INFO=true in the Terminal, before running R + char *print_info = getenv("PRINT_INFO"); // Read from the PATH + Rboolean PRINT_INFO = print_info && ( + strcmp(print_info, "1") == 0 || + strcmp(print_info, "TRUE") == 0 || + strcmp(print_info, "true") == 0 || + strcmp(print_info, "T") == 0 || + strcmp(print_info, "t") == 0 + ); int multiplier = 0; struct timeval start, end; @@ -65,7 +70,7 @@ SEXP CCubes(SEXP tt) { } int *p_tt = INTEGER(tt); - int ttrows = nrows(tt); // number of rows in the data matrix + int ttrows = nrows(tt); // number of rows in the data matrix int ninputs = ncols(tt) - 1; // number of inputs (columns - 1, the last one is the outcome) // calculate the number of positive output rows (the ON set) @@ -112,19 +117,25 @@ SEXP CCubes(SEXP tt) { } } - int value_bit_width = 0; - while (max_value > 0) { - max_value >>= 1; // Shift right until no bits remain - value_bit_width++; + int bits_needed = ceil(log2(max_value)); // Compute the necessary bits + int value_bit_width = 1; + while (value_bit_width < bits_needed) { + value_bit_width *= 2; // Round up to the next power of 2 } + if (value_bit_width > BITS_PER_WORD) { + BITS_PER_WORD = value_bit_width; // Adjust the bits per word + } + + int value_bit_mask = (1U << value_bit_width) - 1U; + // calculate the number of values (biggest number) for each input int nofvalues[ninputs]; int nofpi[ninputs]; for (int i = 0; i < ninputs; i++) { nofvalues[i] = 0; // initialize - nofpi[i] = 0; // initialize + nofpi[i] = 0; // initialize for (int r = 0; r < ttrows; r++) { if (nofvalues[i] < p_tt[i * ttrows + r]) { @@ -146,10 +157,10 @@ SEXP CCubes(SEXP tt) { int estimPI = 250000; // the index of the PIs, in descending order of their number of covered ON-set minterms - int *p_covered = R_Calloc(estimPI, int); + int *p_covered = (int *) R_Calloc(estimPI, int); // many PIs will have the same coverage, but they don't necessarily cover the same minterms - // to employ row dominance when solving the PI chart, we need to compare the coverage of the + // to employ row dominance when solving the coverage matrix, we need to compare the coverage of the // current PI with the coverage of the previous PIs. If this PI survives the comparison, its // coverage has to be added in the p_covered vector, and its order in the p_covered // vector, at the last index of the PI coverage with the same number of minterms @@ -158,9 +169,9 @@ SEXP CCubes(SEXP tt) { // p_pichart = malloc(estimPI * posrows * sizeof(int)); // memset(p_pichart, false, estimPI * posrows * sizeof(int)); int *p_pichart = (int *) R_Calloc(estimPI * posrows, int); - // prefixing (int *) prefills in all values with 0s + // prefixing (int *) before R_Calloc() prefills all values with 0s - int pichart_words = (posrows + BITS_PER_WORD - 1) / BITS_PER_WORD; // Words needed per PI chart columns + int pichart_words = (posrows + BITS_PER_WORD - 1) / BITS_PER_WORD; // Words needed per coverage matrix columns unsigned int *p_pichart_pos = (unsigned int *) R_Calloc(estimPI * pichart_words, unsigned int); int implicant_words = (ninputs + BITS_PER_WORD - 1) / BITS_PER_WORD; // Words needed per PI representation unsigned int *p_implicants_pos = (unsigned int *) R_Calloc(estimPI * implicant_words, unsigned int); @@ -168,15 +179,14 @@ SEXP CCubes(SEXP tt) { int prevfoundPI = 0; // the number of previously found PIs int foundPI = 0; - int prevsolmin = 0; // the minimum number of PIs that solve the PI chart + int prevsolmin = 0; // the previous (level k - 1), minimum number of PIs to solve the coverage matrix int solmin = 0; - // the positions of the PIs solving the PI chart - // a vector which can never be lengthier than the number of ON minterms (posrows) + // the positions of the PIs solving the coverage matrix + // this vector can never be lengthier than the number of ON minterms (posrows) int previndices[posrows]; int indices[posrows]; - for (int i = 0; i < posrows; i++) { previndices[i] = 0; indices[i] = 0; @@ -200,6 +210,11 @@ SEXP CCubes(SEXP tt) { } int n_tasks = nchoosek(ninputs, k); + if (n_tasks == 0) { + // overflow, too many tasks + return (R_NilValue); + } + int n_tasks_batch = 512; for (int task = 0; task < n_tasks; task+=n_tasks_batch) { /* adjust if batch size is larger than total job size */ @@ -277,22 +292,22 @@ SEXP CCubes(SEXP tt) { #pragma omp parallel for schedule(dynamic) #endif - for (int task = 0; task < nchoosek(ninputs, k); task++) { + for (int task = 0; task < n_tasks; task++) { int tempk[k]; int x = 0; - int start_point = task; + int combination = task; - // fill the combination for the current task + // fill the combination for the current task / combination number for (int i = 0; i < k; i++) { - while (nchoosek(ninputs - (x + 1), k - (i + 1)) <= start_point) { - start_point -= nchoosek(ninputs - (x + 1), k - (i + 1)); + while (nchoosek(ninputs - (x + 1), k - (i + 1)) <= combination) { + combination -= nchoosek(ninputs - (x + 1), k - (i + 1)); x++; } tempk[i] = x; x++; } - // allocate vectors of decimal row numbers for the positive and negative rows + // allocate vectors of decimal numbers for the ON-set and OFF-set rows int decpos[posrows]; int decneg[negrows]; @@ -326,7 +341,7 @@ SEXP CCubes(SEXP tt) { int possible_rows[posrows]; Rboolean possible_cover[posrows]; - possible_cover[0] = true; // Rboolean flag, to be set with false if found among the OFF set + possible_cover[0] = true; // boolean flag, to be set with false if found among the OFF set int found = 0; @@ -370,7 +385,6 @@ SEXP CCubes(SEXP tt) { for (int f = 0; f < found; f++) { - // create a temporary vector of length k, containing the values from the initial ON set // plus 1 (because 0 now signals a minimization, it becomes 1, and 1 becomes 2 etc. int tempc[k]; @@ -388,11 +402,11 @@ SEXP CCubes(SEXP tt) { int value = ON_set[tempk[c] * posrows + frows[f]]; tempc[c] = value + 1; - int word_index = tempk[c] / BITS_PER_WORD; - int bit_index = tempk[c] % BITS_PER_WORD; + int word_index = tempk[c] / (BITS_PER_WORD / value_bit_width); + int bit_index = (tempk[c] % (BITS_PER_WORD / value_bit_width)) * value_bit_width; - fixed_bits[word_index] |= 1U << bit_index; - value_bits[word_index] |= (unsigned int)value << (bit_index * value_bit_width); + fixed_bits[word_index] |= (value_bit_mask << bit_index); + value_bits[word_index] |= ((unsigned int)value << bit_index); } // check if the current PI is not redundant @@ -414,13 +428,15 @@ SEXP CCubes(SEXP tt) { for (int w = 0; w < implicant_words; w++) { // If the new PI has values on positions outside the existing PI’s fixed positions, it’s not a subset - if ((fixed_bits[w] & p_implicants_pos[i * implicant_words + w]) != p_implicants_pos[i * implicant_words + w]) { + unsigned int index_mask = p_implicants_pos[i * implicant_words + w]; + + if ((fixed_bits[w] & index_mask) != index_mask) { is_subset = false; break; } // then compare the value bits, if one or more values on those positions are different, it’s not a subset - if ((value_bits[w] & p_implicants_val[i * implicant_words + w]) != p_implicants_val[i * implicant_words + w]) { + if ((value_bits[w] & index_mask) != (p_implicants_val[i * implicant_words + w] & index_mask)) { is_subset = false; break; } @@ -479,7 +495,6 @@ SEXP CCubes(SEXP tt) { #pragma omp critical #endif { - // push the PI information to the global arrays for (int i = foundPI; i > last_index[covsum - 1]; i--) { @@ -497,7 +512,7 @@ SEXP CCubes(SEXP tt) { p_implicants_val[implicant_words * foundPI + w] = value_bits[w]; } - // populate the PI chart + // populate the coverage matrix for (int r = 0; r < posrows; r++) { for (int w = 0; w < pichart_words; w++) { p_pichart_pos[foundPI * pichart_words + w] = pichart_values[w]; @@ -510,21 +525,20 @@ SEXP CCubes(SEXP tt) { // when needed, increase allocated memory if (foundPI / estimPI > 0.9) { - int old_size = estimPI; estimPI += 100000; - p_pichart = R_Realloc(p_pichart, posrows * estimPI, int); - p_pichart_pos = R_Realloc(p_pichart_pos, estimPI, unsigned int); - p_implicants_val = R_Realloc(p_implicants_val, ninputs * estimPI, unsigned int); - p_implicants_pos = R_Realloc(p_implicants_pos, ninputs * estimPI, unsigned int); - p_covered = R_Realloc(p_covered, estimPI, int); + p_pichart = R_Realloc(p_pichart, posrows * estimPI, int); + p_pichart_pos = R_Realloc(p_pichart_pos, pichart_words * estimPI, unsigned int); + p_implicants_val = R_Realloc(p_implicants_val, implicant_words * estimPI, unsigned int); + p_implicants_pos = R_Realloc(p_implicants_pos, implicant_words * estimPI, unsigned int); + p_covered = R_Realloc(p_covered, estimPI, int); - for (unsigned int i = old_size; i < posrows * estimPI; i++) { + for (unsigned int i = foundPI; i < posrows * estimPI; i++) { p_pichart[i] = 0; } - for (unsigned int i = old_size; i < estimPI; i++) { + for (unsigned int i = foundPI; i < pichart_words * estimPI; i++) { p_pichart_pos[i] = 0U; } - for (unsigned int i = old_size; i < ninputs * estimPI; i++) { + for (unsigned int i = foundPI; i < implicant_words * estimPI; i++) { p_implicants_val[i] = 0U; p_implicants_pos[i] = 0U; } @@ -565,7 +579,6 @@ SEXP CCubes(SEXP tt) { if (ON_set_covered) { // Rprintf("posrows: %d; foundPI: %d\n", posrows, foundPI); - if (PRINT_INFO) { gettimeofday(&end, NULL); // End time elapsed_time = (end.tv_sec - start.tv_sec) + (end.tv_usec - start.tv_usec) / 1e6; @@ -573,8 +586,6 @@ SEXP CCubes(SEXP tt) { gettimeofday(&start, NULL); } - - SEXP pic = PROTECT(allocMatrix(INTSXP, posrows, foundPI)); for (long long unsigned int i = 0; i < posrows * foundPI; i++) { INTEGER(pic)[i] = p_pichart[i]; @@ -599,13 +610,12 @@ SEXP CCubes(SEXP tt) { } UNPROTECT(5); + // Rprintf("solution minima: %d\n", solmin); - - if (PRINT_INFO) { gettimeofday(&end, NULL); elapsed_time = (end.tv_sec - start.tv_sec) + (end.tv_usec - start.tv_usec) / 1e6; - Rprintf("Time spent solving the PI chart: %f sec.\n", elapsed_time); + Rprintf("Time spent solving the coverage matrix: %f sec.\n", elapsed_time); gettimeofday(&start, NULL); } @@ -643,7 +653,7 @@ SEXP CCubes(SEXP tt) { } } - // printf("solmin: %d\n", solmin); + // Rprintf("solmin: %d\n", solmin); if (PRINT_INFO) { Rprintf("--- END ---\n"); } @@ -654,11 +664,13 @@ SEXP CCubes(SEXP tt) { for (int c = 0; c < solmin; c++) { for (int r = 0; r < ninputs; r++) { unsigned int value = 0; - int word_index = r / BITS_PER_WORD; // Word index within the implicant - int bit_index = r % BITS_PER_WORD; // Bit position within the word + int word_index = r / (BITS_PER_WORD / value_bit_width); // Word index within the implicant + int bit_index = (r % (BITS_PER_WORD / value_bit_width)) * value_bit_width; // Bit position within the word - if (p_implicants_pos[indices[c] * implicant_words + word_index] & (1U << bit_index)) { - value = 1U + ((p_implicants_val[indices[c] * implicant_words + word_index] >> (bit_index * value_bit_width)) & ((1U << value_bit_width) - 1U)); + // (1U << bit_index) is just as good + if (p_implicants_pos[indices[c] * implicant_words + word_index] & (value_bit_mask << bit_index)) { + value = (p_implicants_val[indices[c] * implicant_words + word_index] >> bit_index) & value_bit_mask; + value++; // 0 indicates a minimization, so we increment the value } p_sol[r * solmin + c] = value; // transposed @@ -677,13 +689,12 @@ SEXP CCubes(SEXP tt) { return (sol); } - long long unsigned int nchoosek( int n, int k ) { + if (k > n) return 0; if (k == 0 || k == n) return 1; - if (k == 1) return n; long long unsigned int result = 1; @@ -692,7 +703,21 @@ long long unsigned int nchoosek( } for (int i = 0; i < k; i++) { - result = result * (n - i) / (i + 1); + // result = result * (n - i) / (i + 1); + + // Check for potential overflow before multiplication + if (result > ULLONG_MAX / (n - i)) { + return 0; // Indicate overflow + } + + result *= (n - i); + + // Check for potential overflow before division + if (result % (i + 1) != 0) { + return 0; // Indicate overflow + } + + result /= (i + 1); } return result; diff --git a/src/ccubes.cl b/src/ccubes.cl index 8d1bf85..1aeb52b 100644 --- a/src/ccubes.cl +++ b/src/ccubes.cl @@ -147,12 +147,12 @@ ccubes_task(int k, int tempk[NINPUTS]; /* max is tempk[ninputs] */ int x = 0; - int start_point = task; + int combination = task; // fill the combination for the current task for (int i = 0; i < k; i++) { - while (nchoosek(NINPUTS - (x + 1), k - (i + 1)) <= start_point) { - start_point -= nchoosek(NINPUTS - (x + 1), k - (i + 1)); + while (nchoosek(NINPUTS - (x + 1), k - (i + 1)) <= combination) { + combination -= nchoosek(NINPUTS - (x + 1), k - (i + 1)); x++; } tempk[i] = x; @@ -255,11 +255,11 @@ ccubes_task(int k, int value = ON_set[tempk[c] * POSROWS + frows[f]]; tempc[c] = value + 1; - int word_index = tempk[c] / BITS_PER_WORD; - int bit_index = tempk[c] % BITS_PER_WORD; + int word_index = tempk[c] / (BITS_PER_WORD / value_bit_width); + int bit_index = (tempk[c] % (BITS_PER_WORD / value_bit_width)) * value_bit_width; - fixed_bits[word_index] |= 1U << bit_index; - value_bits[word_index] |= (unsigned int)value << (bit_index * VALUE_BIT_WIDTH); + fixed_bits[word_index] |= (value_bit_mask << bit_index); + value_bits[word_index] |= ((unsigned int)value << bit_index); } // check if the current PI is not redundant @@ -281,13 +281,15 @@ ccubes_task(int k, for (int w = 0; w < IMPLICANT_WORDS; w++) { // If the new PI has values on positions outside the existing PI’s fixed positions, it’s not a subset - if ((fixed_bits[w] & p_implicants_pos[i * IMPLICANT_WORDS + w]) != p_implicants_pos[i * IMPLICANT_WORDS + w]) { + unsigned int index_mask = p_implicants_pos[i * implicant_words + w]; + + if ((fixed_bits[w] & index_mask) != index_mask) { is_subset = false; break; } // then compare the value bits, if one or more values on those positions are different, it’s not a subset - if ((value_bits[w] & p_implicants_val[i * IMPLICANT_WORDS + w]) != p_implicants_val[i * IMPLICANT_WORDS + w]) { + if ((value_bits[w] & index_mask) != (p_implicants_val[i * implicant_words + w] & index_mask)) { is_subset = false; break; } diff --git a/src/row_dominance.c b/src/row_dominance.c deleted file mode 100644 index 0b2c93a..0000000 --- a/src/row_dominance.c +++ /dev/null @@ -1,100 +0,0 @@ - -#include -#include "row_dominance.h" - - -void row_dominance( - int p_pichart[], - int p_implicants[], - int *p_ck, - int pirows, - unsigned int *foundPI, - int ninputs -) { - - unsigned int picols = *foundPI; - - - // int* survcols = (int *) R_Calloc (picols, int); - - // for (int i = 0; i < picols; i++) { - // survcols[i] = true; - // } - - Rboolean survcols[picols]; - int colsums[picols]; - int sortcol[picols]; - int temp; - - for (unsigned int c = 0; c < picols; c++) { - colsums[c] = 0; - - for (int r = 0; r < pirows; r++) { - colsums[c] += p_pichart[c * pirows + r]; - } - - sortcol[c] = c; - survcols[c] = true; - } - - for (unsigned int c1 = 0; c1 < picols; c1++) { - for (unsigned int c2 = c1 + 1; c2 < picols; c2++) { - if (colsums[sortcol[c1]] < colsums[sortcol[c2]]) { - temp = sortcol[c1]; - sortcol[c1] = sortcol[c2]; - sortcol[c2] = temp; - } - } - } - - for (unsigned int c1 = 0; c1 < picols; c1++) { - if (survcols[sortcol[c1]]) { - for (unsigned int c2 = c1 + 1; c2 < picols; c2++) { - if (survcols[sortcol[c2]]) { - if (colsums[sortcol[c1]] > colsums[sortcol[c2]]) { - - Rboolean itcovers = true; // assume it's covered - int r = 0; - - while (r < pirows && itcovers) { - if (p_pichart[sortcol[c2] * pirows + r]) { - itcovers = p_pichart[sortcol[c1] * pirows + r]; - } - r++; - } - - if (itcovers) { - survcols[sortcol[c2]] = false; - --(*foundPI); - } - } - } - } - } - } - - - if (*foundPI < picols) { - - // move (overwrite) all surviving columns towards the beginning - int s = 0; - for (unsigned int c = 0; c < picols; c++) { - if (survcols[c]) { - for (int r = 0; r < pirows; r++) { - p_pichart[s * pirows + r] = p_pichart[c * pirows + r]; - } - - for (int r = 0; r < ninputs; r++) { - p_implicants[s * ninputs + r] = p_implicants[c * ninputs + r]; - } - - // same with the vector positions - p_ck[s] = p_ck[c]; - - s++; - } - } - } - - return; -} diff --git a/src/row_dominance.h b/src/row_dominance.h deleted file mode 100755 index 0ac39fb..0000000 --- a/src/row_dominance.h +++ /dev/null @@ -1,10 +0,0 @@ -#include - -void row_dominance( - int p_pichart[], - int p_implicants[], - int *p_ck, - int pirows, - unsigned int *foundPI, - int ninputs -);