Problems & Puzzles: Problems
Problem 38.
Sebastián
Martín Ruiz Prime formulas
Sebastián Martín Ruiz, a high school teacher from Chipiona,
Cadiz, Spain, has came up [1] with his own prime formulas p(n)
and nxt(p)
 p(n) calculates the nth prime number. See [2]:
(here
means the largest integer less than or equal to x; in Ubasic this
corresponds to the function "int(x)" or
"floor(x)")
 nxt(p) calculate the following prime to another
given p. See [3]
Sebastian
does not claim that his formulas satisfy one of the desiderata
asked by Hardy and Wright for these mathematical
targets:
"Is
there a simple general formula for the nth prime pn (a formula, that is
to say, by which we can calculate the value of pn for any given n with less labour than by the use of the sieve of
Eratosthenes)?"
because he
recognizes that both functions should work slower than the corresponding
algorithms based in the sieve of
Erathostenes.
(This
happens due to the very nature of one concept at the core of both
formulas: the exhaustive calculation of the number of divisors of a number
n; this calculation is needed to evaluate the socalled Smarandache
Prime function, a binary one, that constitutes  at the end  the
essential calculation)
Nevertheless
 with this feature against  Sebastián hopes indeed that at
least his formulas are better than the previously produced
by others  like Sierpinski, Willans & Gandhi.
See [4]^{ }  in two specific
aspects:
 His
formulas do not suppose any previous knowledge or calculation of any
prime number.
 His
formulas are "faster" than any other formulas alike in order to make the
same calculations.
As a matter
of fact Sebastián claims that both formulas are of a polynomial
order  k*(n log n)^3  while the previous ones are factorial or
exponential order.
Questions:
1) Do you agree with the Sebastian's two before
claims?
In the
Ribenboim's interpretation of the Hardy & Wright desiderata for
the prime formula he says:
"What
is intended... is to find a closed expression for the nth prime p_{n}
in terms of n, by means of functions that are
computable and , if possible, classical" (p.179, op.
cit.)
2) Dou you think that the Sebastian's formulas
are in terms of n, closed, computable
and
classical?
3)
From a computational point of view, are the Sebastian's prime formulas an
improvement or a back step regarding the Liu
Fenghsui' s prime formula? (See Problem 37 in
these pages)
4) Do you devise any way to optimize the Sebastian's prime
formulas?
_________________ Notes:
[1]:
Sebastián Martín Ruiz original publications can be found here: p(n)
and nxt(p). The formulas
were developed in 1999 but published in the SJN at Spring 2000. At
the time of the publications he was considering to be the only author of
the main inner functions [ d(n) and G(n) ]. Only recently he found that
d(n) was already developed, according to "Introducción al la teoría de
números", Ivan Niven & H. S. Zuckerman. Editorial LimusaWiley, S.A.
México, 1969. pp. 888990.
[2] &
[3]: Here are my Ubasiccodes to calculate these functions, and below the
Sebastian's ones for Mathematica.
p(n) 
nxt(p) 
10 input
N:clr time:'pnsmr.ub 20 S1=0:Z=2*(int(N*log(N))+1) 30 for K=1
to Z 40 S2=0 50 for J=2 to K 60 S3=0 70 for L=1 to
J 80 S3=S3+J\L(J1)\L 90 next L 100
S2=S2+1+(2S3)\J 110 next J 120 S1=S1+1S2\N 130 next
K 140 print "The prime";N;"th is=";S1+1; 150 if
S1+1<>prm(N) then print "Error":end 160 print "Verified y
Right",time:end

10 input
P:clr time:'pqsmr.ub 20 S1=0 30 for M=P+1 to 2*P 40
F=1 50 for I=P+1 to M:S2=0 60 for J=1 to
I:S2=S2+I\J(I1)\J:next J 70 S2=2S2:F=F*((S2\I)) 80 next
I 90 S1=S1+F 100 next M 110 print
"nxt(";P;")=";P+1+S1,time:end

