Euler #125

The palindromic number 595 is interesting because it can be written as the sum of consecutive squares: 62 + 72 + 82 + 92 + 102 + 112 + 122...

Euler #64

All square roots are periodic when written as continued fractions and can be written in the form...

Method for calculation of square root continued fraction is taken from  http://web.math.princeton.edu/mathlab/jr02fall/Periodicity/mariusjp.pdf

Brute-forcing Eurler #70

Euler's Totient function, φ(n) [sometimes called the phi function], is used to determine the number of positive numbers less than or equal to n which are relatively prime to n. For example, as 1, 2, 4, 5, 7, and 8, are all less than nine and relatively prime to nine, φ(9)=6....

First I was trying to  brute-force the solution using FT of Euler totient function (http://en.wikipedia.org/wiki/Euler%27s_totient_function#Fourier_transform), than through the definition but all these solutions was too long to wait. In the end finished with product formula.

In the process of solution found quite faster recipe from ActiveState on generating prime numbers (http://code.activestate.com/recipes/366178-a-fast-prime-number-list-generator/).

Also used multiprocessing module to speed up solution.

And after...

$ time ./70.py 
(8319823, 8313928, 1.0007090511248113)

real    235m23.873s
user    1644m0.217s
sys     0m36.206s

Euler #80

It is well known that if the square root of a natural number is not an integer, then it is irrational.

from decimal import *
from math import sqrt

getcontext().prec = 105

reduce(int.__add__, 
       map(lambda a: \
           reduce(lambda x,y: int(x)+int(y), 
                  str(a.sqrt()).replace('.', '')[:100], 
                  0), 
           map(lambda n: Decimal(str(n)), 
               map(lambda k: \
                   0 if int(sqrt(k)) ** 2 == k else k, 
                   range(100))
              )
          ), 
       0)

Surely it's cheating ;)

Euler #74

Determine the number of factorial chains that contain exactly sixty non-repeating terms.

Quite straightforward...

Euler #65

Find the sum of digits in the numerator of the 100th convergent of the continued fraction for e.

Euler #112

Euler #73

Euler #71

#!/usr/bin/env python

from fractions import gcd

N = 1000000

for d in xrange(2, N):
    for n in xrange(3*d/7+1, 1, -1):
        if n*7 < 3*d:
            if gcd(n,d) == 1:
                res = (n, d)
                break;

print(res)

Euler #81

#!/usr/bin/env python

def read_matrix(filename):
    res = []
    f = open(filename, 'r')
    for line in f.readlines():
        res.append(map(lambda x: int(x), line.strip().split(',')))
    return res

matrix = read_matrix('matrix.txt')
matrix_len = len(matrix)

for i in xrange(matrix_len):
    for j in xrange(matrix_len):
        if i > 0 and j > 0:
            s1 = matrix[i][j] + matrix[i - 1][j]
            s2 = matrix[i][j] + matrix[i][j - 1]
            matrix[i][j] = min(s1, s2)
        elif i == 0 and j > 0:
            matrix[i][j] = matrix[i][j - 1] + matrix[i][j]
        elif i > 0 and j == 0:
            matrix[i][j] = matrix[i - 1][j] + matrix[i][j]

print matrix[matrix_len - 1][matrix_len - 1]

Euler #47

Euler #46

#!/usr/bin/env python

from libeuler import is_prime, primes
from math import sqrt

if __name__ == '__main__':
    P = primes(100000)
    n = 1
    while True:
        found = False
        for p in P[1:]:
            a2 = n - 0.5*(p - 1)
            if a2 < 0:
                break
            
            if sqrt(a2) == float(int(sqrt(a2))):
                print "%d -> %d + 2*%d^2" % (2*n + 1, p, sqrt(a2))
                found = True
                n += 1
                break

        if not found and not (2*n+1) in P:
            print 2*n + 1
            break

Project Euler #35

let is_prime n =
  if n < 0 then
    false
  else
    let lim = int_of_float (sqrt (float_of_int n)) in
    let rec check_mod i l =
      if i > l then
        true
      else
        if n mod i = 0 then
          false
        else
          true && (check_mod (succ i) l)
    in
    check_mod 2 lim;;

let primes n =
  let rec inner_primes i =
    if i > n then
      []
    else
      if is_prime i then
        [i] @ inner_primes (succ i)
      else
        inner_primes (succ i)
  in
  inner_primes 2;;

let circulations n = 
  let rec circulate str =
    let str_len = String.length str in
    if str_len = 1 then
      [int_of_string n]
    else
      let variant = (String.sub str 1 (str_len - 1)) ^ (String.sub str 0 1) in
      if variant = n then
        [(int_of_string variant)]
      else
        [(int_of_string variant)] @ (circulate variant)
  in
  circulate n;;

(*
let is_circular n primes = 
  let number_string = string_of_int n in
  let variants = circulations number_string in
  List.fold_left 

    (&&) 

    true 

    (List.map 

      (fun x -> List.exists (fun y -> (x = y)) primes) 

      variants);;
*)

let is_circular n =
  let number_string = string_of_int n in
  let variants = circulations number_string in
  List.fold_left (&&) true (List.map (is_prime) variants);;

let primes_till_1e6 = primes 1000000;;
let test_primes = primes_till_1e6;;

Printf.printf "%d\n" 
    (List.fold_left 
       (+) 
       0 
       (List.map 
         (fun x -> if (is_circular x) then 1 else 0 ) 
         test_primes));;

Strange thing: List.exists is quite slow, so when I switched from using List.exists (commented out) to just plain is_prime in is_circular it became possible to get result in a minute. Or maybe I just can't use List.exists properly.

Project Euler #27

#load "nums.cma"

open Big_int;;

type triple = { a: big_int; 
                b: big_int; 
                len: big_int };;

let is_prime n =
  if lt_big_int n zero_big_int then
    false
  else
    let lim = sqrt_big_int n in
    let rec check_mod i l =
      if gt_big_int i  l then
        true
      else
        if eq_big_int (mod_big_int n i) zero_big_int then
          false
        else
          true && (check_mod (succ_big_int i) l)
    in
    check_mod (big_int_of_int 2) lim;;
      
let primes n =
  let rec inner_primes i =
    if gt_big_int i n then
      []
    else
      if is_prime i then
        [i] @ inner_primes (succ_big_int i)
      else
        inner_primes (succ_big_int i)
  in
  inner_primes (big_int_of_int 2);;

let make_form a b = 
  fun n -> 
    add_big_int 
      (mult_big_int n n) 
      (add_big_int (mult_big_int a n) b);;

let sequence_length f =
  let rec advance_sequence n =
    if not (is_prime (f n)) then
      zero_big_int
    else
      add_big_int 
        unit_big_int (advance_sequence (succ_big_int n))
  in
  advance_sequence (big_int_of_int 0);;

let triple_max a b =
  if gt_big_int a.len b.len then
    a
  else
    b;;

let find_longest_sequence i1 i2 n =
  let rec make_lengths_list i =
    if gt_big_int i i2 then
      []
    else
      [{ a=i; b=n; len=sequence_length (make_form i n); }] @ 
      (make_lengths_list (succ_big_int i))
  in
  List.fold_left (triple_max) 
    { a=zero_big_int; b=zero_big_int; len=zero_big_int; } 
    (make_lengths_list i1);;

let bs = primes (big_int_of_int 1000);;
let a1 = big_int_of_string "-999";;
let a2 = big_int_of_string "999";;

let max_ab = List.fold_left (triple_max) 
    { a=zero_big_int; b=zero_big_int; len=zero_big_int; } 
    (List.map (find_longest_sequence a1 a2) bs);;

Printf.printf "%s\n" 
    (string_of_big_int (mult_big_int max_ab.a max_ab.b));;

Project Euler #29

#load "nums.cma";;

open Big_int;;

let rec generate_combinations a1 a2 b1 b2 =
  let rec generate_single_number a b lim =
    if gt_big_int b lim then
      []
    else
      [power_big_int_positive_big_int a b] @ 
      generate_single_number 
        a (succ_big_int b) lim
  in
  if gt_big_int a1 a2 then
    []
  else
    (generate_single_number 
       a1 b1 b2) @ 
    (generate_combinations 
       (succ_big_int a1) a2 b1 b2);;

let rec uniq lst = match lst with
  [] -> []
| hd :: tl when List.exists (fun x -> eq_big_int x hd) tl -> uniq tl
| hd :: tl -> [hd] @ uniq tl;;

let a1 = big_int_of_int 2;;
let a2 = big_int_of_int 100;;
let b1 = big_int_of_int 2;;
let b2 = big_int_of_int 100;;
Printf.printf "%d\n" 
    (List.length 
       (uniq 
   (generate_combinations a1 a2 b1 b2)));;

Quick to come up but not so fast uniq function.

Project Euler #28

let make_spiral n =
  if n mod 2 = 0 then
    raise (Invalid_argument "Size should be odd");
  let res = Array.make_matrix n n 0 in
  let last_value = n * n in
  let current_value = ref 1 in
  let i = ref (n / 2) in
  let j = ref (n / 2) in
  let shift = ref 1 in
  let delta = ref 1 in
  res.(!i).(!j) <- !current_value;
  while !current_value < last_value do
    for k = 1 to !shift do
      if !current_value < last_value then
        begin
          j := !j + !delta;
          current_value := !current_value + 1;
          res.(!i).(!j) <- !current_value;
        end;
    done;
    for k = 1 to !shift do
      if !current_value < last_value then
        begin
          i := !i + !delta;
          current_value := !current_value + 1;
          res.(!i).(!j) <- !current_value;
        end;
    done;
    delta := - !delta;
    shift := !shift + 1;
  done;
  res;;

let sum_diagonals arr = 
  let size = Array.length arr in
  let sum = ref 0 in
  for i = 0 to size - 1 do
    sum := !sum + arr.(i).(i) + arr.(i).(size - 1 - i)
  done;
  !sum - 1;;
    
Printf.printf "%d\n" (sum_diagonals (make_spiral 1001));;

Very bad OCaml ;)

Project Euler #23

(* generates array with values [1..n] *)
let rec range_array n =
  if n == 1 then
    [| n |]
  else
    Array.append (range_array (n - 1)) [| n |];;

(* calculates sum of proper divisors of given number *)
let sum_of_proper_divisors n =
  let limit = n / 2 in
  let rec sum_divisors k l =
    if k > l then
      0
    else
      if n mod k = 0 then
        k + sum_divisors (k + 1) l
      else
        sum_divisors (k + 1) l
  in
  if n = 1 then
    1
  else
    sum_divisors 1 limit;;

(* checks whether given number is abundant*)
let is_abundant n =
  n < sum_of_proper_divisors n;;

(* generates array of abundant number less than given number *)
let find_abundant_numbers n =
  let rec gen_numbers k =
    if k = n then
      [||]
    else
      if is_abundant k then
        Array.append [| k |] (gen_numbers (k + 1))
      else
        gen_numbers (k + 1)
  in
  gen_numbers 12;;

let binary_search (data:'a array) elem =
  let rec iter a b =
    if a = b then -1
    else
      let median = a + (b - a)/2 in
      match data.(median) with
      | value when value = elem -> median
      | value when value < elem -> iter (median + 1) b
      | _                       -> iter a median
  in
  iter 0 (Array.length data);;

let make_possible_sums a =
  let len = (Array.length a) - 1 in
  let res = Array.make_matrix (len + 1) (len + 1) 0 in
  for i = 0 to len do
    for j = i to len do
      res.(i).(j) <- a.(i) + a.(j)
    done
  done;
  res;;

(* sets all elements of a found in s to 0 *)
let refine_array a s =
  let to_delete = Array.make (Array.length a) false in
  let l1 = Array.length s in
  let l2 = Array.length s.(0) in
  for i = 0 to l1 - 1 do
    for j = 0 to l2 - 1 do
      let idx = binary_search a s.(i).(j) in
      if idx != -1 then
        to_delete.(idx) <- true
    done
  done;
  for i = 0 to (Array.length to_delete) - 1 do
    if to_delete.(i) then
      a.(i) <- 0
  done;;

let all_numbers = range_array 28123;;
let abundant_numbers_sums = make_possible_sums (find_abundant_numbers 28123);;
refine_array all_numbers abundant_numbers_sums;;

Printf.printf "%d\n" (Array.fold_left (+) 0 all_numbers);;

Straightforward and very slow but correct solution ;)

Project Euler #25

Quite obvious

Project Euler #21

let proper_divisors_sum n =
  let limit = n / 2 in
  let rec proper_divisors_list n k =
    if k > limit then
      []
    else
      if n mod k = 0 then
        [k] @ (proper_divisors_list n (k + 1))
      else
        proper_divisors_list n (k + 1)
  in
  List.fold_left (+) 0 (proper_divisors_list n 1);;

let amicable_numbers_sum n =
  let rec amicable_numbers k =
    if k > n then
      []
    else
      let b = proper_divisors_sum k in
      let a = proper_divisors_sum b in
      if (a != b) && (a = k) then
        [k] @ (amicable_numbers (k + 1))
      else
        amicable_numbers (k + 1)
  in
  List.fold_left (+) 0 (amicable_numbers 2);;

Printf.printf "%d\n" (amicable_numbers_sum 10000);;

Quite straightforward and surely using lists ;) Actually we can sum without intermediate lists.

