## Phil Harvey’s Puzzle

### November 15, 2011

According to the binomial theorem, there are 12870 ways to choose 8 items from a set of 16:

```(define (choose n k)   (if (zero? k) 1     (* n (/ k) (choose (- n 1) (- k 1)))))```

```> (choose 16 8) 12870```

A recursive algorithm enumerates the k-subsets of n items. Consider the first element of the list. There are two kinds of subsets: those that include the first element, and those that don’t. Those that do can be enumerated by combining the first element with all (k−1)-subsets of the rest of the elements. Those that don’t can be enumerated by taking all the k-subsets of the rest of the elements. The recursion stops when k is zero or the list is empty. Here is our version of the algorithm:

```(define (combs n xs)   (define (comb x xs)     (if (null? xs) xs       (if (pair? (car xs))           (cons (append (list x) (car xs)) (comb x (cdr xs)))           (cons (list x (car xs)) (comb x (cdr xs))))))   (if (or (zero? n) (null? xs)) (list)     (if (= n 1) xs       (append (comb (car xs) (combs (- n 1) (cdr xs)))               (combs n (cdr xs))))))```

Now it’s easy. The sum of the numbers from 1 to 16 is 136, the sum of their squares is 1496, and the sum of their cubes is 18496, so we enumerate all 12870 subsets and check for half of each of the three sums:

```> (let loop ((xss (combs 8 (range 1 17))))     (let* ((xs (car xss))            (s1 (apply + xs))            (s2 (apply + (map square xs)))            (s3 (apply + (map cube xs))))     (if (and (= s1 68) (= s2 748) (= s3 9248)) xs       (loop (cdr xss))))) (1 4 6 7 10 11 13 16)```

The other partition is thus `(2 3 5 8 9 12 14 15)`.

We used `range` from the Standard Prelude. The definitions of `square` and `cube` are obvious, and are given at http://programmingpraxis.codepad.org/90qQobLz, where you can run the program.

Pages: 1 2

### 18 Responses to “Phil Harvey’s Puzzle”

1. kernelbob said

Here’s my Python interactive session, with one (typo, error) removed.

```~> python
Python 2.7.1+ (r271:86832, Apr 11 2011, 18:13:53)
[GCC 4.5.2] on linux2
>>> a = range(1, 17)
>>> a
[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16]
>>> from itertools import combinations
>>> s1 = sum
>>> s2 = lambda s: sum(i**2 for i in s)
>>> s3 = lambda s: sum(i**3 for i in s)
>>> s3(a)
18496
>>> s2(a)
1496
>>> s1(a)
136
>>> for p in combinations(a, 8):
...  if all(f(p) == t / 2 for f, t in zip((s1, s2, s3), (136, 1496, 18496))):
...   print p
...
(1, 4, 6, 7, 10, 11, 13, 16)
(2, 3, 5, 8, 9, 12, 14, 15)
>>>
```
2. ```
#include <stdio.h>
#include <math.h>

void
process(int s)
{
int asum, asum2, asum3 ;
int bsum, bsum2, bsum3 ;
int i ;

asum = asum2 = asum3 = 0 ;
bsum = bsum2 = bsum3 = 0 ;

for (i=0; i<16; i++) {
int val, val2, val3 ;
val = i ;
val2 = val * val ;
val3 = val2 * val ;
if (s & (1<<i)) {
asum += val ;
asum2 += val2 ;
asum3 += val3 ;
} else {
bsum += val ;
bsum2 += val2 ;
bsum3 += val3 ;
}
}
if (asum == bsum && asum2 == bsum2 && asum3 == bsum3) {
printf("%d %d %d\n", asum, asum2, asum3) ;
printf("A: ") ;
for (i=0; i<16; i++)
if (s & (1<<i)) printf("%d ", i+1) ;
printf("B: ") ;
for (i=0; i<16; i++)
if ((s & (1<<i)) == 0)  printf("%d ", i+1) ;
printf("\n") ;
}
}

main()
{
int s ;

for (s=0; s<(1<<15); s++)
process(s) ;
}
```
3. If you run my version of the program, you’ll see that there are two unique solutions, only one of which divides { 1..16 } into subsets of equal size. Both solutions have the same sums though. Interesting.

