Blum integer: Difference between revisions

Content added Content deleted
(Replace the "blum" sequence with an array.)
m (→‎{{header|Free Pascal}}: changed output format)
Line 420: Line 420:
}
}
const
const
LIMIT = 10*1000*1000;//8802377;//
LIMIT = 10*1000*1000;// >750
type
⚫
tP4n3 = array of Uint32;


function CommaUint(n : Uint64):AnsiString;
function CommaUint(n : Uint64):AnsiString;
Line 454: Line 456:
end;
end;


procedure Sieve4n_3_Primes(Limit:Uint32;var P4n3:tP4n3);
var
var
BlumField : array [0..LIMIT] of boolean;
sieve : array of byte;
BlumPrimes : array[0..332398] of Uint32;
BlPrCnt,idx,n,j: Uint32;
⚫
sieve : array[0..((LIMIT-3) DIV 4)] of boolean;

⚫
EndDigit : array[0..9] of Uint32;
⚫
n,idx,j,k,BlPrCnt : Uint64;
begin
begin
//DIV 3 -> smallest factor of Blum Integer
//sieve primes 4*i+3
n := (LIMIT DIV 3 -3) DIV 4+ 1;
setlength(sieve,n);
setlength(P4n3,n);

BlPrCnt:= 0;
BlPrCnt:= 0;
IDX := 0;
idx := 0;
repeat
repeat
if NOT(sieve[idx]) then
if sieve[idx]= 0 then
begin
begin
n := idx*4+3;
n := idx*4+3;
BlumPrimes[BlPrCnt] := n;
P4n3[BlPrCnt] := n;
inc(BlPrCnt);
inc(BlPrCnt);
j := idx+n;
j := idx+n;
Line 476: Line 479:
while j <= High(sieve) do
while j <= High(sieve) do
begin
begin
sieve[j] := true;
sieve[j] := 1;
inc(j,n);
inc(j,n);
end;
end;
Line 484: Line 487:
//collect the rest
//collect the rest
for idx := idx+1 to High(sieve) do
for idx := idx+1 to High(sieve) do
if NOT(sieve[idx]) then
if sieve[idx]=0 then
Begin
Begin
BlumPrimes[BlPrCnt] := idx*4+3;
P4n3[BlPrCnt] := idx*4+3;
inc(BlPrCnt);
inc(BlPrCnt);
end;
end;
writeln('There ',BlPrCnt,' primes 4*n+3 to Limit ',LIMIT);
setlength(P4n3,BlPrCnt);
dec(BlPrCnt);
setlength(sieve,0);
end;

var
⚫
BlumField : array[0..LIMIT] of boolean;
BlumPrimes : tP4n3;
EndDigit : array[0..9] of Uint32;
k : Uint64;
⚫
n,idx,j,P4n3Cnt : Uint32;
begin
Sieve4n_3_Primes(Limit,BlumPrimes);
P4n3Cnt := length(BlumPrimes);
writeln('There are ',CommaUint(P4n3Cnt),' needed primes 4*n+3 to Limit ',CommaUint(LIMIT));
dec(P4n3Cnt);
writeln;


//generate Blum-Integers
//generate Blum-Integers
For idx := 0 to BlPrCnt do
For idx := 0 to P4n3Cnt do
Begin
Begin
n := BlumPrimes[idx];
n := BlumPrimes[idx];
For j := idx+1 to BlPrCnt do
For j := idx+1 to P4n3Cnt do
Begin
Begin
k := n*BlumPrimes[j];
k := n*BlumPrimes[j];
Line 509: Line 526:
j := 0 ;
j := 0 ;
repeat
repeat
while Not(BlumField[idx]) do
while (idx<LIMIT) AND Not(BlumField[idx]) do
inc(idx);
inc(idx);
⚫
if idx = LIMIT then
BREAK;
if j mod 10 = 0 then
if j mod 10 = 0 then
writeln;
writeln;
Line 520: Line 539:


//count and calc and summate decimal digit
//count and calc and summate decimal digit
writeln(' relative occurence of digit');
writeln(' n.th |Blum-Integer|Digit: 1 3 7 9');
idx :=0;
idx :=0;
j := 0 ;
j := 0 ;
⚫
n := 0;
k := 26828;
k := 26828;
repeat
repeat
while Not(BlumField[idx]) do
while (idx<LIMIT) AND Not(BlumField[idx]) do
inc(idx);
inc(idx);
if idx = LIMIT then
BREAK;
//count last decimal digit
//count last decimal digit
inc(EndDigit[idx MOD 10]);
inc(EndDigit[idx MOD 10]);
Line 531: Line 555:
if j = k then
if j = k then
begin
begin
writeln(CommaUint(j):8,'.th Blum integer is: ',CommaUint(idx):10);
write(CommaUint(j):10,' | ',CommaUint(idx):11,'| ');
⚫
write(EndDigit[1]/j*100:7:3,'% |');
write(EndDigit[3]/j*100:7:3,'% |');
write(EndDigit[7]/j*100:7:3,'% |');
writeln(EndDigit[9]/j*100:7:3,'%');
if k < 100000 then
if k < 100000 then
k := 100000
k := 100000
Line 540: Line 568:
until j >= 400000;
until j >= 400000;
Writeln;
Writeln;

⚫
idx := 0;
For j := 0 to 9 do
inc(idx,EndDigit[j]);
For j := 0 to 9 do
⚫
if EndDigit[j] <> 0 then
⚫
writeln('Digit ',j:2,' rel. ',EndDigit[j]/idx*100:7:3,'%');
end.</syntaxhighlight>
end.</syntaxhighlight>
{{out|@TIO.RUN}}
{{out|@TIO.RUN}}
<pre>
<pre>
There 332398 primes 4*n+3 to Limit 10000000
There are 119,644 needed primes 4*n+3 to Limit 10,000,000

First 50 Blum-Integers
First 50 Blum-Integers


21 33 57 69 77 93 129 133 141 161
21 33 57 69 77 93 129 133 141 161
Line 559: Line 581:
597 633 649 669 681 713 717 721 737 749
597 633 649 669 681 713 717 721 737 749


relative occurence of digit
26,828.th Blum integer is: 524,273
n.th |Blum-Integer| 1 3 7 9
100,000.th Blum integer is: 2,075,217
26,828 | 524,273| 24.933% | 25.071% | 25.030% | 24.966%
200,000.th Blum integer is: 4,275,533
100,000 | 2,075,217| 24.973% | 25.026% | 25.005% | 24.996%
300,000.th Blum integer is: 6,521,629
200,000 | 4,275,533| 24.990% | 24.986% | 25.033% | 24.992%
400,000.th Blum integer is: 8,802,377
300,000 | 6,521,629| 24.982% | 25.014% | 25.033% | 24.970%

400,000 | 8,802,377| 25.001% | 25.017% | 24.997% | 24.985%
Digit 1 rel. 25.001%
Digit 3 rel. 25.017%
Digit 7 rel. 24.997%
Digit 9 rel. 24.985%


Real time: 0.108 s User time: 0.076 s Sys. time: 0.032 s CPU share: 98.90 %
Real time: 0.097 s User time: 0.068 s Sys. time: 0.028 s CPU share: 99.08 %
C-Version //-O3 -marchive=native
C-Version //-O3 -marchive=native
Real time: 1.658 s User time: 1.612 s Sys. time: 0.033 s CPU share: 99.18 %</pre>
Real time: 1.658 s User time: 1.612 s Sys. time: 0.033 s CPU share: 99.18 %</pre>
Retrieved from "https://rosettacode.org/wiki/Blum_integer"