1: package de.luschny.math.factorial;
2:
3: import de.luschny.math.Xmath;
4: import de.luschny.math.arithmetic.Xint;
5: import de.luschny.math.primes.IPrimeIteration;
6: import de.luschny.math.primes.PrimeSieve;
7:
8: public class FactorialPrimeSwingList implements IFactorialFunction
9: {
10:
11: private int[][] listPrime;
12: private int[] listLength;
13: private int[] tower;
14: private int[] bound;
15:
16: public FactorialPrimeSwingList()
17: {
18: }
19:
20: public final String getName()
21: {
22: return "PrimeSwingList ";
23: }
24:
25: public Xint factorial(int n)
26: {
27: // For very small n the 'NaiveFactorial' is OK.
28: if (n < 20) { return Xmath.Factorial(n); }
29:
30: // log2n = floor(log2(n));
31: int log2n = 31 - Integer.numberOfLeadingZeros(n);
32: int j = log2n, hN = n;
33:
34: this.listPrime = new int[log2n][];
35: listLength = new int[log2n];
36: bound = new int[log2n];
37: tower = new int[log2n + 1];
38:
39: while (true)
40: {
41: tower[j] = hN;
42: if (hN == 1)
43: {
44: break;
45: }
46: bound[--j] = hN / 3;
47: listPrime[j] = new int[omegaSwingHighBound(hN)];
48: hN /= 2;
49: }
50: tower[0] = 2;
51:
52: primeFactors(n);
53: return iterQuad().shiftLeft(n - Integer.bitCount(n));
54: }
55:
56: private Xint iterQuad()
57: {
58: int init = listLength[0] == 0 ? 1 : 3;
59: Xint fact = Xint.valueOf(init);
60:
61: int listLen = listPrime.length;
62:
63: for (int i = 1; i < listLen; i++)
64: {
65: fact = (fact.square()).multiply(listPrime[i], listLength[i]);
66: }
67: return fact;
68: }
69:
70: private void primeFactors(int n)
71: {
72: int maxBound = n / 3;
73: int lastList = listPrime.length - 1;
74: int start = tower[1] == 2 ? 1 : 0;
75:
76: PrimeSieve sieve = new PrimeSieve(n);
77:
78: for (int list = start; list < listPrime.length; list++)
79: {
80: IPrimeIteration pIter = sieve.getIteration(tower[list] + 1, tower[list + 1]);
81:
82: for (int prime : pIter)
83: {
84: listPrime[list][listLength[list]++] = prime;
85: if (prime > maxBound)
86: {
87: continue;
88: }
89:
90: int np = n;
91: do
92: {
93: int k = lastList;
94: int q = np /= prime;
95:
96: do
97: {
98: if ((q & 1) == 1)
99: {
100: listPrime[k][listLength[k]++] = prime;
101: }
102: }
103: while (((q /= 2) > 0) && (prime <= bound[--k]));
104:
105: }
106: while (prime <= np);
107: }
108: }
109: }
110:
111: private int omegaSwingHighBound(int n)
112: {
113: return n < 4 ? 6 : (int) (2.0 * (Math.floor(Math.sqrt(n) + n / (Math.log(n) * 1.4426950408889634 - 1.0))));
114: }
115:
116: } // endOfFactorialPrimeSwingList