---
title: "High Performance Benchmarks"
author: "Joseph Wood"
date: "2/8/2020"
output:
  rmarkdown::html_vignette:
    toc: true
    number_sections: false
vignette: >
  %\VignetteIndexEntry{High Performance Benchmarks}
  %\VignetteEngine{knitr::rmarkdown}
  \usepackage[utf8]{inputenc}
---

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
](<https://stackoverflow.com/a/47983855/4408538>). 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](<https://cran.r-project.org/package=arrangements>), [partitions](<https://cran.r-project.org/package=partitions>), and [RcppAlgos](<https://cran.r-project.org/package=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.

```r
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

```r
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

```r
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

```r
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

```r
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

```r
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

```r
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

```r
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

```r
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

```r
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

```r
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.

```r
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

```r

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
```
