Project Euler – Problem 145

The problem: Problem 145

Analysis

Suprisingly enough the brute-force solution works. It takes a couple of minutes but it gets the job done. Suppose we have a function called is_reversible then the following algorithm can be used to count all the reversible numbers below one-billion.

count = 0
for i = 1 upto 1000000000 do
  if (is_reversible(i)) then count++
endfor
output count

Implementation

Uses the OCaml programming language.

(* a left fold over the digits of a positive integer *)
let rec nfoldl f v = function
    0 -> v
  | n -> nfoldl f (f v (n mod 10)) (n / 10)

let reverse = nfoldl (fun n d -> (n * 10) + d) 0

let odd n = n mod 2 = 1
let allodd = nfoldl (fun b d -> b && odd d) true

let is_reversible n =
  n mod 10 <> 0 && allodd (n + reverse n)

let naive_count n =
  let cnt = ref 0 in
    for i = 1 to n-1 do
      if is_reversible i then incr cnt;
    done;
    !cnt

The above code is elegant and it solves the problem but I’m not satisfied with its speed. In that regard I decided to improve on the algorithm a bit.

More Analysis

Notice that since the sum of a number, n, and its reverse, reverse(n) must yield a number with all digits odd, it must be the case that n cannot begin and end with digits of the same parity (i.e. both even or both odd). Hence, a likely candidate must lie in one of these sets A = \{n : n \,\text{begins with}\, 2, 4, 6, 8 \,\text{and}\, n \,\text{ends with}\, 1, 3, 5, 7, 9\} or B = \{n : n \,\text{begins with}\, 1, 3, 5, 7, 9 \,\text{and}\, n \,\text{ends with}\, 2, 4, 6, 8\}. Now, |A| = |B| = (4*5)*10^7 = 200000000, i.e. 200 million. Furthermore, if a number in A is reversible then its reverse is in B and is obviously reversible as well. This means we only need to generate the elements and count those in A that are reversible. Our final answer must then be multiplied by 2 to account for those in B.

A Better Implementation

Uses the Haskell programming language. I switched languages when I realised that the list code I needed for this to work was not implemented using tail-recursion in the OCaml standard library. I did rewrite those pieces using tail recursion but in the end I still decided to switch languages just for the fun of it. Actually, the Haskell version below isn’t much better. Its still slow. But very elegant :D. At least I can use it for testing correctness on small cases since the solution is so simple and easy to read.

nfoldl _ v 0 = v
nfoldl f v n = nfoldl f (f v (n `mod` 10)) (n `div` 10)

nreverse = nfoldl (\n d -> (n * 10) + d) 0
allodd = nfoldl (\b d -> b && odd d) True
isReversible n = n `mod` 10 /= 0 && allodd (n + nreverse n)

countDigits = (+1) . floor . logBase 10 . fromIntegral

genCandidates n =
  let gen m i len = if i > len then
                      if m >= n then [] else [m]
                    else
                      if i == 1 then
                        concatMap (\d -> gen (m*10+d) (i+1) len) [1,3,5,7,9]
                      else if i == len then
                        concatMap (\d -> gen (m*10+d) (i+1) len) [2,4,6,8]
                      else
                        concatMap (\d -> gen (m*10+d) (i+1) len) [0..9]
  in concatMap (gen 0 1) [1..countDigits n]

betterCount = (*2) . length . (filter isReversible) . genCandidates

Finally, I bit my tongue and decided to code up a solution in everyone’s favorite system programming language, C.

#include <stdio.h>
#include <math.h>

#define FALSE 0
#define TRUE  1

#define even(n) ((n)%2==0)

int reverse(int n) {
  int m = 0;
  while (n > 0) { m = m*10 + (n%10); n /= 10; }
  return m;
}

int allodd(int n) {
  while (n > 0) {
    if (even(n%10)) return FALSE;
    n /= 10;
  }
  return TRUE;
}

#define is_reversible(n) ((n)%10>0 && allodd((n) + reverse(n)))
#define count_digits(n) ((int)(floor(log10(n) + 1)))

int __count;

void __better_count(int m, int i, int len, int n) {
  int j;
  if (i > len) {
    if (m < n && is_reversible(m)) __count += 2;
  } else {
    if (i == 1) {
      for (j = 1; j <= 9; j+=2) __better_count(m*10+j, i+1, len, n);
    } else if (i == len) {
      for (j = 2; j <= 8; j+=2) __better_count(m*10+j, i+1, len, n);
    } else {
      for (j = 0; j <= 9; j++) __better_count(m*10+j, i+1, len, n);
    }
  }
}

