jwood000 / RcppAlgos

Tool for Solving Problems in Combinatorics and Computational Mathematics
GNU General Public License v2.0
44 stars 5 forks source link

Multisets under constraint with integer vectors #12

Closed jwood000 closed 4 years ago

jwood000 commented 4 years ago

When trying to replicate the results of this question: Remove repeating numbers from for loop, I found an issue with how version 2.3.5 handles multisets and their respective frequencies when the vector passed isn't sorted. Observe:

set.seed(42)
vals <- sample(-50:50, 20)
vals
 [1]  42  43 -22  31  12  -1  19 -38  11  14  -9  41  33 -28 -10  30
[17]  38 -41 -11  -5

vals <- c(0L, vals)
myfreqs <- c(length(vals) - 1, rep(1, length(vals) - 1))

myfreqs
 [1] 20  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1

combs <- comboGeneral(vals, length(vals), freqs = myfreqs, 
                                            constraintFun = "sum", 
                                            comparisonFun = "==", 
                                            limitConstraints = 170)

dim(combs)
[1] 6561   21

## Brute force verification
allCombs <- comboGeneral(sort(vals), length(vals),
                         freqs = myfreqs[order(vals)], constraintFun = "sum")

verifiedCombs <- allCombs[which(allCombs[,22] == 170), ]

dim(verifiedCombs)
[1] 2135   22

Clearly, something is not right. We obtained 2135 solutions with brute force and 6561 with the constraint algorithm. Let's have a look at the actual results.

head(combs[,1:12])
     [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12]
[1,]  -41  -38  -28  -22  -11  -10   -5    0   11    11    11    11
[2,]  -41  -38  -28  -22  -11  -10   -5   11   11    11    11    11
[3,]  -41  -38  -28  -22  -11  -10   -5   11   11    11    11    11
[4,]  -41  -38  -28  -22  -11  -10   -1   11   11    11    11    11
[5,]  -41  -38  -28  -22  -11  -10    0   11   11    11    11    11
[6,]  -41  -38  -28  -22  -11  -10    0   11   11    11    11    11

head(verifiedCombs)
     [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12]
[1,]  -41  -38  -28  -22  -10   -5    0    0    0     0    11    12
[2,]  -41  -38  -28  -22   -9   -5   -1    0    0     0    11    12
[3,]  -41  -38  -28  -22   -1    0    0    0    0     0     0    11
[4,]  -41  -38  -28  -11  -10   -5    0    0    0     0     0    12
[5,]  -41  -38  -28  -11   -9   -5   -1    0    0     0     0    12
[6,]  -41  -38  -28  -11   -9   -5    0    0    0     0     0    11
     [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21] [,22]
[1,]    14    19    30    31    33    38    41    42    43   170
[2,]    14    19    30    31    33    38    41    42    43   170
[3,]    12    19    30    31    33    38    41    42    43   170
[4,]    14    19    30    31    33    38    41    42    43   170
[5,]    14    19    30    31    33    38    41    42    43   170
[6,]    14    19    30    31    33    38    41    42    43   170

We can see that the brute force solution correctly repeats zero multiple times whereas the constraint algo repeats 11 multiple times. This sorting is taken care of in GetPartitionCase:

https://github.com/jwood000/RcppAlgos/blob/8364c46d275f745c1cde9ad064f9c67f7c86c559/src/ConstraintsUtils.cpp#L397-L408

And we see that it is a templated function that can accept vectors of different types:

https://github.com/jwood000/RcppAlgos/blob/8364c46d275f745c1cde9ad064f9c67f7c86c559/src/ConstraintsUtils.cpp#L377-L382

Now, ConstraintsMaster.cpp calls GetPartitionCase here:

https://github.com/jwood000/RcppAlgos/blob/8364c46d275f745c1cde9ad064f9c67f7c86c559/src/ConstraintsMaster.cpp#L150-L156

And here is the problem. We are only taking care of numeric vectors (vNum and targetVals are vectors of type double). We can verify this by simply wrapping vals above with as.numeric:

combsNum <- comboGeneral(as.numeric(vals), length(vals),
                         freqs = myfreqs, 
                         constraintFun = "sum", 
                         comparisonFun = "==", 
                         limitConstraints = 170)

dim(combsNum)
[1] 2135   21

all.equal(combsNum, verifiedCombs[, 1:21])
[1] TRUE

We simply need to ensure we handle the integer case appropriately.

jwood000 commented 4 years ago

Fixed by https://github.com/jwood000/RcppAlgos/commit/e0b4909050a24c008613877634d7b2e0d0d3a873