77 lines
2.1 KiB
Plaintext
77 lines
2.1 KiB
Plaintext
#include "isprime.bas"
|
|
|
|
Type PrimeHelper
|
|
inc(7) As Integer
|
|
idx As Integer
|
|
End Type
|
|
|
|
Function initPrimeHelper() As PrimeHelper
|
|
Dim helper As 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
|
|
helper.idx = 0
|
|
Return helper
|
|
End Function
|
|
|
|
Function firstPrimeFactor(n As Integer) As Integer
|
|
If n = 1 Then Return 1
|
|
If n Mod 3 = 0 Then Return 3
|
|
If n Mod 5 = 0 Then Return 5
|
|
|
|
Dim helper As PrimeHelper = initPrimeHelper()
|
|
Dim k As Integer = 7
|
|
|
|
While k * k <= n
|
|
If n Mod k = 0 Then Return k
|
|
k += helper.inc(helper.idx)
|
|
helper.idx = (helper.idx + 1) Mod 8
|
|
Wend
|
|
|
|
Return n
|
|
End Function
|
|
|
|
Sub main()
|
|
Dim As Integer blum(49), counts(9)
|
|
Dim As Integer bc = 0, i = 1, p, q
|
|
Dim As Integer j
|
|
|
|
Dim As Double t0 = Timer
|
|
Do
|
|
p = firstPrimeFactor(i)
|
|
If p Mod 4 = 3 Then
|
|
q = i \ p
|
|
If q <> p Andalso q Mod 4 = 3 Andalso isPrime(q) Then
|
|
If bc < 50 Then blum(bc) = i
|
|
counts(i Mod 10) += 1
|
|
bc += 1
|
|
|
|
If bc = 50 Then
|
|
Print "First 50 Blum integers:"
|
|
For j = 0 To 49
|
|
Print Using "####"; blum(j);
|
|
If (j + 1) Mod 10 = 0 Then Print
|
|
Next
|
|
Print
|
|
Elseif bc = 26828 Orelse bc Mod 100000 = 0 Then
|
|
Print Using "The ###,###th Blum integer is: #,###,###"; bc; i
|
|
|
|
If bc = 400000 Then
|
|
Print !"\n% distribution of the first 400,000 Blum integers:"
|
|
For j = 1 To 9 Step 2
|
|
If j <> 5 Then Print Using " ##.###% end in #"; (counts(j)/4000); j
|
|
Next
|
|
Exit Do
|
|
End If
|
|
End If
|
|
End If
|
|
End If
|
|
i += Iif(i Mod 5 = 3, 4, 2)
|
|
Loop
|
|
Print Chr(10); Timer - t0; " sec."
|
|
End Sub
|
|
|
|
main()
|
|
|
|
Sleep
|