# High Performance Benchmarks

#### Joseph Wood

#### 2/8/2020

This document serves as an overview for measuring the performance of `RcppAlgos` against other tools that were determined to be efficient for generating combinations and permutations from the benchmarks found here: Permutations and combinations with/without replacement and for distinct/non-distinct items/multiset. You will note that the examples in that post were relatively small. The benchmarks below will focus on larger examples where performance really matters and for this reason we only consider the packages arrangements, partitions, and RcppAlgos.

## Setup Information

For the benchmarks below, we used a `Macbook Pro i7 16Gb` machine. We also tested on a Windows and Linux machine with similar specs and obtained similar results.

    library(RcppAlgos)
    library(partitions)
    library(arrangements)
    library(microbenchmark)
    pertinent_output <- capture.output(sessionInfo())

    options(digits = 4)
    cat(paste(pertinent_output[1:3], collapse = "\n"))
    #> R version 3.6.2 (2019-12-12)
    #> Platform: x86_64-apple-darwin15.6.0 (64-bit)
    #> Running under: macOS Mojave 10.14.6

    cat(pertinent_output[which(pertinent_output == "other attached packages:") + 1L])
    #> [1] microbenchmark_1.4-7 partitions_1.9-22    arrangements_1.1.8   RcppAlgos_2.3.6

    numThreads <- as.integer(RcppAlgos::stdThreadMax() / 2)
    print(numThreads)
    #> [1] 4

## Combinations

### Combinations - Distinct

    set.seed(13)
    v1 <- sort(sample(100, 30))
    m <- 21
    t1 <- comboGeneral(v1, m, Parallel = T)
    t2 <- combinations(v1, m)
    identical(t1, t2)
    #> [1] TRUE
    dim(t1)
    #> [1] 14307150       21
    rm(t1, t2)
    invisible(gc())
    microbenchmark(cbRcppAlgosPar = comboGeneral(v1, m, nThreads = numThreads),
                   cbRcppAlgosSer = comboGeneral(v1, m),
                   cbArrangements = combinations(v1, m),
                   times = 15, unit = "relative")
    #> Unit: relative
    #>            expr   min    lq  mean median    uq   max neval
    #>  cbRcppAlgosPar 1.000 1.000 1.000  1.000 1.000 1.000    15
    #>  cbRcppAlgosSer 2.403 2.382 2.373  2.389 2.374 2.243    15
    #>  cbArrangements 2.095 2.365 2.343  2.380 2.365 2.237    15

### Combinations - Repetition

    v2 <- v1[1:10]
    m <- 20
    t1 <- comboGeneral(v2, m, repetition = TRUE, nThreads = numThreads)
    t2 <- combinations(v2, m, replace = TRUE)
    identical(t1, t2)
    #> [1] TRUE
    dim(t1)
    #> [1] 10015005       20
    rm(t1, t2)
    invisible(gc())
    microbenchmark(cbRcppAlgosPar = comboGeneral(v2, m, TRUE, nThreads = numThreads),
                   cbRcppAlgosSer = comboGeneral(v2, m, TRUE),
                   cbArrangements = combinations(v2, m, replace = TRUE),
                   times = 15, unit = "relative")
    #> Unit: relative
    #>            expr   min    lq  mean median    uq   max neval
    #>  cbRcppAlgosPar 1.000 1.000 1.000  1.000 1.000 1.000    15
    #>  cbRcppAlgosSer 2.567 2.566 2.562  2.560 2.556 2.560    15
    #>  cbArrangements 2.262 2.562 2.542  2.565 2.552 2.592    15

### Combinations - Multisets

    myFreqs <- c(2, 4, 4, 5, 3, 2, 2, 2, 3, 4, 1, 4, 2, 5)
    v3 <- as.integer(c(1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610))
    t1 <- comboGeneral(v3, 20, freqs = myFreqs, nThreads = numThreads)
    t2 <- combinations(freq = myFreqs, k = 20, x = v3)
    identical(t1, t2)
    #> [1] TRUE
    dim(t1)
    #> [1] 14594082       20
    rm(t1, t2)
    invisible(gc())
    microbenchmark(cbRcppAlgosPar = comboGeneral(v3, 20, freqs = myFreqs, nThreads = numThreads),
                   cbRcppAlgosSer = comboGeneral(v3, 20, freqs = myFreqs),
                   cbArrangements = combinations(freq = myFreqs, k = 20, x = v3),
                   times = 10, unit = "relative")
    #> Unit: relative
    #>            expr   min    lq  mean median    uq   max neval
    #>  cbRcppAlgosPar 1.000 1.000 1.000  1.000 1.000 1.000    10
    #>  cbRcppAlgosSer 3.640 2.737 2.792  2.719 2.719 2.720    10
    #>  cbArrangements 6.461 4.862 4.936  4.817 4.803 4.734    10

