Inhaltsverzeichnis

Development of
Algorithmic Constructions

20:07:04
Deutsch
28.Mar 2024


Prime sieving for the polynomial f(n)=n2+1


0. Abstract

1. Description of the basic algorithm

2. Programming and algorithms written in Mupad and SageMath
a) Programming of the basic algorithm
b) little improvement of the basic algorithm
c) version of the basic algorithm for a better estimation
d) Improved programming for speed I, estimation of the work
e) Improved programming for speed II, second appearance instead the first
f) Improved programming for speed III, Hensel-Lifting
g) Improved programming for speed IV, Hensel-Lifting with a better limit
h) Improved programming for speed V, multiplication instead of division
i) Improved programming for speed VI, Logarithm
j) basic algorithm in the complex field for Gaussian prime elements
k) very simple basic algorithm for an estimation of the irreducible primes

3. Programming and algorithms written in C / C++ for download
a) simple and easy understanding implementation
b) improved and easy understanding implementation
c) Used programm for the investigation in C, 2012
d) Used programm for the investigation in C++ with Heap, 2015

4. Results of the distribution of the primes
a) Table up to 1012 and table up to 241
b) Graphic of the distribution of the primes
c) Graphic of the proportion between the "reducible" and the "irreducible" primes
d) Graphic of the distribution of all primes concerning their huge
e) Table of the distribution of all primes concerning their second appearance
f) Comparison between the amount of primes between their 2. appearance and the "irreducible" primes of the form p=n2+1
g) Table of the distribution of the primes mod 6 and mod 8

5. Runtime of the algorithm
a) Table of the amount of divisions
b) Table of the use of memory
c) Estimation of the runtime

6. Mathematical theory
a) If p | f(n) then also p | f(n+k*p)
b) If p | f(n) then also p | f(-n)
c) If p | f(-n) then also p | f(-n+k*p)
d) If p is a primitiv prime factor then p > n
e) If p is a primitiv prime factor then p > 2n
f) If p is a primitiv prime factor p with p | f(n) then n < p/2
g) Let f°(n+1) ∈ ℕ
which remains after the sieving procedure resp. when f(n+1) / f°(m) with m from 0 up to n
if f°(n+1) > 1 then f°(n+1) is a
primitiv prime factor and a Stormer prime
h) Let f°°(n+1) ∈ ℕ
which remains when f°°(n+1) is divided as often as possible
by the second root of the prime f°°(m) with m from 0 up to n
and f°°(n+1) > 1 then f°°(n+1) is a prime

i) Hensel-lifting
j) if p ∈ ℙ\2 and with p | f(n)=n2+1 then p | f(y)=4y2+1.
k) With the help of the chinese remainder it is possible to calculate a n
where f(n) is either a prime or the product of new primes (Extension to Euclid)
and f(n) is not divisible by the primes found before.

l) if n+I is divisible by a+bI then n-I is divisible by a-bI and vice visa.
m) Every prime p = 1 mod 4 is a factor of f(n)=4n2+1 for a special n ∈ ℕ
n) if p | f(n)=n2+1=(n+I)(n-I) then p | f(n+kp)=(n+kp+I)(n+kp-I)
o) all p=1 mod 4 is "sieved out" for f(n)=n²+1 for n<p/2
p)

7. Sequences of primes and reference to Njas-Sequences
a) Primes of form n2 + 1
b) unsorted list of Primes and 1
c) Primitive prime factors of the sequence k2 + 1 in the order that they are found
d) Prime factors of numbers of the form n2 + 1 which themselves are not of this form
e) Primes congruent to 1 or 2 modulo 4; or, primes of form n2+y2; or, -1 is a square mod p
f) number of primes of the form n2 + 1 < 10n
g) number of distinct prime divisors (taken together) of numbers of the form n2+1 for n<=10n
h) Pythagorean primes: primes of form 4n + 1
i) Numbers n such that n2 + 1 is prime

8. Links
a) A Sieve Method for Factoring Numbers of the Form n2 + 1, by Daniel Shanks, 1959
b) On the Conjecture of Hardy and Littlewood concerning the Number of Primes of the form n2+a, by Daniel Shanks, 1960
c) On the Gaussian Primes on the Line Im(x)=1, by M. C. Wunderlich, 1973
d) On Primes Represented by Quadratic Polynomials, by Stephan Baier and Liangyi Zhao, 2007
e) Search for primes of the form n2 + 1 by Marek Wolf, 2010
f) Batemann Horn Conjecture
g) Gaußsche Zahl
h) Gaussian Integer from MathWorld
i) Gaussian Prime from MathWorld

0. Abstract:


There are several aspects of the following article:

Landau mentioned four "unattackable" problems during the Fifth Congress of Mathematicians in Cambridge in 1912.
One of these is the question, if there are infinitively many "irreducible" primes of the form p=n2+1.

This question is conjoined with the question if there are infinitively many gaussian prime elements of the form p=n+I in the complex plane.

Besides there is still the "unsettled question" "whether the reducible numbers have a definity density" (Shanks)

Daniel Shanks published 1959 a paper in which he elaborated an algorithm with the calculation of the primes up to n=180.000

The following research describes some better algorithms and some numerical results about the distribution of primes resp. of prime ideals concerning the irreducible polynomial f(n)=n²+1.

The main part in this article is the construction of several algorithms. For the illustration the progamming in Mupad was choosen. For the research a program in C with GMP was written.

2012, calculation up to 238 was done in approximately 30 days on Amd Athlon 64 600e Quad-Core processor on one core with 16 Gbyte Ram and 2 Terabyte of disk space.
2015, with some better implementation in C a calucation up to 240 was done on the same computer.
2019, with the help of a more effectiv datastructure calculation up to 241 was done on Amd 6300 Six-Core with 32 Gbyte Ram and 2,5 Terabyte of disk space.
2023, implementing a very easy understanding algorithm in SageMath in order to make an estimation of the amount of the irreducible primes.

I offer three possibilities to understand the algorithms, a descriptive, a programming and a mathematical way.
May the reader choose the best possibility for his own.

1. Description of the basic algorithm


a) A list with the function terms f(n)=n2+1 from n=0 up to n_max is created.

b) The algorithm starts with n=0, f(0)=1, nothing is done, because the function value f(0)=1 is not a prime.

c) The n is increased by one, f(1)=2, 2 is a prime.
    As n+kp and p-n+kp with k element of N, resp. 1+2k and 2-1+2k describe the same values,
    the prime 2 occurs only one time periodically in this sequence as divisor.
    The prime 2 is a singularity, because it divide the discriminant of f(n). All other primes occur twice peridically in the sequence.
    All function values f(1+2k) are divided by 2.

d) The n is increased by one.
    Let f*(n) is denoted as the function value f(n) which remains after dividing by the primes which occurs before.
    if f*(n)=1 nothing is done and the step d) is repeated.

e) If p=f*(n) is a prime, it occurs peridically twice for all function values at f(n+kp) and f(p-n+kp) with k element ℕwith n+kp < n_max and p-n+kp < n_max.
    These function values are divided as often as possible by the prime so that the divisor of p is eleminiated.

f) Go to step d) untill n=n_max

Explication by example

a) Create a list from n=0 to n=list_max=100 with f(n)=n2+1

f(0)= 1
f(1)= 2
f(2)= 5
f(3)= 10 = 2*5
f(4)= 17
f(5)= 26 = 2*13
f(6)= 37
f(7)= 50 = 2*52
f(8)= 65 = 5*13
f(9)= 82 = 2*41
f(10)= 101
f(11)= 122 = 2*61
f(12)= 145 = 5*29
f(13)= 170 = 2*5*17
f(14)= 197
f(15)= 226 = 2*113
f(16)= 257
f(17)= 290 = 2*5*29
f(18)= 325 = 52*13
f(19)= 362 = 2*181
f(20)= 401
f(21)= 442 = 2*13*17
f(22)= 485 = 5*97
f(23)= 530 = 2*5*53
f(24)= 577
f(25)= 626 = 2*313
f(26)= 677
f(27)= 730 = 2*5*73
f(28)= 785 = 5*157
f(29)= 842 = 2*421
f(30)= 901 = 17*53
f(31)= 962 = 2*13*37
f(32)= 1025 = 52*41
f(33)= 1090 = 2*5*109
f(34)= 1157 = 13*89
f(35)= 1226 = 2*613
f(36)= 1297
f(37)= 1370 = 2*5*137
f(38)= 1445 = 52*17
f(39)= 1522 = 2*761
f(40)= 1601
f(41)= 1682 = 2*292
f(42)= 1765 = 5*353
f(43)= 1850 = 2*52*37
f(44)= 1937 = 13*149
f(45)= 2026 = 2*1013
f(46)= 2117 = 29*73
f(47)= 2210 = 2*5*13*17
f(48)= 2305 = 5*461
f(49)= 2402 = 2*1201
f(50)= 2501 = 41*61
f(51)= 2602 = 2*1301
f(52)= 2705 = 5*541
f(53)= 2810 = 2*5*281
f(54)= 2917
f(55)= 3026 = 2*17*89
f(56)= 3137
f(57)= 3250 = 2*53*13
f(58)= 3365 = 5*673
f(59)= 3482 = 2*1741
f(60)= 3601 = 13*277
f(61)= 3722 = 2*1861
f(62)= 3845 = 5*769
f(63)= 3970 = 2*5*397
f(64)= 4097 = 17*241
f(65)= 4226 = 2*2113
f(66)= 4357
f(67)= 4490 = 2*5*449
f(68)= 4625 = 53*37
f(69)= 4762 = 2*2381
f(70)= 4901 = 132*29
f(71)= 5042 = 2*2521
f(72)= 5185 = 5*17*61
f(73)= 5330 = 2*5*13*41
f(74)= 5477
f(75)= 5626 = 2*29*97
f(76)= 5777 = 53*109
f(77)= 5930 = 2*5*593
f(78)= 6085 = 5*1217
f(79)= 6242 = 2*3121
f(80)= 6401 = 37*173
f(81)= 6562 = 2 17 193
f(82)= 6725 = 52*269
f(83)= 6890 = 2*5*13*53
f(84)= 7057
f(85)= 7226 = 2*3613
f(86)= 7397 = 13*569
f(87)= 7570 = 2*5*757
f(88)= 7745 = 5*1549
f(89)= 7922 = 2*17*233
f(90)= 8101
f(91)= 8282 = 2*41*101
f(92)= 8465 = 5*1693
f(93)= 8650 = 2*52*173
f(94)= 8837
f(95)= 9026 = 2*4513
f(96)= 9217 = 13*709
f(97)= 9410 = 2*5*941
f(98)= 9605 = 5*17*113
f(99)= 9802 = 2*132*29
f(100)= 10001 = 73*137

