Counting block combinations I: Project Euler 114

Question

Problem statement (For more detailed description (including pictures), please, visit the link):

A row measuring seven units in length has red blocks with a minimum length of three units placed on it, such that any two red blocks (which are allowed to be different lengths) are separated by at least one grey square. There are exactly seventeen ways of doing this.

How many ways can a row measuring fifty units in length be filled?

Code:

from scipy import special
import time
import itertools


def partitionfunc(n,k,l=3):
    '''n is the integer to partition, k is the 
    length of partitions, l is the min partition element size'''
    if k < 1:
        return 0
    if k == 1:
        if n >= l:
            yield (n,)
        return 0
    for i in range(l,n+1):
        for result in partitionfunc(n-i,k-1,i):
            yield (i,)+result


def valid_partitions(p):
    total =  p
    count = 0 #Max. number of tiles that can be placed on a row.
    while True:
        count += 1
        total -=3
        if total <= 3:
            break
        total-=1
    '''Find all the valid partitions with length [1,count] that can be placed on the row with length p'''
    data = []
    for k in range(1,count+1):
        min_part = k*3 
        for n in range(min_part,p+1):
            Allowed = []
            LIST = list(partitionfunc(n,k))
            for b in LIST:
                if sum(b) + (len(b)-1) <= p: 
                    Allowed.append(b)
            data+= Allowed
    return data


def count_permutations(array):
    '''Counts how many possible permutations are there for the particular partition'''
    get_unique_elements = set(array)
    total_length = len(array)
    lengths = [array.count(x) for x in get_unique_elements]
    answer = 1
    for b in lengths:
        answer*= special.comb(total_length, b)
        total_length-= b
    return answer

def calculate_ways(m,n):
    return special.comb(n-m+1,n-2*m+1)


def final(w):
    total_variations = 0
    data = valid_partitions(w)
    for q in data:
        m = len(q)
        remain = w - sum(q) - m + 1
        n = 2*m -1 + remain
        total_variations+= calculate_ways(m,n)*count_permutations(q)
    return int(total_variations+1)

if __name__ == '__main__':
    start = time.time()
    print('Answer: {}'.format(final(50)))
    print(time.time()-start)

I'll explain the reasoning behind the code using example provided by the Euler:

Step 1

First, we find maximum number of tiles that can be placed on the row with length 7. At most only 2 tiles can be placed, hence number equals 2 (call this number count).

Then we iterate through [1,count] and and find all the valid partitions of the number seven:

In our case:

1 tile: (3), (4), (5), (6), (7) (The reason we omit (1) and (2) is because minimum length of red tile is 3)

2 tiles: (3, 3) (Note, that although (3,4) is a partition of 7 too, but it won't work in our case, because it is specified that there must be at least 1 tile gap between two red tiles.

Step 2.

For each partition obtained in step 1, we calculate number of ways the partition can be placed on the row.

For example:

(5) Represents 1 tile with length 5. There are 3 ways to place such tile.

(3,3) Represents 2 tiles both with length 3. There is 1 way to place them on the row.

When we add up the numbers, we end up with 16. We need to add 1, because (I'm curious why) the row that doesn't contain any red tiles is a valid case too. The final answer is 17.

---

I believe that there are a lot of things that can be improved. I'm glad to hear any suggestions!

P.S The problem I would like to specifically point out is the variable/function names. I believe there is a possibility to make them way more descriptive/clear then they are now.

Answer 1

Firstly, I have to admit to not doing too much in the way of thinking up better variable names; in part because the algorithm was a little hard to follow without solving the problem myself.

The spacing is all over the place: I'm not a PEP 8 purist but you have -= with a space before on one line, with a space after on another, and neither on a third. I fixed all the things my IDE complained about, including a redundant import of itertools.

   LIST = list(partitionfunc(n,k))
   for b in LIST:
      ...

Writing LIST to not shadow the builtin list is evil but it's also unecessary: you don't need to listify the generator to iterate over it. That becomes:

   for b in partition(n, k):
      ...

I changed the partition from tuples to lists, since I thought that was more natural in Python for an homogenous sequence of arbitrary length.

There is a loop in valid_partitions:

count = 0 #Max. number of tiles that can be placed on a row.
while True:
    count += 1
    total -=3
    if total <= 3:
        break
    total-=1

I replaced this with:

max_tiles = (p + 1) // 4

I made valid_partitions a generator, and lost some temporary variables like 'data'. I also changed the range() to start from 1 as the comment suggests.

There is a little logic in the function final that I felt was doing too much, so moved it into calculate_ways()

for q in data:
    m = len(q)
    remain = w - sum(q) - m + 1
    n = 2*m -1 + remain
    total_variations+= calculate_ways(m,n)*count_permutations(q)

final() became:

def total_variations(w):
    return sum(calculate_ways(q, w) * count_permutations(q) for q in valid_partitions(w)) + 1

Mayble the one liner is a little too dense, and one should keep the for loop.

calculate_ways() is now a bit opaque and could probably use tidying but at least it's all in one place now.

I added MIN_SIZE = 3. It's probably overkill to pass 3 around as a parameter but probably worth flagging as a magic number.

from scipy import special
import time

MIN_SIZE = 3


def partition(n, k, min=MIN_SIZE):
    """n is the integer to partition, k is the
    length of partitions, min is the min partition element size"""
    if k == 1:
        if n >= min:
            yield [n]
        return 0
    for i in range(min, n+1):
        for result in partition(n - i, k - 1, i):
            yield [i] + result


def valid_partitions(p):
    """Find all the valid partitions with length [1,max_tiles] that can be placed on the row with length p"""
    max_tiles = (p + 1) // 4
    for k in range(1, max_tiles+1):
        min_part = k * MIN_SIZE
        for n in range(min_part, p+1):
            for b in partition(n, k):
                if sum(b) + (len(b)-1) <= p:
                    yield b


def count_permutations(array):
    """Counts how many possible permutations are there for the particular partition"""
    get_unique_elements = set(array)
    lengths = [array.count(x) for x in get_unique_elements]
    total_length = len(array)
    answer = 1
    for b in lengths:
        answer *= special.comb(total_length, b, exact=True)
        total_length -= b
    return answer


def calculate_ways(q, w):
    m = len(q)
    remain = w - sum(q) - m + 1
    n = 2*m - 1 + remain
    return special.comb(n - m + 1, n - 2*m + 1, exact=True)


def total_variations(w):
    return sum(calculate_ways(q, w) * count_permutations(q) 

if __name__ == '__main__':
    start = time.time()
    print('Answer: {}'.format(total_variations(50)))
    print(time.time()-start)

As always with Project Euler, there may be some deep mathematical insight that can optimise the algorithm. And there still are a lot of single letter variables; it's a start, though.

Answer 2

Here is a better algorithm to solve the problem.

Let \$c(n, m)\$ be the number of ways to cover \$n\$ units with red blocks of minimum length \$m\$. Let \$c_r(n, m)\$ and \$c_b(n, m)\$ be the respective number of covers that ends with a red or black unit. We have $$ c(n, m) = c_r(n, m) + c_b(n, m)\label{f1}\tag{1} $$

Since any length-\$(n-1)\$ cover can be extended with a black unit to yield a length \$n\$ cover, we have $$ c_b(n,m)=c(n-1,m)\label{f2}\tag{2} $$

For any length-\$n\$ cover that ends with a red unit, it either extends a red-ending, length-\$(n-1)\$ cover, or adds a length-\$m\$ red block to a black-ending, length-\$(n-m)\$ cover. Therefore $$ \begin{eqnarray} c_r(n,m) & = & c_r(n-1,m)+c_b(n-m,m) \\ & = & c(n - 1, m) - c_b(n - 1, m) + c_b(n - m, m) \\ & = & c(n - 1, m) - c(n - 2, m) + c(n - m - 1, m) \label{f3}\tag{3} \end{eqnarray} $$

Substituting (\ref{f2}) and (\ref{f3}) into (\ref{f1}) yields: $$ c(n, m) = 2 \cdot c(n - 1, m) - c(n - 2, m) + c(n - m - 1, m) \label{f4} \tag{4} $$

The sequence \$\{c(i, m)\}_{i=0}^n\$ can now be generated from the linear recurrence (\ref{f4}) with initial values \$c(-1,m)=c(0,m)=\ldots=c(m-1,m)=1\$.

While a linear algorithm solves the original problem with ease, it is not efficient enough to solve the extended version where the input \$n\$ can go up to \$10^{18}\$. To speed up computation using the linear recurrence, we define a length-\$(m+1)\$ column vector $$ C_i= \begin{pmatrix} c(i+m,m) & c(i+m-1,m) & \ldots & c(i,m) \end{pmatrix} ^T $$ and a \$(m+1)\times (m+1)\$ coefficient matrix $$ A= \begin{pmatrix} 2 & -1 & 0 & \cdots & 0 & 1 \\ 1 & 0 & 0 & \cdots & 0 & 0 \\ 0 & 1 & 0 & \cdots & 0 & 0 \\ 0 & 0 & 1 & \cdots & 0 & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots & \vdots \\ 0 & 0 & 0 & \cdots & 1 & 0 \\ \end{pmatrix} $$ Then we can see \$C_i=AC_{i-1}=\ldots=A^{i+1}C_{-1}\$, where \$ C_{-1}= \begin{pmatrix} 1 & 1 & \ldots & 1 \end{pmatrix} ^T \$. Therefore \$c(n,m)\$ can be obtained by extracting the first element of $$C_{n-m}=A^{n-m+1}C_{-1}$$ Since the power operation can be computed in \$\Theta(\log(n-m+1))\$ matrix multiplications and the multiplication time complexity is \$O(m^3)\$, the entire algorithm is \$O(m^3\log(n-m))\$ which is sufficient to solve the extended problem.