70 lines
1.6 KiB
Plaintext
70 lines
1.6 KiB
Plaintext
XIncludeFile "isprime.pb"
|
|
|
|
Structure PrimeHelper
|
|
inc.i[8]
|
|
index.i
|
|
EndStructure
|
|
|
|
Procedure.i firstPrimeFactor(n.q)
|
|
If n = 1 : ProcedureReturn 1 : EndIf
|
|
If n % 3 = 0 : ProcedureReturn 3 : EndIf
|
|
If n % 5 = 0 : ProcedureReturn 5 : EndIf
|
|
|
|
Define helper.PrimeHelper
|
|
helper\inc[0] = 4 : helper\inc[1] = 2 : helper\inc[2] = 4
|
|
helper\inc[3] = 2 : helper\inc[4] = 4 : helper\inc[5] = 6
|
|
helper\inc[6] = 2 : helper\inc[7] = 6
|
|
|
|
Define k.q = 7
|
|
While k * k <= n
|
|
If n % k = 0 : ProcedureReturn k : EndIf
|
|
k + helper\inc[helper\index]
|
|
helper\index = (helper\index + 1) % 8
|
|
Wend
|
|
ProcedureReturn n
|
|
EndProcedure
|
|
|
|
OpenConsole()
|
|
Define Dim blum.q(49)
|
|
Define.q bc = 0, i = 1
|
|
Define Dim counts.q(9)
|
|
|
|
Repeat
|
|
Define p.q = firstPrimeFactor(i)
|
|
If p % 4 = 3
|
|
Define q.q = i / p
|
|
If q <> p And q % 4 = 3 And isPrime(q)
|
|
If bc < 50 : blum(bc) = i : EndIf
|
|
counts(i % 10) + 1
|
|
bc + 1
|
|
|
|
If bc = 50
|
|
PrintN("First 50 Blum integers:")
|
|
For j = 0 To 49
|
|
Print(" " + RSet(Str(blum(j)), 3))
|
|
If (j + 1) % 10 = 0 : PrintN("") : EndIf
|
|
Next
|
|
PrintN("")
|
|
ElseIf bc = 26828 Or bc % 100000 = 0
|
|
PrintN("The " + RSet(Str(bc), 6) + "th Blum integer is: " + RSet(Str(i), 7))
|
|
If bc = 400000
|
|
PrintN(#CRLF$ + "% distribution of the first 400,000 Blum integers:")
|
|
For j = 1 To 9 Step 2
|
|
If j <> 5
|
|
PrintN(RSet(StrF(counts(j)/4000, 3), 5) + "% end in " + Str(j))
|
|
EndIf
|
|
Next
|
|
Break
|
|
EndIf
|
|
EndIf
|
|
EndIf
|
|
EndIf
|
|
If i % 5 = 3
|
|
i + 4
|
|
Else
|
|
i + 2
|
|
EndIf
|
|
ForEver
|
|
|
|
PrintN(#CRLF$ + "Press ENTER to exit"): Input()
|