(I added the factorization of the results, because it shows the periodical appearance of the primes.)

b) f(0)=1 nothing is done

c) f(1)=2, divide f(1+2k) / 2 with 1+2k<=n_max. k element N

d) f(2)=5, divide f(2+5k) / 5 as often as there is no factor 5 in the result.
    You have to divide f(7), f(12), f(17), f(22), f(27), f(32), f(37), f(42), f(47), f(52), f(57), f(62) and so on by 5.

    divide f(-2+5k) / 5 as often as there is no factor 5 in the result.
    You have to divide f(3), f(8), f(13), f(18), f(23), f(28), f(33), f(38), f(43), f(48), f(53), f(58) and so on by 5.

e) Go for n from 3 to liste_max
    if p=f*(n) > 1 (f*(n) is the prime which remains after dividing by f(1), f*(2) up to f*(n_max) )
    divide f(n+kp) / p and
    divide f(-n+kp) / p
    as often as there is no factor p in the result.

    You get an unsorted list of primes.

2. a) Programming of the basic algorithm


This is a short implementation for the algorithm concerning the polynom f(n)=n2+1.
Every prime (without the 2) which occurs is sieved out resp. divides periodically at two positions the initialised list. (for the proof see 6. c), d) and e) )
The procedure sieving only divides as often as possible the function values in the list by the sieving primes.

  • liste_max:=1000;

    sieving:=proc (stelle, p)
    // remove all presence of p from elements in liste.
    // Must be called with initial value of stelle pointing to list element divisible by p.
    begin
       while (stelle<=liste_max) do
           erg:=liste[stelle];
           repeat
              erg:=erg /p;
           until (erg mod p>0) end_repeat;
           liste[stelle]:=erg;
           stelle:=stelle+p;
       end_while;
    end_proc;

    for x from 1 to liste_max do
       liste [x]:=x2+1;
    end_for;

    for x from 1 to liste_max do
        p:=liste[x]; 
        if (p>1) then
          print (x, p); 
          sieving (x+p,  p);
          if p > 2 then
            sieving (p-x, p);
          end_if;
        end_if;
    end_for;

You get as result an unsorted list of primes.

2. b) little improvement of the basic algorithm


The initialization of the list is separated in two loops, the prime number 2 is directly sieved out by the initialization.
Therefore one if statement could be deleted.

  • liste_max:=1000;

    sieving:=proc (stelle, p)
    // remove all presence of p from elements in liste.
    // Must be called with initial value of stelle pointing to list element divisible by p.
    begin
       while (stelle<=liste_max) do
           erg:=liste[stelle];
           repeat
              erg:=erg /p;
           until (erg mod p>0) end_repeat;
           liste[stelle]:=erg;
           stelle:=stelle+p;
       end_while;
    end_proc;

    for x from 2 to liste_max step 2 do
       liste [x]:=(x2+1);
    end_for;

    print (1, 2);
    for x from 1 to liste_max step 2 do
       // No Division, only a rightshift is needed
       liste [x]:=(x2+1)/2;
    end_for;

    for x from 2 to liste_max do
        p:=liste[x];  
        if (p>1) then
          print (x, p);
          sieving (x+p,  p);
          sieving (p-x, p);
        end_if;
    end_for;

You get as result an unsorted list of primes.

2. c) version of the basic algorithm for a better estimation


The main loop for x from 1 to liste_max is divided into two loops with the same huge.
In the second loop only the second sieving part is made because the first sieving part is redundant, as the prime p=f(x) > x (for the proof see 6. f) )

  • liste_max:=1000;

    sieving:=proc (stelle, p)
    // remove all presence of p from elements in liste.
    // Must be called with initial value of stelle pointing to list element divisible by p.
    begin
       while (stelle<=liste_max) do
           erg:=liste[stelle];
           repeat
              erg:=erg /p;
           until (erg mod p>0) end_repeat;
           liste[stelle]:=erg;
           stelle:=stelle+p;
       end_while;
    end_proc;

    for x from 2 to liste_max step 2 do
       liste [x]:=(x2+1);
    end_for;

    for x from 1 to liste_max step 2 do
       // No Division, only a rightshift is needed
       liste [x]:=(x2+1)/2;
    end_for;

    for x from 2 to liste_max/2 do
        p:=liste[x];  
        if (p>1) then
          // print (x, p);
          sieving (x+p,  p);
          sieving (p-x, p);
        end_if;
    end_for;

    for x from liste_max/2+1 to liste_max do
        p:=liste[x];  
        if (p>1) then
          // print (x, p);
          sieving (p-x, p);
        end_if;
    end_for;

2. d) Improved programming for speed I


Using that every prime p=f(x) > 2x (for the proof see I. g) ) the loop could be divide in two loops,
where the second loop contains only one sieving.

  • liste_max:=10000;

    sieving:=proc (stelle, p)
    begin
       while (stelle<=liste_max) do
           erg:=liste[stelle];
           repeat
              erg:=erg /p;
              anz_div:=anz_div+2;
           until (erg mod p>0) end_repeat;
           liste[stelle]:=erg;
           stelle:=stelle+p;
       end_while;
    end_proc;

    for x from 2 to liste_max step 2 do
       liste [x]:=(x2+1);
    end_for;

    for x from 1 to liste_max step 2 do
       // No Division, only a right shift is needed
       liste [x]:=(x2+1)/2;
    end_for;

    for x from 2 to trunc (liste_max/3) do
        p:=liste[x];  
        if (p>1) then
          print (x, p);
          sieving (x+p,  p);
          sieving (p-x, p);
        end_if;
    end_for;

    for x from trunc (liste_max/3+1) to liste_max do
        p:=liste[x];  
        if (p>1) then
          print (x, p);
          sieving (p-x, p);
        end_if;
    end_for;

This idea is not used in the following described algorithms.

2. e) Improved programming for speed II


Instead of sieving for the first appearance of the prime,
the sieving is made by the second appearance of the prime.
This is an improvement because less primes have to store and the numbers of divisions is less.

    liste_max:=1000;

    sieving:=proc (stelle, p)
    begin
       while (stelle<=liste_max) do
           erg:=liste[stelle];
           repeat
              erg:=erg /p;
           until (erg mod p>0) end_repeat;
           liste[stelle]:=erg;
           stelle:=stelle+p;
       end_while;
    end_proc;

    for x from 2 to liste_max step 2 do
       liste [x]:=(x2+1);
    end_for;

    for x from 1 to liste_max step 2 do
       // No Division, only a right shift is needed
       liste [x]:=(x2+1)/2;
    end_for;

    for x from 2 to liste_max do
        p:=liste[x];  
        if (p>1) then

          // This decide whether the prime occurs for the second time.
          if x > p-x then
            sieving (x+p,  p);
            sieving (2*p-x, p);
          else

            // The primes which occures for the first time are printed
             print (x, p);
            
          end_if;
        end_if;
    end_for;

2. f) Improved programming for speed III, Hensel-Lifting


Mathematical background : Hensel-Lifting for polynomials

if p(x) is a polynomial and if f | p(x1)
then the Hensel-Lifting calculates the x2, x3, ..., xn so that
f2 | p(x2)
f3 | p(x3)
...
fn | p(xn)

Using the Hensel-Lifting only one division is needed instead of one division and a modulo operation in the sieving function.

  • liste_max:=1000;

    f:=proc (x)
    begin
       return (x2+1);
    end;

    fd:=proc (x)
    begin
       return (2*x);
    end;

    sieving:=proc (s, p, a)
    begin
       while (s<=liste_max) do
          liste[s]:=liste[s]/p;
          s:=s+a;
       end_while;
    end_proc;

    hensel:=proc (x, p)
    begin
         sieving (x+p, p, p);
         sieving (2*p-x, p, p);

         s:=x;
         a:=p;
         inv:=fd(s)^(-1) mod p;
         repeat
            a:=a*p;
            s:=s-f(s)*inv;
            s:=s mod a;
            sieving (s, p, a);
            sieving (a-s, p, a);
         until s>liste_max and (a-s)>liste_max end_repeat;
    end_proc;

    for x from 2 to liste_max step 2 do
       liste [x]:=(x2+1);
    end_for;

    for x from 1 to liste_max step 2 do
       // No Division, only a right shift is needed
       liste [x]:=(x2+1)/2;
    end_for;

    // The prime number 2 is eleminated by the initialisation
    print (1, 2);

    for x from 2 to liste_max-1 do
       p:=liste[x];  

       if (p>1) then
          
          // Second occurance of the prime
          if (x > p-x) then
             hensel (x, p);
          else
          // Primes p at position x
          // print (x, p);
          end_if;
       end_if;
    end_for;


2. g) Improved programming for speed IV, Hensel-Lifting with a better limit


Instead of the calculation the Hensel-lifting for the primes from the second appearance concerning the polynomial
the limit is decreased so that only the Hensel-lifting is made if the primes occurs for the 3. time up to the max limit
of calculation.

  • liste_max:=1000;

    f:=proc (x)
    begin
       return (x2+1);
    end;

    fd:=proc (x)
    begin
       return (2*x);
    end;

    sieving:=proc (s, p, a)
    begin
       while (s<=liste_max) do
          liste[s]:=liste[s]/p;
          s:=s+a;
       end_while;
    end_proc;

    hensel:=proc (x, p)
    begin
         sieving (x+p, p, p);
         sieving (2*p-x, p, p);

         s:=x;
         a:=p;
         inv:=fd(s)^(-1) mod p;
         repeat
            a:=a*p;
            s:=s-f(s)*inv;
            s:=s mod a;
            sieving (s, p, a);
            sieving (a-s, p, a);
         until s>liste_max and (a-s)>liste_max end_repeat;
    end_proc;

    for x from 2 to liste_max step 2 do
       liste [x]:=(x2+1);
    end_for;

    for x from 1 to liste_max step 2 do
       // No Division, only a right shift is needed
       liste [x]:=(x2+1)/2;
    end_for;

    // The prime number 2 is eleminated by the initialisation
    print (1, 2);

    for x from 2 to liste_max-1 do
       p:=liste[x];  

       if (p>1) then
          
          // Second occurance of the prime
          if (x > p-x) then
             // 3. appearance of the prime
            if ((2*p-x)<liste_max) then
               hensel (x, p);
            end_if;
          else
          // Primes p at position x
          // print (x, p);
          end_if;
       end_if;
    end_for;


2. h) Improved programming for speed V, multiplication instead of division


Multiplication is faster as division.
Therefore you can change the initialising and the following sieving.
You need to use the Hensel-lifting in order to determine the power of the factors.

Before : Initialising with f(n) and dividing by the factors.

Better : Initialising with 1 and multiplication by the factors.

2. i) Improved programming for speed V, Logarithm


The divisions of the function value f(n)=n2+1 by the sieving primes is made
by substructions of the logarithm values of the primes from the logarithm of f(n).
By using this replacement you have to pay attention to the error estimation of the logarithm and the maximum amount of substractions / divisions by the primes

  • liste_max:=100;
    log_base:=10;
    // the prime number 2 is eliminated
    count_prime:=1;

    // return the function value of the polynomial f(x)=x2+1
    f:=proc (x)
    begin
       return (x2+1);
    end;

    // return the first derivation / f&#id_180;(x)=2x
    fd:=proc (x)
    begin
       return (2*x);
    end;

    // intialisation of list with the log value base log_base where x is even
    for x from 2 to liste_max step 2 do
       liste [x]:=log (log_base, f(x));
    end_for;

    // intialisation of list with the log value base log_base where x is odd
    for x from 3 to liste_max step 2 do
       liste [x]:=log (log_base, f(x)/2);
    end_for;

    // sieving operation / division of p by substraction of the log values
    sieving:=proc (s, p, a, log_p)
    begin
       while (s<=liste_max) do
          liste[s]:=liste[s]-log_p;
          s:=s+a;
       end_while;
    end_proc;

    // Calculation of the hensel lifting
    hensel:=proc (x, p, log_p)
    begin
         sieving (x+p,   p, p, log_p);
         sieving (2*p-x, p, p, log_p);

         s:=x;
         a:=p;
         inv:=fd(s)^(-1) mod p;
         repeat
            a:=a*p;
            s:=s-f(s)*inv;
            s:=s mod a;
            sieving (  s, p, a, log_p);
            sieving (a-s, p, a, log_p);
         until s>liste_max and (a-s)>liste_max end_repeat;
    end_proc;

    // main loop 
    for x from 2 to liste_max do
       log_p:=liste[x];
       p:=round (log_base^log_p);  
       if (p > 1) then
          if (p > 2*x) then
             // First occurance of the primes / primes p at position x
             print (x, p);
             count_prime:=count_prime+1;
          else
             // Second occurance of the prime
             // 3. appearance of the prime
            if ((2*p-x)<liste_max) then
              hensel (x, p, log_p);
            end_if;
          end_if;
       end_if;
    end_for;

    print (count_prime, liste_max);

2. j) basic algorithm in the complex field for Gaussian prime elements



This algorithm is interesting from the mathematical point of view and the additional results.
The sieving is not made in ℕbut in the complex field ℂ, so that every prime has a relationship
to a complex gaussian number and the conjugated complex number p=(a+bi)*(a-bi)=a2+b2
This algorithm is a possible way to generate gaussian primes, although they are not ordered by huge.

  • liste_max:=100;

    division_complex:=proc (zelle, p)
    begin
       // first complex division
       erg:=zelle/p;
       repeat
          // complex division
          erg_2:=erg/p;
          // if complex result a+bI not a and b element of N return
          if ((domtype (Re (erg_2)) <> DOM_INT) or (domtype (Im (erg_2)) <> DOM_INT)) 
             then return (erg);
          end_if;
          erg:=erg_2;
       until (FALSE) end_repeat;
    end_proc;


    // Sieving by the complex argument p
    sieving:=proc (stelle, p, prim)
    begin
       while (stelle<=liste_max) do
           liste[stelle]:=division_complex (liste[stelle], p);
           stelle:=stelle+prim;
       end_while;
    end_proc;

    for x from 1 to liste_max do
       liste [x]:=x+I;
    end_for;

    for x from 1 to liste_max do
        p:=liste[x]; 
        prim:=p*conjugate (p);
        if (prim>1) then
          print ("x=", x, "prim = ", prim, "=", p,"*", conjugate (p), "=", (Re (p))2, "+", (Im(p))2 ); 
          // First sieving
          sieving (x+prim,  p, prim);
          if (prim>2) then
            // Second sieving by the conjugate argument of p
            sieving (prim-x, conjugate (p), prim);
          end_if;
        end_if;
    end_for;


x= 1 prim =  2 = (1 + I) * (1 - I) = 1 + 1

x= 2 prim =  5 = (2 + I) * (2 - I) = 4 + 1

x= 4 prim =  17 = (4 + I) * (4 - I) = 16 + 1

x= 5 prim =  13 = (3 - 2 I) * (3 + 2 I) = 9 + 4

x= 6 prim =  37 = (6 + I) * (6 - I) = 36 + 1

x= 9 prim =  41 = (5 - 4 I) * (5 + 4 I) = 25 + 16

x= 10 prim =  101 = (10 + I) * (10 - I) = 100 + 1

x= 11 prim =  61 = (6 - 5 I) * (6 + 5 I) = 36 + 25

x= 12 prim =  29 = (5 - 2 I) * (5 + 2 I) = 25 + 4

x= 14 prim =  197 = (14 + I) * (14 - I) = 196 + 1

x= 15 prim =  113 = (8 - 7 I) * (8 + 7 I) = 64 + 49

x= 16 prim =  257 = (16 + I) * (16 - I) = 256 + 1

x= 19 prim =  181 = (10 - 9 I) * (10 + 9 I) = 100 + 81

x= 20 prim =  401 = (20 + I) * (20 - I) = 400 + 1

x= 22 prim =  97 = (9 - 4 I) * (9 + 4 I) = 81 + 16

x= 23 prim =  53 = (7 - 2 I) * (7 + 2 I) = 49 + 4

x= 24 prim =  577 = (24 + I) * (24 - I) = 576 + 1

x= 25 prim =  313 = (13 - 12 I) * (13 + 12 I) = 169 + 144

x= 26 prim =  677 = (26 + I) * (26 - I) = 676 + 1

x= 27 prim =  73 = (3 - 8 I) * (3 + 8 I) = 9 + 64

x= 28 prim =  157 = (11 + 6 I) * (11 - 6 I) = 121 + 36

x= 29 prim =  421 = (15 - 14 I) * (15 + 14 I) = 225 + 196

x= 33 prim =  109 = (10 - 3 I) * (10 + 3 I) = 100 + 9

x= 34 prim =  89 = (8 - 5 I) * (8 + 5 I) = 64 + 25

x= 35 prim =  613 = (18 - 17 I) * (18 + 17 I) = 324 + 289

x= 36 prim =  1297 = (36 + I) * (36 - I) = 1296 + 1

x= 37 prim =  137 = (4 - 11 I) * (4 + 11 I) = 16 + 121

x= 39 prim =  761 = (20 - 19 I) * (20 + 19 I) = 400 + 361

x= 40 prim =  1601 = (40 + I) * (40 - I) = 1600 + 1

x= 42 prim =  353 = (17 - 8 I) * (17 + 8 I) = 289 + 64

x= 44 prim =  149 = (10 + 7 I) * (10 - 7 I) = 100 + 49

x= 45 prim =  1013 = (23 - 22 I) * (23 + 22 I) = 529 + 484

x= 48 prim =  461 = (19 + 10 I) * (19 - 10 I) = 361 + 100

x= 49 prim =  1201 = (25 - 24 I) * (25 + 24 I) = 625 + 576

x= 51 prim =  1301 = (26 - 25 I) * (26 + 25 I) = 676 + 625

x= 52 prim =  541 = (21 - 10 I) * (21 + 10 I) = 441 + 100

x= 53 prim =  281 = (16 - 5 I) * (16 + 5 I) = 256 + 25

x= 54 prim =  2917 = (54 + I) * (54 - I) = 2916 + 1

x= 56 prim =  3137 = (56 + I) * (56 - I) = 3136 + 1

x= 58 prim =  673 = (23 + 12 I) * (23 - 12 I) = 529 + 144

x= 59 prim =  1741 = (30 - 29 I) * (30 + 29 I) = 900 + 841

x= 60 prim =  277 = (14 - 9 I) * (14 + 9 I) = 196 + 81

x= 61 prim =  1861 = (31 - 30 I) * (31 + 30 I) = 961 + 900

x= 62 prim =  769 = (25 - 12 I) * (25 + 12 I) = 625 + 144

x= 63 prim =  397 = (19 - 6 I) * (19 + 6 I) = 361 + 36

x= 64 prim =  241 = (15 + 4 I) * (15 - 4 I) = 225 + 16

x= 65 prim =  2113 = (33 - 32 I) * (33 + 32 I) = 1089 + 1024

x= 66 prim =  4357 = (66 + I) * (66 - I) = 4356 + 1

x= 67 prim =  449 = (7 - 20 I) * (7 + 20 I) = 49 + 400

x= 69 prim =  2381 = (35 - 34 I) * (35 + 34 I) = 1225 + 1156

x= 71 prim =  2521 = (36 - 35 I) * (36 + 35 I) = 1296 + 1225

x= 74 prim =  5477 = (74 + I) * (74 - I) = 5476 + 1

x= 77 prim =  593 = (8 - 23 I) * (8 + 23 I) = 64 + 529

x= 78 prim =  1217 = (31 + 16 I) * (31 - 16 I) = 961 + 256

x= 79 prim =  3121 = (40 - 39 I) * (40 + 39 I) = 1600 + 1521

x= 80 prim =  173 = (13 - 2 I) * (13 + 2 I) = 169 + 4

x= 81 prim =  193 = (12 - 7 I) * (12 + 7 I) = 144 + 49

x= 82 prim =  269 = (10 - 13 I) * (10 + 13 I) = 100 + 169

x= 84 prim =  7057 = (84 + I) * (84 - I) = 7056 + 1

x= 85 prim =  3613 = (43 - 42 I) * (43 + 42 I) = 1849 + 1764

x= 86 prim =  569 = (20 - 13 I) * (20 + 13 I) = 400 + 169

x= 87 prim =  757 = (9 - 26 I) * (9 + 26 I) = 81 + 676

x= 88 prim =  1549 = (35 + 18 I) * (35 - 18 I) = 1225 + 324

x= 89 prim =  233 = (8 - 13 I) * (8 + 13 I) = 64 + 169

x= 90 prim =  8101 = (90 + I) * (90 - I) = 8100 + 1

x= 92 prim =  1693 = (37 - 18 I) * (37 + 18 I) = 1369 + 324

x= 94 prim =  8837 = (94 + I) * (94 - I) = 8836 + 1

x= 95 prim =  4513 = (48 - 47 I) * (48 + 47 I) = 2304 + 2209

x= 96 prim =  709 = (22 + 15 I) * (22 - 15 I) = 484 + 225

x= 97 prim =  941 = (10 - 29 I) * (10 + 29 I) = 100 + 841


2. k) very simple basic algorithm for an estimation of the irreducible primes


The following program is written in SageMath and is very simple.
There is only one sieving function call. Every prime p = 1 mod 4 is a factor of f(n)=n2+1 for a special n ∈ ℕ (Proof)
the sieving is theoretically made on the line n+I in the complex plane
and should give an easier way for an estimation of the amount of the irreducible primes.

liste = []     
# huge of the list
liste_max=129
# initialising the list 
for j in range(liste_max):liste.append (j^2+1)

# devide as often as possible the function value f(j)=j^2+i 
# by the prime p at the position "stelle" in the list
def sieving (stelle, p):
    while (stelle<liste_max):
        while ((liste [stelle])%p==0): (liste [stelle])/=p
        stelle+=p;
    return ()   

for j in range (liste_max):
    p=liste [j]
    print (j,"+ kp + i, prime : ", p)
    if (p>1): sieving (j+p, p) 


k ∈ ℕ, p>1 is the sieved out prime.

0 + kp + i, prime :  1
1 + kp + i, prime :  2
2 + kp + i, prime :  5
3 + kp + i, prime :  5
4 + kp + i, prime :  17
5 + kp + i, prime :  13
6 + kp + i, prime :  37
7 + kp + i, prime :  1
8 + kp + i, prime :  13
9 + kp + i, prime :  41
10 + kp + i, prime :  101
11 + kp + i, prime :  61
12 + kp + i, prime :  29
13 + kp + i, prime :  17
14 + kp + i, prime :  197
15 + kp + i, prime :  113
16 + kp + i, prime :  257
17 + kp + i, prime :  29
18 + kp + i, prime :  1
19 + kp + i, prime :  181
20 + kp + i, prime :  401
21 + kp + i, prime :  1
22 + kp + i, prime :  97
23 + kp + i, prime :  53
24 + kp + i, prime :  577
25 + kp + i, prime :  313
26 + kp + i, prime :  677
27 + kp + i, prime :  73
28 + kp + i, prime :  157
29 + kp + i, prime :  421
30 + kp + i, prime :  53
31 + kp + i, prime :  37
32 + kp + i, prime :  41
33 + kp + i, prime :  109
34 + kp + i, prime :  89
35 + kp + i, prime :  613
36 + kp + i, prime :  1297
37 + kp + i, prime :  137
38 + kp + i, prime :  1
39 + kp + i, prime :  761
40 + kp + i, prime :  1601
41 + kp + i, prime :  1
42 + kp + i, prime :  353
43 + kp + i, prime :  1
44 + kp + i, prime :  149
45 + kp + i, prime :  1013
46 + kp + i, prime :  73
47 + kp + i, prime :  1
48 + kp + i, prime :  461
49 + kp + i, prime :  1201
50 + kp + i, prime :  61
51 + kp + i, prime :  1301
52 + kp + i, prime :  541
53 + kp + i, prime :  281
54 + kp + i, prime :  2917
55 + kp + i, prime :  89
56 + kp + i, prime :  3137
57 + kp + i, prime :  1
58 + kp + i, prime :  673
59 + kp + i, prime :  1741
60 + kp + i, prime :  277
61 + kp + i, prime :  1861
62 + kp + i, prime :  769
63 + kp + i, prime :  397
64 + kp + i, prime :  241
65 + kp + i, prime :  2113
66 + kp + i, prime :  4357
67 + kp + i, prime :  449
68 + kp + i, prime :  1
69 + kp + i, prime :  2381
70 + kp + i, prime :  1
71 + kp + i, prime :  2521
72 + kp + i, prime :  1
73 + kp + i, prime :  1
74 + kp + i, prime :  5477
75 + kp + i, prime :  97
76 + kp + i, prime :  109
77 + kp + i, prime :  593
78 + kp + i, prime :  1217
79 + kp + i, prime :  3121
80 + kp + i, prime :  173
81 + kp + i, prime :  193
82 + kp + i, prime :  269
83 + kp + i, prime :  1
84 + kp + i, prime :  7057
85 + kp + i, prime :  3613
86 + kp + i, prime :  569
87 + kp + i, prime :  757
88 + kp + i, prime :  1549
89 + kp + i, prime :  233
90 + kp + i, prime :  8101
91 + kp + i, prime :  101
92 + kp + i, prime :  1693
93 + kp + i, prime :  173
94 + kp + i, prime :  8837
95 + kp + i, prime :  4513
96 + kp + i, prime :  709
97 + kp + i, prime :  941
98 + kp + i, prime :  113
99 + kp + i, prime :  1
100 + kp + i, prime :  137
101 + kp + i, prime :  5101
102 + kp + i, prime :  2081
103 + kp + i, prime :  1061
104 + kp + i, prime :  373
105 + kp + i, prime :  149
106 + kp + i, prime :  661
107 + kp + i, prime :  229
108 + kp + i, prime :  2333
109 + kp + i, prime :  457
110 + kp + i, prime :  12101
111 + kp + i, prime :  1
112 + kp + i, prime :  193
113 + kp + i, prime :  1277
114 + kp + i, prime :  317
115 + kp + i, prime :  389
116 + kp + i, prime :  13457
117 + kp + i, prime :  1
118 + kp + i, prime :  557
119 + kp + i, prime :  1
120 + kp + i, prime :  14401
121 + kp + i, prime :  7321
122 + kp + i, prime :  229
123 + kp + i, prime :  1
124 + kp + i, prime :  15377
125 + kp + i, prime :  601
126 + kp + i, prime :  15877
127 + kp + i, prime :  1613
128 + kp + i, prime :  1


4. a) Results of the distribution of the primes


Legend of the tables:

Column A = exponent concerning base 10 resp. base 2
Column B = is the interval [0,x]
Column C = "irreducible" primes of the form p=x2+1
Column D = "reducible" primes with p | x2+1 and p < x2+1, counted by the first appearance resp. for the smallest x.
Column E = 1, f(n)=n^2+1 decomposes into all "preceded" primes, the "preceded" primes are found earlier by the described algorithm


Column A = exponent concerning base 10
Column B = is the interval [0,x]
Column C = all primes resp. D+E, the amount of the primes depends of p(x) for the interval [0,x] and is not calculated by primes < x
Column D = "irreducible" primes of the form p=x2+1
Column E = "reducible" primes with p | x2+1 and p < x2+1, counted by the first appearance resp. for the smallest x.

The numbers in the rows F-K are rounded with 6 digits after the decimal point.

A B C D E F G H I J K  
10^n x all Primes P(x)=x^2+1 P(x) | x^2+1 C/B D/B E/B C(n) / C(n-1) D(n) / D(n-1) E(n) / E(n-1)  
1 10 7 5 2 0,700000 0,500000 0,200000        
2 100 70 19 51 0,700000 0,190000 0,510000 10,000000 3,800000 25,500000  
3 1.000 720 112 608 0,720000 0,112000 0,608000 10,285714 5,894737 11,921569  
4 10.000 7.102 841 6.261 0,710200 0,084100 0,626100 9,863889 7,508929 10,297697  
5 100.000 70.780 6.656 64.124 0,707800 0,066560 0,641240 9,966207 7,914388 10,241814  
6 1.000.000 704.537 54.110 650.427 0,704537 0,054110 0,650427 9,953899 8,129507 10,143269  
7 10.000.000 7.026.559 456.362 6.570.197 0,702656 0,045636 0,657020 9,973300 8,433968 10,101360  
8 100.000.000 70.122.424 3.954.181 66.168.243 0,701224 0,039542 0,661682 9,979625 8,664571 10,070968  
9 1.000.000.000 700.184.485 34.900.213 665.284.272 0,700184 0,034900 0,665284 9,985172 8,826155 10,054435  
10 10.000.000.000 6.993.568.566 312.357.934 6.681.210.632 0,699357 0,031236 0,668121 9,988180 8,950029 10,042640  
11 100.000.000.000 69.870.544.960 2.826.683.630 67.043.861.330 0,698705 0,028267 0,670439 9,990686 9,049502 10,034688  
12 1.000.000.000.000 698.175.242.376 25.814.570.672 672.360.671.704 0,698175 0,025815 0,672361 9,992412 9,132458 10,028669  


A B C D E F G H I J K  
2^n x all Primes P(x)=x^2+1 P(x) | x^2+1 C/B D/B E/B C(n) / C(n-1) D(n) / D(n-1) E(n) / E(n-1)  
                       
1 2 2 2 0 1,000000 1,000000 0,000000        
2 4 3 3 0 0,750000 0,750000 0,000000 1,500000 1,500000    
3 8 5 4 1 0,625000 0,500000 0,125000 1,666667 1,333333    
4 16 12 7 5 0,750000 0,437500 0,312500 2,400000 1,750000 5,000000  
5 32 22 10 12 0,687500 0,312500 0,375000 1,833333 1,428571 2,400000  
6 64 46 14 32 0,718750 0,218750 0,500000 2,090909 1,400000 2,666667  
7 128 90 24 66 0,703125 0,187500 0,515625 1,956522 1,714286 2,062500  
8 256 181 43 138 0,707031 0,167969 0,539063 2,011111 1,791667 2,090909  
9 512 364 70 294 0,710938 0,136719 0,574219 2,011050 1,627907 2,130435  
10 1.024 736 114 622 0,718750 0,111328 0,607422 2,021978 1,628571 2,115646  
11 2.048 1.459 212 1.247 0,712402 0,103516 0,608887 1,982337 1,859649 2,004823  
12 4.096 2.920 393 2.527 0,712891 0,095947 0,616943 2,001371 1,853774 2,026464  
13 8.192 5.810 713 5.097 0,709229 0,087036 0,622192 1,989726 1,814249 2,017016  
14 16.384 11.620 1.301 10.319 0,709229 0,079407 0,629822 2,000000 1,824684 2,024524  
15 32.768 23.198 2.459 20.739 0,707947 0,075043 0,632904 1,996386 1,890085 2,009788  
16 65.536 46.339 4.615 41.724 0,707077 0,070419 0,636658 1,997543 1,876779 2,011862  
17 131.072 92.654 8.418 84.236 0,706894 0,064224 0,642670 1,999482 1,824052 2,018886  
18 262.144 185.106 15.867 169.239 0,706123 0,060528 0,645596 1,997820 1,884890 2,009105  
19 524.288 369.794 29.843 339.951 0,705326 0,056921 0,648405 1,997742 1,880822 2,008704  
20 1.048.576 738.766 56.534 682.232 0,704542 0,053915 0,650627 1,997777 1,894381 2,006854  
21 2.097.152 1.476.544 106.787 1.369.757 0,704071 0,050920 0,653151 1,998663 1,888899 2,007758  
22 4.194.304 2.949.935 203.025 2.746.910 0,703319 0,048405 0,654914 1,997865 1,901215 2,005399  
23 8.388.608 5.894.833 387.308 5.507.525 0,702719 0,046171 0,656548 1,998293 1,907686 2,004989  
24 16.777.216 11.782.897 739.671 11.043.226 0,702315 0,044088 0,658228 1,998852 1,909775 2,005116  
25 33.554.432 23.550.908 1.416.635 22.134.273 0,701872 0,042219 0,659653 1,998737 1,915223 2,004330  
26 67.108.864 47.074.786 2.716.922 44.357.864 0,701469 0,040485 0,660984 1,998852 1,917870 2,004035  
27 134.217.728 94.099.962 5.218.926 88.881.036 0,701099 0,038884 0,662215 1,998946 1,920897 2,003727  
28 268.435.456 188.106.701 10.044.585 178.062.116 0,700752 0,037419 0,663333 1,999009 1,924646 2,003376  
29 536.870.912 376.049.307 19.352.155 356.697.152 0,700446 0,036046 0,664400 1,999128 1,926626 2,003218  
30 1.073.741.824 751.784.745 37.339.024 714.445.721 0,700154 0,034775 0,665379 1,999165 1,929450 2,002948  
31 2.147.483.648 1.503.000.381 72.139.395 1.430.860.986 0,699889 0,033593 0,666297 1,999243 1,932011 2,002757  
32 4.294.967.296 3.004.917.702 139.535.723 2.865.381.979 0,699637 0,032488 0,667149 1,999279 1,934251 2,002558  
33 8.589.934.592 6.007.850.623 270.187.320 5.737.663.303 0,699406 0,031454 0,667952 1,999339 1,936331 2,002408  
34 17.179.869.184 12.012.007.361 523.695.185 11.488.312.176 0,699191 0,030483 0,668708 1,999385 1,938267 2,002263  
35 34.359.738.368 24.017.148.553 1.016.029.276 23.001.119.277 0,698991 0,029570 0,669421 1,999428 1,940116 2,002132  
36 68.719.476.736 48.021.305.927 1.973.029.796 46.048.276.131 0,698802 0,028711 0,670091 1,999459 1,941903 2,002002  
37 137.438.953.472 96.018.463.670 3.834.641.365 92.183.822.305 0,698626 0,027901 0,670726 1,999497 1,943529 2,001895  
38 274.877.906.944 191.991.204.563 7.458.662.439 184.532.542.124 0,698460 0,027134 0,671325 1,999524 1,945074 2,001789  
39 549.755.813.888 383.896.376.530 14.518.923.631 369.377.452.899 0,698303 0,026410 0,671894 1,999552 1,946585 2,001693  
40 1.099.511.627.776 767.630.202.986 28.282.415.900 739.347.787.086 0,698156 0,025723 0,672433 1,999577 1,947969 2,001605  
41 2.199.023.255.552 1.534.954.071.943 55.130.064.461 1.479.824.007.482 0,698016 0,025070 0,672946 1,999601 1,949270 2,001526

The quotient between the number of all primes [p=f(x) and p|f(x), counted by their first appearance by increasing x], and x resp. the green column F = C/B goes probably to ln (2)= 0,693147181.

4. b) Graphic of the distribution of the primes


Legend of the graphic:

The green dots / column F are the values for all primes concerning the polynom x2+1
The orange dots / column G are the values for the "irreducible" primes with p=x2+1
The yellow dots / column H are the values for the "reducible" primes by the first occurance with p|x2+1 and p smaller than x2+1
The x-axis is the logarithm to the base of 2 and goes up to 240


