Wednesday, August 28, 2013

Real Root Isolation of Polynomials

Let's first list a few properties about polynomials. Given a polynomial \[    f(x) = \alpha_0 + \alpha_1 x + ... + \alpha_n x^n    \]
the Fundamental Theorem of Algebra tells us that there are $n$ complex roots (counted with multiplicity). A root $r$ has multiplicity $k$ if $f(r) = f'(r) = ... = f^{(k)}(r) = 0 \neq f^{(k+1)}(r)$. A root is simple if it has multiplicity one and a polynomial is square-free if all of its roots are simple. The complex roots must come in conjugate pairs and any rational roots must be of the form (factor of $a_0$) / (factor of $a_n$). Also, the ratio $(-1)^k \alpha_{n-k}/\alpha_n$ gives the sum of the products of all $k$-combinations of roots. Lastly, a result of Cauchy is that the magnitude of any root is bounded by $M = 1 +\max \left| a_i/a_n \right|$.

These properties as well as Descartes' rule of signs (stated below) and synthetic division can all be used to find the roots of a polynomial. Other methods include Newton's, Bisection, Secant, Regula Falsi, etc. Below describes two subdivision methods for isolating all the real roots of a square-free polynomial with real coefficients. They work by determining whether an interval contains zero roots, one root, or multiple roots, and then respectively discarding, storing, or subdividing the interval. The methods can be generalized but not without some effort.

Descartes' Method


Polya Counting cont'd...

In my very first blog post I briefly mentioned Polya counting but I didn't go into the details because I was lazy. Since then I've Tex'd up some notes for my own reference that are accessible in the link below. The proof is only outlined so it is a bit terse but the example problems should be enough to show how the method is used.

Polya Counting Notes


Wednesday, April 3, 2013

Puzzler - Egg Drop

You are given two identical eggs and access to a 100-story building. You want to determine the highest floor from which an egg will not break when dropped. If an egg is dropped and does not break, it is undamaged and can be dropped again. However, once an egg is broken, that’s it for that egg.

If an egg breaks when dropped from floor n, then it would also have broken from any floor above that. If an egg survives a fall, then it will survive any fall shorter than that.

The question is: What strategy should you adopt to minimize the number egg drops it takes to find the solution? (And what is the worst case for the number of drops it will take?)

Solution: Suppose we start by dropping the first egg from floor 50 and it breaks. The only way to determine the correct floor is to then start with floor 1, then 2, 3, ... until our second egg breaks. If the highest floor is 49 then we get 50 total egg drops -- far from optimal. Generalizing, if our first drop is from floor n and the egg breaks, then in the worst case we require n total egg drops. We know that the solution is between 1 and n-1 inclusive. On the other hand, if the egg does not break, we know the solution is between n+1 and 100. Within this range we can use the same technique but instead drop from n-1 floors above floor n. If the egg then breaks, then in the worst case we still total n egg drops. If it doesn't break we continue by taking our third drop to be n-2 floors above the second drop, and so on.. So assuming the first egg does not break, our successive drops are n, 2n-1, 3n-3, 4n-6, ..., n(n+1)/2. -- given by n, n+(n-1), n+(n-1)+(n-2), ..., n+...+1. When the egg finally does break on the kth drop, we drop the second egg incrementally from the floor of the (k-1)th drop up to the floor of the kth drop and in the worst case we get n total drops. Now we simply need to minimize n, which is done by taking the smallest integer n such that n(n+1)/2 > 99. Our solution is thus 14 total drops.


Puzzler - Prisoners and Lightbulbs

There are 100 prisoners in solitary confinement. There's a central living room with one light bulb that is  initially off. No prisoner can see the light bulb from his or her own cell. Everyday, the warden picks a prisoner at random to visit the living room. While there, the prisoner can toggle the bulb if he or she wishes. The prisoner also has the option of asserting that all 100 prisoners have been to the living room by now. If this assertion is false, all 100 prisoners are shot. However, if it is indeed true, all prisoners are set free. So the assertion should only be made if the prisoner is 100% certain of its validity. The prisoners are allowed to get together one night in the courtyard, to discuss a plan. What plan should they agree on, so that eventually, someone will make a correct assertion?

Solution: Here is one sub-optimal solution: Suppose each of the prisoners are assigned a different number from 0-99. Each day is taken modulo 100 and if the prisoner is assigned the same number as the day he turns the light on (or keeps it on) and otherwise he turns the light off (or keeps it off). Whenever a prisoner is called to the room and sees that the light is on, he makes note that the previous days' prisoner has been to the room. Once any of the prisoners has tallied every other prisoner being in the room, he can declare so and they will be set free. For example, say prisoner 50 is called to the room on day 134=34 (mod 100) and sees that the light is on. He concludes that prisoner 33 was in the room the previous day and then proceeds to turn the light off.

