11
\$\begingroup\$

This challenge is partly an algorithms challenge, partly an optimization challenge and partly simply a fastest code challenge.

A T matrix is fully specified by its first row r and first column c. Each remaining element of the matrix is just a copy of the element that is diagonally up and left. That is M[i,j] = M[i-1,j-1]. We will allow T matrices which are not square. We do however always assume that the number of rows is no more than the number of columns. For example, consider the following 3 by 5 T matrix.

10111
11011
11101

We say a matrix has property X if it contains two non-empty sets of columns with non-identical indices which have the same (vector) sum. The vector sum of one or more columns is simply an element-wise summation of their columns. That is the sum of two or more columns containing x elements each is another column containing x elements. The sum of one column is trivially the column itself.

The matrix above trivially has property X as the first and last columns are the same. The identity matrix never has property X.

If we just remove the last column of the matrix above then we get an example which does not have property X and would give a score of 4/3.

1011
1101
1110

The task

The task is to write code to find the highest scoring T matrix with binary entries and which does not have property X. For clarity, a matrix with binary entries has the property that each one of its entries is either 0 or 1.

Score

Your score will be the number columns divided by the number of rows in your best scoring matrix.

Tie Breaker

If two answers have the same score, the one submitted first wins.

In the (very) unlikely event that someone finds a method to get unlimited scores, the first valid proof of such a solution will be accepted. In the even more unlikely event that you can find a proof of optimality of a finite matrix I will of course award the win too.

Hint

All the answers at Find highest scoring matrix without property X are valid here but they are not optimal. There are T matrices without property X which are not cyclic.

For example, there is a 7 by 12 T matrix without property X but no such cyclic matrix.

21/11 would beat all the current answers from this and the previous challenge.

Languages and libraries

You can use any language which has a freely available compiler/interpreter/etc. for Linux and any libraries which are also freely available for Linux.

Bonus The first answer with a score greater than 2 gets an immediate 200 point bounty award. Ton Hospel has now achieved this!


Current leader board

  • C++. Score 31/15 by Ton Hospel
  • Java. Score 36/19 by Peter Taylor
  • Haskell. Score 14/8 by alexander-brett
\$\endgroup\$
  • \$\begingroup\$ By the " two non-empty sets of columns with non-identical indices" you mean two sets of columns that are disjoint? Or, to rephrase this, is {1 ,3}, {1, 5} a valid two subsets of columns? \$\endgroup\$ – pawel.boczarski Apr 26 '15 at 10:52
  • \$\begingroup\$ @pawel.boczarski No not disjoint. Just non-identical. So {1 ,3}, {1, 5} is valid. \$\endgroup\$ – user9206 Apr 26 '15 at 10:53
  • \$\begingroup\$ Ok. What about M[i,1] - is it a "borrow" from the last column of M[i-1] (zero is not a valid matrix column index)? And actually this is "up and left" rather than "up and right". \$\endgroup\$ – pawel.boczarski Apr 26 '15 at 10:58
  • \$\begingroup\$ @pawel.boczarski "right" was a typo. Thanks. The first row and column can be set to anything you like. They define the whole rest of the matrix. Does that answer your question? \$\endgroup\$ – user9206 Apr 26 '15 at 13:13
  • \$\begingroup\$ Ok, I got it. It was my fault to not have read carefully that the first column is also defined. \$\endgroup\$ – pawel.boczarski Apr 26 '15 at 14:28
6
+100
\$\begingroup\$

C++, Score 23/12 25/13 27/14 28/14 31/15

Finally a result with ratio > 2:

rows=15,cols=31
1 1 0 0 1 0 0 0 1 1 1 1 0 1 0 1 1 0 0 1 0 0 0 1 1 1 1 0 1 1 0 
1 1 1 0 0 1 0 0 0 1 1 1 1 0 1 0 1 1 0 0 1 0 0 0 1 1 1 1 0 1 1 
1 1 1 1 0 0 1 0 0 0 1 1 1 1 0 1 0 1 1 0 0 1 0 0 0 1 1 1 1 0 1 
1 1 1 1 1 0 0 1 0 0 0 1 1 1 1 0 1 0 1 1 0 0 1 0 0 0 1 1 1 1 0 
1 1 1 1 1 1 0 0 1 0 0 0 1 1 1 1 0 1 0 1 1 0 0 1 0 0 0 1 1 1 1 
0 1 1 1 1 1 1 0 0 1 0 0 0 1 1 1 1 0 1 0 1 1 0 0 1 0 0 0 1 1 1 
0 0 1 1 1 1 1 1 0 0 1 0 0 0 1 1 1 1 0 1 0 1 1 0 0 1 0 0 0 1 1 
1 0 0 1 1 1 1 1 1 0 0 1 0 0 0 1 1 1 1 0 1 0 1 1 0 0 1 0 0 0 1 
0 1 0 0 1 1 1 1 1 1 0 0 1 0 0 0 1 1 1 1 0 1 0 1 1 0 0 1 0 0 0 
0 0 1 0 0 1 1 1 1 1 1 0 0 1 0 0 0 1 1 1 1 0 1 0 1 1 0 0 1 0 0 
0 0 0 1 0 0 1 1 1 1 1 1 0 0 1 0 0 0 1 1 1 1 0 1 0 1 1 0 0 1 0 
1 0 0 0 1 0 0 1 1 1 1 1 1 0 0 1 0 0 0 1 1 1 1 0 1 0 1 1 0 0 1 
0 1 0 0 0 1 0 0 1 1 1 1 1 1 0 0 1 0 0 0 1 1 1 1 0 1 0 1 1 0 0 
0 0 1 0 0 0 1 0 0 1 1 1 1 1 1 0 0 1 0 0 0 1 1 1 1 0 1 0 1 1 0 
1 0 0 1 0 0 0 1 0 0 1 1 1 1 1 1 0 0 1 0 0 0 1 1 1 1 0 1 0 1 1 

I completely explored 1 to 14 rows. 15 would take too long to completely explore. The results are:

1/1   = 1
2/2   = 1
4/3   = 1.333
5/4   = 1.25
7/5   = 1.4
9/6   = 1.5
12/7  = 1.714
14/8  = 1.75
16/9  = 1.778
18/10 = 1.8
20/11 = 1.818
23/12 = 1.917
25/13 = 1.923
28/14 = 2

The code given below is an older version of the program. The newest version is at https://github.com/thospel/personal-propertyX.

/*
  Compile using something like:
    g++ -Wall -O3 -march=native -fstrict-aliasing -std=c++11 -g propertyX.cpp -lpthread -o propertyX
*/
#include <cstdint>
#include <climits>
#include <ctgmath>
#include <iostream>
#include <vector>
#include <array>
#include <chrono>
#include <mutex>
#include <atomic>
#include <thread>

using namespace std;

const int ELEMENTS = 2;

using uint    = unsigned int;
using Element = uint64_t;
using Column  = array<Element, ELEMENTS>;
using Set     = vector<Column>;
using Sum     = uint8_t;
using Index   = uint32_t;
using sec = chrono::seconds;

int const PERIOD = 5*60;
int const MAX_ROWS = 54;
int const COL_FACTOR = (MAX_ROWS+1) | 1;                // 55
int const ROW_ZERO   = COL_FACTOR/2;                    // 27
int const ROWS_PER_ELEMENT = CHAR_BIT * sizeof(Element) / log2(COL_FACTOR); //11
Element constexpr ELEMENT_FILL(Element v = ROW_ZERO, int n = ROWS_PER_ELEMENT) {
    return n ? ELEMENT_FILL(v, n-1) * COL_FACTOR + v : 0;
}
Element constexpr POWN(Element v, int n) {
    return n ? POWN(v, n-1)*v : 1;
}
Element const ELEMENT_TOP = POWN(COL_FACTOR, ROWS_PER_ELEMENT -1);
int const MAX_COLS = ROWS_PER_ELEMENT * ELEMENTS;       // 22

atomic<Index> col_next;
atomic<uint>  period;
chrono::steady_clock::time_point start;
mutex period_mutex;

uint ratio_row;
uint ratio_col;
mutex ratio_mutex;

auto const nr_threads = thread::hardware_concurrency();
// auto const nr_threads = 1;