## Permutations

### Permutations - Distinct

    v4 <- as.integer(c(2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59))
    t1 <- permuteGeneral(v4, 6, nThreads = numThreads)
    t2 <- permutations(v4, 6)
    identical(t1, t2)
    #> [1] TRUE
    dim(t1)
    #> [1] 8910720       6
    rm(t1, t2)
    invisible(gc())
    microbenchmark(cbRcppAlgosPar = permuteGeneral(v4, 6, nThreads = numThreads),
                   cbRcppAlgosSer = permuteGeneral(v4, 6),
                   cbArrangements = permutations(v4, 6),
                   times = 15, unit = "relative")
    #> Unit: relative
    #>            expr   min    lq  mean median    uq   max neval
    #>  cbRcppAlgosPar 1.000 1.000 1.000  1.000 1.000 1.000    15
    #>  cbRcppAlgosSer 1.828 1.837 1.851  1.863 1.875 2.049    15
    #>  cbArrangements 2.239 2.211 2.115  2.194 2.076 1.888    15

    ## Indexing permutation example with the partitions package
    t1 <- permuteGeneral(11, nThreads = 4)
    t2 <- permutations(11)
    t3 <- perms(11)

    dim(t1)
    #> [1] 39916800       11

    identical(t1, t2)
    #> [1] TRUE
    identical(t1, t(as.matrix(t3)))
    #> [1] TRUE
    rm(t1, t2, t3)
    invisible(gc())

    microbenchmark(cbRcppAlgosPar = permuteGeneral(11, nThreads = 4),
                   cbRcppAlgosSer = permuteGeneral(11),
                   cbArrangements = permutations(11),
                   cbPartitions = perms(11),
                   times = 5, unit = "relative")
    #> Unit: relative
    #>            expr      min       lq     mean   median       uq      max neval
    #>  cbRcppAlgosPar 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000     5
    #>  cbRcppAlgosSer 1.905055 1.904145 1.856310 1.932584 1.943967 1.689898     5
    #>  cbArrangements 2.085293 2.087785 2.118516 2.418879 2.417984 1.769080     5
    #>    cbPartitions 2.572582 2.562715 2.440539 2.737129 2.772614 1.874632     5

### Permutations - Repetition

    v5 <- v3[1:5]
    t1 <- permuteGeneral(v5, 10, repetition = TRUE, nThreads = numThreads)
    t2 <- permutations(v5, 10, replace = TRUE)
    identical(t1, t2)
    #> [1] TRUE
    dim(t1)
    #> [1] 9765625      10
    rm(t1, t2)
    invisible(gc())
    microbenchmark(cbRcppAlgosPar = permuteGeneral(v5, 10, TRUE, nThreads = numThreads),
                   cbRcppAlgosSer = permuteGeneral(v5, 10, TRUE),
                   cbArrangements = permutations(x = v5, k = 10, replace = TRUE),
                   times = 10, unit = "relative")
    #> Unit: relative
    #>            expr   min    lq  mean median    uq   max neval
    #>  cbRcppAlgosPar 1.000 1.000 1.000  1.000 1.000 1.000    10
    #>  cbRcppAlgosSer 1.639 1.590 1.509  1.582 1.380 1.402    10
    #>  cbArrangements 1.530 1.492 1.303  1.477 1.033 1.258    10

### Permutations - Multisets

    v6 <- sort(runif(12))
    t1 <- permuteGeneral(v6, 7, freqs = rep(1:3, 4), nThreads = numThreads)
    t2 <- permutations(freq = rep(1:3, 4), k = 7, x = v6)
    identical(t1, t2)
    #> [1] TRUE
    dim(t1)
    #> [1] 19520760        7
    rm(t1, t2)
    invisible(gc())
    microbenchmark(cbRcppAlgosPar = permuteGeneral(v6, 7, freqs = rep(1:3, 4), nThreads = numThreads),
                   cbRcppAlgosSer = permuteGeneral(v6, 7, freqs = rep(1:3, 4)),
                   cbArrangements = permutations(freq = rep(1:3, 4), k = 7, x = v6),
                   times = 10, unit = "relative")
    #> Unit: relative
    #>            expr   min    lq  mean median    uq   max neval
    #>  cbRcppAlgosPar 1.000 1.000 1.000  1.000 1.000 1.000    10
    #>  cbRcppAlgosSer 3.345 2.392 2.461  2.389 2.393 2.400    10
    #>  cbArrangements 3.455 2.482 2.544  2.477 2.465 2.459    10