int better_count(int n) {
  int len;
  
  __count = 0;
  for (len = 2; len <= count_digits(n); len++) __better_count(0, 1, len, n);
  return __count;
}

int main(void) {
  printf("%d\n", better_count(1000000000));
  return 0;
}

In the end, the problem was easy and I had a lot of fun playing around with the various languages. I learned quite a lot, so it certainly wasn’t a waste of time.

Getting Previous Play Whe Results using Python

Quick Update (17th January, 2012):

  • I have created an extension for Google Chrome that allows you to view the latest Play Whe results. Check it out here: Play Whe Results Viewer
  • You can view even more results and statistics if you follow this link: Play Whe Smarter
  • And finally, you can like the Play Whe Smarter page on facebook to keep up to date with the latest tricks and strategies for winning at Play Whe. Like us here: Play Whe Smarter

Play Whe is a numbers game brought to Trinidad by Chinese immigrants. Its usually played by persons who are influenced by intuition, superstition, dreams and caprice. The rules to the game are pretty simple. You have 36 numbers, 1 through 36. When a number is drawn the person(s) who selected that number wins a sum of money; $24 to every $1 that was played. Everyday (except Sundays and public holidays) there are two drawings, one at 1:00pm and the other at 6:30pm.

The NLCB provides a Play Whe Statistics form that can be used to query the numbers that were drawn on any given day way back to the first draw on July 4th, 1994.

The Play Whe Statistics Form

That is great.

However, the results are not returned and can’t be returned in a program friendly format. json or XML would have been nice. Even if the HTML it returned was properly marked-up I would have been thankful, but that is not the case either. I am still lucky though, because there seems to be some consistency to the data and so it can be extracted using regular expressions.

Let’s assume we want to get a “Monthly Summary” for April, 2011. Well, we will need to figure out how the Play Whe Statistics form is querying the server. If we fire up Live HTTP Headers in Firefox, we will be able to determine the POST request that’s being sent to the server when a query is made for the “Monthly Summary”. [By the way, I know its a POST request due to the fact that a GET request would have passed the arguments in the URL.]

The three highlighted regions in the dialog below gives us all the information we need.

Play Whe Query for April 2011 Monthly Summary

With this information we can now write some Python code to make the query.

import urllib
sel = '/search/pwq/countdateCash.php'
host = 'nlcb.co.tt'
params = urllib.urlencode({
  'month': 'Apr',
  'year': 11
})
data = urllib.urlopen('http://' + host + sel, params).read()
print data

It returns some messed up HTML but with a little regular expression magic we can extract some useful information:

import re
pattern = r"(\d{2})-Apr-11: Draw # (\d+) : (AM|PM)'s draw  was (\d+)"
results = re.findall(pattern, data)
# sorts by draw number and then prints the results
for res in sorted(results, key=lambda r: r[1]):
  print res

Here’s a snippet of the data that’s returned:

('01', '10293', 'AM', '19') # the draw at 1:00pm
('01', '10294', 'PM', '23') # the draw at 6:30pm
('02', '10295', 'AM', '3')
('02', '10296', 'PM', '35')
('04', '10297', 'AM', '1')
('04', '10298', 'PM', '33')
('05', '10299', 'AM', '4')
('05', '10300', 'PM', '32')
('06', '10301', 'AM', '28')
('06', '10302', 'PM', '20')
('07', '10303', 'AM', '29')
('07', '10304', 'PM', '17')
('08', '10305', 'AM', '15')
('08', '10306', 'PM', '13')
('09', '10307', 'AM', '5')
('09', '10308', 'PM', '32')
...

That’s it. You now have a programmatic way for getting your current and past Play Whe results. Enjoy!

Project Euler – Trinidad and Tobago Standings

Below is an image of the current standings in projecteuler.net as of 21/11/2010 for Trinidad and Tobago.

Project Euler standings as of 21/11/2010

As you can see I’m 2nd in the standings. Damn you wallygold. Lol. I need to solve 10 more problems to reach Level 3. Hopefully, I will have time to do it before the year ends.

A note about the languages I use. Currently I’m using racket (a scheme based language) to solve the problems. However, in the past I’ve solved some of the problems using C, Java, Python, Haskell, Scala, LISP and of course the trusty old pencil and paper. I guess I tend to utilize whatever language I am currently into.

Project Euler – Problem 71

The problem: Problem 71

Analysis

Let d be some integer larger than 1, F(d) = \{ \frac{a}{b} : 1 \leq a < b \leq d\} and F_0(d) = \{0\} \cup F(d). Let prev_d be a function from F(d) to F_0(d) defined as follows:

prev_d(t) = s,

