## Nearly-Square Divisors

### August 5, 2016

We have today a fun little problem from number theory:

Given positive integers n, a and b such that n = a · b with ab, find a and b such that the difference ab is as small as possible.

For n square, the solution is just the square root of n; for instance, with n = 36, a = b = 6. Otherwise, a and b will be the two divisors nearest the square root; for instance, with n = 60, a = 10 and b = 6.

Your task is to write a program to find a and b as described above; use your program to find the nearly square divisors of n = 224403121196654400. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.

Pages: 1 2

### 17 Responses to “Nearly-Square Divisors”

1. Ernie said

Another method is to start with the square root of n, rounded down, then decrement by one until a divisor is found. On my computer this took 0.6 sec.
The pitfall with any method is to deal with a prime n.

2. Paul said

In Python. The divisors are scanned for the maximum value less equal isqrt(n).

```def nsd(n):
s = isqrt(n)
a = max(d for d in divisors(n) if d <= s)
return a, n // a
```
3. Matt W said

@Ernie had the same idea! Here’s a solution that does that in Racket.

```(define (near-squares n)
(let loop ((a (floor (sqrt n))))
(let ((b (/ n a)))
(if (and (integer? b)
(= n (* a b)))
(list (inexact->exact a)
(inexact->exact b))
(loop (sub1 a))))))
```
4. Milbrae said

I came up with the same idea as the previous commenters. For large n with a higher number of divisors this may become incredibly time consuming, e.g. n = 100#. So I guess there’s a quicker solution than simply testing all possible numbers from floor(sqrt(n)) downward.

5. Milbrae said

Python:

```import math

def nsd(n):
s = (int)(math.sqrt(n))
while n % s:
s -= 1
return s, n // s

print (nsd(224403121196654400))
```

Output in approx. 0.012s: (473670855, 473753280)

Anyhow, this solution is based on the slow brute-force method. I’m almost certain there’s a much faster number-theoretical way to solve this for huge n.

6. Paul said

@Ernie and @Milbrae. The number N=224403121196654400 has an incredibly large number of divisors and is the worst-case scenario for a solution using divisors and the best-case scenario for the solution walking down from the root. Now try to use both methods in the range [N-10, N+10). Timing is 90 ms for the divisor solution and the other solution takes forever.

7. Ernie said

Paul: That is point I was trying to make in my first post: if N has very few divisors, decrementing from the sqrt may take forever, and if N is prime, well don’t ask.

8. matthew said

224403121196654400 is the 22nd “colossally abundant” number (and mentioned here: https://programmingpraxis.com/2016/04/22/gcd-sum/#comment-59608)

Finding nearly-square divisors seems about as hard as factoring (if we have a method for finding them we can apply it recursively to get a full factorization, and if we have a factorization, finding such divisors is straightforward).

Fermat factorization seems relevant too, but maybe only applies to odd numbers.

9. Milbrae said

@Paul: You may have a point there. Anyhow, consider n having at least 30 distinct prime factors, each appearing once. n will then have 2**30 factors. You wouldn’t want to try almost all of them to find the correct NSD.
There’s a similar problem on Project Euler (https://projecteuler.net/problem=266). I don’t think the people you have already solved this did it through simply factoring.

10. matthew said

@Milbrae: Good point, actually it’s a variant on the knapsack problem – find a maximal subset (or sub-multiset) of prime factors of n, the sum of whose logs is less than log n / 2 (and there are ways of solving the knapsack problem that are better than testing each alternative).

11. Paul said

@Milbrae. Thanks for pointing me to Euler 266. I solved for the product of first 30 primes with divisors and it took about 4 hours. The Euler problem needs the first 42 primes. That will indeed need another approach.

12. Milbrae said

@Matthew: A varint of knapsack? Could you elaborate on this?

13. matthew said

@Milbrae: One form of knapsack problem is: given a set of real numbers S, and a limit k, find the maximal subset sum less (or equal) to k (ie. find a subset with sum less than or equal to k, and greater or equal to the sum of any other such subset). This is exactly what we want for the nearly-square root of a product of single primes (as in the Project Euler problem) – S is the set of logs of the primes, and k is the log of the square root.

Here’s a simple Haskell program that computes the log of the nearly square root of the product of the first 20 primes:

```-- k target s: return largest subset sum of s <= target
k target [] = 0
k target (a:s) =
if a > target then k target s
else max (a + k (target-a) s) (k target s)
primes = [2,3,5,7,11,13,17,19,23,29,
31,37,41,43,47,53,59,61,67,71]
s = map log(reverse primes)
target = (sum s)/2
r = k target s
main = print target >> print (k target s)
```

Scan the tree of subsets, discarding a branch which cannot lead to a solution (here just subsets that would be too large – we could also keep track of the best solution found and discard branches that cannot lead to a better solution).

14. Milbrae said

@Matthew: Wow!! Works fine and extremely fast. I assume correctly the second line of the output is the exponent to e to get the result? I know practically nothing about Haskell.

15. […] a recent exercise, we discussed the problem of computing the nearly square divisors of a composite number n = a · b, […]

16. […] nearly square divisors exercise has generated a considerable amount of interest, and several excellent solutions in the comments. […]

17. akshaya pandey said

package test;

import java.util.*;

public class SquareDivisors {

public static int squareDivisor(int n){
int diff=10000;
if (n < 0){
return diff;
}
List pairs = new ArrayList();
//no need to look beyond sqrt hence check for (i*i )<=n
for(int i=1;(i*i)<=n;i++){
if (n%i==0){
System.out.println("Pairs : "+i +" "+ (n/i));
}
}
for(Pairs pair: pairs){
diff= (diff < (Math.abs(pair.a-pair.b)))?diff:(Math.abs(pair.a-pair.b));
}
return diff;
}

public static void main(String args[]){
System.out.println("Diff : "+squareDivisor(36));

}
}
class Pairs{
public Pairs(int a,int b){
this.a=a;
this.b=b;
}
public int a,b;
}