Re: Pseudo Random Number Generation across PicoLisp implementations

2015-04-28 Thread Christophe Gragnic
On Tue, Apr 28, 2015 at 8:35 AM, Alexander Burger a...@software-lab.de wrote:
 Hi Christophe,
 […]
n = (Seed  32) % n;

 should it be OK?

Yes. It works. Great.

 Thanks for the input!

Thank you for considering it !

The 2)b) of my big email (april, 26th) of this thread was not answered but
I'll copy paste it in the seed thread since they seem related.


chri

-- 

http://profgra.org/lycee/ (site pro)
http://delicious.com/profgraorg (liens, favoris)
https://twitter.com/profgraorg
http://microalg.info
--
UNSUBSCRIBE: mailto:picolisp@software-lab.de?subject=Unsubscribe


Re: Pseudo Random Number Generation across PicoLisp implementations

2015-04-28 Thread Alexander Burger
Hi Christophe,

 Would it be possible to just cast (Seed  32) to a long:
 
 if ((n = evInt(ex.Cdr) + 1 - ((Number)x).Cnt) != 0)
 n = (long)(Seed  32) % n;

Yes, this seems to be a good idea.

The (long) cast should not be needed, as 'Seed' and 'n' already are long
numbers. So if I change this to

   n = (Seed  32) % n;

should it be OK?


 Some naive tests (see here: http://pastebin.com/hnnRLPA0) seem
 to confort my idea, but I may have missed something (again).

Good. I have changed it, and uploaded a new version. Can you perhaps
test it a little more?

Thanks for the input!
♪♫ Alex
-- 
UNSUBSCRIBE: mailto:picolisp@software-lab.de?subject=Unsubscribe


Re: Pseudo Random Number Generation across PicoLisp implementations

2015-04-27 Thread Christophe Gragnic
On Mon, Apr 27, 2015 at 8:01 AM, Alexander Burger a...@software-lab.de wrote:
 […]
 A shift of 32 bits, however, would not suffice, because
 it would put the sign bit of the original 64-bit number into the MSB of
 the result.

OK, I understand.

I have some more ideas to try to make the behaviors of
(rand m n) in pil32 and ersatz the same, but maybe you are
not so much interested. Let me try one.
Would it be possible to just cast (Seed  32) to a long:

if ((n = evInt(ex.Cdr) + 1 - ((Number)x).Cnt) != 0)
n = (long)(Seed  32) % n;
return new Number(((Number)x).Cnt + n);

Some naive tests (see here: http://pastebin.com/hnnRLPA0) seem
to confort my idea, but I may have missed something (again).

If still wrong, I may implement a smaller generator specially for my
language. Like one with just Xn+1 = a Xn, and
a(m−1)2^53 (for ability to use withe JS numbers).
Inspired by Wikipedia and this paper:
http://www.ams.org/journals/mcom/1999-68-225/S0025-5718-99-00996-5/S0025-5718-99-00996-5.pdf
BTW, I wonder why MINSTD is not listed in Table 3.


chri

-- 

http://profgra.org/lycee/ (site pro)
http://delicious.com/profgraorg (liens, favoris)
https://twitter.com/profgraorg
http://microalg.info
--
UNSUBSCRIBE: mailto:picolisp@software-lab.de?subject=Unsubscribe


Re: Pseudo Random Number Generation across PicoLisp implementations

2015-04-27 Thread Alexander Burger
Hi Christophe,

 Alex, if nobody in the mailing list is against modifying the rand of
 Ersatz, may I humbly ask you to use 32 instead of 33 ?

Sorry, I think this is not a good idea. Now I remembered why the shift
is 33 and not 32:

It has to do with the fact that Java doesn't support unsigned 32-it
integers.

But in 'rand', the line

   n = (int)(Seed  33) % n;

requires the result of the shift operation to be a positive number.

Fortunately, Java has with '' a shift operation which inserts zero
bits from the right (i.e. an unsigned shift), thus guaranteeing a
positive result. A shift of 32 bits, however, would not suffice, because
it would put the sign bit of the original 64-bit number into the MSB of
the result.

♪♫ Alex
-- 
UNSUBSCRIBE: mailto:picolisp@software-lab.de?subject=Unsubscribe


Re: Pseudo Random Number Generation across PicoLisp implementations

2015-04-26 Thread Alexander Burger
H Christophe,

 final static long initSeed(Any x) {
 
 Means it gives back a signed 64 bits integer for anything we threw at
 it (called x).

Yes.


 long n; for (n = 0; x instanceof Cell; x = x.Cdr) n += initSeed(x.Car);
 
 Is a recursive trick that walks through x (if walkable, ie Cell) and cumulates
 in n the results of initSeed for anything non walkable that was thrown at it.
 
 
 if (x != Nil) {…
 
 When the non walkable is not a list termination…

Right.


 if (x instanceof Number  ((Number)x).Big == null) n += ((Number)x).Cnt;

Yes. For efficiency, a Number may internally hold a short number (32
bits) in 'Cnt', or a big number in 'Big'.


 If x is a 32 bits integer, just add it to n (not sure about the (Number) cast,
 Number doesn't seem to have children).

The cast is only needed to keep the compiler happy. As 'x' is of type
'Any', it would complain that there is no 'Cnt' property.


In any case, as you observed, we accumulated all that stuff into a long
number 'n'.


 return n=0? n*2 : -n*2+1;
 
 This ensures the return value is non negative (?).
 But why would not the result overflow?

Yes, it may overflow. The above operation is to keep the sign bit in bit
zero of the result. The sign bit is more important than the highest bit
of the number.


 Now if I'm correct until here, initSeed(1) should just return 2.

Correct.


 Then I don't understand why (seed 1) returns -5718471626015965606,
 since (seed x) just multiplies the result of initSeed(x) by 
 6364136223846793005,
 and 2*63641362238467930052^64 (no overflow).

It does overflow. Though 2*6364136223846793005 = 12728272447693586010
fits into 64 bits, it has its most significant bit set. As Java uses
2-complement representation, this is a negative number.

A negative number in 2-complement representation has all bits inverted
plus 1

   :  (inc (x| (dec (** 2 64)) (* 2 6364136223846793005)))
   - 5718471626015965606

♪♫ Alex
-- 
UNSUBSCRIBE: mailto:picolisp@software-lab.de?subject=Unsubscribe


Re: Pseudo Random Number Generation across PicoLisp implementations

2015-04-26 Thread Christophe Gragnic
On Sun, Apr 26, 2015 at 5:48 PM, Alexander Burger a...@software-lab.de wrote:
 Hi Christophe,

 […]

 return n=0? n*2 : -n*2+1;

 This ensures the return value is non negative (?).
 But why would not the result overflow?

 Yes, it may overflow. The above operation is to keep the sign bit in bit
 zero of the result. The sign bit is more important than the highest bit
 of the number.

Interesting. But I don't understand why the sign bit is so important
for the seed. Maybe it's related to what I discovered when working
on (rand), see below.

 […explanation about 2-complement…]

Ah ah ! In fact I realised this when working on the rand function.
It didn't occurred to me that it was also the solution for initSeed.
After sending the email about initSeed, I returned to rand.
I did take more time to test EmuLisp's new rand function,
in particular when sending it two args.

In fact rand worked better than I first thought since it gave
the same results than the Ersatz one. It is not proved, but
when I say the results are the same, in fact the results
matched for ten consecutive calls of (rand 0 1) which
is a quite high probability for «eternal» matching.
The quality of the tests is the same when testing (rand)
and (rand T), with respect to the space of the possible values.

Here I'll discuss
1: (rand 0 1)
2: (rand)
3: (rand T)

1)a) First, for (rand 0 1), my implementation in EmuLisp
was correct when trying to mimic Ersatz.

1)b) Moreover, when removing this line:
shiftedSeed = shiftedSeed.div(2);
which implemented the 33 instead of 32,
the results are the ones of PicoLisp32.

2)a) Then for (rand), my implementation accidentally borrowed
a line of pil32: getting the higher bits instead of the lower ones.
Now that the mistake is fixed it correctly mimics Ersatz.
I had to use the x2^31 ? x : x-2^32
trick (discussed by Alex) to mimic the (int) casting.