4. Ah. That teaches me to code before breakfast. I have a small error in the program above: the line which reads val = i should be “val = i + 1 ;”.

When that change is made, it finds the single, unique solution to the problem. My hint that something was wrong was that the sum of the numbers was 60. The sum of the numbers from 1 to 16 should be 136, and since they are divided into two partitions, they should each sum to 68.

5. Graham said

Since Michi and kernelbob already posted some Python, here’s a Common Lisp solution:

```(defun combinations (xs k)
(cond ((or (zerop k) (null xs)) nil)
((= k 1) (mapcar #'list xs))
(t (append (mapcar #'(lambda (c) (cons (first xs) c))
(combinations (rest xs) (1- k)))
(combinations (rest xs) k)))))

(defun sums (n)
(list (/ (* n (1+ n)) 2)
(/ (* n (1+ n) (1+ (* 2 n))) 6)
(expt (/ (* n (1+ n)) 2) 2)))

(defun sum (seq) (reduce #'+ seq :initial-value 0))
(defun square (x) (* x x))
(defun cube (x) (* x x x))

(defun prop (subseq)
(equal (mapcar #'(lambda (x) (/ x 2)) (sums 16))
(list (sum subseq)
(sum (mapcar #'square subseq))
(sum (mapcar #'cube subseq)))))

(let ((xs (loop for i from 1 to 16 collect i)))
(dolist (c (remove-if-not #'prop (combinations xs 8)))
(print (cons c (list (set-difference xs c))))))
```

```import Data.List ((\\))

combinations :: [Integer] -> Integer -> [[Integer]]
combinations _ 0 = []
combinations [] _ = []
combinations xs 1 = map (:[]) xs
combinations (x:xs) k = map (x:) (combinations xs (k - 1)) ++ combinations xs k

sums :: Integer -> [Integer]
sums n = [(n * (n + 1)) `div` 2,
(n * (n + 1) * (2 * n + 1)) `div` 6,
(n * (n + 1) `div` 2) ^2]

prop :: [Integer] -> Bool
prop xs = map (`div` 2) (sums 16) == map sum [xs, map (^2) xs, map (^3) xs]

main :: IO ()
main = mapM_ printPairs filteredCombs where
printPairs = print . (\ c -> (c, [1..16] \\ c))
filteredCombs =  filter prop \$ combinations [1..16] 8
```

In both, I calculated the sums of the integers, squares of integers, and cubes of integers in {1, …, 16} via formulae for the first three types of figurate numbers.

6. slabounty said

Here’s a ruby version …

```def sum(x)
x.inject(0) { |sum, v| sum += v }
end

def sum_squares(x)
x.inject(0) { |sum, v| sum += v**2 }
end

def sum_cubes(x)
x.inject(0) { |sum, v| sum += v**3 }
end

x = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16]

sum_x = sum(x)
sum_squares_x = sum_squares(x)
sum_cubes_x = sum_cubes(x)

x.permutation(8) do |p|
if (sum(p) == sum_x/2) && (sum_squares(p) == sum_squares_x/2) && (sum_cubes(p) == sum_cubes_x/2)
break;
end
end
```
7. Jan Van lent said

This problem can be formulated as a (binary) integer programming problem.
For example, using GMPL

```var x{i in 1..16}, binary;
s.t. pow{j in 0..3}: sum{i in 1..16} i^j * x[i] = (sum{i in 1..16} i^j)/2;
solve;
display{i in 1..16} x[i];
end;
```
8. kernelbob said

I just found the “snoob” function over here. http://www.steike.com/code/useless/evil-c/ That led to this somewhat golfed C solution.

Basically, snoob(n) returns the next integer after n with the same number of bits set. If we use a bit mask to denote the members of the first partition, we can use snoob() to iterate through the partitions directly.

for (i = 0xFF; i <= 0xFF00; i = snoob(i)) { …}
// 0xFF is the first partition, {1, 2, 3, 4, 5, 6, 7, 8}.
// 0xFF00 is the last partition, {9, 10, 11, 12, 13, 14, 15, 16}.