struct State {
    State(uint cols);
    void process(Index i);
    void extend(uint row);
    void print(uint rows);
    Index nr_columns() const { return static_cast<Index>(1) << cols_; }

    Column last_;
    Element top_;
    int top_offset_;
    uint ratio_row_ = 0;
    uint ratio_col_ = 1;
    uint cols_;
    array<Sum, MAX_ROWS + MAX_COLS -1> side;
    vector<Set> sets;
};

ostream& operator<<(ostream& os, Column const& column) {
    for (int i=0; i<ELEMENTS; ++i) {
        auto v = column[i];
        for (int j=0; j<ROWS_PER_ELEMENT; ++j) {
            auto bit = v / ELEMENT_TOP;
            cout << " " << bit;
            v -= bit * ELEMENT_TOP;
            v *= COL_FACTOR;
        }
    }
    return os;
}

State::State(uint cols) : cols_{cols} {
    sets.resize(MAX_ROWS+2);
    for (int i=0; i<2; ++i) {
        sets[i].resize(2);
        for (int j=0; j < ELEMENTS; ++j) {
            sets[i][0][j] =  ELEMENT_FILL();
            sets[i][1][j] =  static_cast<Element>(-1) - ELEMENT_FILL(1);
        }
    }
    top_ = POWN(COL_FACTOR, (cols_-1) % ROWS_PER_ELEMENT);
    top_offset_ = ELEMENTS - 1 - (cols_-1) / ROWS_PER_ELEMENT;
}

void State::print(uint rows) {
    for (auto c=0U; c<cols_;c++) {
        for (auto r=0U; r<rows;r++) {
            cout << static_cast<int>(side[cols_-c+r-1]) << " ";
        }
        cout << "\n";
    }
    cout << "----------" << endl;
}

void check(uint cols, uint t) {
    State state(cols);

    Index nr_columns = state.nr_columns();
    while (1) {
        Index col = col_next++;
        if (col >= nr_columns) break;
        state.process(col);

        auto now = chrono::steady_clock::now();
        auto elapsed = chrono::duration_cast<sec>(now-start).count();
        if (elapsed >= period) {
            lock_guard<mutex> lock{period_mutex};
            if (elapsed >= period) {
                cout << "col=" << col << "/" << nr_columns << " (" << 100.*col/nr_columns << "% " << elapsed << " s)" << endl;
                period = (elapsed/PERIOD+1)*PERIOD;
            }
        }
    }
}

void State::process(Index col) {
    last_.fill(0);
    for (uint i=0; i<cols_; ++i) {
        Element bit = col >> i & 1;
        side[i] = bit;
        Element carry = 0;
        for (int j=0; j<ELEMENTS; ++j) {
            auto c = last_[j] % COL_FACTOR;
            last_[j] = last_[j] / COL_FACTOR + carry * ELEMENT_TOP;
            if (j == top_offset_ && bit) last_[j] += top_;
            carry = c;
        }
    }
    // cout << "col=" << col << ", value=" << last_ << "\n";
    extend(0);
}

