[R] All the products of common factors

[Stavros Macrakis]

On Mon, Feb 23, 2009 at 9:52 PM, Fox, Gordon <gfox at cas.usf.edu> wrote:
> This is a seemingly simple problem - hopefully someone can help.
> Problem: we have two integers. We want (1) all the common factors, and
> (2) all the possible products of these factors. We know how to get (1),
> but can't figure out a general way to get (2).

Your problem statement is rather complicated, but I believe it reduces
to simply:

   Find the factors of the greatest common divisor of A and B.

Note: *factors*, not *prime factors*.  Finding the prime factors first
is unnecessary.

Euclid's algorithm straightforwardly finds GCDs:

   gcd <- function(a,b) {if (a>b) gcd(b,a) else if (a==0) b else gcd(a,b%%a)}
   gcd(40,80) => 40

Now find all the factors of 40.  Might as well use a simple, naive
algorithm for such small numbers:

   num <- 40
   which(num/1:num == floor(num/1:num))
   [1]  1  2  4  5  8 10 20 40

Or if you want to get really fancy and efficient:

   smallfacs <-  which( (q <- num/1:sqrt(num)) == floor(q) )
   c(smallfacs,rev(num/smallfacs))

You do *not* need to factor the original numbers.  You do *not* need
to iterate through combinations.

Treating your specification as a draft program leads to some rather
complicated and inefficient code...:

system.time({commfac <- c(1,rep(2,20))   # about 1e6
   cfun <- function(i) { combn(commfac,i,prod) }
   L <- lapply(as.list(1:length(commfac)),cfun)
   unique(sort(unlist(L)))})
   user  system elapsed
  24.71    0.01   24.89
> system.time({num <- 2^20; which(num/1:num == floor(num/1:num))})
   user  system elapsed
   0.17    0.00    0.17

> num <- prod(1:16)
> format(num,digits=20)
[1] "20922789888000"
> system.time( { smallfacs <-  which( (q <- num/1:sqrt(num)) == floor(q) ); c(smallfacs,rev(num/smallfacs))})
   user  system elapsed
   0.64    0.06    0.70

Hope this helps,

         -s