There is a simple way to alter this solution to make it way better -- the expected wait goes from ~10,000 to ~1,500 days. Can you think of it?


Saturday, March 16, 2013

Puzzler - Hand Shakes

My wife and I attend a party with four other couples. When we arrive there’s a certain amount of hand-shaking, but no one shakes his own hand and no husband shakes his wife’s hand. I ask the nine other guests how many hands each shook, and I get nine different answers. How many hands did my wife shake?

Solution: There must be exactly one person who shook 0, 1, 2, ..., 8 hands. The person who shook 8 hands shook everyone's except his/her spouses'. So everyone shook atleast one hand except for the spouse, who therefore must have shaken 0 hands. So 8 is paired with 0. Similarly, 7 is paired with 1, 6 with 2, 5 with 3, and 4 is left alone. The only person who didn't have a spouse in the people you asked is your wife, so your wife shook 4 hands.

Friday, March 15, 2013

Discrete Fourier Transform

Suppose we have to multiply two degree polynomials together, $A(x)B(x) = C(x)$ where $A$ and $B$ are given, $C$ is unknown. Let $A,B,C$ be of degree $r,s,r+s$ respectively, and let $n$ be the smallest power of 2 larger than $r+s$.  \[    A(x) = a_0 + a_1x + a_2x^2 + ... + a_{n-1}x^{n-1}     \] \[    B(x) = b_0 + b_1x + b_2x^2 + ... + b_{n-1}x^{n-1}     \] \[    C(x) = c_0 + c_1x + c_2x^2 + ... + c_{2n-2}x^{2n-2}     \]
Note that the coefficients indexed larger than the respective degrees are simply zero. Then we can find our coefficients of $C$ in $\mathcal{O}(n^2)$ by multiplying the coefficients of $A$ and $B$ together in the convoluted sum \[    c_k = \sum_{j=0}^{k}a_j b_{k-j}    \]
But $n=2^{20}$ and $\mathcal{O}(n^2)$ is too costly, so we need a better method. Through the use of the FFT we can find coefficients $c_k$ in $\mathcal{O}(n\log n)$ time by
    1.     (FFT) Evaluate $A$ and $B$ at $n$th roots of unity $\omega^k, \;\;\;\; k=0,...,n-1$
    2.      Find $C(\omega^k) = A(\omega^k) B(\omega^k)$  
    3.      (iFFT) Find $C$ via polynomial interpolation  

Step (1) Suppose we write $A(x) = A_0(x^2) + xA_1(x^2)$ where \[   A_0(x) = a_0 + a_2x + a_4x^2 + ... + a_{n-2}x^{n/2-1} \] \[    A_1(x) = a_1 + a_3x + a_5x^2 + ... + a_{n-1}x^{n/2-1}  \]
To evaluate $A$ at the $n$th roots of unity, we need to evaluate $A_0$ and $A_1$ at the square of each of the roots. But this is the same as evaluating the $(n/2)$th roots of unity at two polynomials of degree $n/2-1$. We immediately get a recursive algorithm (Fast Fourier Transform) with complexity $T(n) = 2T(n/2) + 2n = \mathcal{O}(n\log n).$ We can easily generalize this to work for $n$ as a power of any prime by splitting $A$ into more parts (or any number provided you have its prime factorization).

Step (2) We evaluate $C(\omega^k )$ in $O(n)$ operations.

Step (3) We have formed a linear system $Mc^* = c$, where $M$ is the matrix $(\omega^{kj})$, $c^*$ is the vector of values $ c_k $, and $c$ is the vector of values $C(\omega^{km})$. $M$ is a Vandermonde matrix, so it has an inverse. Thus, we can find $c^{*}$ by applying $M^{-1}$ to both sides, $c^* = M^{-1} c$. It turns out $M^{-1} = (1/n (\omega^{-ij}))$ and applying this matrix to $c$ is the same as step 1, namely the inverse Fourier transform.

Python Code

    from cmath import exp, pi 

    def FFT(X):
        n = len(X)
        w = exp(-2*pi*1j/n)

        if n > 1
            X = FFT(X[::2]) + FFT(X[1::2])

            for k in xrange(n/2):
                xk = X[k]
                X[k] = xk + w**k*X[k+n/2]
                X[k+n/2] = xk - w**k*X[k+n/2]

        return X


Thursday, March 7, 2013

Gravity Toy

Recently found this fun program on the web that lets you create our solar system in 2D. Pretty fun to play with.