void State::extend(uint row) {
    // cout << "Extend row " << row << " " << static_cast<int>(side[cols_+row-1]) << "\n";
    if (row >= MAX_ROWS) throw(range_error("row out of range"));

    // Execute subset sum. The new column is added to set {from} giving {to}
    // {sum} is the other set.
    auto const& set_from = sets[row];
    auto const& set_sum  = sets[row + 1];
    auto      & set_to   = sets[row + 2];
    if (set_to.size() == 0) {
        auto size = 3 * set_from.size() - 2;
        set_to.resize(size);
        for (int j=0; j<ELEMENTS; ++j)
            set_to[size-1][j] = static_cast<Element>(-1) - ELEMENT_FILL(1);
    }

    // Merge sort {set_from - last_} , {set_from} and {set_from + last_}
    auto ptr_sum    = &set_sum[1][0];
    auto ptr_low    = &set_from[0][0];
    auto ptr_middle = &set_from[0][0];
    auto ptr_high   = &set_from[0][0];
    Column col_low, col_high;
    for (int j=0; j<ELEMENTS; ++j) {
        col_low   [j] = *ptr_low++  - last_[j];
        col_high  [j] = *ptr_high++ + last_[j];
    }

    auto ptr_end = &set_to[set_to.size()-1][0];
    auto ptr_to  = &set_to[0][0];
    while (ptr_to < ptr_end) {
        for (int j=0; j<ELEMENTS; ++j) {
            if (col_low[j] < ptr_middle[j]) goto LOW;
            if (col_low[j] > ptr_middle[j]) goto MIDDLE;
        }
        // low == middle
        // cout << "low == middle\n";
        return;

      LOW:
        // cout << "LOW\n";
        for (int j=0; j<ELEMENTS; ++j) {
            if (col_low[j] < col_high[j]) goto LOW0;
            if (col_low[j] > col_high[j]) goto HIGH0;
        }
        // low == high
        // cout << "low == high\n";
        return;

      MIDDLE:
        // cout << "MIDDLE\n";
        for (int j=0; j<ELEMENTS; ++j) {
            if (ptr_middle[j] < col_high[j]) goto MIDDLE0;
            if (ptr_middle[j] > col_high[j]) goto HIGH0;
        }
        // middle == high
        // cout << "middle == high\n";
        return;

      LOW0:
        // cout << "LOW0\n";
        for (int j=0; j<ELEMENTS; ++j) {
            *ptr_to++  = col_low[j];
            col_low[j] = *ptr_low++ - last_[j];
        }
        goto SUM;

      MIDDLE0:
        // cout << "MIDDLE0\n";
        for (int j=0; j<ELEMENTS; ++j)
            *ptr_to++ = *ptr_middle++;
        goto SUM;

      HIGH0:
        // cout << "HIGH0\n";
        for (int j=0; j<ELEMENTS; ++j) {
            *ptr_to++ = col_high[j];
            col_high[j] = *ptr_high++ + last_[j];
        }
        goto SUM;
      SUM:
        for (int j=-ELEMENTS; j<0; ++j) {
            if (ptr_to[j] > ptr_sum[j]) {
                ptr_sum += ELEMENTS;
                goto SUM;
            }
            if (ptr_to[j] < ptr_sum[j]) goto DONE;
        }
        // sum == to
        for (int j=-ELEMENTS; j<0; ++j)
            if (ptr_to[j] != ELEMENT_FILL()) {
                // sum == to and to != 0
                // cout << "sum == to\n";
                // cout << set_sum[(ptr_sum - &set_sum[0][0])/ELEMENTS-1] << "\n";
                return;
            }
      DONE:;
    }
    // cout << "Wee\n";
    auto row1 = row+1;
    if (0)
        for (uint i=0; i<row1+2; ++i) {
            cout << "Set " << i << "\n";
            auto& set = sets[i];
            for (auto& column: set)
                cout << column << "\n";
        }

    if (row1 * ratio_col_ > ratio_row_ * cols_) {
        ratio_row_ = row1;
        ratio_col_ = cols_;
        lock_guard<mutex> lock{ratio_mutex};

        if (ratio_row_ * ratio_col > ratio_row * ratio_col_) {

            auto now = chrono::steady_clock::now();
            auto elapsed = chrono::duration_cast<sec>(now-start).count();
            cout << "cols=" << cols_ << ",rows=" << row1 << " (" << elapsed << " s)\n";
            print(row1);
            ratio_row = ratio_row_;
            ratio_col = ratio_col_;
        }
    }

    auto last = last_;

    Element carry = 0;
    for (int j=0; j<ELEMENTS; ++j) {
        auto c = last_[j] % COL_FACTOR;
        last_[j] = last_[j] / COL_FACTOR + carry * ELEMENT_TOP;
        carry = c;
    }

    side[cols_+row] = 0;
    extend(row1);

    last_[top_offset_] += top_;
    side[cols_+row] = 1;
    extend(row1);

    last_ = last;
}