snoob() should be blazingly fast since it's about 7 arithmetic operations and no memory references.

#include <stdio.h>
#include <string.h>

unsigned snoob(unsigned x) {
unsigned smallest = x & -x;
unsigned ripple = x + smallest;
unsigned ones = x ^ ripple;
ones = (ones >> 2) / smallest;
return ripple | ones;
}

int main()
{
unsigned i, j, k, l;

for (i = 0xFF; i <= 0xFF00; i = snoob(i)) {
unsigned sums = { { 0, 0, 0 }, { 0, 0, 0 } };
for (j = 0; j < 16; j++)
for (k = 0, l = 1; k < 3; k++)
sums[!(i & 1 << j)][k] += l *= j + 1;
if (!memcmp(sums, sums, sizeof sums)) {
for (j = 0; j < 16; j++)
if (i & 1 << j)
printf("%d ", j + 1);
printf("\n");
}
}
return 0;
}

9. kernelbob said

Oops. Formatted.

```#include <stdio.h>
#include <string.h>

unsigned snoob(unsigned x) {
unsigned smallest = x & -x;
unsigned ripple = x + smallest;
unsigned ones = x ^ ripple;
ones = (ones >> 2) / smallest;
return ripple | ones;
}

int main()
{
unsigned i, j, k, l;

for (i = 0xFF; i <= 0xFF00; i = snoob(i)) {
unsigned sums = { { 0, 0, 0 }, { 0, 0, 0 } };
for (j = 0; j < 16; j++)
for (k = 0, l = 1; k < 3; k++)
sums[!(i & 1 << j)][k] += l *= j + 1;
if (!memcmp(sums, sums, sizeof sums)) {
for (j = 0; j < 16; j++)
if (i & 1 << j)
printf("%d ", j + 1);
printf("\n");
}
}
return 0;
}
```
10. neez said

```Prelude Data.List> let ps p = sum.map (^p) in filter (\x -> and [8==length x,68==ps 1 x,748==ps 2 x,9248==ps 3 x]) \$ subsequences [1..16]
[[2,3,5,8,9,12,14,15],[1,4,6,7,10,11,13,16]]
```
11. neez said

And a little sorter solution with list comprehension:

```Prelude Data.List> [x | x <- subsequences [1..16], let ps p = sum.map (^p), 8==length x, 68==ps 1 x, 748==ps 2 x, 9248==ps 3 x]
[[2,3,5,8,9,12,14,15],[1,4,6,7,10,11,13,16]]
```
12. Jan Van lent said

@neez Nice Haskell solution. You can also do it without precalculating the constants.

```Prelude Data.List> let ps p = sum.map (^p) in filter (\x -> and [ps i [1..16]==2*(ps i x) | i <- [0..3]]) \$ subsequences [1..16]
[[2,3,5,8,9,12,14,15],[1,4,6,7,10,11,13,16]]
```
13. neez said

@Jan Van lent
Wow, beautiful soution, especially using the i=0 case to partition into equal sized sets!

14. Jan Van lent said

@neez Well, it is only a small change from your solution. I used the same trick in the GMPL solution. Note that the constraint that the sets be of equal size is actually redundant. You can try with [1..3], [2..3] and even just 3 and they all give the same solution.

15. Yuushi said

Create 2 (random) partitions, and swap between them, making sure the difference is less. Keep doing this until we get a difference of 0.