2)b) I got stuck when trying to mimic the (rand) of pil32.
In doRand, at this line in big.c:
https://code.google.com/p/picolisp/source/browse/src/big.c#1213
return box(hi(Seed));
In Ersatz in rand (x n)
at this line: 
https://code.google.com/p/picolisp/source/browse/ersatz/fun.src#3325
return new Number((int)Seed);
This seemed so similar that I banged my head on the screen showing
the discrepancies between the results of pil32 and Ersatz that I couldn't
make disappear.
Then, I noticed a pattern, that leaded to reversing the sign bit trick you
just explained in initSeed: n=0? n*2 : -n*2+1;
One the one hand, Ersatz takes the lower 32 bits of the seed,
then returns them as a signed integer. But on the other hand,
pil32 takes the higher 32 bits as n (unsigned integer) and returns
n%2==0 ? n/2 : -(n-1)/2.
Correct ?
A question: why did you choose different implementations ?
Is it for a better efficiency in Ersatz ?

3) (rand T) behaves the same across pil32, Ersatz and EmuLisp.

Alex, if nobody in the mailing list is against modifying the rand of
Ersatz, may I humbly ask you to use 32 instead of 33 ?
This change in the official Ersatz will allow me to have reproducibility
across the implementations of PicoLisp I use without asking my users
to use an unofficial version of Ersatz.

Also, it would be nice, but much less important to me since I don't use
calls to rand with no args, to have the same behavior for (rand) with
pil32 and Ersatz.


chri

-- 

http://profgra.org/lycee/ (site pro)
http://delicious.com/profgraorg (liens, favoris)
https://twitter.com/profgraorg
http://microalg.info
--
UNSUBSCRIBE: mailto:picolisp@software-lab.de?subject=Unsubscribe


Re: Pseudo Random Number Generation across PicoLisp implementations

2015-04-25 Thread Christophe Gragnic
Hi,
Still struggling. I'd like some help to find my mistake trying to understand
the Ersatz version of initSeed. The source is here, if you don't have it handy:
https://code.google.com/p/picolisp/source/browse/ersatz/sys.src#234