void my_main(int argc, char** argv) {
    if (!col_next.is_lock_free()) cout << "col_next is not lock free\n";
    if (!period.  is_lock_free()) cout << "period is not lock free\n";

    int min_col = 2;
    int max_col = MAX_COLS;
    if (argc > 1) {
        min_col = atoi(argv[1]);
        if (min_col < 2)
            throw(range_error("Column must be >= 2"));
        if (min_col > MAX_COLS)
            throw(range_error("Column must be <= " + to_string(MAX_COLS)));
    }
    if (argc > 2) {
        max_col = atoi(argv[2]);
        if (max_col < min_col)
            throw(range_error("Column must be >= " + to_string(min_col)));
        if (max_col > MAX_COLS)
            throw(range_error("Column must be <= " + to_string(MAX_COLS)));
    }

    for (int cols = min_col; cols <= max_col; ++cols) {
        cout << "Trying " << cols << " columns" << endl;
        ratio_row = 0;
        ratio_col = 1;
        col_next = 0;
        period = PERIOD;
        start = chrono::steady_clock::now();
        vector<thread> threads;
        for (auto t = 1U; t < nr_threads; t++)
            threads.emplace_back(check, cols, t);
        check(cols, 0);
        for (auto& thread: threads)
            thread.join();
    }
}

int main(int argc, char** argv) {
    try {
        my_main(argc, argv);
    } catch(exception& e) {
        cerr << "Error: " << e.what() << endl;
        exit(EXIT_FAILURE);
    }
    exit(EXIT_SUCCESS);
}
\$\endgroup\$
  • \$\begingroup\$ This is great. The great mystery is if you can ever get a score of 2. \$\endgroup\$ – user9206 Nov 5 '16 at 14:19
  • \$\begingroup\$ By extrapolation, 28/14 ought to exist I think which would be really exciting. But is it just out of reach? \$\endgroup\$ – user9206 Nov 5 '16 at 19:50
  • \$\begingroup\$ n=14 would take about 200 days with my current code on my 8-core CPU. The code can probably be sped up by 30% or so. After that I run out of ideas. And your extrapolation seems rather optimistic anyways... \$\endgroup\$ – Ton Hospel Nov 5 '16 at 19:56
  • \$\begingroup\$ I think the circulant 50 by 25 matrix with first row 01011011100010111101000001100111110011010100011010 might work. This was found by an optimization heuristic which might be useful. \$\endgroup\$ – user9206 Nov 5 '16 at 20:10
  • \$\begingroup\$ 140 hours for exhaustive cover of n=14 is very impressively fast I have to say. \$\endgroup\$ – user9206 Nov 6 '16 at 7:04
2
+50
\$\begingroup\$

Haskell 14/8 = 1.75

1 1 0 0 0 1 0 1 1 0 1 1 0 0
1 1 1 0 0 0 1 0 1 1 0 1 1 0
0 1 1 1 0 0 0 1 0 1 1 0 1 1
1 0 1 1 1 0 0 0 1 0 1 1 0 1
0 1 0 1 1 1 0 0 0 1 0 1 1 0
0 0 1 0 1 1 1 0 0 0 1 0 1 1
0 0 0 1 0 1 1 1 0 0 0 1 0 1
0 0 0 0 1 0 1 1 1 0 0 0 1 0

Previously 9/6 = 1.5

1 0 1 0 1 1 0 0 1
1 1 0 1 0 1 1 0 0
1 1 1 0 1 0 1 1 0
1 1 1 1 0 1 0 1 1
1 1 1 1 1 0 1 0 1
1 1 1 1 1 1 0 1 0

I wrote this, then looked at the answers to the other question and was... discouraged.

import Data.List
import Data.Hashable
import Control.Monad
import Control.Parallel.Strategies
import Control.Parallel
import qualified Data.HashSet as S

matrix§indices = [ matrix!!i | i<-indices ]

powerset :: [a] -> [[a]]
powerset = filterM (const [True, False])

hashNub :: (Hashable a, Eq a) => [a] -> [a]
hashNub l = go S.empty l
    where
      go _ []     = []
      go s (x:xs) = if x `S.member` s
        then go s xs
        else x : go (S.insert x s) xs

getMatrix :: Int -> Int -> [Int] -> [[Int]]
getMatrix width height vector = [ vector § [x..x+width-1] | x<-[0..height-1] ]

hasDuplicate :: (Hashable a, Eq a) => [a] -> Bool
hasDuplicate m = go S.empty m
    where
        go _ [] = False
        go s (x:xs) = if x `S.member` s
            then True
            else go (S.insert x s) xs