4. c) Graphic of the proportion between the "reducible" and the "irreducible" primes



A B C D E F          
2^n x P(x) | x^2+1 P(x)=x^2+1 C/D log (B, 3)          
                     
1 2 0 2 0,0000 0,6309          
2 4 0 3 0,0000 1,2619          
3 8 1 4 0,2500 1,8928          
4 16 5 7 0,7143 2,5237          
5 32 12 10 1,2000 3,1546          
6 64 32 14 2,2857 3,7856          
7 128 66 24 2,7500 4,4165          
8 256 138 43 3,2093 5,0474          
9 512 294 70 4,2000 5,6784          
10 1.024 622 114 5,4561 6,3093          
11 2.048 1.247 212 5,8821 6,9402          
12 4.096 2.527 393 6,4300 7,5712          
13 8.192 5.097 713 7,1487 8,2021          
14 16.384 10.319 1.301 7,9316 8,8330          
15 32.768 20.739 2.459 8,4339 9,4639          
16 65.536 41.724 4.615 9,0410 10,0949          
17 131.072 84.236 8.418 10,0067 10,7258          
18 262.144 169.239 15.867 10,6661 11,3567          
19 524.288 339.951 29.843 11,3913 11,9877          
20 1.048.576 682.232 56.534 12,0676 12,6186          
21 2.097.152 1.369.757 106.787 12,8270 13,2495          
22 4.194.304 2.746.910 203.025 13,5299 13,8805          
23 8.388.608 5.507.525 387.308 14,2200 14,5114          
24 16.777.216 11.043.226 739.671 14,9299 15,1423          
25 33.554.432 22.134.273 1.416.635 15,6245 15,7732          
26 67.108.864 44.357.864 2.716.922 16,3265 16,4042          
27 134.217.728 88.881.036 5.218.926 17,0305 17,0351          
28 268.435.456 178.062.116 10.044.585 17,7272 17,6660          
29 536.870.912 356.697.152 19.352.155 18,4319 18,2970          
30 1.073.741.824 714.445.721 37.339.024 19,1340 18,9279          
31 2.147.483.648 1.430.860.986 72.139.395 19,8347 19,5588          
32 4.294.967.296 2.865.381.979 139.535.723 20,5351 20,1898          
33 8.589.934.592 5.737.663.303 270.187.320 21,2359 20,8207          
34 17.179.869.184 11.488.312.176 523.695.185 21,9370 21,4516          
35 34.359.738.368 23.001.119.277 1.016.029.276 22,6382 22,0825          
36 68.719.476.736 46.048.276.131 1.973.029.796 23,3389 22,7135          
37 137.438.953.472 92.183.822.305 3.834.641.365 24,0398 23,3444          
38 274.877.906.944 184.532.542.124 7.458.662.439 24,7407 23,9753          
39 549.755.813.888 369.377.452.899 14.518.923.631 25,4411 24,6063          
40 1.099.511.627.776 739.347.787.086 28.282.415.900 26,1416 25,2372







4. d) Graphic of the distribution of all primes concerning their huge


2^n n=10 n=20 n=30 n=39
binary logarithm of the prime        
1 0 0 0 0
2 0 0 0 0
3 1 1 1 1
4 2 2 2 2
5 4 4 4 4
6 6 6 6 6
7 10 10 10 10
8 21 21 21 21
9 38 38 38 38
10 66 66 66 66
11 92 127 127 127
12 76 233 233 233
13 87 432 432 432
14 68 805 805 805
15 65 1.511 1.511 1.511
16 64 2.837 2.837 2.837
17 56 5.378 5.378 5.378
18 55 10.186 10.186 10.186
19 23 19.294 19.294 19.294
20 0 36.827 36.827 36.827
21 0 48.789 70.157 70.157
22 0 46.604 133.975 133.975
23 0 44.541 256.852 256.852
24 0 42.639 492.882 492.882
25 0 41.369 946.848 946.848
26 0 39.635 1.823.129 1.823.129
27 0 38.155 3.513.599 3.513.599
28 0 36.694 6.780.412 6.780.412
29 0 35.469 13.103.462 13.103.462
30 0 34.558 25.348.870 25.348.870
31 0 33.258 34.087.506 49.090.415
32 0 32.278 33.039.675 95.167.380
33 0 31.028 32.054.166 184.662.052
34 0 30.181 31.132.063 358.630.671
35 0 28.958 30.244.146 697.097.364
36 0 28.618 29.415.874 1.356.051.342
37 0 25.664 28.635.210 2.639.870.329
38 0 27.183 27.883.394 5.142.823.877
39 0 15.365 27.189.203 10.025.582.289
40 0 0 26.519.512 13.574.819.204
41 0 0 25.874.372 13.247.762.275
42 0 0 25.265.740 12.936.129.149
43 0 0 24.689.787 12.638.529.888
44 0 0 24.125.789 12.354.600.489
45 0 0 23.599.360 12.083.161.608
46 0 0 23.089.555 11.823.171.181
47 0 0 22.607.976 11.574.181.661
48 0 0 22.149.632 11.335.619.027
49 0 0 21.696.082 11.106.498.716
50 0 0 21.251.111 10.886.741.743
51 0 0 20.807.352 10.675.323.639
52 0 0 20.533.816 10.472.016.150
53 0 0 19.952.893 10.276.163.302
54 0 0 19.569.826 10.087.386.117
55 0 0 19.093.140 9.905.692.828
56 0 0 18.892.638 9.730.489.216
57 0 0 17.208.655 9.561.216.748
58 0 0 18.167.763 9.397.807.108
59 0 0 10.460.543 9.239.790.510
60 0 0 0 9.086.871.013
61 0 0 0 8.939.247.717
62 0 0 0 8.796.575.998
63 0 0 0 8.657.577.119
64 0 0 0 8.522.899.197
65 0 0 0 8.392.678.194
66 0 0 0 8.268.245.924
67 0 0 0 8.145.305.721
68 0 0 0 8.021.141.566
69 0 0 0 7.895.311.485
70 0 0 0 7.830.607.024
71 0 0 0 7.645.125.168
72 0 0 0 7.530.870.899
73 0 0 0 7.381.754.182
74 0 0 0 7.335.361.187
75 0 0 0 6.710.801.165
76 0 0 0 7.113.815.439
77 0 0 0 4.113.563.288






4. e) Table of the distribution of all primes concerning their second appearance


A

B

C

D

E

F

G

exponent
=log2 (x)

<=x

number of all primes
by their 1. appearance

number of all primes
by their 2. appearance

C / D

ln (B)

0.7*B/ln(B)

1

2

2

0



2,01

2

4

3

1

3,0000

1,3905

2,01

3

8

5

2

2,5000

2,0858

2,68

4

16

12

3

4,0000

2,7811

4,03

5

32

22

7

3,1429

3,4763

6,44

6

64

46

10

4,6000

4,1716

10,74

7

128

90

19

4,7368

4,8669

18,41

8

256

181

33

5,4848

5,5621

32,22

9

512

364

58

6,2759

6,2574

57,28

10

1.024

736

112

6,5714

6,9527

103,10

11

2.048

1.459

203

7,1872

7,6480

187,45

12

4.096

2.920

355

8,2254

8,3432

343,66

13

8.192

5.810

677

8,5820

9,0385

634,44

14

16.384

11.620

1.255

9,2590

9,7338

1.178,25

15

32.768

23.198

2.333

9,9434

10,4290

2.199,40

16

65.536

46.339

4.375

10,5918

11,1243

4.123,87

17

131.072

92.654

8.261

11,2158

11,8196

7.762,58

18

262.144

185.106

15.446

11,9841

12,5148

14.662,66

19

524.288

369.794

29.212

12,6590

13,2101

27.781,88

20

1.048.576

738.766

55.488

13,3140

13,9054

52.785,58

21

2.097.152

1.476.544

105.123

14,0459

14,6006

100.543,96

22

4.194.304

2.949.935

199.932

14,7547

15,2959

191.947,56

23

8.388.608

5.894.833

382.384

15,4160

15,9912

367.204,02

24

16.777.216

11.782.897

730.623

16,1272

16,6864

703.807,70

25

33.554.432

23.550.908

1.400.694

16,8137

17,3817

1.351.310,79

26

67.108.864

47.074.786

2.687.730

17,5147

18,0770

2.598.674,60

27

134.217.728

94.099.962

5.167.473

18,2101

18,7723

5.004.854,78

28

268.435.456

188.106.701

9.952.863

18,8998

19,4675

9.652.219,94

29

536.870.912

376.049.307

19.193.397

19,5926

20,1628

18.638.769,53

30

1.073.741.824

751.784.745

37.059.261

20,2860

20,8581

36.034.954,43

31

2.147.483.648

1.503.000.381

71.635.400

20,9813

21,5533

69.745.073,09

32

4.294.967.296

3.004.917.702

138.633.544

21,6753

22,2486

135.131.079,12

33

8.589.934.592

6.007.850.623

268.568.123

22,3699

22,9439

262.072.395,86

34

17.179.869.184

12.012.007.361

520.836.610

23,0629

23,6391

508.728.768,44

35

34.359.738.368

24.017.148.553

1.010.957.203

23,7568

24,3344

988.387.321,54

36

68.719.476.736

