72 lines
1.7 KiB
Fortran
72 lines
1.7 KiB
Fortran
!
|
|
! Curzon Numbers
|
|
! tested with Intel ifx (IFX) 2025.2.1 20250806 on Kubuntu 25.04
|
|
! GNU Fortran (Ubuntu 14.2.0-19ubuntu2) 14.2.0 on Kubuntu 25.04
|
|
! VSI Fortran x86-64 V8.6-001 on OpenVMS x86_64 V9.2-3
|
|
! No Non-standard features used, should compile on any fairly recent Fortran.
|
|
! U.B., September 2025
|
|
!
|
|
program Curzon
|
|
|
|
implicit none
|
|
|
|
integer, parameter :: kint = 8 ! 8-byte-integers are big enough, 4 bytes is insufficient.
|
|
|
|
integer (kind=kint) :: base, num, count
|
|
|
|
do base=2,10,2
|
|
print ("( 'Curzon numbers with base ', i2,':')"), base
|
|
count = 0
|
|
num = 1
|
|
do while (count .lt. 50)
|
|
if (isCurzon (num, base)) then
|
|
count = count + 1
|
|
write (*, '(i5)', advance='no') num
|
|
if (mod (count, 10) .eq. 0 ) write (*,*) ! EOL
|
|
end if
|
|
num = num + 1
|
|
end do
|
|
do
|
|
if (isCurzon (num,base)) then
|
|
count = count + 1
|
|
if (count .eq. 1000) exit
|
|
endif
|
|
num = num + 1
|
|
end do
|
|
print ("('1000th Curzon number with base ', i0, ': ', i0/)") , base, num
|
|
|
|
end do
|
|
|
|
contains
|
|
|
|
|
|
function isCurzon (n, k) result (YN)
|
|
|
|
integer(kind=kint), intent(in) :: n,k
|
|
integer(kind=kint) :: m, p
|
|
logical :: YN
|
|
|
|
m = k*n+1
|
|
p = ModularExp (k,n,m) + 1
|
|
YN = m .eq. p
|
|
|
|
end function isCurzon
|
|
|
|
function ModularExp(a, b, n) result(d)
|
|
|
|
integer(kind=kint) , intent(in) :: a,b,n
|
|
integer (kind=kint) :: d
|
|
integer (kind=kint) :: i, k
|
|
d = 1
|
|
k = 0
|
|
do while (ishft (b, -k) .gt. 0)
|
|
k = k + 1
|
|
end do
|
|
do i = k - 1, 0, -1
|
|
d = mod (d*d, n)
|
|
if (iand (ishft (b, -i) , 1_kint) .gt. 0) d = mod (d*a ,n);
|
|
end do
|
|
|
|
end function ModularExp
|
|
|
|
end program Curzon
|