Problem 69 Analytic - 100th post

So, problem 69 took a long time to brute force, but as is the case with most totient problems as I am starting to discover, its much easier to think about properties of the totient function than have to actually calculate it. Turns out, it is rather simple to find what the answer to problem 69 must be - it has to be primorial, the product of the first n primes for some n. Anyway, this post, which frees up J yet again, is my 100th on this blog. yay.

Problem 70 - Perl

I hadn't touched Perl except in a very, very brief rushed learning session at the end of the hackathon where I solved the first 20 problems. 50 problems later, I have finally returned to it. And though it is capable and all, boy is it a weird language. My solution to this problem is definitely smarter than a brute force (brute forcing the totient of all numbers below 10^7 would not be fun/fast), but still not as smart as some of the solutions on the Euler problem. My solution works using the fact that the solution must have very few prime factors, but I do an "exhaustive" search over numbers with "few" prime factors, instead of choosing a much nicer, more restrictive set of conditions to search that other people found. Anyway, the solution runs in about 30s, and this is essentially my birthday present to me, today is my birthday, and 70 problems are now done.

sub isPrime {
    my (@ops) = @_;
    my $x = $ops[0];
    if ($x < 2 || $x % 2 == 0) {
        return x == 2
    }
    for (my $i = 3; $i * $i <= $x; $i += 2) {
        if ($x % $i == 0) {
            return 0
        }
    }
    return 1
}

sub nextPrime {
    my (@ops) = @_;
    my $x = $ops[0];
    my $y = $x + 1;
    while (not(isPrime($y))) {
        $y += 1;
    }
    return $y;
}

sub lstprod {
    my (@ops) = @_;
    my $r = 1;
    for my $x (@ops) {
        $r *= $x;
    }
    return $r
}

sub totientIsPermut {
    my (@ops) = @_;
    my $n = lstprod(@ops);
    my $prev = 0;
    my $t = $n;
    #print "-$n, @ops\n";
    for my $p (@ops) {
        if ($p != $prev) {
            $t = $t / $p;
            $t = ($p - 1) * $t;
            $prev = $p
        }
    }
    my @tarr = sort split('', "$t");
    my @narr = sort split('', "$n");
    my $ts   = join('',@tarr);
    my $ns   = join('',@narr);
    if ($ts eq $ns) {
        $rat = $n / $t;
    }
    return ($ts eq $ns) 
}

sub update {
    my (@ops) = @_;
    my $size = @ops;
    for (my $x = 0; $x < $size; $x += 1) {
        my $next = nextPrime($ops[$x]);
        my @slice = @ops[($x+1).. ($size-1)];
        my $nprod = $next * lstprod(@slice);
        if ($nprod < 10000000 && $nprod != 0) {
            @slice = @ops[($x) .. ($size-1)];
            $slice[$x] = $next;
            return @slice
        }
    }
    @slice = @ops[1 .. ($size-1)];
    return update(@slice)
}
my @nums = (2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2);
my $done = 0;
my $minrat = 2;
$rat = 2;
my $ans  = 0;
while (not($done)) {
    if (totientIsPermut(@nums)) {
        if ($rat < $minrat) {
            $minrat = $rat;
            $n = lstprod(@nums);
            if ($n > 10000000) {
                $done = 1
            } else {
                $ans = lstprod(@nums)
            }
        }
    }
    @nums = update(@nums)
}
print "$ans\n"

Problem 68 - Fantom

Fantom is yet another language that runs on the JVM, though the designers of Fantom don't stop there - it can also compile to code for the .NET framework or javascript. Anyway, its a fairly simple language with a common Scala/Java/Gosu/C/C#/Squirrel/Ceylon/How many curly-braced languages are there again? syntax. I decided to use it on problem 68 as the problem is not all that hard using a nice object oriented approach and ordering the different kinds of rings with lexicographic permutations (I think I have written that NthLP method in at least 4 different languages by now). Runs in about 1.8s on my machine.

class Ring {
    Int[] elems
    
    static Int fac(Int x) {
        if (x < 2) {
            return 1
        } else {
            return x * fac(x-1)
        }
    }

    Int[] NthLP(Int N, Int[] p, Int[] rem) {
        if (rem.size == 1) {
            p.addAll(rem)
            return p
        }
        Int k := fac(rem.size - 1)
        Int timesDivided := N / k
        p.add(rem[timesDivided])
        rem.removeAt(timesDivided)
        return NthLP(N % k, p, rem) 
    }   
    
    Bool isMagic() {
        Int target := elems[4] + elems[9] + elems[5]
        for (i := 0; i < 4; ++i) {
            if (elems[i] + elems[5+i] + elems[6+i] != target) {
                return false
            }
        }
        return true
    }
    Str description() {
        Int minIndex := 0
        Str desc := ""
        5.times |i| { if (elems[i] < elems[minIndex]) { 
                          minIndex = i} }
        5.times |i| { 
            Int j := (i + minIndex) % 5
            [elems[j], elems[5+j], 
        elems[(j == 4 ? 5 : 6 + j)]].each |d| { desc += "$d" }
        }
        return desc
    }
    new make(Int N) {
        //construct a Ring with elements
        //in the formation of the Nth Lexicographic permutation
        this.elems = NthLP(N, [,], [10,1,2,3,4,5,6,7,8,9])
    }
}

class E68 {

    static Void main() {
        Str maxstr := ""
        for (i := 0; i < Ring.fac(10); ++i) {
            Ring r := Ring(i)
            Str s := ""
            if (r.isMagic()) {
                s = r.description()
                if (s.size() > 16) {
                    break
                } else if (s.compare(maxstr) > 0) {
                    maxstr = s
                }
            }
        }
        echo(maxstr)
    }
}
        

Problem 19 redo - F#

F# is a language based on OCaml but which also does other stuff and such...I didn't have a very exciting experience with it for this problem, but in my delaying solving problem 68, I decided to free up some of my more capable language. Even though I have only used Perl for about 10 minutes in the past, it seemed like a language I might want to use again, so I ported my solution to problem 19 from Perl to F#. Solution runs in about 75ms on my machine:

let rec answer day sundays mo y =
    if y = 101  then sundays else
    if mo = 12  then (answer day sundays 0 (y + 1)) else
        let newsun = (if (day = 0) && not (y = 0) then 
                      sundays + 1 else sundays)
        if (mo = 0 || mo = 2 || mo = 4 || mo = 6 || mo = 7
         || mo = 9 || mo = 11) then
            (answer ((day + 3) % 7) newsun (mo + 1) y) else
        if (mo = 1) then
            (answer (if ((y % 4 = 0) && not (y = 0)) then 
            ((day + 1) % 7) else day) newsun (mo + 1) y)
        else
            (answer ((day + 2) % 7) newsun (mo + 1) y)

printfn "%d" (answer 1 0 0 0)

Problem 69 - J

I decided to skip problem 68, at least for a little bit - it involves magic and polygons and all sorts of weird stuff. So, instead I decided to solve problem 69. This problem is rather simple in J, where a reasonably efficient implementation of Euler's totient function is part of the language: (5 p: y) is the totient of y. So, with that, writing a solution was not hard at all. Now, I am sure there are better solutions, as mine runs in about 20 seconds, but what works works.

f =: 3 : 'y % (5 p: y)'
lst =: f i.1000000
lst i.{.\:~lst

Problem 35 - ChucK

ChucK is a programming language that was never intended to be used for Project Euler. It is a"Strongly-timed, concurrent, and On-The-Fly Music Programming Language." Obviously, I did not have to make music to solve this problem, so I used few of ChucK's features. However, it is a programming language, so it is fair game for this challenge. My solution in ChucK to problem 35 has, unfortunately, 25 times as many lines as my previous solution, which was in J (I wanted to use J for problem 69, as will be clear in my next post). Nonetheless, ChucK was a reasonable enough language, even if the chuck operator => is rather odd. Solution runs in about 25s on my machine...not fast, but not Euler-rules-breakingly slow either.

fun int isPrime(int x) {
    if (x < 2) {
        return 0;
    } else if (x % 2 == 0) {
        return x == 2;
    }

    for (3 => int i; i * i <= x; 2+=>i) {
        if (x % i == 0) {
            return 0;
        }
    }
    return 1;
}

fun int mag(int x) {
    if (x < 10) {
        return 1;
    } else {
        return 10 * mag(x/10);
    }
}

fun int rotate(int x) {
    x % 10 => int d;
    return (d * mag(x)) + (x / 10);
}

fun int isCircPrime(int x) {
    if (!isPrime(x)) {
        return 0;
    }
    x => rotate => int r;
    while (x != r) {
        if (!isPrime(r)) {
            return 0;
        }
        r => rotate => r;
    }
    return 1;
}

//main
0 => int ans;
for (0 => int i; i < 1000000; ++i) {
    if (isCircPrime(i)) {
        ++ans;
    }
}
<<<ans>>>;

Problem 67 - Prolog

Problem 67 is really identical to problem 18, other than that the data set is larger. Luckily, my solution in Prolog to problem 18 was plenty fast for it to scale well to this problem. My solution runs in about half a second, having to perform 34,211 inferences in Prolog. And yeah...maybe I should find a way to load the input from a text file, but for now, the entire data table is explicitly in the program...therefore, below is an abbreviated version, the full solution, with the 100 line triangle description, is linked to from the Progress table.

getRow(R, Row) :-
nth0(R,[[59],
[73,41],
[52,40,9],
[26,53,6,34],
...,59,63,35]], Row).

max(X, Y, Y) :- Y >= X.
max(X, Y, X) :- X >= Y.

prepend([], Z, [Z]).
prepend([X|Y], Z, [Z,X|Y]).

combineRows([X, Y], [Z], Row) :-
    XZ is X + Z,
    YZ is Y + Z,
    max(XZ, YZ, N),
    prepend([], N, Row).

combineRows([X,Y|T1], [Z|T2], Row) :-
    XZ is X + Z,
    YZ is Y + Z,
    max(XZ, YZ, N),
    combineRows([Y|T1], T2, Tail),
    prepend(Tail, N, Row).

ans(99, W) :- getRow(99, W).
ans(N, W)  :-
    N1 is N + 1,
    getRow(N, Top),
    ans(N1, Bottom),
    combineRows(Bottom, Top, W).

%% Command to Find the Answer:
e67(W) :- ans(0, List),
    nth0(0, List, W).

Problem 18 - WIND

The WIND is another fictional machine introduced in my CS237 class, with an Intel-like instruction set. I re-solved problem 18 in WIND, in order to re-use my original solution to problem 18 (in Prolog) for problem 67. While WIND instructions don't exactly run super quickly (the only interpreter is one written in Java, the fastest of languages), but luckily this program was fast enough that it still ran in about .5 seconds (taking about 14000 assembly instructions)

jmp main
.origin 85
.data 100,101,103,106,110,115,121,128,136,145,155,166,178,191,205
.data 75,95,64,17,47,82,18,35,87,10,20,4,82,47,65,19,1,23,75,3,34,
88,2,77,73,7,63,67,99,65,4,28,6,16,70,92,41,41,26,56,83,40,80,70,
33,41,48,72,33,47,32,37,16,94,29,53,71,44,65,25,43,91,52,97,51,14,
70,11,33,28,77,73,17,78,39,68,17,57,91,71,52,38,17,14,91,43,58,50,
27,29,48,63,66,4,68,89,53,67,30,73,16,69,87,40,31,4,62,98,27,23,
9,70,98,73,93,38,53,60,4,23

main:
    mov     $85,r3
    mov     $13,r1     
oloop:
    cmp     $0,r1
    jl      overoloop
    mov     $0,r2
    mov     r1,r0
    add     $1,r0
    mov     0(r3,r0),r4     ;; get ptr into row below
iloop:
    cmp     r2,r1
    jl      overiloop
    mov     0(r4,r2),r5     ;; elem below
    mov     r2,r0
    add     $1,r0
    mov     0(r4,r0),r15    ;; elem below and over
    cmp     r15,r5
    jl      r15big
    mov     r5,r15
r15big:
    mov     0(r3,r1),r14     ;; ptr into cur row
    add     r15,0(r14,r2)
    add     $1,r2
    jmp     iloop
overiloop:
    sub     $1,r1
    jmp     oloop
overoloop:
    mov     $100,r0
    mov     (r0),r0
    trap    $SysPutNum
    mov     $10,r0
    trap    $SysPutChar
    trap    $SysHalt

Progress

This table records all of the languages I have used on this challenge, along with which problems have been solved in which languages. The statement of the problems can all be found on Project Euler, and my solutions are described in this blog, as well as being linked to from the table below. This table keeps growing, and keeps getting pushed lower and lower in the blog, so I keep reprinting it. The links all lead to postings of the code on github: so if you wanted all of the code, you could clone the whole repository (https://github.com/drsam94/PELC). A * denotes a solution that has been replaced by a solution in another language. If a language occurs more than once in the table, that is only because it would make the table too wide to list all solutions in one row.

BrainFuck	1*
Haskell	2*,37
J	3*,24*,35*,69*
Ruby	4*,30*,56
Whenever	2,5*
Ook!	6*
Racket	7*,26*,57
K	8*,30
Javascript	9*,31*,60
sML	10*,42
PHP	11
Go	12*,32
cLisp	13*,62
C	48*,50
Fortran95	3,15*
Python	16*,24*,39*
Python	65*,66
ELM	17*,48
Scala	18*,44
Perl	19*, 70
Java	20
FALSE	4
Squirrel	21*,52*,55
D	22

Ceylon	23
Postscript	7
Befunge	27
Boo	29
Frink	16
Forth	9
Shakespeare	12
Bash	33
Batch	34
Whitespace	36
LOLCODE	10
Lua	38
Erlang	24
Rust	40
Ada	41
Clojure	14
Coffeescript	43
Prolog	18*,67
R	45
Julia	46
x86-64	47
INTERCAL	17
OCaml	49
COBOL	21

Cobra	51
C++	14*,53
APL	52
Chapel	54
C#	26*, 59
WARM	58
Kotlin	26
Pascal	61
ALGOL68	31
Smalltalk	13
Dart	63*,65
Tcl	64
Gosu	39
Rebol	63
WIND	18
ChucK	35
F#	19
Fantom	68
Analytic	1,5,6,8
Analytic	15,25,28,69

Problem 66 - Python

So, problem 66...It involves solving Diophantine equations, which is arguably the scariest math word to appear in a Project Euler problem up to this point (yes, it just means we are looking for solutions over the integers, but Diophantine...big word). Trying to answer this problem by brute force gets one nowhere...many of the minimal solutions are absolutely enormous, and would take forever to find by brute force (especially as one can't exactly binary search for a minimal solution...). Instead, one must find out by some means (my source: Wikipedia), that the equation in question is Pell's equation, which has nice solutions in terms of the convergents of the continued fraction expansion of the square root of D. Now, Euler did a good job of preparing most users for this problem: Problem 64 was about the continued fraction expansion of square roots, and problem 65 was about finding convergents of continued fractions. So if you are most people, by the time you reached problem 66, you should have had most of the code written already. However, as I must solve all my problems in different languages, this was not the case for me. Now, I decided to actually reuse my code for problem 65, but that meant resolving problem 65, as explained in the posts below. I did however have to port my problem 64 solution in Tcl into python, which was not too difficult. All in all, the combination of code reuse and wikipedia reading made this problem not terribly hard, but it would certainly be hard in an arbitrary language. My solution takes advantage of Python's lists, dictionaries, classes, arbitrary precision integer arithmetic, and potential for functions to have multiple return values. Resolving this problem in another language would certainly not be an easy task, and is not one I plan on doing any time too soon - this problem makes for the 5th Python solution, so its probably not a bad thing, either that I will probably not recover Python for a while. Anyway, now that this problem is done, which required writing three solutions, onwards to problem 67: which is essentially identical to problem 18, explaining why 5x as many people have solved that problem as had solved 66. Below is my solution, which runs in about 1.3s on my machine.

import math

def gcd(a, b):
    return a if b == 0 else gcd(b, a%b)

class rational:
    def __add__(self, other):
        return rational(self.num*other.denom+other.num*self.denom,
                        self.denom*other.denom)

    def inverse(self):
        return rational(self.denom, self.num)

    def __init__(self, n, d):
        g = gcd(n, d)
        self.num   = n / g
        self.denom = d / g

coefDict = {}

def coeff(n,D):
    if n == 1:
        return int(math.sqrt(D))
    else:
        lst = coefDict[D]
        return lst[(n - 2)%len(lst)]

def nextDepth(n,N,D):
    if n == N:
        return rational(coeff(n,D),1)
    else:
        return rational(coeff(n,D),1) + nextDepth(n+1,N,D).inverse()

def nthConv(n,D):
    return nextDepth(1,n,D)

def mag(x):
    return 100 if x < 10 else 10 * mag(x / 10)

def nextTerm(N, m, b):
    d  = (N - b ** 2) / m
    k  = int(math.ceil(2.0 * b / d))
    bp = k*d - b
    if bp ** 2 > N:
        k -= 1
        bp = k*d - b
    return k,d,bp

def computeExpansion(D):
    b = int(math.sqrt(D))
    coefDict[D] = []
    m = 1
    computations = []
    while True:
        a, m, b = nextTerm(D, m, b)
        c = mag(b)*m + b
        if c not in computations:
            coefDict[D].append(a)
            computations.append(c)
        else:
            break

def minx(D):
    computeExpansion(D)
    n = 1
    while True:
        test = nthConv(n,D)
        if test.num ** 2 - D * (test.denom ** 2) == 1:
            return test.num
        n += 1

ans = 0
maxsol = 0
for D in range(2,1000):
    if int(math.sqrt(D)) ** 2 != D:
        x = minx(D)
        if x > maxsol:
            ans = D
            maxsol = x
print ans

65 redone in Dart

After I made a hole for problem 65 by using Python in 66, I was in a bit of a pickle thinking about what language to use, as this problem is fairly hard in a language that doesn't either have support for arbitrary-integer precision rational numbers or can be made to have such capabilities. Thus, I used Dart, which has arbitrarily large integers and nice classes and such. Though I liked Dart enough the last time I used it, on problem 63, using it on problem 65 showed me why people think it would be a better replacement for Javascript...classes and objects in Javascript can be a real pain to deal with, but in Dart they are very reasonable to work with, it was not hard at all to make a class, do operator overloading, and all that sort of stuff. (also, I discovered that Dart has a succinct syntax for single line functions, which is always nice...and a ternary operator, because every good language needs one). Solution below runs in about 80ms on my machine:

class Rational {
    int n, d;
    static int gcd(a, b) =>
        b == 0 ? a : gcd(b, a%b);

    Rational(int top, int bottom) :
        n = top ~/ gcd(top, bottom), d = bottom ~/ gcd(top, bottom) {}

    operator +(Rational other) => 
        new Rational(this.n*other.d + other.n*this.d, other.d*this.d);

    Rational inverse() => new Rational(this.d, this.n);
}

int digisum(x) => x < 10 ? x : x%10 + digisum(x~/10);

int coeff(n) => n == 1 ? 2 : (n % 3 == 0 ? 2*(n~/3) : 1);

Rational nextDepth(n,N) {
    if (n == N) {
        return new Rational(1,coeff(n));
    } else {
        return new Rational(coeff(n),1) + nextDepth(n+1,N).inverse();
    }
}

Rational nthConv(n) => nextDepth(1,n); 

main() {
    print(digisum(nthConv(100).n));
}

Problem 63 - Rebol (a long sequence of redos)

This and the next few posts are actually going to be in the reverse order that I solved a set of problems, in order to keep with the technical rules of this challenge. Problem 66 involves solving Pell's Equation, which as some quick Wikipedia ready is enough to find out, is solved relatively easily using the convergent's of continued fractions.  So, I solved problem 66 in Python, reusing some of the code from Problem 65, then resolved problem 65 using Dart, as my solution was most portable to a nice language like dart with arbitrary precision integer arithmetic and classes, and then to fill the last hole, I solved problem 63 in Rebol.

Rebol is an interesting language...it's advertising scheme revolves around the fact that it is a very fast, very small (the compiler is <1 MB), and yet very powerful language...and though I didn't need to use too much of the arguable 'power' (most of which involves its capabilities for large-scale programs, including GUI stuff), it seems like a nice enough language, once you get past some of the odd syntax (everything is in square brackets, "either" replaces if-else, functions can be declared in multiple ways, though I still only understand the 'func' identifier (which is in some ways different from 'does')). Anyway, here is a solution, runs in about 20ms on my machine:

log-10: func [x] [return (log-e x) / (log-e 10)]

ans: 0
j: 1
b: 1
while [b < 10] [
    i: b
    while [i < 10] [
        digs: 1 + (to integer! (j * log-10 (i)))
        either digs = j [
            ans: ans + 1
        ] [
            b: i + 1
        ]
        i: i + 1
    ]
    j: j + 1
]
print ans

Problem 65 - Python, again

After using Gosu this morning, Python was free again, so I decided to use it. This might have been a little brash, as problem 65 wasn't too too hard, but it was especially easy given arbitrary precision integer arithmetic and classes. This is the fourth problem that I have solved in Python (yes, I kind of like Python), and there may even be more to come in the future.

def gcd(a, b):
    return a if b == 0 else gcd(b, a%b)

class rational:
    def __add__(self, other):
        return rational(self.num*other.denom+other.num*self.denom,
                        self.denom*other.denom)

    def inverse(self):
        return rational(self.denom, self.num)

    def __init__(self, n, d):
        g = gcd(n, d)
        self.num   = n / g
        self.denom = d / g

def digisum(x):
    if x < 10:
        return x
    else:
        return x % 10 + digisum(x / 10)

def coeff(n):
    if n == 1:
        return 2
    elif n % 3 == 0:
        return 2 * (n / 3)
    else:
        return 1

def nextDepth(n,N):
    if n == N:
        return rational(1,coeff(n))
    else:
        return rational(coeff(n),1) + nextDepth(n+1,N).inverse()

def nthConv(n):
    return nextDepth(1,n)

print digisum(nthConv(100).num)

Problem 39 - Redo in Gosu

Gosu is a fairly fun language based that runs on the JVM. As it runs on the JVM, it has access to Java classes and stuff, so it has access to collections like maps, and all that fun stuff. As such, I decided to use Gosu in order to resolve problem 39, which I had previously solved in python by leveraging the Python dicts in order to make things simple: in Gosu, I just used Java HashMaps (though with a bit less of the clunky Java syntax). The hardest part of this solution was finding how to use Math.sqrt in gosu...oh well. Solution runs in about 2.5 seconds on my machine.

uses java.lang.*
    
var sols = { 0 -> 0 }

for (a in 2..500) {
    for (b in a..500) {
        var c2 = a*a + b*b
        var c  = Math.sqrt(c2) as int
        if (c2 == c * c) {
            var p = a + b + c
            if (p <= 1000) {
                if (!sols.containsKey(p)) {
                    sols.put( p, 1)
                } else {
                    sols.put( p, sols.get(p) + 1)
                }
            }
        }
    }
}

var ans = 0
var maxsols = 0
for (p in sols.keySet()) {
    if (maxsols < sols.get(p)) {
        ans = p
        maxsols = sols.get(p)
    }
}
print(ans)
     

Problem 64 - Tcl

Tcl is a weird little language. It's not all that bad once you get used to it, but things like how everything is actually a list of strings, and calculations need to be wrapped in an expr lead to some weird stuff. But, after plenty of difficulty, I managed to solve problem 64 using Tcl. Part of the added difficulty was my lack of prior exposure to continued fractions, which resulted in me having to derive the method entirely from the single given example, which caused me to state the method incorrectly more than once. Also, I messed up my strange method of storing a pair of numbers in a list a couple times (to check for repeats, I checked for repeats of what I call in the code m and b, by storing m * (10^(1 + the order of magnitude of b) + b in a list...its an effective method, despite any oddness. (and me forgetting the +1 at first, causing some false positives). The solution is a tad slow, taking about 3 seconds, but I feel like that is probably mostly Tcl's fault...as I said, everything is a list of strings apparently, which gives me no confidence in the speed of code written in the language.

proc nextTerm {N m b} {
    set d [expr {$N - ($b ** 2)}]
    set d [expr {int($d/ $m)}]
    set k [expr {ceil((2 * $b) / $d)}]
    set bp [expr {($k * $d) - $b}]
    if [expr {$bp ** 2 > $N}] {
        set k [expr {$k-1}]
        set bp [expr {($k * $d) - $b}]
    }
    return [list $d $bp]
}

proc mag {x} {
    if [expr {$x < 10}] {
        return 100
    } else {
        set b [expr {$x / 10}]
        set g [mag $b]
        return [expr {10 * $g}]
    }
}

set ans 0
for {set N 2} {$N <= 10000 } {incr N} {
    set a0 [expr {int(sqrt($N))}]
    if [expr {!($a0 * $a0 == $N)}] {
        set done -1
        set b $a0
        set m 1
        set coeffs {}
        while {$done == -1} {
            set data [nextTerm $N $m $b]
            set m    [lindex $data 0]
            set b    [lindex $data 1]
            set mg   [mag $b]
            set v    [expr {($m * $mg) + $b}]
            set done [lsearch $coeffs $v]
            if [expr {$done == -1 }] {
                lappend coeffs $v
            }
        }
        set l [llength $coeffs]
        if [expr {$l % 2 == 1}] {
            incr ans
        }
    }
}
puts $ans

Problem 63 - Dart

Dart is a programming language that plans on replacing Javascript...maybe it will be able to that, it has a lot of engineers from Google trying to make that happen. But either way, it makes a nice language to solve a Project Euler problem in. I discovered the language through http://helloworldquiz.com/, which a friend showed me, and which had surprisingly many languages which I had yet to consider for this challenge. Solution runs in about 44ms on my machine.

import "dart:math";

main() {
    num ans = 0;
    num j = 1;
    num b = 1;
    while (b < 10) {
        for (num i = b; i <= 9; ++i) {
            num pdigs = (log(pow(i,j)) / LN10 + 1).floor();
            if (pdigs == j) {
                ++ans;
            } else {
                b = (i + 1);
            }
        }
        ++j;
    }
    print(ans);
}

Problem 62 - Lisp

Lisp is a nice language, I personally prefer the Scheme dialect to Common Lisp, but I have no major qualms with it. This problem wasn't too bad...however, my solution right now is pretty slow (~2-3 min). The last time I had this issue (the other day of Prob. 60), I found a way to make the solution much faster the next day...so I will take the same approach here, and consider myself done, but I will look tomorrow to see if I can tweak it or find a better solution.

(defun digits (x)
  (if (< x 10) (cons x '()) 
      (cons (mod x 10) (digits (floor (/ x 10))))))

(defun cubes (size)
  (defun help (elem acc)
    (if (< (log (expt elem 3) 10) (- size 1)) acc 
        (help (- elem 1) (cons (expt elem 3) acc))))
  (help (floor (expt 10 (/ size 3))) '()))

(defun lsum (lst)
  (defun help (l tot)
    (if (equal l '()) tot (help (cdr l) (+ tot (car l)))))
  (help lst 0))

(defun ans (lst)
  (if (equal lst '()) 0
    (let ((digs (sort (digits (car lst)) #'>)))
    (let ((a (lsum (map 'list #'(lambda (x) 
        (if (equal digs (sort (digits x) #'>)) 1 0)) (cdr lst)))))
        (if (= a 4) (car lst) (ans (cdr lst)))))))

(defun e62 (size)
  (let ((a (ans (cubes size)))) (if (= a 0) (e62 (+ size 1)) a)))

(e62 5)

Problem 13 re-do - Smalltalk cop-out

I have never been a big fan of Smalltalk - I was first introduced to it in Principles of Programming Languages in the context of Squeak, a very strange IDE that only helped to exaggerate the strangeness of the message-sending object-oriented language. However, it was just about the only language left out there which I had actually used at some point in my life, yet I had not used in this challenge yet. Luckily, turns out Smalltalk has arbitrary precision integer arithmetic, which allowed me to use it to provide a very simple solution to problem 13, freeing up Common Lisp for the first time since the initial 24-hour hackathon during which I solved the first 20 problems. My wonderfully unexciting program is below.

((37107287533902102798797998220837590246510135740250 + 
46376937677490009712648124896970078050417018260538 + 
74324986199524741059474233309513058123726617309629 + 
91942213363574161572522430563301811072406154908250 + 
23067588207539346171171980310421047513778063246676 + 
89261670696623633820136378418383684178734361726757 + 
28112879812849979408065481931592621691275889832738 + 
44274228917432520321923589422876796487670272189318 + 
47451445736001306439091167216856844588711603153276 + 
70386486105843025439939619828917593665686757934951 + 
62176457141856560629502157223196586755079324193331 + 
64906352462741904929101432445813822663347944758178 + 
92575867718337217661963751590579239728245598838407 + 
58203565325359399008402633568948830189458628227828 + 
80181199384826282014278194139940567587151170094390 + 
35398664372827112653829987240784473053190104293586 + 
86515506006295864861532075273371959191420517255829 + 
71693888707715466499115593487603532921714970056938 + 
54370070576826684624621495650076471787294438377604 + 
53282654108756828443191190634694037855217779295145 + 
36123272525000296071075082563815656710885258350721 + 
45876576172410976447339110607218265236877223636045 + 
17423706905851860660448207621209813287860733969412 + 
81142660418086830619328460811191061556940512689692 + 
51934325451728388641918047049293215058642563049483 + 
62467221648435076201727918039944693004732956340691 + 
15732444386908125794514089057706229429197107928209 + 
55037687525678773091862540744969844508330393682126 + 
18336384825330154686196124348767681297534375946515 + 
80386287592878490201521685554828717201219257766954 + 
78182833757993103614740356856449095527097864797581 + 
16726320100436897842553539920931837441497806860984 + 
48403098129077791799088218795327364475675590848030 + 
87086987551392711854517078544161852424320693150332 + 
59959406895756536782107074926966537676326235447210 + 
69793950679652694742597709739166693763042633987085 + 
41052684708299085211399427365734116182760315001271 + 
65378607361501080857009149939512557028198746004375 + 
35829035317434717326932123578154982629742552737307 + 
94953759765105305946966067683156574377167401875275 + 
88902802571733229619176668713819931811048770190271 + 
25267680276078003013678680992525463401061632866526 + 
36270218540497705585629946580636237993140746255962 + 
24074486908231174977792365466257246923322810917141 + 
91430288197103288597806669760892938638285025333403 + 
34413065578016127815921815005561868836468420090470 + 
23053081172816430487623791969842487255036638784583 + 
11487696932154902810424020138335124462181441773470 + 
63783299490636259666498587618221225225512486764533 + 
67720186971698544312419572409913959008952310058822 + 
95548255300263520781532296796249481641953868218774 + 
76085327132285723110424803456124867697064507995236 + 
37774242535411291684276865538926205024910326572967 + 
23701913275725675285653248258265463092207058596522 + 
29798860272258331913126375147341994889534765745501 + 
18495701454879288984856827726077713721403798879715 + 
38298203783031473527721580348144513491373226651381 + 
34829543829199918180278916522431027392251122869539 + 
40957953066405232632538044100059654939159879593635 + 
29746152185502371307642255121183693803580388584903 + 
41698116222072977186158236678424689157993532961922 + 
62467957194401269043877107275048102390895523597457 + 
23189706772547915061505504953922979530901129967519 + 
86188088225875314529584099251203829009407770775672 + 
11306739708304724483816533873502340845647058077308 + 
82959174767140363198008187129011875491310547126581 + 
97623331044818386269515456334926366572897563400500 + 
42846280183517070527831839425882145521227251250327 + 
55121603546981200581762165212827652751691296897789 + 
32238195734329339946437501907836945765883352399886 + 
75506164965184775180738168837861091527357929701337 + 
62177842752192623401942399639168044983993173312731 + 
32924185707147349566916674687634660915035914677504 + 
99518671430235219628894890102423325116913619626622 + 
73267460800591547471830798392868535206946944540724 + 
76841822524674417161514036427982273348055556214818 + 
97142617910342598647204516893989422179826088076852 + 
87783646182799346313767754307809363333018982642090 + 
10848802521674670883215120185883543223812876952786 + 
71329612474782464538636993009049310363619763878039 + 
62184073572399794223406235393808339651327408011116 + 
66627891981488087797941876876144230030984490851411 + 
60661826293682836764744779239180335110989069790714 + 
85786944089552990653640447425576083659976645795096 + 
66024396409905389607120198219976047599490197230297 + 
64913982680032973156037120041377903785566085089252 + 
16730939319872750275468906903707539413042652315011 + 
94809377245048795150954100921645863754710598436791 + 
78639167021187492431995700641917969777599028300699 + 
15368713711936614952811305876380278410754449733078 + 
40789923115535562561142322423255033685442488917353 + 
44889911501440648020369068063960672322193204149535 + 
41503128880339536053299340368006977710650566631954 + 
81234880673210146739058568557934581403627822703280 + 
82616570773948327592232845941706525094512325230608 + 
22918802058777319719839450180888072429661980811197 + 
77158542502016545090413245809786882778948721859617 + 
72107838435069186155435662884062257473692284509516 + 
20849603980134001723930671666823555245252804609722 + 
53503534226472524250874054075591789781264330331690 ) 
asString copyFrom: 1 to: 10) printNl

Problem 60 - Javascript

Problem 60 was a bit annoying - its a problem that I had never solved in my past life of doing Euler problems without language restrictions, but everything turned out O.K. My first solution was excessively slow (2-5min), but I noticed something really dumb I was doing that caused me to search way more combinations than necessary. The code below runs in <1s. And now I am formally up in another multiple of 10 - the 60's are upon us.

function isPrime(x) {
    if (x < 2) {
        return false
    } else if (x % 2 == 0) {
        return x == 2
    }

    for (var i = 3; i * i <= x; ++i) {
        if (x % i == 0) {
            return false;
        }
    }
    return true;
}
function mag(x) {
    return x < 10 ? 10 : 10 * mag(x / 10)
}
function concat(x, y) {
    return x * mag(y) + y
}

function getPrimes(N) {
    p = new Array()
    for (var i = 2; i <= N; ++i) {
        if (isPrime(i)) {
            p.push(i)
        }
    }
    return p
}

function satisfies(p,s) {
    for (var i = 0; i < s.length; ++i) {
        if (!(isPrime(concat(s[i],p)) && isPrime(concat(p,s[i])))) {
            return false
        }
    }
    return true
}
function setSum(set) {
    sum = 0;
    for (var i = 0; i < set.length; ++i) {
        sum += set[i]
    }
    return sum 
}
function scopy(set) {
    return set.slice(0)
}

primes = getPrimes(10000)
function finishSet(set, target, start) {
    if (set.length == 5) {
        test = setSum(set)
        return test < target ? test : target
    } else {
        for (var i = start; i < primes.length; ++i) {
            p = primes[i]
            if (satisfies(p,set)) {
                cset = scopy(set)
                cset.push(p)
                test = finishSet(cset,target,i)
                if (test < target) {
                    target = test
                    break
                }
            }
        }
    }
    return target
}
/** this is 'main' */
function euler60() {
    ans = finishSet(new Array(),10000000,0)
    console.log(ans)
    alert(ans)
}

Progress

This table records all of the languages I have used on this challenge, along with which problems have been solved in which languages. The statement of the problems can all be found on Project Euler, and my solutions are described in this blog, as well as being linked to from the table below. This table keeps growing, and keeps getting pushed lower and lower in the blog, so I keep reprinting it. The links all lead to postings of the code on github: so if you wanted all of the code, you could clone the whole repository (https://github.com/drsam94/PELC). A * denotes a solution that has been replaced by a solution in another language.

BrainFuck	1*
Haskell	2*,37
J	3*,24*,35
Ruby	4*,30*,56
Whenever	2,5*
Ook!	6*
Racket	7*,26*,57
K	8*,30
Javascript	9*,31*,60
sML	10*,42
PHP	11
Go	12*,32
cLisp	13*,62
C	48*,50
Fortran95	3,15*
Python	16*,24*,39*, 65
ELM	17*,48
Scala	18*,44
Perl	19
Java	20
FALSE	4
Squirrel	21*,52*,55
D	22
Ceylon	23

Postscript	7
Befunge	27
Boo	29
Frink	16
Forth	9
Shakespeare	12
Bash	33
Batch	34
Whitespace	36
LOLCODE	10
Lua	38
Erlang	24
Rust	40
Ada	41
Clojure	14
Coffeescript	43
Prolog	18
R	45
Julia	46
x86-64	47
INTERCAL	17
OCaml	49
COBOL	21
Cobra	51

C++	14*,53
APL	52
Chapel	54
C#	26*, 59
WARM	58
Kotlin	26
Pascal	61
ALGOL68	31
Smalltalk	13
Dart	63
Tcl	64
Gosu	39
Analytic	1,5,6,8
Analytic	15,25,28

ALGOL 68 to re-do problem 31

Continuing with the trend of using more historical languages after using Pascal yesterday, I decided to use ALGOL 68 today, in order to re-solve problem 31 (so that I have a language I know as well as Javascript in order to solve problem 60, which seems a little annoying / I haven't ever solved it even in my prior Euler experience, so yeah, Javascript will be nice). Anyway, ALGOL 68 wasn't that bad to actually write code in...and a nice paper by Tannenbaum basically gave all the description I needed (http://maben.homeip.net/static/S100/research%20machines/software/algol/ACMSurveys0676-Tanenbaum.pdf), but I ran into a lot of issues not related to the actual writing of code...rather, I ran into issues with compilers. After taking a little while to find a compiler at all, I install ALGOL 68 Genie, which purports to be one of the better/more full compilers of ALGOL. Turns out it actually is a reasonably good compiler, but there was one large difference between the compiler and the Tannenbaum paper: Tanenbaum uses all lowercase words, and a68g requires all reserved word to be written in all caps. Once I discovered this, it wasn't much of an issue...but it took a long time to discover this. It wasn't documented particularly well anywhere, and the compiler errors related to using lower case words were utterly unhelpful (most looked like missing-semicolon errors)...but anyway, after getting past all of that, things were not too bad. Below is my solution which takes 27s to run...I imagine that is more due to the compiler not being the most efficient than because my code is actually slow.

BEGIN

[1:7] INT a := (1,2,5,10,20,50,100);
INT ways := 1;

PROC helper = (INT amountindex, inheritedamount) VOID:
BEGIN
    INT current;
    BOOL done := FALSE;   
    FOR i FROM 0 BY 1 TO 200 WHILE NOT done
    DO
        current := inheritedamount + (i * a[amountindex]);
        IF current > 200
            THEN done := TRUE
            ELIF current = 200
            THEN ways +:= 1;
                 done := TRUE
            ELIF amountindex < 7
            THEN helper(amountindex + 1, current)
        FI
    OD
END;

helper(1,0);
write((ways, new line))          
END

Problem 61 - Pascal

Pascal is one of those languages that comes up all the time as an older language that used to be hip and cool, but nobody uses anymore. So, that makes it a perfect target for use in this blog. My experience seemed to validate Pascal's current status: its a fine language and all, but its a bit stiff and has a lot of weird features that have been smoothed out in later languages (variable declarations need to go in a certain place, function parameters are grouped by type, semicolons are necessary in weird places, and must be omitted in even weirder places [defining weird as unnatural to my sensibilities, though the rules do themselves make some sense]). Anyway, most of my difficulty with this problem came not in difficulties with Pascal, but rather difficulties with a corner case of the problem itself: it is necessary to find a cycle of numbers where one is triangular, one square, one pentagonal, one hexagonal, one heptagonal, and one octagonal. I originally decided to make the beginning of my cycle an octagonal number, as there are the fewest of them, but I realized this to be foolhardy due to an issue involving all hexagonal numbers also being triangular. This was solved by beginning cycles with triangular numbers. Anyway, my Pascal solution is below, and it runs fairly quickly (~ 20ms on my machine). Also, if anyone noticed, yes, I skipped problem 60. I'll get to it soon enough.

program Euler61;

var
    temp : integer;
    ans  : integer;
    t    : integer;
    i    : integer;

function nthTri(x : integer) : integer;
begin
    nthTri := (x * (x + 1)) div  2
end;
function isTri(x : integer) : Boolean;
begin
    temp  := round(sqrt(2*x + 1/4) + 1/2);
    isTri := temp * (temp + 1) / 2 = x
end;

function isSquare(x : integer) : Boolean;
begin
    temp     := round(sqrt(x));
    isSquare := temp * temp = x
end;

function isPent(x : integer) : Boolean;
begin
    temp   := round(sqrt(2/3 * x + 1 / 36) + 1/6);
    isPent := temp * (3 * temp - 1) / 2 = x
end;

function isHex(x : integer) : Boolean;
begin
    temp  := round(sqrt(x / 2 + 1 / 16) + 1/4);
    isHex := temp * (2 * temp - 1) = x
end;

function isHept(x : integer) : Boolean;
begin
    temp   := round(sqrt(2 / 5 * x + 9 / 100) + 3 / 10);
    isHept := temp * (5 * temp - 3) / 2 = x
end;

function isOct(x : integer) : Boolean;
begin
    temp  := round(sqrt(x / 3 + 1/9) + 1/3);
    isOct := temp * (3 * temp - 2) = x
end;

function shape(x,found : integer) : integer;
begin
    if isHept(x) then
        shape := 5
    else
        if isHex(x) then
            shape := 4
        else
            if isPent(x) then
                shape := 3
            else
                if isSquare(x) then
                    shape := 2
                else
                    if isOct(x) then 
                        shape := 1
                    else
                        shape := 0;
end;

function Next(x, found, sum, origin : integer) : integer;
var
    base : integer;
    s : integer;
    r : integer;
    test : integer;
    tfound : integer;
    tsum : integer;
begin
    r := 10;
    base := (x mod 100) * 100;
    while r <= 99 do
    begin
        s := 1 shl shape(base + r, found);
        if (found and s) = 0 then
        begin
            tfound := found + s;
            tsum   := sum + (base + r);
            if tfound = 63 then
            begin
                if r = origin div 100 then
                    Next := tsum
                else
                    Next := 0;
                exit;
            end

            else
            begin
                test := Next(base + r, tfound, tsum, origin);
                if test > 0 then
                begin
                    Next := test;
                    exit
                end;
            end;
        end;
        r := r + 1;
    end;
    Next := 0
end;

begin
    ans := 0;
    i := 45;
    while ans = 0 do
    begin
        t   := nthTri(i);
        if t mod 100 > 10 then
            ans := Next(t,1,t,t);
        i   := i + 1;
        if t > 9999 then
            exit;
    end;
    writeln(ans)
end.

Problem 59 - C# again

Despite the fact that this is the second C# solution in about 5 solutions, I really am not a huge fan of the language - It is a nice language with C-like syntax though which I had yet to use for a hard enough problem that it was attractive to continue uprooting it. This problem was mildly annoying though, require things like file i/o that can sometimes be frustrating or near-impossible with more esoteric languages. This problem involves XOR encryption of a file, and you must brute-force find a key for the message. Though I remembered the message and key from the last time I solved this problem, I didn't take any advantage of that knowledge in my solution, as indicated by the fact that my test for correctness is simply containing " the " somewhere in the message. Anyway, this should put C# to rest for a while...solution runs in 2 seconds, but I could probably speed it up if I wanted to by not re-reading the file on every key. Also, for anyone checking my solutions, I reformatted the input file to have one number per line.

using System.IO;
public class Euler59 {
    
    public static char[] decode(char[] key) {
        StreamReader reader = new StreamReader("cipher1.txt");
        char[] message = new char[1201];
        int i = 0;
        do {
            char c = (char) System.Convert.ToInt16(reader.ReadLine());
            message[i] = (char)(((int)c) ^ ((int)key[i%3]));
            i++;
        } while (reader.Peek() != -1);
        return message;
    }


    public static int Main(string[] args) {
        for (char a = 'a'; a <= 'z'; ++a) {
            for (char b = 'a'; b <= 'z'; ++b) {
                for (char c = 'b'; c <= 'z'; ++c) {
                    char[] key = new char[3] {a,b,c};
                    char[] asciis = decode(key);
                    string message = new string(asciis);
                    if (message.Contains(" the ")) {
                        int ans = 0;
                        for (int i = 0; i < asciis.Length; ++i) {
                            ans += (int) asciis[i];
                        }
                        System.Console.WriteLine(ans);
                        return 0;
                    }
                }
            }
        }
        return 0;
    }
}

Kotlin on problem 26

Kotlin is a language that compiles to both the JVM and javascript. It seems like a fairly reasonable language, using some features that seem similar to Scala, but overall retaining a very Java-like syntax. I had some difficulty installing Kotlin (aka my hard drive does not have enough free space for it), so I only used the free online in-browser version...as such, Kotlin would not have worked very well for the problem I am up to, problem 59, which requires file I/O. So instead, I used Kotlin to rewrite my solution to problem 26 that was in C#, and C# will be used for a solution to 59.

import java.util.Collections.*
fun biggestCycle(d: Int): Int {
  val rems: Array< int > = Array< int >(1000){i -> 0}
  var r: Int = 1
  var i: Int = 0
  while (r != 0) {
    r = (r * 10) % d
    for (j in 0..i) {
      if (rems[j] == r) {
        return i
      }
    }
    rems[i++] = r
  }
  return i
}
  
  
fun main(args : Array< string >) {
  var ans = 0
  for (d in 1..1000) {
    if (ans < biggestCycle(d)) {
      ans = d
    } 
  }
  println(ans)
}

After a long hiatus, WARM for Problem 58

I have not gotten around to solving too many project Euler problems as of late, as unfortunately I have had ample supplies of classwork to do. However, after finishing a midterm earlier today, I decided that I finally had enough time to solve a problem. I couldn't entirely leave my studies though, so I decided to solve this problem in WARM, a Risc assembly language based on ARM created for Williams CS237 classes. The solution is rather slow due to the lack of a direct implementation of WARM, it has to run in a virtual machine. But nonetheless, it produces the correct answer, bringing me up to 58 problems solved.

    b       main
isPrime:
    cmp     r0,#2
    blt     finif
    beq     finit
    tst     r0,#1
    beq     finif
    mov     r1,#3
loop:
    mul     r2,r1,r1
    cmp     r2,r0
    bgt     finit
    div     r2,r0,r1
    mul     r2,r2,r1
    cmp     r2,r0
    beq     finif
    add     r1,r1,#2
    b       loop
finif:
    mov     r0,#0
    mov     pc, lr
finit:
    mov     r1,#0
    mov     pc, lr

main:
    mov     v1,#1  ;; level
    mov     v3,#0  ;; primes
loop2:
    mov     v2,#1  ;; index in level
loop3:
    cmp     v2,#4
    bge     overloop3
    mul     r0,v1,#2
    sub     r0,r0,#1
    mul     r0,r0,r0
    mul     v4,v1,v2, lsl #1
    add     r0,r0,v4 ;; finished computation
    bl      isPrime
    cmp     r0,#0
    addne   v3,v3,#1
    add     v2,v2,#1
    b       loop3
overloop3:
    mul     v4,v1,v2
    add     v4,v4,#1
    div     r0,v4,v3
    cmp     r0,#10
    bge     overloop2
    add     v1,v1,#1
    b       loop2
overloop2:
    mov     r0,v1, lsl #1 
    add     r0,r0,#1
    swi     #SysPutNum
    mov     r0,#10
    swi     #SysPutChar
    swi     #SysHalt