final static long initSeed(Any x) {

Means it gives back a signed 64 bits integer for anything we threw at
it (called x).


long n; for (n = 0; x instanceof Cell; x = x.Cdr) n += initSeed(x.Car);

Is a recursive trick that walks through x (if walkable, ie Cell) and cumulates
in n the results of initSeed for anything non walkable that was thrown at it.


if (x != Nil) {…

When the non walkable is not a list termination…


if (x instanceof Number  ((Number)x).Big == null) n += ((Number)x).Cnt;

If x is a 32 bits integer, just add it to n (not sure about the (Number) cast,
Number doesn't seem to have children).


else {
byte b[] = x instanceof Symbol? x.name().getBytes() :
((Number)x).Big.toByteArray();
for (int i = 0; i  b.length; ++i)
n += b[i];
}

If x is not a 32 bits int, get a bytes array (via its name if it is a
sym, or via its
array of bytes if it is a Java BigInteger, and cumulates the bytes in n.


return n=0? n*2 : -n*2+1;

This ensures the return value is non negative (?).
But why would not the result overflow?


Now if I'm correct until here, initSeed(1) should just return 2.
Then I don't understand why (seed 1) returns -5718471626015965606,
since (seed x) just multiplies the result of initSeed(x) by 6364136223846793005,
and 2*63641362238467930052^64 (no overflow).

See here for the source:
https://code.google.com/p/picolisp/source/browse/ersatz/fun.src#3308
only two lines:
n = initSeed(ex.Cdr.Car.eval()) * 6364136223846793005L;
return new Number(Seed = n);

Thanks reading so far !


chri

-- 

http://profgra.org/lycee/ (site pro)
http://delicious.com/profgraorg (liens, favoris)
https://twitter.com/profgraorg
http://microalg.info
--
UNSUBSCRIBE: mailto:picolisp@software-lab.de?subject=Unsubscribe


Re: Pseudo Random Number Generation across PicoLisp implementations

2015-04-24 Thread Christophe Gragnic
On Thu, Apr 23, 2015 at 11:12 AM, Jon Kleiser jon.klei...@fsat.no wrote:

 Feel free to give it a try. As EmuLisp uses JavaScript numbers,
 6364136223846793005 is too many digits, but you may find a way around it.

Progress is being made.

I use this JS lib for arbitrary big integers:
https://github.com/defunctzombie/int

The changes in the EmuLisp code (inspired by Ersatz, as always) are here:
https://github.com/Microalg/microalg/commit/09e1a435c670c16047b67a3349a0599d1ef06bb4

It quite works but I must refine the shiftings (or maybe something else)
to try to get exactly the same values than PicoLisp or Ersatz.


chri

-- 

http://profgra.org/lycee/ (site pro)
http://delicious.com/profgraorg (liens, favoris)
https://twitter.com/profgraorg
http://microalg.info
-- 
UNSUBSCRIBE: mailto:picolisp@software-lab.de?subject=Unsubscribe


Re: Pseudo Random Number Generation across PicoLisp implementations

2015-04-23 Thread Christophe Gragnic
Hi,
After some more time to investigate, I found out that the LCG used in
PicoLisp was not the Haynes one, but the one referred as «Newlib» here:
http://en.wikipedia.org/wiki/Linear_congruential_generator#Parameters_in_common_use
More auto-answers and other questions below.

On Wed, Apr 22, 2015 at 10:50 PM, Christophe Gragnic
christophegrag...@gmail.com wrote:
 Now the question. I found discrepancies between PicoLisp:
 : (rand 1 1000)
 - 1
 : (rand 1 1000)
 - 934
 : (bye)

 And Ersatz
 : (rand 1 1000)
 - 1
 : (rand 1 1000)
 - 967
 : (bye)

 I tried to inspect the source for getting some explanation with no luck.
 PicoLisp
 https://code.google.com/p/picolisp/source/browse/src/big.c#1213
 Ersatz:
 https://code.google.com/p/picolisp/source/browse/ersatz/fun.src#3325
 Is there a reason for those different behaviours ?

Maybe I was too tired. The difference is 33, and we see in the Ersatz src:
n = (int)(Seed  33) % n;
But it is a coincidence since a bit shift should not have this effect
(confirmed when trying with the args 0 and 1000):
PicoLisp:
: (rand 0 1000)
- 0
: (rand 0 1000)
- 648
: (rand 0 1000)
- 760
Ersatz:
: (rand 0 1000)
- 0
: (rand 0 1000)
- 824
: (rand 0 1000)
- 880

Then I had to find the shift in the C code, which is hidden in the `hi` macro
(for high bits):
https://code.google.com/p/picolisp/source/browse/src/pico.h#161
#define hi(w)   num((w)BITS)
Where:
#define WORD ((int)sizeof(long))
#define BITS (8*WORD)

What I don't understand is why Ersatz shifts by 33, and PicoLisp by BITS,
(which should IMHO be 32). This connects to my next remarks/questions:
* The Wikipedia page mentions taking «bits 63...32». I guess that
  they start at bit 0, thus should the shift be 32 ?
* Is 33 used because of the sign bit ?
* In the PicoLisp src, why use LL and not ULL to define
6364136223846793005?


chri

-- 

http://profgra.org/lycee/ (site pro)
http://delicious.com/profgraorg (liens, favoris)
https://twitter.com/profgraorg
http://microalg.info
--
UNSUBSCRIBE: mailto:picolisp@software-lab.de?subject=Unsubscribe


Re: Pseudo Random Number Generation across PicoLisp implementations

2015-04-23 Thread Alexander Burger
Hi Christophe,

our mails overlapped :)


 PicoLisp was not the Haynes one, but the one referred as «Newlib» here:
 http://en.wikipedia.org/wiki/Linear_congruential_generator#Parameters_in_common_use

Right. I took it originally from the book by Donald Knuth.


 Maybe I was too tired. The difference is 33, and we see in the Ersatz src:
 n = (int)(Seed  33) % n;

Right.


 What I don't understand is why Ersatz shifts by 33, and PicoLisp by BITS,
 (which should IMHO be 32).

True. In fact, as I wrote, I don't remember that either. Perhaps a sign
issue? Would it make sense to change it? Should be well tested though.



This connects to my next remarks/questions:

 * The Wikipedia page mentions taking «bits 63...32». I guess that
   they start at bit 0, thus should the shift be 32 ?

Knuth writes that the highest bits are more random than the lower bits
in a linear congruential generator.


 * Is 33 used because of the sign bit ?

Yes, perhaps.


 * In the PicoLisp src, why use LL and not ULL to define
 6364136223846793005?

   : (hex 6364136223846793005)
   - 5851F42D4C957F2D

As you see, this number doesn't have its sign bit set in 64-bit
two-complement representation (it would be 8 or higher in its most
siginificant digit). So unsigned makes no difference.

♪♫ Alex
-- 
UNSUBSCRIBE: mailto:picolisp@software-lab.de?subject=Unsubscribe


Re: Pseudo Random Number Generation across PicoLisp implementations

2015-04-23 Thread Alexander Burger
Hi Christophe,

 I use PicoLisp to implement a pedagogical language, directly
 embedded in PicoLisp. In an exercise about dices, I noticed a
 kind of similarity among the results of my students.
 ...

Nice! This sounds really interesting.


 Now the question. I found discrepancies between PicoLisp:
 ...
 - 934

 And Ersatz
 ...
 - 967

Right. In fact, three of the four implementations (pil64, pil32, mini
and ersatz) will give different results. Only mini and pil32 behave
identically.


 Is there a reason for those different behaviours ?
 
 Last thing (mostly for Jon), I read that the PicoLisp generator is
 a 64-bit linear congruential random generator (both pil64 and pil32)
 with the multiplier 6364136223846793005.

That's correct. All four basically do

   Seed = Seed * 6364136223846793005 + 1

The differences arise from how the Seed is reduced to a 32-bit number
after that.

pil64 uses the middle 64 bits of the 128-bit result of the
multiplication, while pil32 and mini use the highest 32 bits of the
64-bit result. Ersatz does similar to pil32 (i.e. a 64-bit result), but
shifts by 33 (I don't know if there was a special reason).

♪♫ Alex
-- 
UNSUBSCRIBE: mailto:picolisp@software-lab.de?subject=Unsubscribe


Re: Pseudo Random Number Generation across PicoLisp implementations

2015-04-23 Thread Jon Kleiser
Hi Chri,

On 22. Apr, 2015, at 22:50, Christophe Gragnic christophegrag...@gmail.com 
wrote:

 . .
 
 Now the question. I found discrepancies between PicoLisp:
 : (rand 1 1000)
 - 1
 : (rand 1 1000)
 - 934
 : (bye)
 
 And Ersatz
 : (rand 1 1000)
 - 1
 : (rand 1 1000)
 - 967
 : (bye)
 
 I tried to inspect the source for getting some explanation with no luck.
 PicoLisp
 https://code.google.com/p/picolisp/source/browse/src/big.c#1213
 Ersatz:
 https://code.google.com/p/picolisp/source/browse/ersatz/fun.src#3325
 Is there a reason for those different behaviours ?
 
 Last thing (mostly for Jon), I read that the PicoLisp generator is
 a 64-bit linear congruential random generator (both pil64 and pil32)
 with the multiplier 6364136223846793005.
 I'd like to implement this generator for EmuLisp. Any advice ?
 Any objection ?
 
 
 chri

Feel free to give it a try. As EmuLisp uses JavaScript numbers, 
6364136223846793005 is too many digits, but you may find a way around it.

/Jon--
UNSUBSCRIBE: mailto:picolisp@software-lab.de?subject=Unsubscribe