where s is the largest fraction in F_0(d) that is less than t. As the example in the problem illustrates, prev_8(\frac{3}{7}) = \frac{2}{5}. We need to find prev_{1000000}(\frac{3}{7}).

Let t \in F(d) and n be an integer between 1 and d inclusive. We now define p_d(t,n) as follows:

p_d(t,n) = s,

where s is the largest fraction with a denominator of n in F_0(d) that is less than t.

Clearly then,

prev_d(t) = max_{2 \leq n \leq d}\{ p_d(t,n) \}.

Thus it suffices to show how to compute p_d(t,n). The implementation which will be given below uses a modified binary search to compute p_d(t,n). The running time is O(lg\,n). Therefore, it takes

lg\,2 + lg\,3 + \cdots + lg\,d = lg\,d! = O(d \cdot lg\,d)

time to compute prev_d(t) for any t.

Implementation of p_d(t,n)

Uses the racket programming language.

;; finds the largest fraction with a denominator of n
;; that is less than t
(define (p t n) ;; assumes n <= d
  (let ([a (numerator t)]
        [b (denominator t)])
    (/ (let loop ([l 0] [h (- n 1)])
         (if (> l h)
             0
             (let ([m (quotient (+ l h) 2)])
               (if (< (* m b) (* a n))
                   (max m (loop (+ m 1) h))
                   (loop l (- m 1))))))
       n)))

An example showing how to use p to compute prev_8(\frac{3}{7}).

> (p (/ 3 7) 2)
0
> (p (/ 3 7) 3)
1/3
> (p (/ 3 7) 4)
1/4
> (p (/ 3 7) 5)
2/5
> (p (/ 3 7) 6)
1/3
> (p (/ 3 7) 7)
2/7
> (p (/ 3 7) 8)
3/8

Since \frac{2}{5} is the maximum of the values that were computed, we can conclude that prev_8(\frac{3}{7}) = \frac{2}{5}.

Finally, we use the same idea to get prev_{1000000}(\frac{3}{7}). The full implementation takes about 7 seconds to compute the answer.

SICP – A Proof for Exercise 1.13

Exercise 1.13 in Section 1.2.2 asks for a proof that Fib(n) is the closest integer to {\phi^n}/{\sqrt{5}}, where \phi=(1+\sqrt{5})/2.

Proof Idea: To show that Fib(n) is the closest integer to {\phi^n}/{\sqrt{5}} we simply need to show that |Fib(n) - \phi^n/\sqrt{5}| < 1/2 for all n \geq 0. We will establish the result by proving two lemmas.

Lemma 1: Fib(n) = (\phi^n - \psi^n)/\sqrt{5} for all n \geq 0.

Lemma 2: |\psi^n|/\sqrt{5} < 1/2 for all n \geq 0.

Notice that if Lemma 1 and Lemma 2 are true it will follow that |Fib(n) - \phi^n/\sqrt{5}| = |(\phi^n - \psi^n)/\sqrt{5} - \phi^n/\sqrt{5}| = |\psi^n|/\sqrt{5} < 1/2 for all n \geq 0.

It remains to be shown that Lemma 1 and Lemma 2 are indeed true.

Proof of Lemma 1: We will prove the result by Induction on n. When n = 0 and n = 1, the result follows trivially. Let n \geq 0 and suppose the result holds for n and n+1. Consider, Fib(n+2). By definition, Fib(n+2) = Fib(n+1) + Fib(n). Then by using the Induction Hypothesis, it follows that

Fib(n+2) = \frac{\phi^{n+1} - \psi^{n+1}}{\sqrt{5}} + \frac{\phi^n - \psi^n}{\sqrt{5}},

which can be further reduced to

\frac{\phi^n(\phi + 1) - \psi^n(\psi + 1)}{\sqrt{5}}.

Since \phi^2 = \phi + 1 and \psi^2 = \psi + 1, the result holds for Fib(n+2). Therefore, by Induction, Fib(n) = (\phi^n - \psi^n)/\sqrt{5} for all n \geq 0.

Q.E.D.

Proof of Lemma 2: \sqrt{5} > 2 implies that \frac{1}{\sqrt{5}} < \frac{1}{2}, i.e. \frac{|\psi^0|}{\sqrt{5}} < \frac{1}{2}. Also, \sqrt{5} < 3 implies that \sqrt{5} - 1 < 2 which implies that \frac{\sqrt{5} - 1}{2} < 1. So, |\psi^n| \cdot \frac{\sqrt{5} - 1}{2} < |\psi^n| for any n \geq 0. But, |\psi^n| \cdot \frac{\sqrt{5} - 1}{2} = |\psi^{n+1}|. Hence, |\psi^{n+1}| < |\psi^n| for all n \geq 0 and it trivially follows that \frac{|\psi^{n+1}|}{\sqrt{5}} < \frac{|\psi^n|}{\sqrt{5}} for all n \geq 0. Since \frac{|\psi^0|}{\sqrt{5}} < \frac{1}{2} and \frac{|\psi^n|}{\sqrt{5}} decreases as n increases, we can conclude that \frac{|\psi^n|}{\sqrt{5}} < \frac{1}{2} for all n \geq 0.