## Partitions

### Partitions - Distinct

#### All Distinct Partitions

    t1 <- comboGeneral(0:140, freqs=c(140, rep(1, 140)),
                       constraintFun = "sum", comparisonFun = "==",
                       limitConstraints = 140)
    t2 <- partitions(140, distinct = TRUE)
    t3 <- diffparts(140)

    # Each package has different output formats... we only examine dimensions
    # and that each result is a partition of 140
    identical(dim(t1), dim(t2))
    #> [1] TRUE
    identical(dim(t1), dim(t(t3)))
    #> [1] TRUE
    all(rowSums(t1) == 140)
    #> [1] TRUE
    all(rowSums(t2) == 140)
    #> [1] TRUE
    all(colSums(t3) == 140)
    #> [1] TRUE

    dim(t1)
    #> [1] 9617150      16
    rm(t1, t2, t3)
    invisible(gc())
    microbenchmark(cbRcppAlgosSer = comboGeneral(0:140, freqs=c(140, rep(1, 140)),
                                                 constraintFun = "sum",
                                                 comparisonFun = "==",
                                                 limitConstraints = 140),
                   cbArrangements = partitions(140, distinct = TRUE),
                   cbPartitions = diffparts(140),
                   times = 10, unit = "relative")
    #> Unit: relative
    #>            expr   min    lq  mean median    uq   max neval
    #>  cbRcppAlgosSer 1.241 1.237 1.224  1.203 1.178 1.378    10
    #>  cbArrangements 1.000 1.000 1.000  1.000 1.000 1.000    10
    #>    cbPartitions 5.807 5.795 5.658  5.771 5.670 5.205    10

#### Restricted Distinct Partitions

    t1 <- comboGeneral(160, 10,
                       constraintFun = "sum", comparisonFun = "==",
                       limitConstraints = 160)
    t2 <- partitions(160, 10, distinct = TRUE)
    identical(t1, t2)
    #> [1] TRUE
    dim(t1)
    #> [1] 8942920      10
    rm(t1, t2)
    invisible(gc())
    microbenchmark(cbRcppAlgosSer = comboGeneral(160, 10,
                                                 constraintFun = "sum",
                                                 comparisonFun = "==",
                                                 limitConstraints = 160),
                   cbArrangements = partitions(160, 10, distinct = TRUE),
                   times = 10, unit = "relative")
    #> Unit: relative
    #>            expr   min    lq  mean median    uq   max neval
    #>  cbRcppAlgosSer 1.000 1.000 1.000  1.000 1.000 1.000    10
    #>  cbArrangements 1.084 1.081 1.208  1.087 1.467 1.486    10

### Partitions - Repetition

#### All Partitions

    t1 <- comboGeneral(0:65, repetition = TRUE, constraintFun = "sum",
                       comparisonFun = "==", limitConstraints = 65)
    t2 <- partitions(65)
    t3 <- parts(65)

    # Each package has different output formats... we only dimensions
    # and that each result is a partition of 65
    identical(dim(t1), dim(t2))
    #> [1] TRUE
    identical(dim(t1), dim(t(t3)))
    #> [1] TRUE
    all(rowSums(t1) == 65)
    #> [1] TRUE
    all(rowSums(t2) == 65)
    #> [1] TRUE
    all(colSums(t3) == 65)
    #> [1] TRUE
    dim(t1)
    #> [1] 2012558      65
    rm(t1, t2, t3)
    invisible(gc())
    microbenchmark(cbRcppAlgosSer = comboGeneral(0:65, repetition = TRUE,
                                                 constraintFun = "sum",
                                                 comparisonFun = "==",
                                                 limitConstraints = 65),
                   cbArrangements = partitions(65),
                   cbPartitions = parts(65),
                   times = 20, unit = "relative")
    #> Unit: relative
    #>            expr   min    lq  mean median    uq    max neval
    #>  cbRcppAlgosSer 1.000 1.000 1.000  1.000 1.000 1.0000    20
    #>  cbArrangements 1.042 1.049 1.015  1.041 1.014 0.8623    20
    #>    cbPartitions 1.723 1.710 1.618  1.695 1.596 1.3554    20