Project Euler #24

(* swaps elements of indexes i and j in array a*)
let swap a i j =
  let t = a.(i) in
  a.(i) <- a.(j); a.(j) <- t;;

(* reverses part of array a from index i to the end of array*)
let reverse a i =
  let l = (Array.length a) - 1 in
  for k = 0 to (l - i)/2 do
    swap a (i + k) (l - k)
  done;;

(* generates next permutation of array *)
let make_next_permutation a =
  let j = ref 0 in
  let l = ref 0 in
  for i = 0 to (Array.length a) - 2 do
    if a.(i) < a.(i+1) then
      j.contents <- i
  done;
  for i = 0 to (Array.length a) - 1 do
    if a.(!j) < a.(i) then
      l.contents <- i
  done;
  swap a !j !l; 
  reverse a (!j+1);;

let make_permutations a =
  let t = ref 1 in
  while !t < 1000000 do
    make_next_permutation a; 
    t := (!t + 1); 
  done;;

let t = [| 0; 1; 2; 3; 4; 5; 6; 7; 8; 9 |];;
make_permutations t;;
Array.iter (Printf.printf "%d") t;;
Printf.printf "\n";;

http://en.wikipedia.org/wiki/Permutation#Systematic_generation_of_all_permutations
Весьма императивное решение, с использованием ссылок.