Q.E.D.

As planned, this completes the proof that Fib(n) is the closest integer to {\phi^n}/{\sqrt{5}}.

Project SICP

SICP book cover

Structure and Interpretation of Computer Programs” is by far the greatest programming book I have ever read. The majority of other programming books don’t teach you about the big ideas in programming or about managing complexity in large software systems. These books spend most of the time teaching you about the peculiarities of the syntax of a particular programming language. And when they do go beyond that the problems they tackle are so fucking mundane and mind numbing. It’s no wonder a majority of the people who program don’t see the beauty in it.

What’s so awesome about this book?

Well firstly, it uses Scheme as the programming language for discussing the various ideas. It never really teaches Scheme because there isn’t much syntactic complexity to the language, you just learn it as you progress through the book. This is great because it means that Scheme doesn’t get in your way when the big ideas are being taught. And mind you this book is all about the big ideas. Take chess as an example. Anyone can be taught in 5 minutes (if so much) the rules of chess and how each of the individual pieces move on the board. But the complexity of chess comes from the interaction of the pieces and the various ways in which the moves can be combined. So to become a master of chess you don’t spend your time studying the rules, you spend your time studying tactics, strategies, openings and end game scenarios. These are the big picture ideas in chess. Similarly, anyone can be taught in 1 hour (if so much) the rules of Scheme and how most of its primitive constructs and special forms work. But the complexity of Scheme comes from the means it provides for forming combinations and making abstractions. So to become a master programmer (not Scheme programmer but programmer in general) it is clear to see that we should spend most of our time understanding how to combine elements (procedures and data) to form compound elements and more importantly understanding techniques for abstraction which allow us to control complexity in both small and large software systems. SICP devotes each printed page to get this knowledge across and because of this the complexity of the programs you begin studying from the start are orders of magnitude more interesting and intellectually enriching.

Secondly, the book is well written and the authors don’t treat you like a “Dummy” or an “Idiot”. This book expands your mind and changes the way you think about programming. They are not afraid to tackle complex issues. Maybe you won’t understand all of the ideas in a 1st, 2nd or 3rd reading but because you’ve been exposed to it your mind is forever changed. They also assume you are willing to spend some time digesting and understanding this material. Its not a become a master programmer in 21 days type of book. Who the fuck can master anything in 21 days?

Anyway here are some of the topics covered in the book: functional programming, object oriented programming, logic programming, programming with constraints, lazy programming, nondeterministic computing, interpreters and compilers, data-directed programming.

And here are some of the interesting algorithms, data structures and programs that are used as examples: algorithms for numerical computing, interval arithmetic, a picture manipulation language, symbolic differentiation, stacks, queues, sets, trees, symbolic algebra, a simulator for digital circuits, a scheme interpreter, an interpreter with lazy evaluation, an interpreter for nondeterministic computing, an interpreter for logic programming, a scheme compiler.

Seems like a lot but its what any good book on programming should be about.

Now onto the reason I started this post.

Project SICP is my attempt to program solutions to every problem in the book. I have already completed a vast majority of them but I want to organize them in a structured way and solve the remaining ones. And yes. I am completely aware of the successful attempts made by other people (for example: Eli Bendersky) to do exactly this but I learn by doing, not by seeing how other people did it.

Check out my SICP repository on github if you want to view my solutions to these problems.

A ray tracer in Racket

Racket is a very cool programming language based on the scheme/lisp family of languages. For the past couple of months I’ve been having lots of fun implementing and re-implementing various programs in Racket. My latest project is a ray tracer. A ray tracer is a program for generating realistic images. The code I’m writing for the ray tracer is based on C++ code taken from the book “Realistic Ray Tracing (2nd edition) by Peter Shirley and R. Keith Morley“. I’m now on the first section of the third chapter which is the part where they show you how to implement a single sample ray tracer with an orthographic camera. Its just beautiful to see how nicely the ideas in the book seamlessly translate into readable Racket code.

By the way, I plan to put the source code for this project up on github. So keep checking there for the code if you’re interested.

In the mean time, here’s my first image ever generated by the program:

And, yes I know its not looking realistic but that will change as I add more features.

Update: I added the repository to github. Check it out at http://github.com/dwayne/rt.