hasProperty :: [[Int]] -> Bool
hasProperty matrix =
    let
        base = replicate (length (matrix !! 0)) 0::[Int]
    in
        if elem base matrix then
            False
        else
            if hasDuplicate matrix then
                False
            else
                if hasDuplicate (map (foldl (zipWith (+)) base) (powerset matrix)) then
                    False
                else
                    True


pmap = parMap rseq

matricesWithProperty :: Int -> Int -> [[[Int]]]
matricesWithProperty n m =
    let
        base = replicate n 0::[Int]
    in
    filter (hasProperty) $
    map (getMatrix n m) $
    sequence [ [0,1] | x<-[0..n+m-1] ]

firstMatrixWithProperty :: Int -> Int -> [[Int]]
firstMatrixWithProperty n m = head $ matricesWithProperty n m

main = mapM (putStrLn. show) $ map (firstMatrixWithProperty 8) [1..]
\$\endgroup\$
  • \$\begingroup\$ Thank you! A Haskell answer is always very welcome. I think the first interesting case is 12/7. Can you get that? \$\endgroup\$ – user9206 Apr 28 '15 at 13:44
  • \$\begingroup\$ I was running it on a dual-core 2009 laptop, so no :) I'll try again on a faster machine though \$\endgroup\$ – alexander-brett Apr 28 '15 at 14:03
  • \$\begingroup\$ Very nice. I just added a comment that 21/11 would beat all previous answers. \$\endgroup\$ – user9206 Apr 28 '15 at 19:15
  • \$\begingroup\$ Could you explain what your code outputs exactly please? \$\endgroup\$ – user9206 Apr 29 '15 at 11:27
  • \$\begingroup\$ In main, the number (in this case 8) is the height of the matrix, and the program iterates through widths [1..]. For each height/width combination, it prints out an array of columns of a valid matrix. \$\endgroup\$ – alexander-brett Apr 30 '15 at 6:56
1
\$\begingroup\$

C

Here's an answer which works and is significantly more memory-efficient than the Haskell. Unfortunately, my computer is still too slow to achieve better than 14/8 in any reasonable length of time.

Try compiling with gcc -std=c99 -O2 -fopenmp -o matrices.exe matrices.c and running with matrices.exe width height or similar. The output is an integer, the bits of which form the basis for the matrix in question, for instance:

$ matrices.exe 8 14
...
valid i: 1650223

Then, since 1650223 = 0b110010010111000101111, the matrix in question is:

0 0 1 1 1 0 1 0 0 1 0 0 1 1
0 ...
1 ...
0
1
1
1
1

If someone with 8 cores and time on their hands wants to run this for a while, I think some good things could result :)


#include <stdio.h>
#include <stdlib.h>
#include <unistd.h>

/*
 * BEGIN WIKIPEDIA CODE
 */
const long long m1  = 0x5555555555555555; //binary: 0101...
const long long m2  = 0x3333333333333333; //binary: 00110011..
const long long m4  = 0x0f0f0f0f0f0f0f0f; //binary:  4 zeros,  4 ones ...
const long long m8  = 0x00ff00ff00ff00ff; //binary:  8 zeros,  8 ones ...
const long long m16 = 0x0000ffff0000ffff; //binary: 16 zeros, 16 ones ...
const long long m32 = 0x00000000ffffffff; //binary: 32 zeros, 32 ones
const long long hff = 0xffffffffffffffff; //binary: all ones
const long long h01 = 0x0101010101010101; //the sum of 256 to the power of 0,1,2,3...
//This uses fewer arithmetic operations than any other known
//implementation on machines with fast multiplication.
//It uses 12 arithmetic operations, one of which is a multiply.
long long hamming(long long x) {
    x -= (x >> 1) & m1;             //put count of each 2 bits into those 2 bits
    x = (x & m2) + ((x >> 2) & m2); //put count of each 4 bits into those 4 bits
    x = (x + (x >> 4)) & m4;        //put count of each 8 bits into those 8 bits
    return (x * h01)>>56;  //returns left 8 bits of x + (x<<8) + (x<<16) + (x<<24) + ...
}
/*
 * END WIKIPEDIA CODE
 */

