Find all integers between m and n whose sum of squared divisors is itself a square

Question

Problem Question

Divisors of 42 are : 1, 2, 3, 6, 7, 14, 21, 42. These divisors squared are: 1, 4, 9, 36, 49, 196, 441, 1764. The sum of the squared divisors is 2500 which is 50 * 50, a square!

Given two integers m, n (1 <= m <= n) we want to find all integers between m and n whose sum of squared divisors is itself a square. 42 is such a number.

The result will be an array of arrays, each subarray having two elements, first the number whose squared divisors is a square and then the sum of the squared divisors.

Code below

How can I make this specific program run faster? My current code times out after n > 9999.

#returns the divisors of each number in an array of arrays

r = (m..n).to_a.map { |z| (1..z).select { |x| z % x == 0} }

#this finds all integers between m and n whose sum of squared divisors is itself a square

squarenumbers = r.map { |x| x.map { |c| c**2 }.inject(:+) }.select { |x| Math.sqrt(x) % 1 == 0 }

#returns an array of booleans. 

booleans = r.map { |x| x.map { |c| c**2 }.inject(:+) }.map { |x| Math.sqrt(x) % 1 == 0 }

#returns the index of each of the true values in booleans as an array

indexer = booleans.map.with_index{|x, i| i if x == true }.compact

#returns the numbers whose squared divisors is a square in an array

unsqr = indexer.map { |x| (m..n).to_a[x] }

#merges the two arrays together, element for element and creates an array of arrays

unsqr.zip(squarenumbers) 

# for m = 1 and n = 1000 the result would be

# [[1, 1], [42, 2500], [246, 84100], [287, 84100], [728, 722500]]

Answer 1

Brute-force calculatioins of factors

You begin by calculating:

m, n = 40, 42
r = (m..n).to_a.map { |z| (1..z).select { |x| z % x == 0} }
  #=> [[1, 2, 4, 5, 8, 10, 20, 40], [1, 41], [1, 2, 3, 6, 7, 14, 21, 42]]

That's OK, but you don't need .to_a:

r = (m..n).map { |z| (1..z).select { |x| z % x == 0} }
  #=> [[1, 2, 4, 5, 8, 10, 20, 40], [1, 41], [1, 2, 3, 6, 7, 14, 21, 42]]

This avoids an extra step, which is the creation of the temporary array^1,2:

(m..n).to_a #=> [40, 41, 42]

Structure of a solution

Let's work backwards to come up with our code. First, concentrate on determining, for any given number q, if the sum of squares of the factors of q is itself a perfect square. Suppose we construct a method magic_number? which takes q as its only argument and returns true if q satisfies the required property and false otherwise. Then we will compute:

(m..n).select { |q| magic_number?(q) }

to return an array of all numbers between m and n that satisfy the property. magic_number? can be written like this:

def magic_number?(q)
  return true if q == 1
  s = sum_of_squared_factors(q)
  s == Math.sqrt(s).round**2
end

Calculating sum of squared factors

So now we are left with writing the method sum_of_squared_factors. We can use your code to obtain the factors:

def factors(q)
  (1..q).select { |x| q % x == 0 }
end

factors(40) #=> [1, 2, 4, 5, 8, 10, 20, 40] 
factors(41) #=> [1, 41] 
factors(42) #=> [1, 2, 3, 6, 7, 14, 21, 42]

and then write:

def sum_of_squared_factors(q)
  factors(q).reduce(0) { |t,i| t + i*i }
end

sum_of_squared_factors(40) #=> 2210 
sum_of_squared_factors(41) #=> 1682 
sum_of_squared_factors(42) #=> 2500 

Speeding the calculation of factors

There's something more we can do to speed up the calculation of factors. If f is a factor of n, f and n/f, are both factors of n. (For example, since 3 is a factor of 42, so is 42/3 #=> 14). We therefore need only obtain the smaller of each pair.

There is one exception to this rule. If n is a perfect square and f == n**0.5, then f = n/f, so we only include f among the factors of n (not n/f as well).

If turns out that if f is the smaller of the pair, f <=(n**0.5).round³. We therefore need only check to see which of the numbers (1..(n**0.5).round) are factors and include their complements (unless n is a perfect square, in which case we do not double-count (n**0.5).round):

q = 42
arr = (1..Math.sqrt(q).round).select { |x| q % x == 0 }
  #=> [1, 2, 3, 6] 
arr = arr.flat_map { |n| [n, q/n] }
  #=> [1, 42, 2, 21, 3, 14, 6, 7] 
arr.pop if a[-2] == a[-1]
arr
  #=> [1, 42, 2, 21, 3, 14, 6, 7]

