Anyone interested in helping with a new random number generator?

[Nitorami]

3. ... please explain what I'm doing wrong.

Chiming in here... When talking cryptography, we think of secure algorithms for widespread use which should resist determined attack even from state actors. The outcome of wars may depend on them. It is essential that such algorithms are PROVABLY secure. Therefore they are constructed of cryptographic "primitives" whose characteristics are accessible to mathematical analysis. Their strength can be calculated and expressed numerically. Brute forcing it is always equivalent to solving one of the known difficult mathematical problems like factoring large numbers. And then, the algorithm needs to go through processes: expert reviews, independent assessment, crypto experts on the entire world trying to take it apart and find weaknesses, often over decades.

None of that applies to amateur attempts which are mostly based on intuition, complicated number shuffling, obfuscation and empirical testing. That can be fun, and may certainly be good enough for personal needs, which is all you wanted. So considering your goal, you probably did not do anything "wrong"; but the term "cryptography" is just related to much higher ambitions.

[Thaddy]

Indeed. Try to understand the chacha20 code I gave you: that is proven secure and because it is relatively simple about the best place to start. In your case you need just key and nonce, because you want it reversible. My Pascal implementation is a direct translation of the standard. The three of us do understand what you want, but you keep using simple lcg's and you can't do serious encryption with those. It may be enough for your purpose, though.

[ad1mt]

I would like to ask for help from anyone who is running Linux or BSD.
Please could you compile the following pascal program.
Then place it along with the following bash shell script in their own folder.
Then run the bash script, which then runs the program many times, to produce an output file "Test_microsecs.log".
Then send the "Test_microsecs.log" file back to me.
I would do this for myself, but I only have access to Linux on a Raspberry Pi, and that hardware, I've discovered, is a bit flakey.
Many thanks.
The Pascal program:

Code: Pascal  [Select][+][-]

1. {$MODE OBJFPC}
2. program Test_microsecs;
3. uses    sysutils;
4. const   MICROSECS_MOD = 1000;
5. type    Int64u = QWord;
6. var
7. EMULATED_MICROSECS      :int32 = -1;
8. MILLISECS                       :int32;
9. DATETIME                        :TTimeStamp;
10. (*------------------------------*)
11. procedure EMULATE_MICROSECS_FROM_TIMER;
12. var
13. COUNTER         :uint32;
14. DATETIMENEXT:TTimeStamp;
15. begin
16. DATETIME:= DateTimeToTimeStamp(Now);
17. COUNTER:= 0;
18. repeat
19.         Inc(COUNTER);
20.         DATETIMENEXT:= DateTimeToTimeStamp(Now);
21. until   ((DATETIMENEXT.Time MOD 12) = 0)
22. and             (DATETIMENEXT.Time <> DATETIME.Time)
23. and             (COUNTER > 1)
24. ;
25. EMULATED_MICROSECS:= (COUNTER MOD MICROSECS_MOD);
26. end;
27. (*------------------------------*)
28. begin
29. EMULATE_MICROSECS_FROM_TIMER;
30. writeln(EMULATED_MICROSECS);
31. end.
32.  

The Bash shell script:

Code: Bash  [Select][+][-]

1. N=0
2. while [ $N -lt 300 ]
3. do
4.         ./Test_microsecs >Test_microsecs_$N.log &
5.         N=$(expr $N + 1)
6. done
7. sleep 10
8. cat Test_microsecs_[0-9].log Test_microsecs_[0-9][0-9].log Test_microsecs_[0-9][0-9][0-9].log > Test_microsecs.log
9. rm -f Test_microsecs_[0-9].log Test_microsecs_[0-9][0-9].log Test_microsecs_[0-9][0-9][0-9].log
10.  

In case anyone is interested, the Pascal code is attempting to emulate a micoseconds value from the OS timer, where this is not available natively. It does this by counting the number of times round a small loop until the next change of milliseconds.
Suggestions welcome from anyone who has better ideas.

[Thaddy]

In case anyone is interested, the Pascal code is attempting to emulate a micoseconds value from the OS timer, where this is not available natively. It does this by counting the number of times round a small loop until the next change of milliseconds.

That is not true. The timer resolution is uSecs or better.

Code: Pascal  [Select][+][-]

1. program HighResTimerExample;// save in lowercase
2.  
3. {$mode objfpc}{$H+}{$inline on}
4.  
5. uses
6.   SysUtils, BaseUnix, UnixType, Linux;
7.  
8. function GetHighResTime: Int64;
9. var
10.   tp: TTimespec;
11. begin
12.   // Use CLOCK_MONOTONIC for consistent timing
13.   if clock_gettime(CLOCK_MONOTONIC, @tp) = 0 then
14.     Result := Int64(tp.tv_sec) * 1000000000 + Int64(tp.tv_nsec)
15.   else
16.     Result := 0;
17. end;
18.  
19. procedure SleepNanoseconds(ns: Int64);
20. var
21.   req, rem: TTimespec;
22. begin
23.   req.tv_sec := ns div 1000000000;
24.   req.tv_nsec := ns mod 1000000000;
25.   rem.tv_sec := 0;
26.   rem.tv_nsec := 0;
27.  
28.   fpNanoSleep(@req, @rem);
29. end;
30.  
31. var
32.   StartTime, EndTime, Elapsed: Int64;
33.   i: Integer;
34.   Sum: Double;
35. begin
36.   WriteLn('High Resolution Timer Example');
37.   WriteLn('=============================');
38.  
39.   // Measure the time for a simple operation
40.   StartTime := GetHighResTime;
41.  
42.   // Do some work
43.   Sum := 0;
44.   for i := 1 to 1000000 do
45.     Sum := Sum + Sin(i * 0.001);
46.  
47.   EndTime := GetHighResTime;
48.   Elapsed := EndTime - StartTime;
49.  
50.   WriteLn('Operation took: ', Elapsed, ' nanoseconds');
51.   WriteLn('              = ', Elapsed / 1000000:0:3, ' milliseconds');
52.   WriteLn('              = ', Elapsed / 1000000000:0:6, ' seconds');
53.   WriteLn('Dummy result: ', Sum:0:4);
54.  
55.   WriteLn;
56.   WriteLn('Testing high resolution sleep:');
57.  
58.   // Test high resolution sleep
59.   StartTime := GetHighResTime;
60.   SleepNanoseconds(100000000); // 100ms
61.   EndTime := GetHighResTime;
62.   Elapsed := EndTime - StartTime;
63.  
64.   WriteLn('Requested sleep: 100 ms (100,000,000 ns)');
65.   WriteLn('Actual sleep:    ', Elapsed, ' ns');
66.   WriteLn('Difference:      ', Abs(Elapsed - 100000000), ' ns');
67. end.

run as

Code: Bash  [Select][+][-]

1. ./highrestimerexample

[ad1mt]

It is essential that such algorithms are PROVABLY secure. Therefore they are constructed of cryptographic "primitives" whose characteristics are accessible to mathematical analysis. Their strength can be calculated and expressed numerically.

I agree with everything you've said. I just get annoyed when I'm told that I am not allowed to play around with programming ideas as a hobby. Especially when I've invited constructive criticism. Its a bit like signing up for a woodworking class, only to be told on the first day "go away, you are not clever enough to do woodworking".
That said, however....
There are some things in mathematics that cannot be proved/solved.
For some of those things we waiting for a proof (sometimes for 100's years), but for others, it is proven that they cannot be proven/solved.
And some of the ideas in cryptography rely on the fact that some things cannot yet be solved/proved, although they might be one day.
I think I've got all that correct.

[Thaddy]

The terms you are looking for are theorum and conjecture.
Btw look at my timer code above. I am surprised you did not find how to use the high resolution timer on Linux.
https://www.freepascal.org/docs-html/rtl/linux/clock_gettime.html

[ad1mt]

That is not true. The timer resolution is uSecs or better.

Thanks.
I didn't know that... I looked in the documentation but couldn't see a way of doing it.
Is it possible that it's only available on Linux?  I'm mainly Windows only, except when playing with my Raspberry Pi.

[Thaddy]

The windows equivalents are QueryPerformanceCounter and QueryPerformanceFrequency.
The resolution is usually the same. There are many examples for those.

(Btw, the code comes partly from a timer I needed on my raspberry pies.)

Btw, it is not about being clever, but about having followed CS classes related to security.
(Or calling e.g. MathMan if you are stuck...  )

[Nitorami]

Hi ad1mt
here is a proposal:

I would use MWC instead of LCG generators. The MWC is just as simple, but instead of a fixed increment uses the high bits of the previous call. The simple lag-1 version uses 2x32 bits as state and 1x32 bit as multiplier. It has similar speed as the LCG, needs no expensive mod operation and is statistically much better. No lattice structure, no weak bits, period around 2^63. It does have its flaws and I think does not pass Big Crush, but I am fairly confident that two xored generators will.

Prof. Marsaglia has provided the theoretical fundament, the rules to select multipliers and to calculate the period. A bigint library is necessary to calculate suitable multipliers, and it will take a few second to calculate a couple of them so they should be precalculated. With these multipliers, the period is known and linearly related to the multiplier. So I would store say 256 multipliers in a table; with four MWC generators in parallel, the possibilities are then broadly 256 x 255 x 254 x 253 x 2^63^4 > 1e80, way more than necessary. The period of a single generator is guaranteed, you'd need no repeat detection, and with 2^63 the period is also long enough to be sure it will never be run through on any normal PC.

[Thaddy]

Good proposal, but I still think chacha20 is better and probably faster even if the 4 mwc's run in parallel.
Because it passes all tests...and is part of TLS.(which in itself means it is good enough for a real-time handshake)
It is a good thing you addressed the modulo operations, though: that is a no, no.
I will put together a demo along your lines and compare to chacha20.

Alternatively one could run chacha8 over, say,  four parallel operations: that is proven to be as good as chacha20 on its own and would be faster.

One other approach, which I hinted on before, is using algorithms that promote diffusion together with algorithms that promote entropy.
I think that may be the most interesting approach when designing something new.
This is currently a recent development in some area's. Will post a link to some discussions.
The usual names show up...

[ad1mt]

I would use MWC instead of LCG generators.

I looked for a Pascal implemenation I could test, and the first 3 hits were, wikipedia, a haskell implementation, and this...
   https://forum.lazarus.freepascal.org/index.php?topic=35918.0
I copied the code, but was unsure if all the problems were fixed.
If you have an updated version of the code, please could you post it? Thanks.
I'm asking if the code is correct because the wikipedia C implementation looks very different, with a 4096 entry state table.

[ad1mt]

@Thaddy
I later found this...
   https://wiki.freepascal.org/Marsaglia%27s_pseudo_random_number_generators
I think you might have mentioned this earlier.
I've copied the MWC code, and am going to do some testing.
I've got a question about your conversion from 32bit int to single/float... to do this you multiply by 4.656613e-10
I can't figure out why. Could you explain please?  Thanks.

[Nitorami]

I tried to send you an example code plus a test routine for safe primes, but I can currently only logon via mobile. Very odd. From 3 available PCs, I get the message that cookies must be enabled, even when just trying to read a post. But cookies ARE definitely enabled. Anybody knows what is going on?

Edit: same from every browser.

[srvaldez]

I asked Grok to translate the C code from Wikipedia and this is the result https://grok.com/share/c2hhcmQtMg%3D%3D_386e63c4-3c5a-49d9-ae77-726f1b61702e
however if I set seed=123 the results differ

Code: [Select]

C          Random CMWC: 2300145674
freePascal Random CMWC: 2674568288
<edit> the result is likely to be different unless the C rand() and the FP Random functions produce the same result
<edit2> indeed, if I change the Rand32 function to return a certain constant then the results from C and FP are the same, so the translation looks good

[Nitorami]

Still can only access the forum from my mobile, so let's try it via mobile

Thaddy: Agreed that chacha20 would be the professional solution. But I think the OP has fun in designing his own algorithm. I sympathise with that, it is less productive than using a library but forces one to think and better understand the problems and pitfalls.

@ad1mt: Not the CMWC4096 "mother of all", just a basic lag1 MWC. Below a simple example, fixed seed, and for real use I would pack the generator into an advanced record.
I also attached  a function to test or suitable multipliers that satisfy Prof. Marsaglias condition. It needs the GMP library which comes packed with FPC. I think it is correct.

Code: Pascal  [Select][+][-]

1. uses gmp;  //library is required for function MPZIsSafePrime() only
2. {$R-,Q-}
3.  
4. function MWC32r1: dword;
5. {Basic MWC lag-1 generator
6. High quality random numbers based on Multiply With Carry, by Prof. Marsaglia
7. x(n)=a*x(n-1)+carry mod 2^32.  Period is a*2^31-1. Multiplier (any one will do):
8. 1791398085 1683268614 1965537969 1675393560 1967773755 1517746329 1447497129 1655692410 1606218150
9. 2051013963 1075433238 1557985959 1781943330 1893513180 1631296680 2131995753 2083801278 1873196400 1554115554
10. or any a for which both a*2^32-1 and a*2^31-1 are prime
11. }
12. var t: qword ;
13. const A: dword = 23;
14.   B: dword = 12318;  //Mind the conditions; A and B shall not be both zero, and not equal multiplier-1;
15.  //not sure -> consult wikipedia
16. begin
17.   t :=  qword (2131995753)*A+B;
18.   A := lo(t);
19.   B := hi(t);
20.   result := A;
21. end;
22.  
23.  
24. function MPZIsSafePrime (const a, b,r: dword): boolean;
25. {decides whether for a given multiplier a the terms p=a*b^r-1 and (p-1)/2 = (a*b^r/2-1) are both prime.
26. In other words, if p = a*(2^32)-1 is a "safe prime".
27.  
28. b is in bits e.g. 32, r is the lag}
29.  
30. var a1,b1,p,one,two: mpz_t;
31. begin
32.   mpz_init (a1);
33.   mpz_init (b1);
34.   mpz_init (p);
35.   mpz_set_ui(a1,a);
36.   mpz_ui_pow_ui(b1,2,b*r); //b1 = 2^(b*r)
37.   mpz_mul (p, a1,b1);
38.   mpz_sub_ui (p, p, 1);
39.   result := mpz_probab_prime_p (p,15) > 0;
40.   //nun testen wir noch (p-1)/2 = a*b^r/2-1.
41.   mpz_ui_pow_ui(b1,2,b*r-1);
42.   mpz_mul (p, a1,b1); // p = a1*b1^r
43.   mpz_sub_ui (p, p, 1);
44.   result := result and (mpz_probab_prime_p (p,15) >0);
45.   mpz_clear (p);
46.   mpz_clear (b1);
47.   mpz_clear (a1);
48. end;
49.  
50.  
51. var a,i: dword;
52.  
53. begin
54.   writeln ('output 10 random values from MWC');
55.   for i := 1 to 10 do writeln (MWC32r1);
56.   writeln;
57.   writeln ('calculate 10 suitable multipliers for 32 bit lag 1 MWC');
58.   i := 0;
59.   repeat
60.     a := high (dword) div 4 * 3 + random (high (dword) div 4); //try some reasonably large number
61.     if MPZIsSafePrime (a, 32, 1) then begin
62.       writeln (a);
63.       inc (i);
64.     end;
65.   until i>10;
66. end.
67.  

Edit: Thaddys version is from 16bit days, essentially it combines two 16-bit MCWs with different multiplier to one 32-bit generator. It is accordingly slow in comparison to a native 32bit Version.

But it is interesting to see that Prof Marsaglia combined two MCWs. That means as long as the multipliers are different, the streams are uncorrelated, otherwise the combined generator would be bad. In turn we may conclude that paralleling several 32bit MWCs is also ok.
On a side note - Marsaglias choice of multiplier (18000) was bad. While the 18000 satisfies his conditions, it causes correlations between subsequent words which can even be visualised in a randogram, which means it is really quite bad. I tested that earlier, and it is true for all even multipliers. So when calculating multipliers, avoid the even ones !