int main ( int argc, char *argv[] ) {
    int height;
    int width;

    sscanf(argv[1], "%d", &height);
    sscanf(argv[2], "%d", &width);

    #pragma omp parallel for
    for (
        /*
         * We know that there are 2^(h+w-1) T-matrices, defined by the entries
         * in the first row and first column. We'll let the long long i
         * represent these entries, with 1s represented by set bits.
         *
         * The first (0) and last (1) matrix we will ignore.
         */
        long long i = 1;
        i < (1 << (height+width-1))-1;
        i++
    ) {
        // Flag for keeping track as we go along.
        int isvalid = 1;

        /*
         * Start by representing the matrix as an array of columns, with each
         * non-zero matrix entry as a bit. This allows us to construct them and
         * check equality very quickly.
         */
        long *cols = malloc(sizeof(long)*width);
        long colmask = (1 << height)-1;
        for (int j = 0; j < width; j++) {
            cols[j] = (i >> j) & colmask;
            if (cols[j] == 0) {
                //check no zero rows
                isvalid = 0;
            } else {
                //check no duplicate rows
                for (int k = 0; k < j; k++) {
                    if (cols[j] == cols[k]) {
                        isvalid = 0;
                    }
                }
            }
        }

        if (isvalid == 1) {
            /*
             * We'll also represent the matrix as an array of rows, in a
             * similar manner.
             */
            long *rows = malloc(sizeof(long)*height);
            long rowmask = (1 << width)-1;
            for (int j = 0; j < height; j++) {
                rows[j] = (i >> j) & rowmask;
            }

            int *sums[(1 << width)];
            for (long j = 0; j < 1<<width; j++) {
                sums[j] = (int*)malloc(sizeof(int)*height);
            }

            for (
                /*
                 * The powerset of columns has size 2^width. Again with the
                 * long; this time each bit represents whether the
                 * corresponding row is a member of the subset. The nice thing
                 * about this is we can xor the permutation with each row,
                 * then take the hamming number of the resulting number to get
                 * the sum.
                 */
                long permutation = 1;
                (isvalid == 1) && (permutation < (1 << width)-1);
                permutation ++
            ) {
                for (int j = 0; j < height; j++) {
                    sums[permutation][j] = hamming( rows[j] & permutation);
                }
                for (int j = permutation-1; (isvalid == 1) && (j > -1); j--) {
                    if (memcmp(sums[j], sums[permutation], sizeof(int)*height) == 0) {
                        isvalid = 0;
                    }
                }
            }

            for (long j = 0; j < 1<<width; j++) {
                free(sums[j]);
            }

            free(rows);

        }

        if (isvalid == 1) {
            printf ("valid i: %ld\n", i);
        }

        free(cols);
    }

    return 0;
}
\$\endgroup\$
  • \$\begingroup\$ I get alexander-brett.c: In function ‘main’: alexander-brett.c:107:21: warning: implicit declaration of function ‘memcmp’ [-Wimplicit-function-declaration] if (memcmp(sums[j], sums[permutation], sizeof(int)*height) == 0) { ^ alexander-brett.c:122:13: warning: format ‘%ld’ expects argument of type ‘long int’, but argument 2 has type ‘long long int’ [-Wformat=] printf ("valid i: %ld\n", i); \$\endgroup\$ – user9206 May 1 '15 at 16:50
  • \$\begingroup\$ How long does ./alexander-brett 8 14 take for you? \$\endgroup\$ – user9206 May 1 '15 at 18:24
  • \$\begingroup\$ Hi Lembik, 8 14 got 5 answers in an hour for me on a quad core machine. I managed to compile with those headers on windows, it'd be wierd if memcmp went missing... \$\endgroup\$ – alexander-brett May 1 '15 at 18:26
  • \$\begingroup\$ stackoverflow.com/questions/7132039/… \$\endgroup\$ – user9206 May 1 '15 at 18:29
  • \$\begingroup\$ I tried your code with 7 12 and one of the answers it outputs is valid i: 7481. In python bin(7481) is 0b1110100111001 which isn't long enough. Any idea what is going on? \$\endgroup\$ – user9206 May 19 '15 at 19:36

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy