Archive of articles classified as Uncategorized

Certificate Verification for Costas Arrays

2014.04.25

If you never heard of Costas arrays before you may want to watch this video.

First of all, I want to be clear with the way that I’m saying “Costas arrays” but I’m not really working with arrays. It’s more common to think of them as permutations.

This is the definition that I found the easiest to understand:

Definition Let f\colon [n] \rightarrow [n]. f has the Costas property iff

    \[f(i+k)-f(i) = f(j+k) - f(j) \implies i = j\]

for all i, j, k \in [n].

This is all we need. Here’s code.

import itertools

def is_costas(f):
    n = len(f)
    return all(f[i+k]-f[i] != f[j+k]-f[j]
               for i, j in itertools.combinations(range(n), 2)
               for k in range(1, n-max(i, j)))

I was a bit more judicious instead of checking for every i, j and k like the definition says. There’s symmetry here and I tried to explore it by assigning i and j to elements of 2-combinations.

The definition leaves the restriction of k implicit, we need to make it explicit thus we keep k limited with max(i, j).

And lastly, you need to be careful with the fact that, in contrast with the definition, Python’s indexing is zero-based.

The decision problem of determining if there is at least one permutation of a given length n satisfying the Costas property is sometimes called CAP (Costas Array Problem). I just showed you an O(n3) certificate verification algorithm for this problem, so that means that CAP is in NP.

As far as I know it’s an open question if it’s NP-complete.

You can use this code to search for Costas arrays of a given n by checking permutations. In 6 minutes I could find all Costas arrays up to n=10. Using the same idea in a constraint programming solver I found all Costas arrays up to n=11 in 1 minute.

Computing the Binomial Coefficient

2014.04.22

Looking at the recursive definition looks promising.

    \[ \binom nk = \begin{cases} 1 & k = 0 \text{ or } k = n\\ \binom{n-1}k + \binom{n-1}{k-1} & 0 < k < n \end{cases} \]

It’s easy to code an almost direct version of it.

import functools

def binomial_coefficient(n, k):
    return _bc(n, min(k, n-k))

@functools.lru_cache(maxsize=None)
def _bc(n, k):
    if k in (0, n):
        return 1
    return _bc(n-1, k) + _bc(n-1, k-1)

As an example this image shows all the values that are called to calculate binomial_coefficient(7, 4).

There’s two nontrivial things in this implementation that isn’t in the original definition. First I’m taking advantage of the symmetry:

    \[ \binom nk = \binom n{n-k} \]

So this is why I start computing the binomial coefficient with the second argument being the minimum among k and n-k.

The second thing is that I made _bc a memoized function. For this case python3’s functools.lru_cache decorator works great.

Want more buzzwords? Of course you do. This recursive implementation with the use of memoization characterizes dynamic programming with a top-down approach.

We can push a little further with the symmetry if we just keep using min(k, n-k) at every step.

@functools.lru_cache(maxsize=None)
def binomial_coefficient(n, k):
    min_k = min(k, n-k)
    if min_k == 0:
        return 1
    return (binomial_coefficient(n-1, min_k) + 
            binomial_coefficient(n-1, min_k-1))

Then we won’t need to ever check if k is n. Also no need for a separate _bc function, the idea of having it just to wrap that min was kind of silly anyway.

But hey, we can do considerably better if we change our approach a bit. The problem here is that our recursive formula gives us only a partial order. It’d be nice if we had a total order so that we can have an expression in tail position.

Conveniently the binomial coefficient delivers with this identity:

    \[ \binom nk = \frac{n}{k} \binom {n-1}{k-1} \]

Now let’s revise our definition:

    \[ \binom nk = \begin{cases} 1 & k = 0 \text{ or } k = n\\ \frac{n}{k} \binom {n-1}{k-1} & 0 < k < n \end{cases} \]

The best part is that with this we’ll be able to eliminate a lot of calculations.

Notice that this has nothing to do with top-down versus bottom-up approach. The important thing is that now we can define a tail-recursive function. Or in other words, because we’re coding in python, just a simple accumulator loop will do.

def binomial_coefficient(n, k):
    min_k = min(k, n-k)
    acc = 1
    for i in range(1, min_k+1):
        acc = (acc*(n-i+1))//i
    return acc 

With a little attention to keep the accumulator as an integer, we have in essence implemented the following:

    \[ \binom nk = \prod_{i=1}^k \frac{n-i+1}{i} \]

Which is the way I like best to define the binomial coefficient.