Open in CoCalc with one click!
In [1]:
import numpy as np from itertools import chain
In [2]:
def memoize(f): memo = {} def helper(x): if x not in memo: memo[x] = f(x) return memo[x] return helper

The following is a terrible function being recursive. Running through the n element subsets requires rebuilding all the m element subsets for m < n. A dynamical algorithm or memoization needs to be used, but at the expense of a huge memory cost.

A better, bounded search algorithm should be looked for here.

In [3]:
def subSet(n, m): if m == 1: return ([i] for i in range(n)) else: S = subSet(n, m - 1) S = (a + [i] for a in S if max(a) < n - 1 for i in range(max(a) + 1, n)) return S
In [4]:
A = np.array([[1 + 2*j + 2*i for j in range(6)] for i in range(6)])
In [5]:
LEFT = np.array([5,7,8,12,15,15]) TOP = np.array([6,7,11,12,14,17])
In [6]:
def err(B, left, top): leftError = left - np.sum(B,1)/np.sum(B != 0,1) topError = top - np.sum(B,0)/np.sum(B != 0,0) err = sum(leftError**2) + sum(topError**2) return err
In [7]:
def sub(B,S): C = B.flatten() C[S] = 0 C = C.reshape(B.shape) return C
In [8]:
for i in range(1,36): g = subSet(36,i) while (True): try: a = next(g) if err(sub(A,a),LEFT,TOP) == 0: print(sub(A,a)) break except: print("No %s element solution"%i) break
No 1 element solution No 2 element solution No 3 element solution No 4 element solution No 5 element solution
/usr/local/lib/python3.6/dist-packages/ipykernel/ RuntimeWarning: invalid value encountered in true_divide from ipykernel import kernelapp as app /usr/local/lib/python3.6/dist-packages/ipykernel/ RuntimeWarning: invalid value encountered in true_divide app.launch_new_instance()
No 6 element solution No 7 element solution No 8 element solution
In [ ]: