Factorizations to get a very special number
Mathematical recreation is fun. I came across a nice Mind your Decisions video and thought it would be cool to use Julia
to find that special number. Have a look at the video here before we continue
using Primes
function getBigestFactors{T}(a::T,b::T)
c=T()
for k in keys(a)
c[k]=max(a[k],b[k])
end
for k in keys(b)
c[k]=max(a[k],b[k])
end
c
end
superfactors=reduce(getBigestFactors,map(factor,BigInt.(1:100)))
n=prodfactors(superfactors)
Julia
can be amazingly expressive. In line 14 are getting all the prime factors of all the numbers from 1 to 100 using the map
function. We are using the type BigInt
since will be generating pretty large numbers. The standard Int
will simply overflow. We then use reduce
to get all the common factors for the numbers from 1 to 100.
If your are slight confused lets take a bit slower. First we generate a one dimensional array of BigInt
s from 1 to 100.
BigInt.(1:100)
This makes use of Julia
‘s very powerful dot factorization syntax. We then map
that array into another array of Primes.Factorization{BigInt}
through the factor
function in the Primes
package. The Primes.Factorization
type is subtype of Associative
. I learned that from reading the implementation.
The cool thing about Associative
types that the can be very conveniently accessed. Lets have a closer look.
a=factor(24) println(a)
yields
Primes.Factorization(2=>3,3=>1)
That is \(2^3 3\)
You can then do the following
print(a[2])
yielding
3
This means that the number 24 has a 2 raised to the third power as a factor.
Or, be more adventurous and try
print(a[5])
which will print a
0
This makes perfect sense as \(5^0=1\) and 1
is always a factor of any number. We will always get 0 exponent for any number that is not a prime factor of a given number. The Primes.jl
package authors have done great use of Julia
‘s very powerful type system.
To solve this problem we need to look at all the prime factors for all the numbers from 1 to 100 find the highest exponent of the of each of those prime factors. To do that we implement the getBigestFactor
function when does that for any two prime factorizations. We plug that into reduce
and et voilà!
The superfactors in line 14 are
2^6 ⋅ 3^4 ⋅ 5^2 ⋅ 7^2 ⋅ 11 ⋅ 13 ⋅ 17 ⋅ 19 ⋅ 23 ⋅ 29 ⋅ 31 ⋅ 37 ⋅ 41 ⋅ 43 ⋅ 47 ⋅ 53 ⋅ 59 ⋅ 61 ⋅ 67 ⋅ 71 ⋅ 73 ⋅ 79 ⋅ 83 ⋅ 89 ⋅ 97
Finally we multiply them to get our special number at line 16
n=69720375229712477164533808935312303556800
Are we done yet?
As the saying goes “there is more than one way to skin a cat”. Some are more elegant than others.
Here is a more elegant and computationally more efficient way.
Using Primes
function getAllFactors{T}(n::T)
p_factors=filter(isprime,T.(1:n))
exponents=[floor(T,log(p,n)) for p in p_factors]
zip(p_factors,exponents) |> Dict |> Primes.Factorization
end
Here we just get the prime numbers in the range from 1 to n
(line 3). Using those primes we then get the highest exponents that will not yield a number outside the range (line 4). Finally, we package everything as a Primes.Factorization
(line 5).
In Julia
REPL we get
julia> getAllFactors(big"100") 2^6 ⋅ 3^4 ⋅ 5^2 ⋅ 7^2 ⋅ 11 ⋅ 13 ⋅ 17 ⋅ 19 ⋅ 23 ⋅ 29 ⋅ 31 ⋅ 37 ⋅ 41 ⋅ 43 ⋅ 47 ⋅ 53 ⋅ 59 ⋅ 61 ⋅ 67 ⋅ 71 ⋅ 73 ⋅ 79 ⋅ 83 ⋅ 89 ⋅ 97 julia> prodfactors(ans) 69720375229712477164533808935312303556800
At that point, you might ask, why stop at a 100. Julia
is lightening fast and powerful. Just for fun we can do the same for numbers from 1 to 1000.
julia> getAllFactors(big"1000") 2^9 ⋅ 3^6 ⋅ 5^4 ⋅ 7^3 ⋅ 11^2 ⋅ 13^2 ⋅ 17^2 ⋅ 19^2 ⋅ 23^2 ⋅ 29^2 ⋅ 31^2 ⋅ 37 ⋅ 41 ⋅ 43 ⋅ 47 ⋅ 53 ⋅ 59 ⋅ 61 ⋅ 67 ⋅ 71 ⋅ 73 ⋅ 79 ⋅ 83 ⋅ 89 ⋅ 97 ⋅ 101 ⋅ 103 ⋅ 107 ⋅ 109 ⋅ 113 ⋅ 127 ⋅ 131 ⋅ 137 ⋅ 139 ⋅ 149 ⋅ 151 ⋅ 157 ⋅ 163 ⋅ 167 ⋅ 173 ⋅ 179 ⋅ 181 ⋅ 191 ⋅ 193 ⋅ 197 ⋅ 199 ⋅ 211 ⋅ 223 ⋅ 227 ⋅ 229 ⋅ 233 ⋅ 239 ⋅ 241 ⋅ 251 ⋅ 257 ⋅ 263 ⋅ 269 ⋅ 271 ⋅ 277 ⋅ 281 ⋅ 283 ⋅ 293 ⋅ 307 ⋅ 311 ⋅ 313 ⋅ 317 ⋅ 331 ⋅ 337 ⋅ 347 ⋅ 349 ⋅ 353 ⋅ 359 ⋅ 367 ⋅ 373 ⋅ 379 ⋅ 383 ⋅ 389 ⋅ 397 ⋅ 401 ⋅ 409 ⋅ 419 ⋅ 421 ⋅ 431 ⋅ 433 ⋅ 439 ⋅ 443 ⋅ 449 ⋅ 457 ⋅ 461 ⋅ 463 ⋅ 467 ⋅ 479 ⋅ 487 ⋅ 491 ⋅ 499 ⋅ 503 ⋅ 509 ⋅ 521 ⋅ 523 ⋅ 541 ⋅ 547 ⋅ 557 ⋅ 563 ⋅ 569 ⋅ 571 ⋅ 577 ⋅ 587 ⋅ 593 ⋅ 599 ⋅ 601 ⋅ 607 ⋅ 613 ⋅ 617 ⋅ 619 ⋅ 631 ⋅ 641 ⋅ 643 ⋅ 647 ⋅ 653 ⋅ 659 ⋅ 661 ⋅ 673 ⋅ 677 ⋅ 683 ⋅ 691 ⋅ 701 ⋅ 709 ⋅ 719 ⋅ 727 ⋅ 733 ⋅ 739 ⋅ 743 ⋅ 751 ⋅ 757 ⋅ 761 ⋅ 769 ⋅ 773 ⋅ 787 ⋅ 797 ⋅ 809 ⋅ 811 ⋅ 821 ⋅ 823 ⋅ 827 ⋅ 829 ⋅ 839 ⋅ 853 ⋅ 857 ⋅ 859 ⋅ 863 ⋅ 877 ⋅ 881 ⋅ 883 ⋅ 887 ⋅ 907 ⋅ 911 ⋅ 919 ⋅ 929 ⋅ 937 ⋅ 941 ⋅ 947 ⋅ 953 ⋅ 967 ⋅ 971 ⋅ 977 ⋅ 983 ⋅ 991 ⋅ 997 julia> prodfactors(ans) 7128865274665093053166384155714272920668358861885893040452001991154324087581111499476444151913871586911717817019575256512980264067621009251465871004305131072686268143200196609974862745937188343705015434452523739745298963145674982128236956232823794011068809262317708861979540791247754558049326475737829923352751796735248042463638051137034331214781746850878453485678021888075373249921995672056932029099390891687487672697950931603520000
Now we get super special number that is 443 digits long. Lets finish by a Julia
one liner to get the number of digits
julia> getAllFactors(big"1000") |> prodfactors|> digits |>length 433