48.021.305.927

1.964.003.096

24,4507

25,0297

1.921.864.236,32

37

137.438.953.472

96.018.463.670

3.818.594.577

25,1450

25,7249

3.739.843.919,33

38

274.877.906.944

191.991.204.563

7.430.384.395

25,8387

26,4202

7.282.853.948,16

39

549.755.813.888

383.896.376.530

14.468.881.009

26,5326

27,1155

14.192.228.206,68

4. f) Comparison between the amount of primes between their 2. appearance and the "irreducible" primes of the form p=n2+1



A B C D E F G
log (n, 2) n amount of "irred." primes amount of "irred." primes amount of "red." primes amount of "red." primes C/F


for n = first root of p for n = sec. root of p for n = first root of p for n = sec. root of p
2 4 3 1 0 0
3 8 4 1 1 1 4
4 16 7 2 5 1 7
5 32 10 3 12 4 2,5
6 64 14 3 32 7 2
7 128 24 4 66 15 1,6
8 256 43 6 138 27 1,59259259259259
9 512 70 7 294 51 1,37254901960784
10 1.024 114 9 622 103 1,10679611650485
11 2.048 212 11 1.247 192 1,10416666666667
12 4.096 393 13 2.527 342 1,14912280701754
13 8.192 713 17 5.097 660 1,08030303030303
14 16.384 1.301 23 10.319 1.232 1,05600649350649
15 32.768 2.459 32 20.739 2.301 1,06866579747936
16 65.536 4.615 42 41.724 4.333 1,06508192937918
17 131.072 8.418 53 84.236 8.208 1,02558479532164
18 262.144 15.867 69 169.239 15.377 1,03186577355791
19 524.288 29.843 90 339.951 29.122 1,02475791497837
20 1.048.576 56.534 113 682.232 55.375 1,02093002257336
21 2.097.152 106.787 157 1.369.757 104.966 1,0173484747442
22 4.194.304 203.024 211 2.746.909 199.721 1,0165380706085
23 8.388.608 387.307 292 5.507.525 382.092 1,01364854537651
24 16.777.216 739.670 392 11.043.226 730.230 1,01292743382222
25 33.554.432 1.416.634 538 22.134.273 1.400.156 1,01176868863184
26 67.108.864 2.716.921 712 44.357.863 2.687.018 1,01112869359267
27 134.217.728 5.218.925 956 88.881.036 5.166.517 1,0101437777133
28 268.435.456 10.044.584 1.300 178.062.115 9.951.563 1,00934737588457
29 536.870.912 19.352.154 1.791 356.697.152 19.191.605 1,00836558484817
30 1.073.741.824 37.339.023 2.458 714.445.720 37.056.803 1,00761587555192
31 2.147.483.648 72.139.394 3.377 1.430.860.986 71.632.023 1,00708301928036
32 4.294.967.296 139.535.722 4.614 2.865.381.978 138.628.930 1,00654114548817
33 8.589.934.592 270.187.319 6.232 5.737.663.303 268.561.891 1,00605234046405
34 17.179.869.184 523.695.184 8.417 11.488.312.175 520.828.193 1,00550467704808
35 34.359.738.368 1.016.029.275 11.539 23.001.119.277 1.010.945.664 1,00502856996279
36 68.719.476.736 1.973.029.795 15.866 46.048.276.130 1.963.987.230 1,00460418726857
37 137.438.953.472 3.834.641.364 21.728 92.183.822.305 3.818.572.849 1,0042079896431
38 274.877.906.944 7.458.662.438 29.842 184.532.542.123 7.430.354.553 1,00380976234688
39 549.755.813.888 14.518.923.630 41.168 369.377.452.899 14.468.839.841 1,00346149308102
40 1.099.511.627.776 28.282.415.899 56.533 739.347.787.085 28.194.567.830 1,00311577994491

4. g) Table of the distribution of the primes mod 6 and mod 8


ABCDEFGHI
exponent =log2 (x)<=xnumber of primes with p=f(x) number of primes with p=f(x) and p%6=1 number of primes with p=f(x) and p%6=5 number of primes with p=f(x) and p%8=1 number of primes with p=f(x) and p%8=3 number of primes with p=f(x) and p%8=5 number of primes with p=f(x) and p%8=7
243021010
384121020
4167152040
53210274050
66414497060
712824914110120
8256431527210210
9512702445360330
101.0241144172600530
112.0482127014110901020
124.09639312726519401980
138.19271323148136003520
1416.3841.30143686465106490
1532.7682.4598431.6151.24401.2140
1665.5364.6151.5473.0672.33002.2840
17131.0728.4182.8545.5634.24604.1710
18262.14415.8675.32310.5437.93807.9280
19524.28829.84310.01319.82914.955014.8870
201.048.57656.53418.83937.69428.171028.3620
212.097.152106.78735.39271.39453.338053.4480
224.194.304203.02567.653135.371101.1190101.9050
238.388.608387.308128.939258.368193.2350194.0720
2416.777.216739.671246.669493.001369.2730370.3970
2533.554.4321.416.635472.296944.338708.1670708.4670
2667.108.8642.716.922905.5071.811.4141.358.36901.358.5520
27134.217.7285.218.9261.739.9163.479.0092.608.73602.610.1890
28268.435.45610.044.5853.348.6446.695.9405.023.04305.021.5410
29536.870.91219.352.1556.451.98812.900.1669.675.25309.676.9010
301.073.741.82437.339.02412.447.47124.891.55218.667.885018.671.1380
312.147.483.64872.139.39524.049.06348.090.33136.067.539036.071.8550
324.294.967.296139.535.72346.514.86293.020.86069.768.975069.766.7470
338.589.934.592270.187.32090.058.651180.128.668135.092.2840135.095.0350
3417.179.869.184523.695.184174.558.447349.136.736261.845.6380261.849.5450
3534.359.738.3681.016.029.274338.653.157677.376.116508.003.5360508.025.7370
3668.719.476.7361.973.029.782657.637.4591.315.392.322986.505.2530986.524.5280

ABCDEFGHI
exponent =log2 (x)<=xnumber of primes with p|f(x) number of primes with p=f(x) and p%6=1 number of primes with p=f(x) and p%6=5 number of primes with p=f(x) and p%8=1 number of primes with p=f(x) and p%8=3 number of primes with p=f(x) and p%8=5 number of primes with p=f(x) and p%8=7
240000000
381100010
4165232030
53212845070
664321913140180
7128663630280380
82561387365630750
951229414614814101530
101.02462232829431103110
112.0481.24766458363106160
124.0962.5271.3111.2161.27101.2560
138.1925.0972.6462.4512.54802.5490
1416.38410.3195.3414.9785.16105.1580
1532.76820.73910.7699.97010.366010.3730
1665.53641.72421.57920.14520.764020.9600
17131.07284.23643.45740.77942.123042.1130
18262.144169.23986.95582.28484.619084.6200
19524.288339.951174.611165.340169.8210170.1300
201.048.576682.232350.420331.812341.1260341.1060
212.097.1521.369.757702.725667.032684.6940685.0630
224.194.3042.746.9101.408.1461.338.7641.373.50401.373.4060
238.388.6085.507.5252.817.9102.689.6152.753.63502.753.8900
2416.777.21611.043.2265.645.0235.398.2035.520.15705.523.0690
2533.554.43222.134.27311.305.69010.828.58311.067.464011.066.8090
2667.108.86444.357.86422.637.75021.720.11422.179.283022.178.5810
27134.217.72888.881.03645.319.44943.561.58744.442.896044.438.1400
28268.435.456178.062.11690.718.10587.344.01189.028.573089.033.5430
29536.870.912356.697.152181.617.979175.079.173178.349.4040178.347.7480
301.073.741.824714.445.721363.493.073350.952.648357.229.4330357.216.2880
312.147.483.6481.430.860.986727.560.774703.300.212715.429.4920715.431.4940
324.294.967.2962.865.381.9791.456.149.2361.409.232.7431.432.663.67101.432.718.3080
338.589.934.5925.737.663.3032.914.225.6142.823.437.6892.868.806.51602.868.856.7870
3417.179.869.18411.488.312.1755.832.117.9655.656.194.2105.744.168.11905.744.144.0560
3534.359.738.36823.001.119.27211.671.204.53711.329.914.73511.500.598.218011.500.521.0540
3668.719.476.73646.048.276.10123.355.396.02622.692.880.07523.024.208.905023.024.067.1960

5. a) Table of the amount of divisions


The values from 2 up to 224 are not listed because the used algorithm had a first level cache segment in order to speed up the algorithm and
therefore the resulting values were not right.

A B C D E
2^n x amount of divisons B*log( log (B; 3); 3) C/B
25 33.554.432 83.524.162 84.246.004 2,49
26 67.108.864 161.330.817 170.887.811 2,40
27 134.217.728 332.155.442 346.386.366 2,47
28 268.435.456 662.288.896 701.658.821 2,47
29 536.870.912 1.348.539.023 1.420.466.101 2,51
30 1.073.741.824 2.728.193.442 2.874.066.288 2,54
31 2.147.483.648 5.532.461.967 5.812.227.625 2,58
32 4.294.967.296 11.199.801.794 11.748.575.126 2,61
33 8.589.934.592 22.682.503.676 23.737.750.633 2,64
34 17.179.869.184 45.891.918.142 47.942.335.660 2,67
35 34.359.738.368 92.828.721.324 96.791.273.397 2,70
36 68.719.476.736 187.668.985.522 195.344.667.920 2,73
37 137.438.953.472 379.269.618.308 394.117.011.261 2,76
38 274.877.906.944 766.168.931.008 794.906.542.041 2,79
39 549.755.813.888 1.547.202.551.607 1.602.811.459.228 2,81
40 1.099.511.627.776 3.123.330.927.170 3.230.961.451.946 2,84
41 2.199.023.255.552 6.303.055.060.506 6.511.348.556.497 2,87





