```import random

gen = random.Random()

def diff(s, t):
return sum(abs(sum((k**x for k in s)) - sum((i**x for i in t))) for x \
in range(1,4))

def swap(s, t):
start = diff(s,t)
end = start + 1
while end > start:
ind1 = gen.randint(0, 7)
ind2 = gen.randint(0, 7)
s[ind1], t[ind2] = t[ind2], s[ind1]
end = diff(s,t)

def main():
x = [i for i in range(1,17)]
s = random.sample(x, 8)
t = [k for k in x if k not in s]

while diff(s,t) != 0:
swap(s,t)

s.sort()
t.sort()

print(s)
print(t)

if __name__ == '__main__':
main()
```
16. Simone said
```public class HarveyPuzzle {
public static final int ASIZE=16;

public HarveyPuzzle(){
puzzleStart();
}

public static void main(String[] args) {
new HarveyPuzzle();
}

private void puzzleStart() {
boolean [] pick = new boolean[ASIZE];
recursiveFinder(pick, 0);
}

private void recursiveFinder(boolean [] pick, int k){
if(k==pick.length-1){
isSolution(pick);
pick[k]=!pick[k];
isSolution(pick);
} else {
recursiveFinder(pick, k+1);
pick[k]=!pick[k];
recursiveFinder(pick, k+1);
}
}

private void isSolution(boolean [] pick){
if(sum(pick, true)==sum(pick, false)
&&
squaresSum(pick, true)==squaresSum(pick, false)
&&
cubesSum(pick, true)==cubesSum(pick, false)){

System.out.println("Solution found!");
printSolution(pick);
}
}

private int sum(boolean [] pick, boolean w) {
int somma=0;

for(int i=0;i<pick.length;i++){
if(pick[i]==w)
somma+=i+1;
}

return somma;
}

private int squaresSum(boolean [] pick, boolean w) {
int somma=0;

for(int i=0;i<pick.length;i++){
if(pick[i]==w)
somma+=Math.pow(i+1, 2);
}

return somma;
}

private int cubesSum(boolean [] pick, boolean w) {
int somma=0;

for(int i=0;i<pick.length;i++){
if(pick[i]==w)
somma+=Math.pow(i+1, 3);
}

return somma;
}

private void printSolution(boolean[] pick) {
System.out.print("[");

for(int i=0;i<pick.length;i++){
if(pick[i]){
System.out.print((i+1)+" ");
}
}

System.out.print("], [");

for(int i=0;i<pick.length;i++){
if(!pick[i]){
System.out.print((i+1)+" ");
}
}

System.out.println("]\n");

System.err.println(sum(pick, true) + "==" + sum(pick, false));
System.err.println(squaresSum(pick, true) + "==" + squaresSum(pick, false));
System.err.println(cubesSum(pick, true) + "==" + cubesSum(pick, false)+"\n\n");
}

}
```
17. Vaibhav Bansal said

package com.algorithm;

import java.util.Arrays;

public class SameSumSquareCube {

private static int[] set1;
private static int[] set2;

public static void main(String[] args) {

set1 = new int[] {1, 2, 3, 4, 5, 6, 7, 8};
set2 = new int[] {9, 10, 11, 12, 13, 14, 15, 16};

recursive(7);
System.out.println(“completed”);
}

private static void recursive(int count) {
for (int i=count; i>=0; i–) {
for (int j=7; j>=0; j–) {
int temp = set1[i];
set1[i] = set2[j];
set2[j] = temp;
if (isCubeSumEqual(set1, set2) && isSquareSumEqual(set1, set2) && isSumEqual(set1, set2)) {
System.out.println(Arrays.toString(set1));
System.out.println(Arrays.toString(set2));
}
}
recursive(–count);
}
}

private static boolean isCubeSumEqual(int[] set1, int[] set2) {

int set1CubeSum = 0;
for (int number : set1) {
set1CubeSum += number * number * number;
}

int set2CubeSum = 0;
for (int number : set2) {
set2CubeSum += number * number * number;
}
return set1CubeSum == set2CubeSum;
}

private static boolean isSquareSumEqual(int[] set1, int[] set2) {

int set1SquareSum = 0;
for (int number : set1) {
set1SquareSum += number * number;
}

int set2SquareSum = 0;
for (int number : set2) {
set2SquareSum += number * number;
}
return set1SquareSum == set2SquareSum;
}

private static boolean isSumEqual(int[] set1, int[] set2) {

int set1Sum = 0;
for (int number : set1) {
set1Sum += number;
}

int set2Sum = 0;
for (int number : set2) {
set2Sum += number;
}
return set1Sum == set2Sum;
}
}

18. […] the weekend, I was amusing myself with this problem from the always entertaining Programming Praxis. The problem is to partition the set {1, 2, […]