Skip to main content
3 of 7
added 61 characters in body
primo
  • 33.5k
  • 5
  • 61
  • 139

#RPython 5.4.1, n ≈ 33 (48 seconds)

from rpython.rlib.rtime import time
from rpython.rlib.rarithmetic import r_int, r_uint
from rpython.rlib.rrandom import Random
from rpython.rlib.rposix import pipe, close, read, write, fork, waitpid


bitcounts = bytearray([0])
for i in range(16):
  b = bytearray([j+1 for j in bitcounts])
  bitcounts += b


def bitcount(n):
  bits = 0
  while n:
    bits += bitcounts[n & 65535]
    n >>= 16
  return bits


def main(argv):
  if len(argv) < 2:
    write(2, 'Usage: %s NUM_THREADS [N]'%argv[0])
    return 1
  threads = int(argv[1])

  if len(argv) > 2:
    n = int(argv[2])
    rnd = Random(r_uint(time()*1000))
    m = []
    for i in range(n):
      row = []
      for j in range(n):
        row.append(1 - r_int(rnd.genrand32() & 2))
      m.append(row)
  else:
    m = []
    strm = ""
    while True:
      buf = read(0, 4096)
      if len(buf) == 0:
        break
      strm += buf
    rows = strm.split("\n")
    for row in rows:
      r = []
      for val in row.split(' '):
        r.append(int(val))
      m.append(r)
    n = len(m)

  a = []
  for row in m:
    val = 0
    for v in row:
      val = (val << 1) | -(v >> 1)
    a.append(val)

  pids = []
  handles = []
  total = 0
  for i in xrange(threads):
    r, w = pipe()
    pid = fork()
    if pid:
      close(w)
      pids.append(pid)
      handles.append(r)
    else:
      close(r)
      total = run(n, a, i, threads)
      write(w, str(total))
      close(w)
      return 0

  for pid in pids:
    waitpid(pid, 0)

  for handle in handles:
    strval = read(handle, 256)
    total += int(strval)
    close(handle)

  print total >> n-1

  return 0


def run(n, a, mynum, threads):
  start = (1 << n-1) * mynum / threads
  end = (1 << n-1) * (mynum+1) / threads

  dtotal = 0
  for delta in xrange(start, end):
    pdelta = 1 - ((bitcount(delta) & 1) << 1)
    for v in a:
      pdelta *= n - (bitcount(delta ^ v) << 1)
    dtotal += pdelta

  return dtotal


def target(*args):
  return main

To compile, download the most recent PyPy source, and execute the following:

pypy /path/to/pypy-src/rpython/bin/rpython matrix-permanent.py

The resulting executable will be named matrix-permanent-c or similiar in the current working directory.

As of PyPy 5.0, RPython's threading primitives are a lot less primitive than they used to be. Newly spawned threads require the GIL, which is more or less useless for parallel computations. I've used fork instead, so it may not work as expected on Windows, although I haven't tested fails to compile (unresolved external symbol _fork).

The executable accepts up to two command line parameters. The first is the number of threads, the second optional parameter is n. If it is provided, a random matrix will be generated, otherwise it will be read from stdin. Each row must be newline separated (without a trailing newline), and each value space separated. The third example input would be given as:

1 -1 1 -1 -1 1 1 1 -1 -1 -1 -1 1 1 1 1 -1 1 1 -1
1 -1 1 1 1 1 1 -1 1 -1 -1 1 1 1 -1 -1 1 1 1 -1
-1 -1 1 1 1 -1 -1 -1 -1 1 -1 1 1 1 -1 -1 -1 1 -1 -1
-1 -1 -1 1 1 -1 1 1 1 1 1 1 -1 -1 -1 -1 -1 -1 1 -1
-1 1 1 1 -1 1 1 1 -1 -1 -1 1 -1 1 -1 1 1 1 1 1
1 -1 1 1 -1 -1 1 -1 1 1 1 1 -1 1 1 -1 1 -1 -1 -1
1 -1 -1 1 -1 -1 -1 1 -1 1 1 1 1 -1 -1 -1 1 1 1 -1
1 -1 -1 1 -1 1 1 -1 1 1 1 -1 1 -1 1 1 1 -1 1 1
1 -1 -1 -1 -1 -1 1 1 1 -1 -1 -1 -1 -1 1 1 -1 1 1 -1
-1 -1 1 -1 1 -1 1 1 -1 1 -1 1 1 1 1 1 1 -1 1 1
-1 -1 -1 -1 -1 -1 -1 1 -1 -1 -1 -1 1 1 1 1 -1 -1 -1 -1
1 1 -1 -1 -1 1 1 -1 -1 1 -1 1 1 -1 1 1 1 1 1 1
-1 1 1 -1 -1 -1 -1 -1 1 1 1 1 -1 -1 -1 -1 -1 1 -1 1
1 1 -1 -1 -1 1 -1 1 -1 -1 -1 -1 1 -1 1 1 -1 1 -1 1
1 1 1 1 1 -1 -1 -1 1 1 1 -1 1 -1 1 1 1 -1 1 1
1 -1 -1 1 -1 -1 -1 -1 1 -1 -1 1 1 -1 1 -1 -1 -1 -1 -1
-1 1 1 1 -1 1 1 -1 -1 1 1 1 -1 -1 1 1 -1 -1 1 1
1 1 -1 -1 1 1 -1 1 1 -1 1 1 1 -1 1 1 -1 1 -1 1
1 1 1 -1 -1 -1 1 -1 -1 1 1 -1 -1 -1 1 -1 -1 -1 -1 1
-1 1 1 1 -1 -1 -1 -1 -1 -1 -1 1 1 -1 1 1 -1 1 -1 -1

Sample Usage

$ time ./matrix-permanent-c 8 30
-15770583040

real    0m5.922s
user    0m47.376s
sys     0m0.000s

Method

I've used the Balasubramanian-Bax/Franklin-Glynn formula, with a runtime complexity of O(2nn). However, instead of iterating the δ in grey code order, I've instead replaced vector-row multiplication with a single xor operation (mapping (1, -1) → (0, 1)). The vector sum can likewise be found in a single operation, by taking n minus twice the popcount.

primo
  • 33.5k
  • 5
  • 61
  • 139