5. b) Table of the use of memory


The first column gives the amount of memory if the primes are stored for the first appearance,
the second column gives the amount of memory if the primes are stored for the second appearance.

n kbyte kbyte B(n+1)/B(n) C(n+1)/C(n)
218 1472 1400

220 2000 1740 1,36 1,24
222 5520 2940 2,76 1,69
224 12800 9196 2,32 3,13
226 37000 25946 2,89 2,82
228 131000 83940 3,54 3,24
230 481000 300000 3,67 3,57
232 1800000 1125000 3,74 3,75

5. c) Estimation of the runtime and efficiency


An upper limit for the runtime is O(N*ln(ln(N)) where N is the huge of the field.
The algorithm finds approximately 0.7*N primes, the irreducible primes in order, the reducible primes in a disorder by huge.

For comparison to the sieve of Eratosthenes:
The runtime of the sieve of Eratosthenes is also O(N*ln(ln(N)), by finding approximately N/ln(N) primes in order by huge.

For comparison to the sieve of Atkin:
The runtime of the sieve of Atkin is O(N/log log (N)), by finding approximately N/ln(N) primes in order by huge.

6. Mathematical theory


Let f(n) = n2 + 1 with n ∈ ℕ

The following proofs explain the mathematical background which is used for the described algorithms.

a) Proof: If p ∈ ℙ and 
             p | f(n) then also p | f(n+p)

   p | f(n)   <=> n2 + 1 = 0                  mod p
   p | f(n+p) <=> (n+p)2 + 1 = 0              mod p
              <=> n2 + 2np + p2 + 1 = 0       mod p
              <=> n2 + 1  = 0                 mod p

   Thus if p | f(n) then p | f(n+p)

   Example: f(4)=42+1=17 so 17 | f(4+17)=(4+17)2+1 = 442 = 2*13*17
   
   Therefore : If p | f(n) then p | f(n+k*p)) for k ∈ ℤ
   

b) Proof: If p | f(n) then also p | f(-n)

   p | f(n)    <=> n2 + 1 = 0                 mod p
   p | f(-n)   <=> (-n)2 + 1 = 0              mod p
               <=> n2 + 1 = 0                 mod p

   Thus if p | f(n) then p | f(-n)


c) Proof: If p | f(-n) then also p | f(-n+p)

   This is a simple conclusion of b) and c) and means that
   if p is a divisor of f(n) then p appears periodically concerning the function values of f(n)=n2+1
   at f(n) and f(-n) with the period length p

   Example: p(-4)=(-4)2+1=17 so 17 | p(-4+17)=(-4+17)2+1 = 170 =2*5*17
   
   Therefore : If p | f(-n) then p | f(-n+k*p)) for k ∈ ℤ


d) Lower estimation: If p is a primitiv prime factor 
                     resp. if p | f(n) with the smallest n > 0
                     then p > n

   supposing that p <= n then p would be found earlier by the described algorithm and would be sieved out.
   
   Assumption of the contrary : if p <= n then p | f(m) for at least one f(m) with m=0 up to m=n which implies a contradiction.
   
   Therefore if p | f(n) then p > n


e) Lower estimation: If p > 2 is a primitiv prime factor 
                     resp. if p | f(n) with the smallest n > 0
                     then p > 2n  

   As f(n)=f(-n).
   
   Assumption of the contrary : if p <= 2n then p | f(m) for at least one f(m) with m=-n up to m=n which implies a contradiction.
   
   Proof d) shows that this is not possible.


f) Upper estimation: If p is a primitiv prime factor p with p | f(n) then n < p/2 
Since every prime p=1 mod 4 has exact two 4.roots with r^4=1 mod p and r²<>1 mod p and r1 + r2 =0 mod p with r1 g) Proof: Let f°(n+1) ∈ ℕ which remains when f(n+1) is divided as often as possible by f°(m) with m from 0 up to n and f°(m)| f°(n+1). If f°(n+1) > 1 then f°(n+1) is a primitiv prime factor and a Stormer prime. Supposing f°(n+1) is not a prime, f°(n+1)=p*q with p, q ∈ ℕ greater than 1 p > 2(n+1), q > 2(n+1)) (Proof e)) f(n)=n2+1 > f°(n+1) = p*q > (2n)(2n) = 4n2 which is a contradiction. h) Proof: If f°°(n+1) ∈ ℕ which remains when f°°(n+1) is divided as often as possible by the second appearance resp. by the 2. root of the prime f°°(m) with m from 0 up to n and f°°(m)| f°°(n+1). If f°°(n+1) > 1 then f°°(n+1) is a prime. proof by induction: for n=0 f(0)=1 n=1 f(1)=2 n=2 f(2)=5 n=3 f(3)=10=2*5 (second appearance for p=5 resp. 2. root of p) Supposing f°°(n) is not equal 1 or not a prime after dividing all prime factors f°°(n) with n from 0 up to n-1 As f(n)=n2+1 is symmetrical to the y-axis, let there be the following construction : f°(-n) is the natural number which remains when f(-n) is divided as often as possible by the first appearance of the prime f°(-n) with n from 0 up to -n+1 f°°(n) is the natural number which remains when f(n) is divided as often as possible by the second appearance of the prime f°°(n) with n from 0 up to n-1 p | f°(-n) <=> p | f°°(p-n) (see Proof c) ) if f°°(p-n) is not equal 1 or a prime, then f°(-n) would also not be equal 1 or a prime, that would be a contradiction to Proof f) i) The Hensel-lifting explains for every p|f(n) where p2|f(y) p3|f(z) ... pn|f(n) can be found. j) Proof: if p ∈ ℕ\2 and with p | f(n)=n2+1 then p | f(y)=4y2+1. The substitution with n=2y gives this result. Therefore all primes without the 2 which are divisor of the f(n)=n2+1 are a divisor of f(y)=4y2+1. Nevertheless the sequence of the primes found by the described algorithm is not the same and the proof e) is not true for the function f(y)=4y2+1 k) With the help of the chinese remainder it is possible to calculate a n where f(n) is either a prime or the product of new primes and f(n) is not divisible by the primes found before. Example: The first prime values after the sieving process for the polynom fn)=n2+1 are f(0)=1 f(1)=2 f(2)=5 f(3)=1 f(4)=17 Supposing that the sequence f(5)-f(oo)=1 p1=2 p2=5 p3=17 are the first primes choose a special n1, n2 and n3 for the primes n1=2 n2=4 n3=5 (Other values are also possible, f(ni) should not be divisible by pi) With the help of the chinese remainder theorem you could solve the equitations : n=0 mod 2 f(0)=1 not divisible by 2 n=4 mod 5 f(4)=17 not divisible by 5 n=5 mod 17 f(5)=26 not divisible by 17 The smallest solution calculated by the chinese remainder theorem is n=124 mod 170 f(124) is not divisible by the primes 2, 5, 17 f(124)=15377 which is by chance a primes. By this way you could calculate a special n taking the first primes of the sequence calculated by the algorithm and you get as result a special f(n) which is either a prime or consists of primes which are not yet in the sequence. By the way the amount of reducible primes / 2^m for m=2 up to oo is monoton increasing. l) Proof: if n+I is divisible by a+bI (a!=0 and b!=0) then n-I is divisible by a-bI and vice visa. Proof: (n+I) / (a+bI) = (n+I)(a-bI)/(a2+b2) = (an+b+aI-bnI)/(a2+b2) = (an+b)/(a2+b2) + (a-bn)I / (a2+b2) = c + dI with a, b, c ,d ∈ ℤ\0, n ∈ ℤ (n-I) / (a-bI) = (n-I)(a+bI)/(a2+b2) = (an+b-aI+bnI)/(a2+b2) = (an+b)/(a2+b2) - (a-bn)I / (a2+b2) = c - dI if (a+bI) / (c+dI) = e+fI then (a-bI) / (c-dI) = e-fI the proof is analogue to the sentence before. m) Statement: Every prime p = 1 mod 4 is a factor of f(n)=n2+1 for a special n ∈ ℕ Proof: if p = 1 mod 4 then -1 is a square (according the quadratic reciprocity law, 1. supplement) in Fp* and hence there is some a element Fp* such that a2 = -1. Put n´:=a element Fp* and choose some n element N such that n´= n mod p then n2 = n´2 = a2 = -1 mod p. Hence p divides f(n) = n2+1. n) Proof: if p | f(n)=n2+1=(n+I)(n-I) then p | f(n+kp)=(n+kp+I+lpI)(n+kp-I+lpI) because n+kp+I+lpI = n+I mod p. o) Let a1, a2, ..., an be primes with ai | f(n)=n^2+1 with ai