about summary refs log tree commit diff
path: root/vendor/modernc.org/mathutil/primes.go
diff options
context:
space:
mode:
Diffstat (limited to 'vendor/modernc.org/mathutil/primes.go')
-rw-r--r--vendor/modernc.org/mathutil/primes.go331
1 files changed, 331 insertions, 0 deletions
diff --git a/vendor/modernc.org/mathutil/primes.go b/vendor/modernc.org/mathutil/primes.go
new file mode 100644
index 0000000..ec1e5f5
--- /dev/null
+++ b/vendor/modernc.org/mathutil/primes.go
@@ -0,0 +1,331 @@
+// Copyright (c) 2014 The mathutil Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package mathutil // import "modernc.org/mathutil"
+
+import (
+	"math"
+)
+
+// IsPrimeUint16 returns true if n is prime. Typical run time is few ns.
+func IsPrimeUint16(n uint16) bool {
+	return n > 0 && primes16[n-1] == 1
+}
+
+// NextPrimeUint16 returns first prime > n and true if successful or an
+// undefined value and false if there is no next prime in the uint16 limits.
+// Typical run time is few ns.
+func NextPrimeUint16(n uint16) (p uint16, ok bool) {
+	return n + uint16(primes16[n]), n < 65521
+}
+
+// IsPrime returns true if n is prime. Typical run time is about 100 ns.
+func IsPrime(n uint32) bool {
+	switch {
+	case n&1 == 0:
+		return n == 2
+	case n%3 == 0:
+		return n == 3
+	case n%5 == 0:
+		return n == 5
+	case n%7 == 0:
+		return n == 7
+	case n%11 == 0:
+		return n == 11
+	case n%13 == 0:
+		return n == 13
+	case n%17 == 0:
+		return n == 17
+	case n%19 == 0:
+		return n == 19
+	case n%23 == 0:
+		return n == 23
+	case n%29 == 0:
+		return n == 29
+	case n%31 == 0:
+		return n == 31
+	case n%37 == 0:
+		return n == 37
+	case n%41 == 0:
+		return n == 41
+	case n%43 == 0:
+		return n == 43
+	case n%47 == 0:
+		return n == 47
+	case n%53 == 0:
+		return n == 53 // Benchmarked optimum
+	case n < 65536:
+		// use table data
+		return IsPrimeUint16(uint16(n))
+	default:
+		mod := ModPowUint32(2, (n+1)/2, n)
+		if mod != 2 && mod != n-2 {
+			return false
+		}
+		blk := &lohi[n>>24]
+		lo, hi := blk.lo, blk.hi
+		for lo <= hi {
+			index := (lo + hi) >> 1
+			liar := liars[index]
+			switch {
+			case n > liar:
+				lo = index + 1
+			case n < liar:
+				hi = index - 1
+			default:
+				return false
+			}
+		}
+		return true
+	}
+}
+
+// IsPrimeUint64 returns true if n is prime. Typical run time is few tens of µs.
+//
+// SPRP bases: http://miller-rabin.appspot.com
+func IsPrimeUint64(n uint64) bool {
+	switch {
+	case n%2 == 0:
+		return n == 2
+	case n%3 == 0:
+		return n == 3
+	case n%5 == 0:
+		return n == 5
+	case n%7 == 0:
+		return n == 7
+	case n%11 == 0:
+		return n == 11
+	case n%13 == 0:
+		return n == 13
+	case n%17 == 0:
+		return n == 17
+	case n%19 == 0:
+		return n == 19
+	case n%23 == 0:
+		return n == 23
+	case n%29 == 0:
+		return n == 29
+	case n%31 == 0:
+		return n == 31
+	case n%37 == 0:
+		return n == 37
+	case n%41 == 0:
+		return n == 41
+	case n%43 == 0:
+		return n == 43
+	case n%47 == 0:
+		return n == 47
+	case n%53 == 0:
+		return n == 53
+	case n%59 == 0:
+		return n == 59
+	case n%61 == 0:
+		return n == 61
+	case n%67 == 0:
+		return n == 67
+	case n%71 == 0:
+		return n == 71
+	case n%73 == 0:
+		return n == 73
+	case n%79 == 0:
+		return n == 79
+	case n%83 == 0:
+		return n == 83
+	case n%89 == 0:
+		return n == 89 // Benchmarked optimum
+	case n <= math.MaxUint16:
+		return IsPrimeUint16(uint16(n))
+	case n <= math.MaxUint32:
+		return ProbablyPrimeUint32(uint32(n), 11000544) &&
+			ProbablyPrimeUint32(uint32(n), 31481107)
+	case n < 105936894253:
+		return ProbablyPrimeUint64_32(n, 2) &&
+			ProbablyPrimeUint64_32(n, 1005905886) &&
+			ProbablyPrimeUint64_32(n, 1340600841)
+	case n < 31858317218647:
+		return ProbablyPrimeUint64_32(n, 2) &&
+			ProbablyPrimeUint64_32(n, 642735) &&
+			ProbablyPrimeUint64_32(n, 553174392) &&
+			ProbablyPrimeUint64_32(n, 3046413974)
+	case n < 3071837692357849:
+		return ProbablyPrimeUint64_32(n, 2) &&
+			ProbablyPrimeUint64_32(n, 75088) &&
+			ProbablyPrimeUint64_32(n, 642735) &&
+			ProbablyPrimeUint64_32(n, 203659041) &&
+			ProbablyPrimeUint64_32(n, 3613982119)
+	default:
+		return ProbablyPrimeUint64_32(n, 2) &&
+			ProbablyPrimeUint64_32(n, 325) &&
+			ProbablyPrimeUint64_32(n, 9375) &&
+			ProbablyPrimeUint64_32(n, 28178) &&
+			ProbablyPrimeUint64_32(n, 450775) &&
+			ProbablyPrimeUint64_32(n, 9780504) &&
+			ProbablyPrimeUint64_32(n, 1795265022)
+	}
+}
+
+// NextPrime returns first prime > n and true if successful or an undefined value and false if there
+// is no next prime in the uint32 limits. Typical run time is about 2 µs.
+func NextPrime(n uint32) (p uint32, ok bool) {
+	switch {
+	case n < 65521:
+		p16, _ := NextPrimeUint16(uint16(n))
+		return uint32(p16), true
+	case n >= math.MaxUint32-4:
+		return
+	}
+
+	n++
+	var d0, d uint32
+	switch mod := n % 6; mod {
+	case 0:
+		d0, d = 1, 4
+	case 1:
+		d = 4
+	case 2, 3, 4:
+		d0, d = 5-mod, 2
+	case 5:
+		d = 2
+	}
+
+	p = n + d0
+	if p < n { // overflow
+		return
+	}
+
+	for {
+		if IsPrime(p) {
+			return p, true
+		}
+
+		p0 := p
+		p += d
+		if p < p0 { // overflow
+			break
+		}
+
+		d ^= 6
+	}
+	return
+}
+
+// NextPrimeUint64 returns first prime > n and true if successful or an undefined value and false if there
+// is no next prime in the uint64 limits. Typical run time is in hundreds of µs.
+func NextPrimeUint64(n uint64) (p uint64, ok bool) {
+	switch {
+	case n < 65521:
+		p16, _ := NextPrimeUint16(uint16(n))
+		return uint64(p16), true
+	case n >= 18446744073709551557: // last uint64 prime
+		return
+	}
+
+	n++
+	var d0, d uint64
+	switch mod := n % 6; mod {
+	case 0:
+		d0, d = 1, 4
+	case 1:
+		d = 4
+	case 2, 3, 4:
+		d0, d = 5-mod, 2
+	case 5:
+		d = 2
+	}
+
+	p = n + d0
+	if p < n { // overflow
+		return
+	}
+
+	for {
+		if ok = IsPrimeUint64(p); ok {
+			break
+		}
+
+		p0 := p
+		p += d
+		if p < p0 { // overflow
+			break
+		}
+
+		d ^= 6
+	}
+	return
+}
+
+// FactorTerm is one term of an integer factorization.
+type FactorTerm struct {
+	Prime uint32 // The divisor
+	Power uint32 // Term == Prime^Power
+}
+
+// FactorTerms represent a factorization of an integer
+type FactorTerms []FactorTerm
+
+// FactorInt returns prime factorization of n > 1 or nil otherwise.
+// Resulting factors are ordered by Prime. Typical run time is few µs.
+func FactorInt(n uint32) (f FactorTerms) {
+	switch {
+	case n < 2:
+		return
+	case IsPrime(n):
+		return []FactorTerm{{n, 1}}
+	}
+
+	f, w := make([]FactorTerm, 9), 0
+	for p := 2; p < len(primes16); p += int(primes16[p]) {
+		if uint(p*p) > uint(n) {
+			break
+		}
+
+		power := uint32(0)
+		for n%uint32(p) == 0 {
+			n /= uint32(p)
+			power++
+		}
+		if power != 0 {
+			f[w] = FactorTerm{uint32(p), power}
+			w++
+		}
+		if n == 1 {
+			break
+		}
+	}
+	if n != 1 {
+		f[w] = FactorTerm{n, 1}
+		w++
+	}
+	return f[:w]
+}
+
+// PrimorialProductsUint32 returns a slice of numbers in [lo, hi] which are a
+// product of max 'max' primorials. The slice is not sorted.
+//
+// See also: http://en.wikipedia.org/wiki/Primorial
+func PrimorialProductsUint32(lo, hi, max uint32) (r []uint32) {
+	lo64, hi64 := int64(lo), int64(hi)
+	if max > 31 { // N/A
+		max = 31
+	}
+
+	var f func(int64, int64, uint32)
+	f = func(n, p int64, emax uint32) {
+		e := uint32(1)
+		for n <= hi64 && e <= emax {
+			n *= p
+			if n >= lo64 && n <= hi64 {
+				r = append(r, uint32(n))
+			}
+			if n < hi64 {
+				p, _ := NextPrime(uint32(p))
+				f(n, int64(p), e)
+			}
+			e++
+		}
+	}
+
+	f(1, 2, max)
+	return
+}
generated by cgit-pink 1.4.1
rss feed of @hanemile@chaos.social mastodon lieu webring search engine XXIIVV webring activitypub