Help:=proc(): print(`PNT(n) , YAP(L), EucSeq(n) `): 
print(`TPL(n) , piJ(n), thetaJ(x) `):
end:

with(numtheory):

HelpOld:=proc(): print(`SE(n), IsDiv(n,L), OneStep(L)`) :
print(`PrimeList(N), PrimeListM(N) `):
end :


#SE(n): Sieve of Era...
SE:=proc(n) local i,j:
{seq(i,i=1..n)} minus 
{seq(seq(j*i,i=2..n/j),j=2..trunc(sqrt(n)))} minus {1}:

end:


#IsDiv(n,L): inputs a pos. integer n and a list
#of pos. integers L, and outputs true if
#n is divisible by any of the members of L

IsDiv:=proc(n,L) local i:
for i from 1 to nops(L) do 
 if n mod L[i]=0 then
  RETURN(true):
 fi:
od:
false:

end:


#OneStep(L): inputs the list of the first primes
#and outputs the same list with one more prime that
#comes right after
OneStep:=proc(L) local i,n:
#n is the largest prime so-far
n:=L[nops(L)]:

 for i from n+1 while IsDiv(i,L) do od:

[op(L),i]:

end:

#The list of the first N prime numbers, our way
PrimeList:=proc(N) local i,L:
L:=[2]:
for i from 2 to N do
 L:=OneStep(L):
od:
L:
end:

#The list of the first N prime numbers, the Maple way
PrimeListM:=proc(N) local i:
 [seq(ithprime(i),i=1..N)]:
end:

#PNT(n): The ratio of #primes <=n and n/log(n)
#(it should go to 1, if the the Prime Number Theorem is
#correct)
PNT:=proc(n):evalf( nops(SE(n))/(n/log(n))): end:

#inputs a list of primes and outputs another prime
#not in L (an algorithmic implementation of Euclid's
#"proof"
YAP:=proc(L) local B,p:

B:=mul(p,p in L)+1:

max(op(factorset(B))):

end:

#The first n terms in the sequence of primes
#implementing Euclid's "proof" that there are infinitely
#many primes
EucSeq:=proc(n) local L,i:
L:=[2]:

for i from 1 to n-1 do
 L:= [op(L),YAP(L)]:
od:

L:

end:

#The first n pairs of twin-primes
TPL:=proc(n) local L,a:
a:=2:
L:=[]:

while nops(L)<n do

  if nextprime(a)-a=2 then
    L:=[op(L),a]:
  fi:
a:=nextprime(a):

od:

L:

end:

#piJ(x): the number of primes less than x
piJ:=proc(x) local i,s:
s:=0:

for i from 1 to x do

 if isprime(i) then
    s:=s+1:
 fi:

od:

s:

end:



#thetaJ(x): the number of primes less than x
thetaJ:=proc(x) local i,s:
s:=0.:

for i from 1 to x do

 if isprime(i) then
    s:=evalf(s+log(i)):
 fi:

od:

s:

end: