#include #include typedef unsigned int num; #define NUMFMT "%u" typedef struct numptr { num n; struct numptr *next; } numptr; typedef struct primeent { num p; numptr *numlist; } primeent; #define NHASHSIZE 300007 num nhashtable[NHASHSIZE]; #define PHASHSIZE 20011 primeent phashtable[PHASHSIZE]; #define SEQLIM 200000 num sequence[SEQLIM]; int nseq = 0; static int isprime (num x) { num d; if (x < 2) return 0; if (x == 2) return 1; if (x % 2 == 0) return 0; for (d=3; d*d <= x; d += 2) { if (x % d == 0) return 0; } return 1; } static num maxprimediv (num x) { num p; if (x == 1) return 1; while (x % 2 == 0) x /= 2; if (x == 1) return 2; for (p = 3; p * p <= x; p += 2) { while (x % p == 0) x /= p; if (x == 1) return p; } return x; } static num gcd (num a, num b) { if (a > b) return gcd(b,a); if (a == 0) return b; return gcd(b%a, a); } static numptr *numptr_alloc (void) { static numptr *freelist; static int nfree = 0; if (!nfree) { nfree = 1024*1024-8; freelist = (numptr *) malloc (nfree * sizeof (numptr)); if (!freelist) { fprintf (stderr, "Out of memory\n"); exit (1); } } nfree--; return freelist++; } static int find_num (num n) { unsigned int hashloc = n % NHASHSIZE; for (;;) { if (nhashtable[hashloc] == 0) return 0; if (nhashtable[hashloc] == n) return 1; hashloc++; if (hashloc >= NHASHSIZE) hashloc = 0; } } static void register_num (num n) { unsigned int hashloc = n % NHASHSIZE; for (;;) { if (nhashtable[hashloc] == 0) { nhashtable[hashloc] = n; return; } hashloc++; if (hashloc >= NHASHSIZE) hashloc = 0; } } static numptr *find_pmults (num p) { unsigned int hashloc = p % PHASHSIZE; for (;;) { if (phashtable[hashloc].p == 0) return (numptr *) NULL; if (phashtable[hashloc].p == p) return phashtable[hashloc].numlist; hashloc++; if (hashloc >= PHASHSIZE) hashloc = 0; } } static void add_pmults (num p, num n) { unsigned int hashloc = p % PHASHSIZE; numptr *s = numptr_alloc (); for (;;) { if (phashtable[hashloc].p == 0) { phashtable[hashloc].p = p; phashtable[hashloc].numlist = s; s->n = n; s->next = (numptr *) NULL; return; } if (phashtable[hashloc].p == p) { s->n = n; s->next = phashtable[hashloc].numlist; phashtable[hashloc].numlist = s; return; } hashloc++; if (hashloc >= PHASHSIZE) hashloc = 0; } } static void add_all_pmults (num n) { num x = n; num p; if (x % 2 == 0) add_pmults ((num) 2, n); while (x % 2 == 0) x /= 2; for (p = 3; p * p <= x; p += 2) { if (x % p == 0) add_pmults (p, n); while (x % p == 0) x /= p; } if (x != 1) add_pmults (x, n); } static int product_covers (num n, num a, num d, num *p_rem) { num g; g = gcd (a, n); a /= g; n /= g; g = gcd (d, n); d /= g; n /= g; if (n > 1) return 0; *p_rem = a * d; return 1; } static int have_num (num n) { if (n == 1) return 1; if (isprime (n)) return 1; if (find_num (n)) return 1; return 0; } static int find_combo (num n) { num bigp; num bigp2; numptr *s; numptr *s2; num rem; bigp = maxprimediv (n); bigp2 = maxprimediv (n/bigp); if (bigp2 == 1) { fprintf (stderr, "find_combo with n prime\n"); exit (1); } if (product_covers (n, bigp, bigp2, &rem) && have_num (rem)) { return 1; } for (s2 = find_pmults (bigp2); s2; s2 = s2->next) { if (product_covers (n, bigp, s2->n, &rem) && have_num (rem)) { return 1; } } for (s = find_pmults (bigp); s; s = s->next) { bigp2 = maxprimediv (n/gcd (n, s->n)); if (product_covers (n, s->n, bigp2, &rem) && have_num (rem)) { return 1; } for (s2 = find_pmults (bigp2); s2; s2 = s2->next) { if (product_covers (n, s->n, s2->n, &rem) && have_num (rem)) { return 1; } } } return 0; } static void print_factorization (num n) { num p; int pow; int first = 1; if (n % 2 == 0) { if (!first) printf ("*"); first = 0; printf ("("); printf (NUMFMT, (num) 2); printf (")"); pow = 0; while (n % 2 == 0) { pow++; n /= 2; } if (pow > 1) printf ("^%d", pow); } for (p = 3; p * p <= n; p += 2) { if (n % p == 0) { if (!first) printf ("*"); first = 0; printf ("("); printf (NUMFMT, p); printf (")"); pow = 0; while (n % p == 0) { pow++; n /= p; } if (pow > 1) printf ("^%d", pow); } } if (n > 1) { if (!first) printf ("*"); first = 0; printf ("("); printf (NUMFMT, n); printf (")"); } } int main (int ac, char **av) { num n; printf ("1, (1)\n"); nseq = 0; for (n=2; nseq < SEQLIM; n++) { if (isprime (n)) continue; if (!find_combo (n)) { sequence[nseq] = n; register_num (n); add_all_pmults (n); nseq++; printf (NUMFMT, n); printf (", "); print_factorization (n); printf ("\n"); fflush (stdout); } } return 0; }