SCALA |
C# |
GO |
SWING SCALA |
SWING C# |
SWING GO |
class Swing(n: Int) {
// OEIS: A056040 Swing numbers
def name() = "Swing"
private var primes: Primes = null
private var factors: Array[Long] = null
if (n >= Small.oddSwing.length)
{
primes = new Primes(n)
factors = new Array[Long](n)
}
def value(m: Int = n):BigInt =
{
return oddSwing(m) << MathFun.bitCount(m)
}
def oddSwing(k: Int): BigInt =
{
if (k > n)
{
throw new IllegalArgumentException(
"k has to be >= n, but was " + k)
}
if (k < Small.oddSwing.length)
return BigInt(Small.oddSwing(k));
val rootK = MathFun.floorSqrt(k)
var i = 0
primes.iterator(3, rootK, p =>
{
var q = k / p
while (q > 0)
{
if ((q & 1) == 1) {factors(i) = p; i += 1}
q /= p
}
})
primes.iterator(rootK + 1, k / 3, p =>
{
if (((k / p) & 1) == 1) {factors(i) = p; i += 1}
})
primes.iterator(k / 2 + 1, k, p =>
{
factors(i) = p; i += 1
})
return MathFun.product(factors,0,i)
}
}
|
public class Swing
{
public string Name
{
get { return "Swing"; }
}
private uint n;
private Primes primes;
private Factors factors;
// OEIS: A056040 Swinging factorial
public Swing(uint n)
{
if (n >= Small.OddSwing.Length)
{
primes = new Primes(n);
factors = new Factors(n);
}
this.n = n;
}
public BigInt Value()
{
return OddSwing(n) << (int)MathFun.BitCount(n);
}
public BigInt OddSwing(uint k)
{
if (k > n)
{
throw new ArgumentOutOfRangeException(
"k has to be >= n, but was " + k);
}
if (k < Small.OddSwing.Length) return Small.OddSwing[k];
uint rootN = MathFun.FloorSqrt(k);
factors.Init();
factors.SetMax(rootN);
primes.Iterator(3, rootN, p =>
{
uint q = k;
while ((q /= p) > 0)
if ((q & 1) == 1) { factors.Add(p); }
});
factors.SetMax(k / 3);
primes.Iterator(rootN + 1, k / 3, p =>
{
if (((k / p) & 1) == 1) { factors.Add(p); }
});
factors.SetMax(k);
primes.Iterator(k / 2 + 1, k, p =>
{
factors.Add(p);
});
return factors.Product();
}
}
|
type swing struct {
primes *primes
factors []uint64
}
func makeSwing(n uint) (ps *swing) {
ps = new(swing)
ps.primes = makePrimes(uint64(n))
if n >= uint(len(smallOddSwing)) {
ps.factors = make([]uint64, n)
}
return
}
func (ps *swing) swing(m uint) *big.Int {
if uint64(m) > ps.primes.sieveLen {
return nil
}
r := ps.oddSwing(m)
return r.Lsh(r, bitCount(uint64(m)))
}
func (ps *swing) oddSwing(k uint) *big.Int {
if k < uint(len(smallOddSwing)) {
return big.NewInt(smallOddSwing[k])
}
rootK := floorSqrt(k)
var i int
ps.primes.iterator(3, uint64(rootK), func(p uint64) {
q := uint64(k) / p
for q > 0 {
if q & 1 == 1 {
ps.factors[i] = p
i++
}
q /= p
}
})
ps.primes.iterator(uint64(rootK+1),uint64(k/3),func(p uint64){
if (uint64(k) / p & 1) == 1 {
ps.factors[i] = p
i++
}
})
ps.primes.iterator(uint64(k/2 + 1), uint64(k), func(p uint64){
ps.factors[i] = p
i++
})
return product(ps.factors[0:i])
}
|
PRIME FACTORIAL SCALA |
PRIME FACTORIAL C# |
PRIME FACTORIAL GO |
class PrimeFactorial(n: Int) {
def name() = "PrimeFactorial"
var swing: Swing = new Swing(n)
def value(k : Int = n): BigInt =
{
if (k > n)
{
throw new IllegalArgumentException(
"k has to be >= n, but was " + k)
}
val exp2 = k - MathFun.bitCount(k)
if (k < Small.oddFactorial.length)
return BigInt(Small.oddFactorial(k)) << exp2
return oddFactorial(k) << exp2
}
private def oddFactorial(n: Int): BigInt =
{
if (n < Small.oddFactorial.length)
{
return BigInt(Small.oddFactorial(n))
}
val of = oddFactorial(n / 2)
return (of * of) * swing.oddSwing(n)
}
}
|
public class PrimeFactorial
{
public string Name
{
get { return "PrimeFactorial"; }
}
public PrimeFactorial() { }
private Swing swing;
// OEIS: A000142 Factorial
public BigInt Value(int n)
{
if (n < 0)
{
throw new ArgumentOutOfRangeException(
"Factorial: n has to be >= 0, but was " + n);
}
int exp2 = (int)(n - MathFun.BitCount((uint)n));
if (n < Small.OddFactorial.Length)
return new BigInt(Small.OddFactorial[n]) << exp2;
swing = new Swing((uint)n);
return OddFactorial((uint)n) << exp2;
}
private BigInt OddFactorial(uint n)
{
if (n < Small.OddFactorial.Length)
{
return Small.OddFactorial[n];
}
return BigInt.Pow(OddFactorial(n / 2), 2) * swing.OddSwing(n);
}
}
|
func (ps *swing) primeFactorial(k uint) (r *big.Int) {
if uint64(k) > ps.primes.sieveLen {
return nil
}
r = ps.oddFactorial(k)
return r.Lsh(r, k - bitCount(uint64(k)))
}
func (ps *swing) oddFactorial(n uint) *big.Int {
if n < uint(len(smallOddFactorial)) {
return big.NewInt(int64(smallOddFactorial[n]))
}
of := ps.oddFactorial(n / 2)
return of.Mul(of.Mul(of, of), ps.oddSwing(n))
}
|
val start2 = System.nanoTime
var ps = new PrimeFactorial(n)
val v = ps.value(n)
val stop2 = System.nanoTime
val sec2 = DoubleFormat.fix2((stop2 - start2)/1e9) + " sec"
n = 234567 -> 68,5 sec
|
var watch = new System.Diagnostics.Stopwatch();
watch.Reset(); watch.Start();
var pFact = new PrimeFactorial();
pFact.Value(n);
watch.Stop(); Console.WriteLine(pFact.Name + " " + watch.Elapsed);
n = 234567 -> PrimeFactorial 00:00:38.6056398
|
func main() {
const n = 234567
start := time.Nanoseconds()
fs := makeSwing(n)
a := fs.primeFactorial(n)
stop := time.Nanoseconds()
fmt.Printf("n = %d -> PrimeFactorial %.2f sec\n",
n, float64(stop - start) / 1e9)
dtrunc := int64(float64(a.BitLen()) * .30103) - 10
var first, rest big.Int
rest.Exp(first.SetInt64(10), rest.SetInt64(dtrunc), nil)
first.Quo(a, &rest)
fstr := first.String()
fmt.Printf("%d! begins %s... and has %d digits.\n",
n, fstr, int64(len(fstr)) + dtrunc)
}
|
SWING FACTORIAL SCALA |
SWING FACTORIAL C# |
SWING FACTORIAL GO |
class SwingFactorial(n: Int) {
def name() = "SwingFactorial"
def value(): BigInt =
{
if (n < 0)
{
throw new IllegalArgumentException(
"Factorial: n has to be >= 0, but was " + n)
}
ndiv2OddFact = BigInt(1)
ndiv4OddFact = ndiv2OddFact
return oddFactorial(n) << (n - MathFun.bitCount(n))
}
private def oddFactorial(n: Int): BigInt =
{
val oddFact =
if (n < Small.oddFactorial.length)
{
BigInt(Small.oddFactorial(n))
}
else
{
val of = oddFactorial(n / 2)
(of * of) * oddSwing(n)
}
ndiv4OddFact = ndiv2OddFact
ndiv2OddFact = oddFact
return oddFact
}
private def oddSwing(n: Int): BigInt =
{
if (n < Small.oddSwing.length)
{
return BigInt(Small.oddSwing(n))
}
val len = if ((n % 4) != 2) (n - 1) / 4 + 1 else (n - 1) / 4
val high = n - ((n + 1) & 1)
val ndiv4 = n / 4
val oddFact = if (ndiv4 < Small.oddFactorial.length)
BigInt(Small.oddFactorial(ndiv4)) else ndiv4OddFact
return product(high, len) / oddFact
}
private def product(m: Int, len: Int): BigInt =
{
if (len == 1) return BigInt(m)
if (len == 2) {val M = m.toLong; return BigInt(M * (M - 2))}
val hlen = len >>> 1
return product(m - hlen * 2, len - hlen) * product(m, hlen)
}
private var ndiv4OddFact = BigInt(1)
private var ndiv2OddFact = BigInt(1)
}
|
public class SwingFactorial
{
public string Name
{
get { return "SwingFactorial"; }
}
public SwingFactorial() { }
private BigInt oddFactNdiv4, oddFactNdiv2;
private int SMALLSWING = Small.OddSwing.Length;
private int SMALLFACT = Small.OddFactorial.Length;
public BigInt Value(int n)
{
if (n < 0)
{
throw new ArgumentOutOfRangeException(
"Factorial: n has to be >= 0, but was " + n);
}
oddFactNdiv4 = oddFactNdiv2 = BigInt.One;
return OddFactorial(n) << (n-(int)MathFun.BitCount((uint)n));
}
private BigInt OddFactorial(int n)
{
BigInt oddFact;
if (n < SMALLFACT)
{
oddFact = Small.OddFactorial[n];
}
else
{
BigInt sqrOddFact = OddFactorial(n / 2);
oddFact = BigInt.Pow(sqrOddFact, 2) * OddSwing(n);
}
oddFactNdiv4 = oddFactNdiv2;
oddFactNdiv2 = oddFact;
return oddFact;
}
private BigInt OddSwing(int n)
{
if (n < SMALLSWING) return Small.OddSwing[n];
int len = (n - 1) / 4;
if ((n % 4) != 2) len++;
int high = n - ((n + 1) & 1);
int ndiv4 = n / 4;
BigInt oddFact = ndiv4 < SMALLFACT ?
Small.OddFactorial[ndiv4] : oddFactNdiv4;
return Product(high, len) / oddFact;
}
private static BigInt Product(int m, int len)
{
if (len == 1) return new BigInt(m);
if (len == 2) return new BigInt((long)m * (m - 2));
int hlen = len >> 1;
return Product(m - hlen * 2, len - hlen) * Product(m, hlen);
}
}
|
func swingFactorial(n uint) (r *big.Int) {
var oddFactNDiv2, oddFactNDiv4 big.Int
// closes on oddFactNDiv2, oddFactNDiv4
oddSwing := func(n uint) (r *big.Int) {
if n < uint(len(smallOddSwing)) {
return big.NewInt(smallOddSwing[n])
}
length := (n - 1) / 4
if n%4 != 2 {
length++
}
high := n - (n+1)&1
ndiv4 := n / 4
var oddFact big.Int
if ndiv4 < uint(len(smallOddFactorial)) {
oddFact.SetInt64(smallOddFactorial[ndiv4])
r = &oddFact
} else {
r = &oddFactNDiv4
}
return oddFact.Quo(oddProduct(high, length), r)
}
var oddFactorial func(uint) *big.Int
// closes on oddFactNDiv2, oddFactNDiv4, oddSwing, and itself
oddFactorial = func(n uint) (oddFact *big.Int) {
if n < uint(len(smallOddFactorial)) {
oddFact = big.NewInt(smallOddFactorial[n])
} else {
oddFact = oddFactorial(n / 2)
oddFact.Mul(oddFact.Mul(oddFact, oddFact), oddSwing(n))
}
oddFactNDiv4.Set(&oddFactNDiv2)
oddFactNDiv2.Set(oddFact)
return oddFact
}
oddFactNDiv2.SetInt64(1)
oddFactNDiv4.SetInt64(1)
r = oddFactorial(n)
return r.Lsh(r, n-bitCount(uint64(n)))
}
func oddProduct(m, length uint) *big.Int {
switch length {
case 1:
return big.NewInt(int64(m))
case 2:
var mb big.Int
mb.SetInt64(int64(m))
mb2 := big.NewInt(int64(m - 2))
return mb2.Mul(&mb, mb2)
}
hlen := length / 2
h := oddProduct(m-hlen*2, length-hlen)
return h.Mul(h, oddProduct(m, hlen))
}
|
val start = System.nanoTime
var fs = new SwingFactorial(n)
val a = fs.value()
val stop = System.nanoTime
val sec = DoubleFormat.fix2((stop - start)/1e9) + " sec"
n = 234567 -> 82,3 sec
|
watch.Reset(); watch.Start();
var sFact = new SwingFactorial();
sFact.Value(n);
watch.Stop();
Console.WriteLine(sFact.Name + " " + watch.Elapsed);
n = 234567 -> SwingFactorial 00:00:49.5329715
|
func testSwingFactorial() {
const n = 234567
start := time.Nanoseconds()
a := swingFactorial(n)
stop := time.Nanoseconds()
fmt.Printf("n = %d -> SwingFactorial %.2f sec\n",
n, float64(stop-start)/1e9)
dtrunc := int64(float64(a.BitLen())*.30103) - 10
var first, rest big.Int
rest.Exp(first.SetInt64(10), rest.SetInt64(dtrunc), nil)
first.Quo(a, &rest)
fstr := first.String()
fmt.Printf("%d! begins %s... and has %d digits.\n",
n, fstr, int64(len(fstr))+dtrunc)
}
|
BINOMIAL SCALA |
BINOMIAL C# |
BINOMIAL GO |
class Binomial(n: Int) {
val primes = new Primes(n)
var factors: Array[Long] = new Array[Long](n)
val rootN = MathFun.floorSqrt(n)
def value(j: Int): BigInt =
{
var k = j; var i = 0
if (0 > k || k > n) return BigInt(0)
if (k > n / 2) {k = n - k; }
if (k < 3)
{
if (k == 0) return BigInt(1)
if (k == 1) return BigInt(n)
val N = n.toLong
if (k == 2) return BigInt((N * (N - 1)) / 2)
}
primes.iterator(2, rootN, p =>
{
var r = 0; var N = n; var K = k
while (N > 0)
{
r = if ((N % p) < (K % p + r)) 1 else 0
if (r == 1) { factors(i) = p; i += 1 }
N /= p; K /= p
}
})
primes.iterator(rootN + 1, n / 2, p =>
{
if (n % p < k % p)
{
factors(i) = p; i += 1
}
})
primes.iterator(n - k + 1, n, p =>
{
factors(i) = p; i += 1
})
return MathFun.product(factors,0,i)
}
}
|
public class Binomial
{
Primes prime;
Factors factors;
uint rootN;
uint n;
public Binomial(uint n)
{
prime = new Primes(n);
factors = new Factors(n);
rootN = MathFun.FloorSqrt(n);
this.n = n;
}
// OEIS: A007318 Pascal's triangle
public BigInt Value(uint k)
{
if (0 > k || k > n) return BigInt.Zero;
if (k > n / 2) { k = n - k; }
if (k < 3)
{
if (k == 0) return BigInt.One;
if (k == 1) return new BigInt(n);
if (k == 2) return new BigInt(((ulong)n * (n - 1)) / 2);
}
factors.Init();
factors.SetMax(rootN);
prime.Iterator(2, rootN, p =>
{
uint r = 0, N = n, K = k;
while (N > 0)
{
r = (N % p) < (K % p + r) ? 1u : 0u;
if (r == 1)
{
factors.Add(p);
}
N /= p; K /= p;
}
});
factors.SetMax(n / 2);
prime.Iterator(rootN + 1, n / 2, p =>
{
if (n % p < k % p)
{
factors.Add(p);
}
});
factors.SetMax(n);
prime.Iterator(n - k + 1, n, p =>
{
factors.Add(p);
});
return factors.Product();
}
}
|
func (p *primes) binomial(n, k uint) *big.Int {
if uint64(n) > p.sieveLen {
return nil
}
var r big.Int
if k > n {
return &r
}
if k > n/2 {
k = n - k
}
if k < 3 {
switch k {
case 0:
return r.SetInt64(1)
case 1:
return r.SetInt64(int64(n))
case 2:
var n1 big.Int
return r.Rsh(r.Mul(r.SetInt64(int64(n)), n1.SetInt64(int64(n-1))), 1)
}
}
var i int
rootN := uint64(floorSqrt(n))
factors := make([]uint64, n)
p.iterator(2, rootN, func(p uint64) {
var r, nn, kk uint64 = 0, uint64(n), uint64(k)
for nn > 0 {
if nn%p < kk%p+r {
r = 1
factors[i] = p
i++
} else {
r = 0
}
nn /= p
kk /= p
}
})
p.iterator(rootN+1, uint64(n/2), func(p uint64) {
if uint64(n)%p < uint64(k)%p {
factors[i] = p
i++
}
})
p.iterator(uint64(n-k+1), uint64(n), func(p uint64) {
factors[i] = p
i++
})
return product(factors[0:i])
}
|
val start = System.nanoTime
val binomial = new Binomial(n)
for (k <- 33 until n by 10000)
val expected = binomial.value(k)
val stop = System.nanoTime
val sec = DoubleFormat.fix2((stop - start)/1e9) + " sec "
n = 234567 -> 5,45 sec
|
watch.Start();
var binomial = new Binomial(234567);
for (uint k = 33; k < n; k += 10000)
var expected = binomial.Value(k);
watch.Stop();
n = 234567 -> 00:00:03.45
|
func testBinomial() {
const n = 234567
start := time.Nanoseconds()
p := makePrimes(n)
var a *big.Int
for k := uint(33); k <= n; k += 10000 {
b := p.binomial(n, k)
if k == 120033 {
a = b
}
}
stop := time.Nanoseconds()
fmt.Printf("n = %d -> %d Binomials %.2f sec\n",
n, (n-32)/10000, float64(stop-start)/1e9)
dtrunc := int64(float64(a.BitLen())*.30103) - 10
var first, rest big.Int
rest.Exp(first.SetInt64(10), rest.SetInt64(dtrunc), nil)
first.Quo(a, &rest)
fstr := first.String()
fmt.Printf("C(%d, 120033) begins %s... and has %d digits.\n",
n, fstr, int64(len(fstr))+dtrunc)
}
|
DOUBLE-FACTORIAL SCALA |
DOUBLE-FACTORIAL C# |
DOUBLE-FACTORIAL GO |
class DoubleFactorial(n: Int) {
def name() = "DoubleFactorial"
val N = if ((n & 1) == 0) n / 2 else n + 1
val swing: Swing = new Swing(N)
// OEIS: A006882 Double factorials n!!: a(n) = n*a(n-2)
def value: BigInt =
{
var dblFact = if (n < Small.oddDoubleFactorial.length)
BigInt(Small.oddDoubleFactorial(n)) else oddFactorial(N, n)
if ((n & 1) == 0)
{
val exp2 = n - MathFun.bitCount(n / 2)
dblFact = dblFact << exp2
}
return dblFact
}
def oddFactorial(n: Int, m: Int): BigInt =
{
if (n < Small.oddFactorial.length)
{
return BigInt(Small.oddFactorial(n))
}
var oddFact = oddFactorial(n / 2, m)
if (n < m) { oddFact *= oddFact }
return oddFact * swing.oddSwing(n)
}
}
|
public class DoubleFactorial
{
private Swing swing;
public DoubleFactorial() { }
// OEIS: A006882 Double factorials n!!: a(n) = n*a(n-2).
public BigInt Value(uint n)
{
if (n < 0)
{
throw new ArgumentOutOfRangeException(
"DoubleFactorial: n has to be >= 0, but was " + n);
}
BigInt dblFact;
uint N = ((n & 1) == 0) ? n / 2 : n + 1;
if (n < Small.OddDoubleFactorial.Length)
{
dblFact = (BigInt)Small.OddDoubleFactorial[n];
}
else
{
if (N >= Small.OddSwing.Length)
{
swing = new Swing(N);
}
dblFact = OddFactorial(N, n);
}
if ((n & 1) == 0)
{
int exp2 = (int)(n - MathFun.BitCount(n / 2));
dblFact = dblFact << exp2;
}
return dblFact;
}
private BigInt OddFactorial(uint n, uint m)
{
if (n < Small.OddFactorial.Length)
{
return Small.OddFactorial[n];
}
var oddSwing = n < Small.OddSwing.Length
? Small.OddSwing[n] : swing.OddSwing(n);
var oddFact = OddFactorial(n / 2, m);
if (n < m) { oddFact *= oddFact; }
return oddFact * oddSwing;
}
}
|
// returns nil if sieve not big enough
func (p *primes) doubleFactorial(n uint) (r *big.Int) {
nEven := n&1 == 0
if n < uint(len(smallOddDoubleFactorial)) {
r = big.NewInt(smallOddDoubleFactorial[n])
} else {
var nn uint
if nEven {
nn = n / 2
} else {
nn = n + 1
}
if uint64(nn) > p.sieveLen && nn > uint(len(smallOddSwing)){
return nil
}
r = p.oddDoubleFactorial(nn, n)
}
if nEven {
r.Lsh(r, n - bitCount(uint64(n/2)))
}
return
}
func (p *primes) oddDoubleFactorial(n, m uint) *big.Int {
if n < uint(len(smallOddFactorial)) {
return big.NewInt(smallOddFactorial[n])
}
of := p.oddDoubleFactorial(n/2, m)
if n < m {
of.Mul(of, of)
}
return of.Mul(of, p.oddSwing(n))
}
|
val start = System.nanoTime
val df = new DoubleFactorial(n)
val v = df.value
val stop = System.nanoTime
val sec = DoubleFormat.fix2((stop - start)/1e9) + " sec "
n = 234566 -> 14,8 sec
n = 234567 -> 20,6 sec
|
watch.Start();
var df = new DoubleFactorial();
var b = df.Value(n);
watch.Stop();
Console.Write("primedouble " + watch.Elapsed + " ");
n = 234566 -> 00:00:09.28
n = 234566 -> 00:00:13.19
|
func testDoubleFactorial() {
const nEven = 234566
const nOdd = 234567
start := time.Nanoseconds()
p := makePrimes(nOdd+1)
sieveDone := time.Nanoseconds()
aEven := p.doubleFactorial(nEven)
evenDone := time.Nanoseconds()
aOdd := p.doubleFactorial(nOdd)
oddDone := time.Nanoseconds()
fmt.Printf("sieve to %d -> %.2f sec\n",
nOdd+1, float64(sieveDone-start)/1e9)
fmt.Printf(" %d!! -> %.2f sec\n",
nEven, float64(evenDone-sieveDone)/1e9)
fmt.Printf(" %d!! -> %.2f sec\n",
nOdd, float64(oddDone-evenDone)/1e9)
dtrunc := int64(float64(aEven.BitLen())*.30103) - 10
var first, rest big.Int
rest.Exp(first.SetInt64(10), rest.SetInt64(dtrunc), nil)
first.Quo(aEven, &rest)
fstr := first.String()
fmt.Printf("%d!! begins %s... and has %d digits.\n",
nEven, fstr, int64(len(fstr))+dtrunc)
dtrunc = int64(float64(aOdd.BitLen())*.30103) - 10
rest.Exp(first.SetInt64(10), rest.SetInt64(dtrunc), nil)
first.Quo(aOdd, &rest)
fstr = first.String()
fmt.Printf("%d!! begins %s... and has %d digits.\n",
nOdd, fstr, int64(len(fstr))+dtrunc)
}
|
ERATOSTHENES SCALA |
ERATOSTHENES C# |
ERATOSTHENES GO |
// Prime number sieve, Eratosthenes (276-194 b.t. )
// This implementation considers only multiples of primes
// greater than 3, so the smallest value has to be mapped to 5.
// Note: There is no multiplication operation in this function
// and no call to a sqrt function.
class Primes(val n: Int) {
private val primesOnBits =
Array[Long](1762821248, 848611808, 3299549660L, 2510511646L)
private val bitsPerInt: Int = 32
private val mask: Int = bitsPerInt - 1
private val log2Int: Int = 5
private var isComposite: Array[Long] = null
if (n < 386) {isComposite = primesOnBits; }
isComposite = new Array[Long]((n / (3 * bitsPerInt)) + 1: Int)
var d1 = 8; var d2 = 8; var p1 = 3; var p2 = 7; var s = 7
var s2 = 3; var l = 0; var c = 1; var max = n / 3; var inc = 0
var toggle: Boolean = false
while (s < max) // -- scan the sieve
{
// -- if a prime is found ...
if ((isComposite(l >> log2Int) & (1 << (l & mask))) == 0)
{
inc = p1 + p2 // -- ... cancel its multiples
// for (c <- s until max by inc) // slow!
c = s; while (c < max) // -- set c as composite
{
isComposite(c >> log2Int) |= 1 << (c & mask)
c += inc
}
c = s + s2; while (c < max)
{
isComposite(c >> log2Int) |= 1 << (c & mask)
c += inc
}
}
l = l + 1
toggle = !toggle
if (toggle) {s += d2; d1 += 16; p1 += 2; p2 += 2; s2 = p2}
else {s += d1; d2 += 8; p1 += 2; p2 += 6; s2 = p1}
}
def iterator(min: Int, max: Int, visitor: (Int => Unit)): Unit =
{
// isComposite[0] ... isComposite[n] includes
// 5 <= prime numbers <= 96*(n+1) + 1
if (max > n) throw new IllegalArgumentException(
"Max not in sieve.")
val min2 = if (min < 2) 2 else min
if (min2 == 2) visitor(2)
if (min2 <= 3) visitor(3)
val absPos = (min2 + (min2 + 1) % 2) / 3 - 1
var index = absPos / bitsPerInt
var bitPos = absPos % bitsPerInt
var prime = 5+3*(bitsPerInt*index+bitPos)-(bitPos & 1)
var toggle: Boolean = (bitPos & 1) == 1
while (prime <= max)
{
var bitField = isComposite(index) >> bitPos
index += 1
while (bitPos < bitsPerInt)
{
if ((bitField & 1) == 0)
{
visitor(prime)
}
toggle = !toggle
prime = if (toggle) prime + 2 else prime + 4
if (prime > max) return
bitField >>= 1
bitPos += 1
}
bitPos = 0
}
}
}
|
/// Prime number sieve, Eratosthenes (276-194 b.t. )
/// This implementation considers only multiples of primes
/// greater than 3, so the smallest value has to be mapped to 5.
/// Note: There is no multiplication operation in this function
/// and no call to a sqrt function.
public class Primes
{
const int bitsPerInt = 32;
const int mask = bitsPerInt - 1;
const int log2Int = 5;
private uint sieveLen = 0;
private static uint[] PrimesOnBits = {
1762821248u, 848611808u, 3299549660u, 2510511646u };
private uint[] isComposite;
public delegate void Visitor(uint x);
public Primes(uint n)
{
sieveLen = n;
if (n < 386) { isComposite = PrimesOnBits; return; }
isComposite = new uint[(n / (3 * bitsPerInt)) + 1];
int d1 = 8, d2 = 8, p1 = 3, p2 = 7, s = 7, s2 = 3;
int l = 0, c = 1, max = (int)n / 3, inc;
bool toggle = false;
while (s < max) // -- scan the sieve
{
// -- if a prime is found ...
if ((isComposite[l >> log2Int] & (1u << (l++ & mask))) == 0)
{
inc = p1 + p2; // -- ... cancel its multiples
for (c = s; c < max; c += inc)
{ // -- ... set c as composite
isComposite[c >> log2Int] |= 1u << (c & mask);
}
for (c = s + s2; c < max; c += inc)
{
isComposite[c >> log2Int] |= 1u << (c & mask);
}
}
if (toggle = !toggle)
{ s += d2; d1 += 16; p1 += 2; p2 += 2; s2 = p2; }
else { s += d1; d2 += 8; p1 += 2; p2 += 6; s2 = p1; }
}
}
public void Iterator(uint min, uint max, Visitor visitor)
{
// isComposite[0] ... isComposite[n] includes
// 5 <= primes numbers <= 96*(n+1) + 1
if (min < 2) min = 2;
if (max > sieveLen) throw
new ArgumentOutOfRangeException("Max larger than sieve.");
if (min <= 2) visitor(2);
if (min <= 3) visitor(3);
int absPos = (int)((min + (min + 1) % 2) / 3 - 1);
int index = absPos / bitsPerInt;
int bitPos = absPos % bitsPerInt;
bool toggle = (bitPos & 1) == 1;
uint prime = (uint)(5+3*(bitsPerInt*index+bitPos)-(bitPos & 1));
while (prime <= max)
{
uint bitField = isComposite[index++] >> bitPos;
for (; bitPos < bitsPerInt; bitPos++)
{
if ((bitField & 1) == 0)
{
visitor(prime);
}
prime += (toggle = !toggle) ? 2u : 4u;
if (prime > max) return;
bitField >>= 1;
}
bitPos = 0;
}
}
}
|
// Prime number sieve, Eratosthenes (276-194 b.t. )
// This implementation considers only multiples of primes
// greater than 3, so the smallest value has to be mapped to 5.
// word size dependent constants
const (
bitsPerInt = 32
mask = bitsPerInt - 1
log2Int = 5
)
// holds completed sieve
type primes struct {
sieveLen uint
isComposite []uint32
}
// constructor, completes sieve.
func makePrimes(n uint) (ps *primes) {
ps = new(primes)
ps.sieveLen = n
if n < 386 {
ps.isComposite = []uint32{1762821248, 848611808,
3299549660, 2510511646}
return
}
ps.isComposite = make([]uint32, (n/(3*bitsPerInt))+1)
var (
d1, d2, p1, p2, s, s2 uint = 8, 8, 3, 7, 7, 3
l, c, max, inc uint = 0, 1, n / 3, 0
toggle bool
)
for s < max { // -- scan the sieve
// -- if a prime is found ...
if (ps.isComposite[l>>log2Int] & (1 << (l & mask))) == 0 {
inc = p1 + p2 // -- ... cancel its multiples
// -- ... set c as composite
for c = s; c < max; c += inc {
ps.isComposite[c>>log2Int] |= 1 << (c & mask)
}
for c = s + s2; c < max; c += inc {
ps.isComposite[c>>log2Int] |= 1 << (c & mask)
}
}
l++
toggle = !toggle
if toggle {
s += d2; d1 += 16; p1 += 2; p2 += 2; s2 = p2
} else {
s += d1; d2 += 8; p1 += 2; p2 += 6; s2 = p1
}
}
return
}
func (ps *primes) iterator(min, max uint, visitor func(uint)) bool {
// isComposite[0] ... isComposite[n] includes
// 5 <= primes numbers <= 96*(n+1) + 1
if min < 2 {
min = 2
}
if max > ps.sieveLen {
return false // Max larger than sieve
}
if min <= 2 {
visitor(2)
}
if min <= 3 {
visitor(3)
}
absPos := uint((min+(min+1)%2)/3 - 1)
index := absPos / bitsPerInt
bitPos := absPos % bitsPerInt
toggle := (bitPos & 1) == 1
prime := uint(5 + 3*(bitsPerInt*index+bitPos) - (bitPos & 1))
for prime <= max {
bitField := ps.isComposite[index] >> bitPos
index++
for ; bitPos < bitsPerInt; bitPos++ {
if (bitField & 1) == 0 {
visitor(prime)
}
toggle = !toggle
if toggle {
prime += 2
} else {
prime += 4
}
if prime > max {
return true
}
bitField >>= 1
}
bitPos = 0
}
return true
}
|
var count = 0
def visitor(prime: Int): Unit = {
count = count + 1
}
def sieve(): Unit = {
val n = 100000000 // Correct answer: 5,761,455
val start = System.nanoTime
val p = new Primes(n)
val stop = System.nanoTime
println("Time taken to sieve = %s sec".format((stop-start)/1e9))
val start2 = System.nanoTime
p.iterator(1, n, visitor)
val stop2 = System.nanoTime
println("Time taken to visit = %s sec".format((stop2-start2)/1e9))
println("Number of primes: " + count)
}
Time taken to sieve = 1.58 sec
Time taken to visit = 0.44 sec
Number of primes: 5761455
|
static int count = 0;
static void visitor(uint prime) {
count = count + 1;
}
static void PrimeSieveTest() {
var watch = new System.Diagnostics.Stopwatch();
var n = 100000000u; // Correct answer: 5,761,455
watch.Reset(); watch.Start();
var p = new Primes(n);
watch.Stop();
Console.WriteLine("Time taken to sieve " + watch.Elapsed);
watch.Reset(); watch.Start();
p.Iterator(1, n, visitor);
watch.Stop();
Console.WriteLine("Time taken to visit " + watch.Elapsed);
Console.WriteLine("Number of primes: " + count);
}
Time taken to sieve 00:00:01.10
Time taken to visit 00:00:00.49
Number of primes: 5761455
|
package main
import (
"fmt"
"time"
)
const n = 100000000 // Correct answer: 5,761,455
func main() {
var count uint
visitor := func(_ uint) {
count = count + 1
}
start := time.Nanoseconds()
p := makePrimes(n)
stop := time.Nanoseconds()
fmt.Printf("Time taken to sieve %.2f sec\n",
float64(stop-start)/1e9)
start = time.Nanoseconds()
p.iterator(1, n, visitor)
stop = time.Nanoseconds()
fmt.Printf("Time taken to visit %.2f sec\n",
float64(stop-start)/1e9)
fmt.Println("Number of primes: ", count)
}
|
MathFun SCALA |
MathFun C# |
MathFun GO |
object MathFun {
// Calibrate the treshhold
private val THRESHOLD_PRODUCT_SERIAL: Int = 28
def product(seq: Array[Long], start: Int, len: Int): BigInt =
{
if (len <= THRESHOLD_PRODUCT_SERIAL)
{
var sprod = BigInt(seq(start))
// for (i <- start + 1 until start + len) // slow!
val stop = start+len; var i = start+1; while(i < stop)
{
sprod *= BigInt(seq(i))
i += 1
}
return sprod
}
val halfLen = len >> 1
val lprod = product(seq, start, halfLen)
var rprod = product(seq, start + halfLen, len - halfLen)
return lprod * rprod
}
def bitCount(v: Int): Int =
{
var w = v - ((0xaaaaaaaa & v) >> 1)
w = (w & 0x33333333) + ((w >> 2) & 0x33333333)
w = w + (w >> 4) & 0x0f0f0f0f
w += w >> 8
w += w >> 16
return w & 0xff
}
def floorSqrt(n: Int): Int =
{
var a: Int = n
var b: Int = n
do {
a = b
b = (n / a + a) / 2
} while (b < a)
return a
}
}
|
public class MathFun
{
// Calibrate the treshholds
private const int THRESHOLD_PRODUCT_SERIAL = 24;
private const int THRESHOLD_PRODUCT_PARALLEL = 8000;
static public BigInt Product(long[] seq, int start, int len)
{
BigInt rprod, lprod;
int halfLen = len >> 1;
if (len > THRESHOLD_PRODUCT_PARALLEL)
{
rprod = BigInt.Zero; lprod = BigInt.Zero;
Parallel.Invoke(
() => { rprod = Product(seq, 0, halfLen); },
() => { lprod = Product(seq, halfLen, len - halfLen); }
);
return lprod * rprod;
}
if (len <= THRESHOLD_PRODUCT_SERIAL)
{
var sprod = new BigInt(seq[start]);
for (int i = start + 1; i < start + len; i++)
{
sprod *= seq[i];
}
return sprod;
}
lprod = Product(seq, start, halfLen);
rprod = Product(seq, start + halfLen, len - halfLen);
return lprod * rprod;
}
/// Bit count, loopfree!
public static uint BitCount(uint w)
{
w = w - ((0xaaaaaaaa & w) >> 1);
w = (w & 0x33333333) + ((w >> 2) & 0x33333333);
w = w + (w >> 4) & 0x0f0f0f0f;
w += w >> 8;
w += w >> 16;
return w & 0xff;
}
public static uint FloorSqrt(uint n)
{
uint a, b;
a = b = n;
do
{
a = b;
b = (n / a + a) / 2;
} while (b < a);
return a;
}
private MathFun() { }
}
}
|
const productSerialThreshold = 24
func product(seq []uint64) *big.Int {
if len(seq) <= productSerialThreshold {
var b big.Int
sprod := big.NewInt(int64(seq[0]))
for _, s := range seq[1:] {
b.SetInt64(int64(s))
sprod.Mul(sprod, &b)
}
return sprod
}
halfLen := len(seq) / 2
lprod := product(seq[0:halfLen])
rprod := product(seq[halfLen:])
return lprod.Mul(lprod, rprod)
}
func bitCount(w uint64) uint { // loopfree!
const (
ff = 1 << 64 - 1
mask1 = ff / 3
mask3 = ff / 5
maskf = ff / 17
maskp = maskf >> 3 & maskf
)
w -= w >> 1 & mask1
w = w & mask3 + w >> 2 & mask3
w = (w + w >> 4) & maskf
return uint(w * maskp >> 56)
}
func floorSqrt(n uint) uint {
for b := n; ; {
a := b
b = (n / a + a) / 2
if b >= a {
return a
}
}
return 0 // unreachable. required by current compiler.
}
|
~ |
~ |
~ |
object Small {
val oddSwing = Array[Long](
1,1,1,3,3,15,5,35,35,315,63,693,231,3003,429,6435,6435,
109395,12155,230945,46189,969969,88179,2028117,676039,
16900975,1300075,35102025,5014575,145422675,9694845,
300540195,300540195,9917826435L,583401555,20419054425L,
2268783825L,83945001525L,4418157975L,172308161025L,
34461632205L,1412926920405L,67282234305L,2893136075115L,
263012370465L,11835556670925L,514589420475L,24185702762325L,
8061900920775L,395033145117975L,15801325804719L,
805867616040669L,61989816618513L,3285460280781189L,
121683714103007L,6692604275665385L,956086325095055L,
54496920530418135L,1879204156221315L,110873045217057585L,
7391536347803839L,450883717216034179L,14544636039226909L,
916312070471295267L,916312070471295267L)
val oddFactorial = Array[Long](
1,1,1,3,3,15,45,315,315,2835,14175,155925,467775,
6081075,42567525,638512875,638512875,10854718875L,
97692469875L,1856156927625L,9280784638125L,194896477400625L,
2143861251406875L,49308808782358125L,147926426347074375L,
3698160658676859375L)
val oddDoubleFactorial = Array[Long](
1,1,1,3,1,15,3,105,3,945,15,10395,45,135135,315,2027025,315,
34459425,2835,654729075,14175,13749310575L,155925,
316234143225L,467775,7905853580625L,6081075,213458046676875L,
42567525,6190283353629375L,638512875,191898783962510625L,
638512875,6332659870762850625L,10854718875L)
}
|
/// Note: This is an attempt to optimize by accumulating
/// small products. Uses a trick I learned from Marco Bodrato.
class Factors
{
private long[] factors;
private long maxProd, prod;
private int count;
public Factors(uint n)
{
int mem = (int)(0.63 * n / System.Math.Log(n));
factors = new long[mem];
}
public void Init()
{
maxProd = 1;
prod = 1;
count = 0;
}
public void SetMax(long max)
{
maxProd = long.MaxValue / max;
if (prod >= maxProd)
{
factors[count++] = prod;
prod = 1;
}
}
public void Add(long prime)
{
if (prod < maxProd)
{ // the Bodrato trick:
prod *= prime; // accumulate small products first.
}
else
{
factors[count++] = prod;
prod = prime;
}
}
public BigInt Product()
{
factors[count++] = prod;
return MathFun.Product(factors, 0, count);
}
}
|
var smallOddSwing []int64 = []int64{1, 1, 1, 3, 3, 15, 5,
35, 35, 315, 63, 693, 231, 3003, 429, 6435, 6435,
109395, 12155, 230945, 46189, 969969, 88179, 2028117, 676039,
16900975, 1300075, 35102025, 5014575, 145422675, 9694845,
300540195, 300540195, 9917826435, 583401555, 20419054425,
2268783825, 83945001525, 4418157975, 172308161025,
34461632205, 1412926920405, 67282234305, 2893136075115,
263012370465, 11835556670925, 514589420475, 24185702762325,
8061900920775, 395033145117975, 15801325804719,
805867616040669, 61989816618513, 3285460280781189,
121683714103007, 6692604275665385, 956086325095055,
54496920530418135, 1879204156221315, 110873045217057585,
7391536347803839, 450883717216034179, 14544636039226909,
916312070471295267, 916312070471295267}
var smallOddFactorial []int64 = []int64{1, 1, 1, 3, 3,
15, 45, 315, 315, 2835, 14175, 155925, 467775,
6081075, 42567525, 638512875, 638512875, 10854718875,
97692469875, 1856156927625, 9280784638125, 194896477400625,
2143861251406875, 49308808782358125, 147926426347074375,
3698160658676859375}
var smallOddDoubleFactorial []int64 = []int64{1, 1, 1, 3, 1,
15, 3, 105, 3, 945, 15, 10395, 45, 135135, 315, 2027025, 315,
34459425, 2835, 654729075, 14175, 13749310575, 155925,
316234143225, 467775, 7905853580625, 6081075, 213458046676875,
42567525, 6190283353629375, 638512875, 191898783962510625,
638512875, 6332659870762850625, 10854718875}
|