Re: [S] Prime divisors

Paul A Tukey (paul@bellcore.com)
Wed, 16 Sep 1998 18:27:15 -0400 (EDT)


This discussion has been fun.

Seems to me we've been gradually heading toward
computing the prime factorization of a number.
That is, a collection of prime numbers (with possible
duplication) whose product is the given number.

A recursive layer on top of Frank Harrell's factorize() does
it -- but the code below uses a slightly shortened version
of factorize() that only returns the smallest divisor > 1.

fac <- function(n) {
p <- n/(z <- 1:floor(sqrt(n)))
z <- z[trunc(p) == p]
c(z, rev(n/z))[2]
}

pfac <- function(n, nn = 0)
{
if(nn == 0)
pfac(n, fac(n))
else if(n <= nn)
nn
else c(nn, pfac(n/nn))
}

Note that prod(pfac(n)) == n.

Now I'm sure someone can write a more elegant version.
Also, recursion is probably neither memory-efficient
nor CPU-efficient in Splus.

-- Paul Tukey
Bellcore

-----------------------------------------------------------------------
This message was distributed by s-news@wubios.wustl.edu. To unsubscribe
send e-mail to s-news-request@wubios.wustl.edu with the BODY of the
message: unsubscribe s-news