DD[i_]:=Sum[Floor[i/k]Floor[(i1)/k],{k,1,i}]
FG[n_]:=Sum[1+Floor[(DD[j]2)/j],{j,2,n}]
VIC[n_]:=1+Sum[1Floor[FG[k]/n],
{k,1,2*(Floor[n*N[Log[n],15]]+1)}] Do[Print[VIC[n],"
",Prime[n]],{n,1,25}]

p=Input["Input a prime number:"]
DD[i_]:=Sum[Floor[i/j]Floor[(i1)/j],{j,1,i}] G[i_]:=Floor[(2DD[i])/i]
F[m_]:=Product[G[i],{i,p+1,m}]
S[n_]:=Sum[F[m],{m,n+1,2*n}]
Print["nxt(",p,")=",p+1+S[p]]

[4] P.
Ribenboim, The New Book of Prime Number Records, pp.
179186
Solution:
Jud McCranie has
suggested to eliminate the 'floor' or the 'int' functions in
the codes that I originally wrote in Ubasic, using  instead  the
integer division operator (x\y). I have followed his suggestion and
the savings in time of computation has been of .... 33.3%... nothing
bad!
Now the codes shown in the
table above reflect that idea.
***
The Jud's idea also can
be applied to the formulas writing, eliminating the clumsy floor brackets
so hard to write in a basic word processor... And I would like if that
idea can be implemented also in the codes for Mathematica, with the
same good results in time savings that we got in Ubasic
***
I (CR) have devised
another important improvement in the codes based in the fact that the
calculation of the number of divisors of a number n, exhaustive via, needs
run only from 1 to isqrt(n) and then multiply for 2 the result obtained.
The changes needed in my codes are:
Function 
Changes in
red 
p(n) 
70 for L=1 to isqrt(J) 100
S2=S2+1+(22*S3)\J 
nxt(p) 
60 for J=1 to isqrt(I):S2=S2+I\J(I1)\J:next J 70 S2=22*S2:F=F*((S2\I))

The times spent in calculation
drop dramatically!
BTW nxt(p): while nxt(p) is the
following prime to the number p, p itself needs not to be
prime!
***
The Sebastian's
formulas introducing the ideas given by Jud and me (CR), and
homogenizing the subindexes in both formulas, look the following
way:
and
in these formulas all the divisions are
integerdivisions
***
The new programs for
Mathematica, sent by Sebastián, are:
Function 
Program 
p(n) 
DD[i_]:=Sum[Quotient[i,k]Quotient[(i1),k],{k,1,Floor[Sqrt[i]]}]
FG[n_]:=Sum[1+Quotient[(22*DD[j]),j],{j,2,n}]
VIC[n_]:=1+Sum[1Quotient[FG[k],n],{k,1,2*(Floor[n*Log[n]]+1)}]
Do[Print[VIC[n],"
",Prime[n]],{n,1,50}] 
nxt(p) 
p=Input["Input a prime
number:"] DD[i_]:=Sum[Quotient[i,j]Quotient[(i1),j],
{j,1,Floor[Sqrt[i]]}]
G[i_]:=Quotient[(22*DD[i]),i]
F[m_]:=Product[G[i],{i,p+1,m}]
S[n_]:=Sum[F[m],{m,n+1,2*n}]
Print["nxt(",p,")=",p+1+S[p]]

***
Patrick De Geest sent the
following ideas to improve the programs a bit more:
1. Using C like notation. The variables thus need to
be read only once instead of twice. (in large loops the gain in speed
can become important) Therefore I changed:
S3 = S3 + J\L(J1)\L into S3 += J\L(J1)\L S2 =
S2 + 1 + (22*S3)\J into S2 += 1 +(22*S3)\J S1 = S1 + 1  S2\N into
S1 += 1  S2\N
2. Using shortest possible variable
names. (in large loops the gain can become noticable) Therefore I
changed S1 to S, S2 to T and S3 to U.
3. Cramming all the fornext lines into one long
line saves a little time when dealing with large loops. Of course that
clarity in structure is lost but maybe that is not so important in small
routines as p(n).
Combining all these 'shortcuts' gave me the above
mentioned runtime reduction. My final version of p(n) looks like
:
10 input N:clr time:'pnsmr.ub 20
S=0:Z=2*(int(N*log(N))+1) 30 for K=1 to Z:T=0:for J=2 to K:U=0:for
L=1 to isqrt(J):U+=J\L(J1)\L:next L:T+=1+(22*U)\J:next
J:S+=1T\N:next K 40 print "The prime";N;"th is =";S+1; 50 if
S+1<>prm(N) then print "Error":end 60 print "Verified y
Right",time:end
The final version of nxt(p) looks like
: (gain in speed is less pronounced here. I give them anyway as
these optimizations might come handy when developing other future time
critical algorithms)
10 input P:clr time:'pqsmr.ub 20 S=0:for M=P+1
to 2*P:F=1:for I=P+1 to M:T=0:for J=1 to isqrt(I):T=T+I\J(I1)\J:next
J:T=22*T:F*=(T\I):next I:S+=F:next M 30 print
"nxt(";P;")=";P+1+S,time:end
***
Sebastian Ruiz has made a very
interesting discovery to his prime formula p(n). You may calculate p(n)
not starting the summa for k since 1 but since any number m less or equal
to p(n)
Before 
After 
p(n)=1 + summa from k=1
to... 
p(n)=m + summa from k=m
to... 
But this means that of m is the previous
prime to p(n), then this formula is also a formula to calculate
recursively prime by prime starting in any of these. This also means that
this formula used this way is a competitor formula for the nxt(p) formula
of Sebastian.
***
As a matter of fact, according to my
(C.R.) implementation of the Ubasic code of the p(n) formula
modified to calculate the series of primes, p(n) is better (faster) than
nxt(p)...!!!.
Jud McCranie was asked to confirm
the Rivera's result and this is what he aswered:
"Using p(n), it takes 0.4 seconds
to calculate p(75), but it takes 8.4 seconds to calculate p(1) ...p(75).
Using NxtPrm, it takes 0.8 seconds to calculate p(75) from p(74), and
16.4seconds to calculate p(1) then p(2) ... p(75). So for these tests,
NxtPrm takes twice as long."
***
Sebastian announces
that his nxt(p) function can be formulated in order to obtain,
recursively, only the primes of an arithmetic progression a+d.n,
n=>0, gcd(a, d)=1. If you want to know the expression of that formula
and/or the corresponding code in Mathematica, please request them to Sebastian.
***
One more communication came
from S. M. Ruiz (31/10/2002):
We can transform the sums of
products in the formula for NXT(p) in a nested
product:
a+ab+abc=a(1+b(1+c))
For a nested product NP[n]
with n levels, we have NP[n]=O(n). Since d[i]=O(i^(1/2)) with Rivera's
improvement, we get NestedNXT[n]=NP[d[i],{i,1,n}]=O(n^(3/2)). And since
pn=O(n log n), we get NestedNXT[pn]=O((n log n)^(3/2)) which is
better than the explicit formula's #pn=O((n log n)^(5/2)). (By
Jonathan Sondow and Sebastián Martín
Ruiz)
Here are the corresponding
codes in Mathematica and Ubasic, sent by himself:
Mathematica 
Ubasic 
DD[i_] :=
Sum[Quotient[i, j]  Quotient[i  1, j], {j, 1, Sqrt[i]}] G[i_]
:= Quotient[2  2*DD[i], i] (* n=Bound for Pi[p] prime counting
function *) (* A=Rosser and Schoenfeld inequalities
*) (* B=Nested Product O(nlogn)^(3/2) *) p =
114 114 While[p < 115, {n = Floor[N[1.045 +
NIntegrate[1/Log[x], {x, 2, p}], 15]]; A = Floor[(n + 1)*Log[n +
1] + (n + 1)*(Log[Log[n + 1]]  0.9385)]; T = G[A]; m = A  1; B
= Timing[While[m > p, (T = (T + 1)*G[m]; m)]]; k = p; p = p
+ 1 + T; Print["NXT(", k, ")=", p, " ", PrimeQ[p], " ",
B]}]

10 input P:clr
time:'pnsmr.ub 15 N=1.39*P/log(P) 17
A=int((N+1)*log(N+1)+(N+1)*(log(log(N+1))0.9385) 18 S3=0 20
for L=1 to isqrt(A) 30 S3=S3+A\L(A1)\L 40 next 50 T=
((22*S3)\A) 120 Z=A1 130 for M=Z to P+1 Step 1 160
S3=0 170 for L=1 to isqrt(M) 180 S3=S3+M\L(M1)\L 190 next
L 200 G=((22*S3)\M) 210 T=(T+1)*G 230 next M 240 print
"The next prime is=";P+1+T; 250 print time:end 
We can also reduce the order of the explicit
formula p(n) from O(n log n)^(5/2) to O(n log
n)^(3/2) compute Pi(k) recursively as:
pi(k)=pi(k1)+F(k) where F(k) is the
characteristic function of prime numbers. There is
more improvements in the article than we have written Jonathan Sondow and me jointly. You can see
the version PDF here : http://arxiv.org/PS_cache/math/pdf/0210/0210312.pdf . (This article has been send to The American Mathematical
Monthly). We thank P. Sebah for discussions on the
bounds.
***
SMR reported (May, 2003) the following
interesting old reference:
An explicit Formula for the kth prime number,
Sthephen Regimbal, Mathematics Magazine, Volume 48, Issue 4 (Sep.
1975), pp. 230232.
***
New improvements to the nextprm function came from
William Bouris. When Sebastian revised the Bouris improvement he
made one more improvement to the Bouris work. Here are the successive
codes and their outputs, as ran by Sebastian (who asks where does come
from the limit for the line 15?)
Original
code 
Bouris improvement
to the original code 
Sebastian
improvement to the Bouris improvement 
10 input P:clr
time:'pqsmr.ub 20 S1=0 30
for M=P+1 to 2*P 40 F=1 50
for I=P+1 to M:S2=0 60 for J=1 to
isqrt(I):S2=S2+I\J(I1)\J:next J 70
S2=22*S2:F=F*((S2\I)) 80 next
I 90 S1=S1+F 100 next
M 110 print
"nxt(";P;")=";P+1+S1,time:end 
5 cls 10
C=0 12 input N:clr
time:'pqsmr.ub 15
D=int(N+sqrt(N)+9) 20
S=0 30 F=1 40 for I=(N+1)
to D 50 T=0 60 for J=1 to
int(sqrt(I)) 65
T=T+int(I/J)int((I1)/J) 67
C=C+1 70 next J 75
T=2(2*T) 80
F=F*(int(T/I)) 85
S=S+F 90 if T=F then I=D 100
next I 110 print "The next prime
after";N;"is";N+1+S;"found in";C;"steps",time:end 120
goto 12 
5 cls 10
C=0 12 input N:clr
time:'pqsmr.ub 15
D=int(N+sqrt(N)+9) 30
F=1 40 for I=(N+1) to D 50
T=0 60 for J=1 to
int(sqrt(I)) 65
T=T+(I\J)((I1)\J) 67
C=C+1 70 next J 75
T=2(2*T) 80
F=F*((T\I)) 90 if T=F then goto
110 100 next I 110 print "The next
prime after";N;"is";I;"found
in";C;"steps",time:end 120 goto
10 
nxtp(100)=101 time =0 sec
nxtp(500)=503 time=4
sec
nxtp(1000)=1009
time=26 sec

nxtp(100)=101
time=0 sec
nxtp(1000)=1009 time=0
sec
nxtp(10000)=10007 time=0
sec
nxtp(100000)=100003 time=1
sec
nxtp(1000000)=1000003 time=2
sec
nxtp(10000000)=10000019 time=25
sec

nxtp(100)=101
time=0 sec
.
.
.
nxtp(10000)=10007
time=0 sec
nxtp(100000)=100003
time=0 sec
nxtp(1000000)=1000003
time=0 sec
nxtp(10000000)=10000019 time=0
sec
nxtp(100000000)=100000007
time=0 sec
nxtp(1000000000)=1000000007 time=0 sec
nxtp(10000000000)=10000000019
time=4 sec
nxtp(100000000000)=100000000003 time=3
sec
nxtp(1000000000000)=1000000000039 time: 98
sec

***
SMR has
developed (March 2004) another formula for the p(n) function, now based in
a totally distinct
characteristic function for the prime numbers recalling the basic
arithmetic function LCM:
His new p(n) formula is this
one:
which he improved as a
calculating tool using the Rosen and Schoenfeld
bound:
obtaining:
This last formula (for the
prime p(200) ) is approximately 8.5 times faster than the previous one,
and 65 times faster than the faster p(n) function not using the LCM
characteristic function.
SMR now is investigating the
following question "Which is the order of this
algorithm?"
"***
