diff --git a/run-tests.py b/run-tests.py index b0ab93f..aa6ad37 100644 --- a/run-tests.py +++ b/run-tests.py @@ -7,9 +7,10 @@ import time dimension=2 real=np.ones((dimension, dimension)) -imag=np.ones((dimension, dimension)) +imag=np.ones((dimension, dimension))*0 submatrix=real+1j*imag -print submatrix.dtype +submatrix[0,0]*=2 +print submatrix p=lib.permanent(submatrix) print p diff --git a/src/permanent.c b/src/permanent.c index db5c9e6..eaf4bd0 100644 --- a/src/permanent.c +++ b/src/permanent.c @@ -10,8 +10,8 @@ #define SM_shape(x0) (int) PyArray_DIM(submatrix, x0) // Complex numbers -static const npy_complex128 complex_zero = {.real=0, .imag=0}; static const npy_complex128 complex_one = {.real=1, .imag=0}; +static const npy_complex128 complex_zero = {.real=0, .imag=0}; static void complex_inc(npy_complex128 a, npy_complex128 b) { a.real+=b.real; a.imag+=b.imag;} static npy_complex128 complex_prod(npy_complex128 a, npy_complex128 b) { npy_complex128 x; @@ -39,17 +39,23 @@ static npy_complex128 perm_ryser(PyArrayObject *submatrix) { npy_complex128 perm = complex_zero; npy_complex128 sum = complex_zero; int exp = 1 << n; + printf("2^%d = %d\n", n, exp); // Iterate over exponentially many index strings for (i=0; i 0) { complex_inc(sum, SM(z,y)); } + printf("sum %.3f %.3f \n", sum.real, sum.imag); } prod = complex_prod(prod, sum); } + printf("prod %.3f %.3f \n", prod.real, prod.imag); if (i%2 == 0) {prod.real*=-1; prod.imag*=-1;} complex_inc(perm, prod); }