safety0ff:

> I think we need a solid integer math module more than anything on this list.
> I find myself often reimplementing ilog2, isqrt, isPowerOf2, etc, etc.

I implemented many of those in my old D1 nonstandard library. I still port some of them to D2 now and then. So I agree they are needed.


> I think we should be careful about adding an isPrime method, I think adding an isProbablePrime plus the prime sieve should cover most use cases.

This is an interesting opinion. Can you explain why an isPrime is less useful  than a isProbablePrime (like this: http://en.wikipedia.org/wiki/Miller-Rabin_primality_test ) plus a sieve? Can't you have both a isPrime and a isProbablePrime in Phobos?

Bye,
bearophile

On Sunday, 30 March 2014 at 23:40:14 UTC, bearophile wrote:
>> I think we should be careful about adding an isPrime method, I think adding an isProbablePrime plus the prime sieve should cover most use cases.
>
> This is an interesting opinion. Can you explain why an isPrime is less useful  than a isProbablePrime (like this: http://en.wikipedia.org/wiki/Miller-Rabin_primality_test ) plus a sieve? Can't you have both a isPrime and a isProbablePrime in Phobos?
>
> Bye,
> bearophile

You can have both, it's not less useful.
I just wonder if there are real-world use-cases worth supporting in Phobos where you want need a one-off primality test where a definite answer is require.
I.E. Are there non-toy uses to support inclusion?

safety0ff:

> I.E. Are there non-toy uses to support inclusion?

I have several non-toy uses for the GCD, binomial, and the permutations/combinations/pairwise ranges.

The other little discrete numerical functions like the integer square root, etc, are useful sufficiently often.

The other two (isPrime and a sieve) are less practical, so their inclusion is more debatable. Sometimes you want some small primes for hashing, or in the unittests of a Rational type, and few other situations. But D is for mathematics usage too. Given the commonality of prime numbers in discrete mathematics (example: if you want to study certain graphs), I think such two/three functions/generators are acceptable, if they aren't too much complex.

Bye,
bearophile

Am Sun, 30 Mar 2014 23:05:55 +0000
schrieb "safety0ff" <safety0ff.dev@gmail.com>:

> I think we need a solid integer math module more than anything on
> this list.
> I find myself often reimplementing ilog2, isqrt, isPowerOf2, etc,
> etc.

More or less real world stuff:

minMax(x, low, high)
  Returns x clamped the low and high boundary so that
  low <= x <= high.
setMax(x, high)
  Same but only clamps high values.
setMin(x, low)
  Same but only clamps low values.
isWithin(x, low, high)
  Checks if x is within a given range.
isPowerOfTwo(x)
  Basically checks if only one bit is set.
nextPowerOfTwo(x)
  Returns x rounded up to be a power of two.
bitMask(T)(T bitNum)
  Returns cast(T) 1 << x. The result is always the type of the
  parameter.
bitMaskRange(firstBitNum, bitCount)
  Generates a mask of consecutive set bits starting with
  'firstBitNum' and counting 'bitCount' bits.
log2(x)
  Returns the highest power of 2 in 'x'.
roundUp(x, multiple)
  Rounds 'x' up to the next multiple of 'multiple'.
  A version that works with pointers is also useful.
roundDown(x, multiple)
  Rounds 'x' down to the next multiple of 'multiple'.
roundUpDiv(x, divisor)
  Normal integer division rounds down. This version rounds up.
sqr(x)
  x*x
KiB(x)
  x*1024
MiB(x)
  x*1024*1024
isEven(x)
  !(x & 1)

From Project Euler solutions:

sum_1_to_n(n)
  Returns the sum of all integers [1 .. n] in O(1) time
  complexity.
sum_1_to_n_dividable_by(n, d)
  Same as above but counts only values dividable by 'd'.
generalized_pentagonal_number(n)
  Don't know what this is, but it is O(1) as well.
gcd(u, v)
  Greatest common divisor. Some copied&pasted algorithm, that
  works by comparing set bits in 'u' and 'v'. Very fast for
  being implemented in standard C without assembly hacks.
are_coprime(u, v)
  Returns true, iff the greatest common divisor of 'u' and 'v'
  is 1.
pascal_row(n)
  Returns the n-th row of Pascal's triangle. This gives you
  "n choose k" for a fixed 'n' and all k in 0<=k<=n.
  (This function allocates an integer array.)

Also there is a Pascal's triangle struct, that allows marching through the triangle quickly using methods like 'rightUp' or 'left'. Pascal's triangle operations are useful in combinatorics.

----------8<-------------

/**
 * The ultimate tool to move around in Pascal's triangle. Pascal(0, 0) is the
 * top of the triangle with the value of 1.
 */
struct Pascal
{
	this(ℕ n, ℕ k) {
		// internally we offset these by +1, to simplify the math
		++n; ++k;
		_n = n;
		if (k > n / 2) {
			_k = _n;
			left(_n - k);
		} else {
			right(k - 1);
		}
	}

	ℕ left(ℕ amount) {
		while (amount--) {
			_value *= --_k;
			_value /= _n - _k;
		}
		return _value;
	}

	ℕ right(ℕ amount) {
		while (amount--) {
			_value *= _n - _k;
			_value /= _k++;
		}
		return _value;
	}

	ℕ rightDown(ℕ amount) {
		while (amount--) {
			_value *= _n++;
			_value /= _k++;
		}
		return _value;
	}

	ℕ leftDown(ℕ amount) {
		while (amount--) {
			_value *= _n++;
			_value /= _n - _k;
		}
		return _value;
	}

	ℕ leftUp(ℕ amount) {
		while (amount--) {
			_value *= --_k;
			_value /= --_n;
		}
		return _value;
	}

	ℕ rightUp(ℕ amount) {
		while (amount--) {
			_value *= _n - _k;
			_value /= --_n;
		}
		return _value;
	}

	@property
	ℕ value() { return _value; }

	alias value this;

	private:

	ℕ _n = 1, _k = 1;
	ℕ _value = 1;
}

-- 
Marco

I think we should get started with a Phobos integer module proposal/implementation.

What is the best way to move forward?
A community editable document with a repo to merge our contributions to?

I am not familiar with using dub.

>>      if (p ^^ 2 > b)
>>            si = si.max(p ^^ 2 - b);

Can't that be written as:

auto t = p ^^ 2 - b;
if (t > 0 && t > si) si = t;


--
Ziad


On Tue, Mar 18, 2014 at 7:23 AM, bearophile <bearophileHUGS@lycos.com>wrote:

> There is a efficient Sieve implementation in C++ here:
>
> http://code.activestate.com/recipes/577966-even-faster- prime-generator/?in=lang-cpp
>
> There are of course far faster implementations, but its performance is not bad, while being still simple and quite short.
>
> A D implementation (it's not a final version because the global enums should be removed and replaced with run-time variables or template arguments. And something different from just counting must be added):
>
>
> import std.stdio, std.algorithm, std.typecons;
>
> enum uint MAX_N = 1_000_000_000U;
> enum uint RT_MAX_N = 32_000; // square of max prime under this limit
> should exceed MAX_N.
> enum uint B_SIZE = 20_000;   // not sure what optimal value for this is;
>                              // currently must avoid overflow when squared.
>
> // mod 30, (odd) primes have remainders 1,7,11,13,17,19,23,29.
> // e.g. start with mark[B_SIZE/30]
> // and offset[] = {1, 7, 11, ...}
> // then mark[i] corresponds to 30 * (i / 8) + b - 1 + offset[i % 8].
> Tuple!(uint, size_t, uint) calcPrimes() pure nothrow {
>     // Assumes p, b odd and p * p won't overflow.
>     static void crossOut(in uint p, in uint b, bool[] mark)
>     pure nothrow {
>         uint si = (p - (b % p)) % p;
>         if (si & 1)
>             si += p;
>         if (p ^^ 2 > b)
>             si = si.max(p ^^ 2 - b);
>
>         for (uint i = si / 2; i < B_SIZE / 2; i += p)
>             mark[i] = true;
>     }
>
>     uint pCount = 1; uint lastP = 2;
>     // Do something with first prime (2) here.
>
>     uint[] smallP = [2];
>
>     bool[B_SIZE / 2] mark = false;
>     foreach (immutable uint i; 1 .. B_SIZE / 2) {
>         if (!mark[i]) {
>             pCount++; lastP = 2 * i + 1;
>             // Do something with lastP here.
>
>             smallP ~= lastP;
>             if (lastP ^^ 2 <= B_SIZE)
>                 crossOut(2 * i + 1, 1, mark);
>         } else
>             mark[i] = false;
>     }
>
>     for (uint b = 1 + B_SIZE; b < MAX_N; b += B_SIZE) {
>         for (uint i = 1; smallP[i] ^^ 2 < b + B_SIZE; ++i)
>             crossOut(smallP[i], b, mark);
>         foreach (immutable uint i; 0 .. B_SIZE / 2) {
>             if (!mark[i]) {
>                 pCount++; lastP = 2 * i + b;
>                 // Do something with lastP here.
>
>                 if (lastP <= RT_MAX_N)
>                     smallP ~= lastP;
>             } else
>                 mark[i] = false;
>         }
>     }
>
>     return tuple(pCount, smallP.length, lastP);
> }
>
> void main() {
>     immutable result = calcPrimes;
>
>     writeln("Found ", result[0], " primes in total.");
>     writeln("Recorded ", result[1], " small primes, up to ", RT_MAX_N);
>     writeln("Last prime found was ", result[2]);
> }
>
>
> Is it a good idea to add a simple but reasonably fast Sieve implementation to Phobos? I have needed a prime numbers lazy range, and a isPrime() function numerous times. (And as for std.numeric.fft, people that need even more performance will use code outside Phobos).
>
> Bye,
> bearophile
>

On Sunday, 30 March 2014 at 23:05:57 UTC, safety0ff wrote:

> I find myself often reimplementing ilog2, isqrt, isPowerOf2,
isPowerOf2 { x => (x && !(x & (x-1))) }

On Monday, 31 March 2014 at 06:10:08 UTC, Dominikus Dittes Scherkl wrote:
> On Sunday, 30 March 2014 at 23:05:57 UTC, safety0ff wrote:
>
>> I find myself often reimplementing ilog2, isqrt, isPowerOf2,
> isPowerOf2 { x => (x && !(x & (x-1))) }

The point is that it should be in the standard library instead of copy pasted every time I have a short 1-file program that requires it.

On Mon, 2014-03-31 at 00:55 +0000, bearophile wrote:
> safety0ff:
> 
> > I.E. Are there non-toy uses to support inclusion?
> 
> I have several non-toy uses for the GCD, binomial, and the permutations/combinations/pairwise ranges.

Finance, bioinformatics, signal processing, are some of the areas I know of where all this sort of stuff is useful. Python, Mathematica, Julia all provide these either fast using hardware numbers or to arbitrary accuracy using software numbers.

> The other little discrete numerical functions like the integer square root, etc, are useful sufficiently often.
> 
> The other two (isPrime and a sieve) are less practical, so their inclusion is more debatable. Sometimes you want some small primes for hashing, or in the unittests of a Rational type, and few other situations. But D is for mathematics usage too. Given the commonality of prime numbers in discrete mathematics (example: if you want to study certain graphs), I think such two/three functions/generators are acceptable, if they aren't too much complex.

Or for force cracking of encodings? ;-)

-- 
Russel. ============================================================================= Dr Russel Winder      t: +44 20 7585 2200   voip: sip:russel.winder@ekiga.net 41 Buckmaster Road    m: +44 7770 465 077   xmpp: russel@winder.org.uk London SW11 1EN, UK   w: www.russel.org.uk  skype: russel_winder