RosettaCodeData/Task/Benfords-law/FreeBASIC/benfords-law.basic
2023-07-01 13:44:08 -04:00

90 lines
2.1 KiB
Text

' version 27-10-2016
' compile with: fbc -s console
#Define max 1000 ' total number of Fibonacci numbers
#Define max_sieve 15485863 ' should give 1,000,000
#Include Once "gmp.bi" ' uses the GMP libary
Dim As ZString Ptr z_str
Dim As ULong n, d
ReDim As ULong digit(1 To 9)
Dim As Double expect, found
Dim As mpz_ptr fib1, fib2
fib1 = Allocate(Len(__mpz_struct)) : Mpz_init_set_ui(fib1, 0)
fib2 = Allocate(Len(__mpz_struct)) : Mpz_init_set_ui(fib2, 1)
digit(1) = 1 ' fib2
For n = 2 To max
Swap fib1, fib2 ' fib1 = 1, fib2 = 0
mpz_add(fib2, fib1, fib2) ' fib1 = 1, fib2 = 1 (fib1 + fib2)
z_str = mpz_get_str(0, 10, fib2)
d = Val(Left(*z_str, 1)) ' strip the 1 digit on the left off
digit(d) = digit(d) +1
Next
mpz_clear(fib1) : DeAllocate(fib1)
mpz_clear(fib2) : DeAllocate(fib2)
Print
Print "First 1000 Fibonacci numbers"
Print "nr: total found expected difference"
For d = 1 To 9
n = digit(d)
found = n / 10
expect = (Log(1 + 1 / d) / Log(10)) * 100
Print Using " ## ##### ###.## % ###.## % ##.### %"; _
d; n ; found; expect; expect - found
Next
ReDim digit(1 To 9)
ReDim As UByte sieve(max_sieve)
'For d = 4 To max_sieve Step 2
' sieve(d) = 1
'Next
Print : Print "start sieve"
For d = 3 To sqr(max_sieve)
If sieve(d) = 0 Then
For n = d * d To max_sieve Step d * 2
sieve(n) = 1
Next
End If
Next
digit(2) = 1 ' 2
Print "start collecting first digits"
For n = 3 To max_sieve Step 2
If sieve(n) = 0 Then
d = Val(Left(Trim(Str(n)), 1))
digit(d) = digit(d) +1
End If
Next
Dim As ulong total
For n = 1 To 9
total = total + digit(n)
Next
Print
Print "First";total; " primes"
Print "nr: total found expected difference"
For d = 1 To 9
n = digit(d)
found = n / total * 100
expect = (Log(1 + 1 / d) / Log(10)) * 100
Print Using " ## ######## ###.## % ###.## % ###.### %"; _
d; n ; found; expect; expect - found
Next
' empty keyboard buffer
While InKey <> "" : Wend
Print : Print "hit any key to end program"
Sleep
End