January 24, 2014 9:00 AM
Today’s exercise sounds, from the title, like another exercise involving prime numbers and integer factorization, and it is, but it’s really a puzzle from the realm of recreational mathematics: Given a positive integer n, compute the prime factorization, including multiplicities, of n! = 1 · 2 · … · n. You should be able to handle very large n, which means that you should not compute the factorial before computing the factors, as the intermediate result will be extremely large.
Your task is to write the function described above. 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.
Posted by programmingpraxis
Categories: Exercises
Tags:
Mobile Site | Full Site
Get a free blog at WordPress.com Theme: WordPress Mobile Edition by Alex King.
#include
#include
// reducing numbers from biggest to 2
// 16 -> 2*8 … 8 -> 2*4 … 4 -> 2*2 …
void main(int argc, char **argv) {
int n = atoi(argv[1]);
int *p = (int *) malloc(sizeof(int) * (n + 1));
int i, j, d;
for(i = 0; i 1; i–)
if(p[i]) {
for(j = i + i, d = 2; j <= n; j += i, d++) {
if(p[j]) {
p[i] += p[j];
p[d] += p[j];
p[j] = 0;
}
}
}
printf("1");
for(i = 2; i <= n; i++)
if(p[i])
printf(" * %i^%i", i, p[i]);
printf("\n");
}
By josejuan (@__josejuan__) on January 24, 2014 at 9:54 AM
// (sorry...) #include <stdio.h> #include <stdlib.h> // reducing numbers from biggest to 2 // 16 -> 2*8 ... 8 -> 2*4 ... 4 -> 2*2 ... void main(int argc, char **argv) { int n = atoi(argv[1]); int *p = (int *) malloc(sizeof(int) * (n + 1)); int i, j, d; for(i = 0; i <= n; i++) p[i] = 1; for(i = n; i > 1; i--) if(p[i]) { for(j = i + i, d = 2; j <= n; j += i, d++) { if(p[j]) { p[i] += p[j]; p[d] += p[j]; p[j] = 0; } } } printf("1"); for(i = 2; i <= n; i++) if(p[i]) printf(" * %i^%i", i, p[i]); printf("\n"); }By josejuan (@__josejuan__) on January 24, 2014 at 9:56 AM
Recursivelly adding factorizations
import Data.List import Data.Monoid import Data.Numbers.Primes facfac 1 = Factors [(1, 1)] facfac n = factorization n `mappend` facfac (n - 1) -- Factorizations Monoid data Factorization a = Factors [(a, a)] deriving Show factorization = Factors . map (\ps -> (head ps, length ps)) . group . primeFactors instance Integral a => Monoid (Factorization a) where mempty = Factors [] mappend (Factors []) fs = fs mappend fs (Factors []) = fs mappend xfs@(Factors (xf@(x, a):xs)) yfs@(Factors (yf@(y, b):ys)) = case x `compare` y of EQ -> Factors ((x, a + b):ms) LT -> Factors (xf:ns) GT -> Factors (yf:os) where Factors ms = mappend (Factors xs) (Factors ys) Factors ns = mappend (Factors xs) yfs Factors os = mappend xfs (Factors ys)By josejuan (@__josejuan__) on January 24, 2014 at 10:40 AM
(Excuse me, I can’t modify previous post)
In C version, we can reduce main loop from “i = n” to “i = n >> 1”.
That version has O(n log log n) and are not needed primes calculation nor mul nor div operations (space is O(n)).
By josejuan (@__josejuan__) on January 24, 2014 at 9:33 PM
Did in Go using my Project Euler library. Since it’s designed to solve PE problems everything is done with type uint64.
package main import ( . "fmt" . "github.com/dgottner/euler" "os" "strconv" ) func main() { var n, p uint64 if len(os.Args) < 2 { Println("Usage: factfact n") return } else { val, err := strconv.ParseInt(os.Args[1], 10, 64) if err != nil { Println(err.Error()) return } n = uint64(val) } Printf("%d! = 1", n) for p = 2; p <= n; p = NextPrime(p) { exp := n / p for fac := p * p; fac <= n; fac *= p { exp += n / fac } Printf(" * %d", p) if exp != 1 { Printf("^%d", exp) } } Println() }By David on January 24, 2014 at 10:32 PM
That’s a fun little problem. Here is my solution. It is quick and dirty, not efficient.
from itertools import count, takewhile def primes(): yield 2 primes = [2] n = 3 while True: if all(n % p for p in primes): yield n primes.append(n) n += 2 def pff(n): """prime factorization of n! returns a list of pairs (prime, power). """ up_to_n = lambda iter: takewhile(lambda p: p <= n, iter) return [(prime, sum(n/w for w in up_to_n(prime**power for power in count(1)))) for prime in up_to_n(primes())] from math import factorial from operator import mul print factorial(100) == reduce(mul, (a**b for (a, b) in pff(100)))By kernelbob on January 25, 2014 at 8:05 AM
Assuming `primes` is already implemented, returning a stream of prime numbers in order, in Haskell,
By Will Ness on January 26, 2014 at 11:08 AM
Another Python version. The optimization suggested by Will Ness did not work (in Python).
def facfac(n): for p in primes(n): nfacs = 0 f = n // p while f: nfacs += f f //= p yield p, nfacsBy Paul on January 27, 2014 at 12:30 PM