q = 36
arr = (1..Math.sqrt(q).round).select { |x| q % x == 0 }
  #=> [1, 2, 3, 4, 6] 
arr = arr.flat_map { |n| [n, q/n] }
  #=> [1, 36, 2, 18, 3, 12, 4, 9, 6, 6] 
arr.pop if a[-2] == a[-1]
  #=> 6 
arr
  #=> [1, 36, 2, 18, 3, 12, 4, 9, 6] 

so we can write:

def factors(q)
  arr = (1..Math.sqrt(q)).select { |x| q % x == 0 }
  arr = arr.flat_map { |n| [n, q/n] }
  arr.pop if arr[-2] == arr[-1]
  arr
end

Substituting out arr ("chaining" expressions), we obtain a typical Ruby expression:

def factors(q)
  (1..Math.sqrt(q)).select { |x| q % x == 0 }.
    flat_map { |n| [n, q/n] }.
    tap { |a| a.pop if a[-2] == a[-1] }
end

factors(42)
  #=> [1, 42, 2, 21, 3, 14, 6, 7] 
factors(36)
  #=> [1, 36, 2, 18, 3, 12, 4, 9, 6] 

See Enumerable#flat_map and Object#tap. (There's no need for this array to be sorted. In applications where it needs to be sorted, just tack .sort onto the end of flat_maps block.)

Wrapping up

In sum, we are left with the following:

def magic_number?(q)
  return true if q == 1
  s = sum_of_squared_factors(q)
  s == Math.sqrt(s).round**2
end

def sum_of_squared_factors(q)
  factors(q).reduce(0) { |t,i| t + i*i }
end

def factors(q)
  (1..Math.sqrt(q)).select { |x| q % x == 0 }.
    flat_map { |n| [n, q/n] }.
    tap { |a| a.pop if a[-2] == a[-1] }
end

m, n = 1, 1000
(m..n).select { |q| magic_number?(q) }
  #=> `[1, 42, 246, 287, 728]

This calculation was completed in a blink of an eye.

Compute primes to further speed calculation of factors

Lastly, let me describe an even faster way to compute the factors of a number, using the method Prime::prime_division. That method decomposes any number into its prime components. Consider, for example, n = 360.

require 'prime'

Prime.prime_division(360)
  #=> [[2, 3], [3, 2], [5, 1]]

This tells us that:

360 == 2**3 * 3**2 * 5**1
  #=> true

It also tells us that every factor of 360 is the product of between 0 and 3 2's, multiplied by between 0 and 2 3's, multiplied by 0 or 1 5's. Therefore:

 def factors(n)
   Prime.prime_division(n).reduce([1]) do |a,(prime,pow)|
     a.product((0..pow).map { |po| prime**po }).map { |x,y| x*y }
   end
 end

 a = factors(360).sort
   #=> [ 1,  2,  3,  4,  5,  6,  8,  9, 10,  12,  15,  18,
   #    20, 24, 30, 36, 40, 45, 60, 72, 90, 120, 180, 360]

We can check that:

 a == (1..360).select { |n| (360 % n).zero? }
   #=> true

One other check:

 factors(40).sort
   #=> [1, 2, 4, 5, 8, 10, 20, 40]            

^{1. You could instead write that [*m..n] #=> [40, 41, 42].2. Why is it not necessary to convert the range to an array? Enumerable#map, being an instance method of the module Enumerable, is available for use by every class that includes Enumerable. Array is one, but (m..n).class #=> Range is another. (See the second paragraph at Range).3. Suppose f is smaller than n/f and f > n**0.5, then n/f < n/(n**0.5) = n**0.5 < f, a contradiction.}

Answer 2

I don't know Ruby but the problem lies with the algorithm used in finding the divisors of a number (which is not specific to the language used, i.e. Ruby in this case).

r = (m..n).to_a.map { |z| (1..z).select { |x| z % x == 0} }

To find the divisors of an integer n you are dividing n by all positive integers unto n - 1 which means the loop runs n - 1 times. However, it is enough to divide upto sort(n) to calculate the divisors. In pseudocode this looks like below:

for i = 1 to i <= sqrt(n)
    r = n % i
    if r == 0 then
        i is a divisor
        if n / i != i then
            n / i is another divisor

For example:

sqrt_42 = 6.48074069840786
i = 1 => 1 and 42 are two divisors
i = 2 => 2 and 21
i = 3 => 3 and 14
i = 4 => no divisor
i = 5 => no divisor
i = 6 => 6 and 7

And thats all.

This will improve the performance a lot since now the loop runs only sort(n) times instead of n - 1 times which is a big difference for large n.