#### Restricted Partitions

    t1 <- comboGeneral(100, 15, TRUE, constraintFun = "sum",
                       comparisonFun = "==", limitConstraints = 100)
    t2 <- partitions(100, 15)
    identical(t1, t2)
    #> [1] TRUE
    dim(t1)
    #> [1] 9921212      15
    rm(t1, t2)

    # This takes a really long time... not because of restrictedparts,
    # but because apply is not that fast. This transformation is
    # needed for proper comparisons. As a result, we will compare
    # a smaller example
    # t3 <- t(apply(as.matrix(restrictedparts(100, 15, include.zero = F)), 2, sort))
    t3 <- t(apply(as.matrix(restrictedparts(50, 15, include.zero = F)), 2, sort))
    identical(partitions(50, 15), t3)
    #> [1] TRUE
    rm(t3)
    invisible(gc())
    microbenchmark(cbRcppAlgosSer = comboGeneral(100, 15, TRUE,
                                                 constraintFun = "sum",
                                                 comparisonFun = "==",
                                                 limitConstraints = 100),
                           cbArrangements = partitions(100, 15),
                           cbPartitions = restrictedparts(100, 15, include.zero=FALSE),
                           times = 10, unit = "relative")
    #> Unit: relative
    #>            expr   min    lq  mean median     uq   max neval
    #>  cbRcppAlgosSer 1.092 1.087 1.048  1.009 0.9924 1.075    10
    #>  cbArrangements 1.000 1.000 1.000  1.000 1.0000 1.000    10
    #>    cbPartitions 3.256 3.240 2.910  2.853 2.6818 2.646    10

### Partitions - Multisets

Currenlty, `RcppAlgos` is the only package capable of efficiently generating partitions of multisets. Therefore, we will only time `RcppAlgos` and use this as a reference for future improvements.

    t1 <- comboGeneral(120, 10, freqs=rep(1:8, 15),
                       constraintFun = "sum", comparisonFun = "==",
                       limitConstraints = 120)
    dim(t1)
    #> [1] 7340225      10
    all(rowSums(t1) == 120)
    #> [1] TRUE
    microbenchmark(cbRcppAlgos = comboGeneral(120, 10, freqs=rep(1:8, 15),
                                              constraintFun = "sum", comparisonFun = "==",
                                              limitConstraints = 120), times = 10)
    #> Unit: seconds
    #>         expr   min  lq  mean median    uq   max neval
    #>  cbRcppAlgos 1.196 1.2 1.215  1.214 1.231 1.236    10

## Extreme Test

    set.seed(2187)
    testVector7 <- sort(sample(10^7, 10^3 + 100))

    for (i in 1:1) {
        gc()
        cat("\n************ ", prettyNum(comboCount(1100, 3), big.mark = ","),
        " COMBINATIONS **************\n")
        cat("\nRcppAlgos::comboGeneral with ", numThreads, " threads:\n")
        print(system.time(RcppAlgos::comboGeneral(testVector7, 3, nThreads = numThreads)))
        gc()
        cat("\nRcppAlgos::comboGeneral single-threaded:\n")
        print(system.time(RcppAlgos::comboGeneral(testVector7, 3)))
        gc()
        cat("\narrangements::combinations:\n")
        print(system.time(arrangements::combinations(x = testVector7, k = 3)))
        gc()
        cat("\n\n************** ", prettyNum(permuteCount(600, 3), big.mark = ","),
        " PERMUTATIONS **************\n")
        cat("\nRcppAlgos::permuteGeneral with ", numThreads, " threads:\n")
        print(system.time(RcppAlgos::permuteGeneral(testVector7[1:600], 3, nThreads = numThreads)))
        gc()
        cat("\nRcppAlgos::permuteGeneral single-threaded:\n")
        print(system.time(RcppAlgos::permuteGeneral(testVector7[1:600], 3)))
        gc()
        cat("\narrangements::permutations:\n")
        print(system.time(arrangements::permutations(x = testVector7[1:600], k = 3)))
        gc()
    }
    #>
    #> ************  221,228,700  COMBINATIONS **************
    #>
    #> RcppAlgos::comboGeneral with  4  threads:
    #>    user  system elapsed
    #>   1.271   0.736   0.534
    #>
    #> RcppAlgos::comboGeneral single-threaded:
    #>    user  system elapsed
    #>   1.091   0.474   1.566
    #>
    #> arrangements::combinations:
    #>    user  system elapsed
    #>   1.356   0.471   1.827
    #>
    #>
    #> **************  214,921,200  PERMUTATIONS **************
    #>
    #> RcppAlgos::permuteGeneral with  4  threads:
    #>    user  system elapsed
    #>   0.987   1.692   0.780
    #>
    #> RcppAlgos::permuteGeneral single-threaded:
    #>    user  system elapsed
    #>   0.930   0.481   1.412
    #>
    #> arrangements::permutations:
    #>    user  system elapsed
    #>   1.645   0.460   2.105