0
$\begingroup$

I want a simple method of visualizing the prime numbers that reveals the apparent tendency of certain quadratic polynomials to generate unusually large numbers of primes.

I was able to display the prime numbers using this

Table[Prime[n], {n, 41538}] 

I am new to Mathematica and I would like a guide on how I can do this in mathematica. I am working on something similar and i wanted to see some of the different done ulam spiral.. is there is any. Thanks in advance.

$\endgroup$

closed as unclear what you're asking by Dr. belisarius, Sjoerd C. de Vries, dr.blochwave, MarcoB, Bob Hanlon Jun 12 '15 at 0:22

Please clarify your specific problem or add additional details to highlight exactly what you need. As it's currently written, it’s hard to tell exactly what you're asking. See the How to Ask page for help clarifying this question. If this question can be reworded to fit the rules in the help center, please edit the question.

  • $\begingroup$ How would you like us to help? You don't tell us what the pattern is, nor what you want to achieve... Table[Prime[n], {n, 41538}] will generate all the primes up to 500,000. Maybe you can start from there and explain more about what you would like to do? $\endgroup$ – MarcoB Jun 11 '15 at 16:52
  • $\begingroup$ Thanks @MarcoB am actually new .. i didnt know where to start now i have an idea. Thanks for this . I believe once i try this i can come up with what i want to see. $\endgroup$ – Madona Syombua Jun 11 '15 at 16:55
  • $\begingroup$ Alternatively, you can use Prime /@ Range[41538] or since Prime is Listable: Prime@Range[41538] $\endgroup$ – Bob Hanlon Jun 11 '15 at 18:33
  • $\begingroup$ I managed to edit my question, i don't know if i can find help on the edited one? $\endgroup$ – Madona Syombua Jun 14 '15 at 18:15
5
$\begingroup$

The fastest way is probably a compiled sieve of Eratosthenes.

PrimesUpTo = Compile[{{n, _Integer}},
  Block[{S = Range[2, n]},
    Do[
      If[S[[i]] != 0,
        Do[
          S[[k]] = 0,
          {k, 2i+1, n-1, i+1}
        ]
      ],
      {i, Sqrt[n]}
    ];
    Select[S, Positive]
  ],
  CompilationTarget -> "C",
  RuntimeOptions -> "Speed"
];

It's fast for your input:

PrimesUpTo[500000] // Length // AbsoluteTiming
{0.009409, 41538}

It's fast(ish) for larger inputs too, but if you go any higher than the following examples, you'll probably want to use a segmented sieve.

PrimesUpTo[10^7] // Length // AbsoluteTiming
{0.337327, 664579}
PrimesUpTo[10^8] // Length // AbsoluteTiming
{3.546054, 5761455}

To find all primes less than $n$ in a segmented fashion, you'll need to first find all primes $\leq \sqrt n$.

SegmentedPrimeSieve = Compile[{{l, _Integer}, {hi, _Integer}, {primes, _Integer, 1}},
  Module[{lo = Max[l, 2], S, sqrt = Round[Sqrt[hi]], p, st},
    S = Range[lo, hi];
      Do[
        p = primes[[i]];

        If[p > sqrt, Break[]];
        st = lo + Mod[-lo, p];

        If[st > hi, Continue[]];
        If[st == p, st += p];
        Do[
          S[[k]] = 0,
          {k, st - lo + 1, hi - lo + 1, p}
        ],
        {i, 1, Length[primes]}
      ];

      Select[S, Positive]
    ],
    CompilationTarget -> "C",
    RuntimeOptions -> "Speed"
];

Now let's count all primes less than 1 billion in batches of length 10^4. This takes ~17.5 sec on my computer.

n = 10^9;
ps = PrimesUpTo[Ceiling[Sqrt[n]]];
Δ = 10^4;

Total @ Table[
  Length[SegmentedPrimeSieve[1 + i*Δ, (i + 1)Δ, ps]], 
  {i, 0, n/Δ - 1}] // AbsoluteTiming
{17.582762, 50847534}

If you have a lot of cores you can speed things up with ParallelTable. With 12 cores, I get about a 7.25 times speedup.

LaunchKernels[];
DistributeDefinitions[SegmentedPrimeSieve];

Total @ ParallelTable[
  Length[SegmentedPrimeSieve[1 + i*Δ, (i + 1)Δ, ps]], 
  {i, 0, n/Δ - 1}] // AbsoluteTiming
{2.413624, 50847534}

Let's try 10 billion in batches of length 10^5.

n = 10^10;
ps = PrimesUpTo[Ceiling[Sqrt[n]]];
Δ = 10^5;

Total @ ParallelTable[
  Length[SegmentedPrimeSieve[1 + i*Δ, (i + 1)Δ, ps]], 
  {i, 0, n/Δ - 1}] // AbsoluteTiming
{30.046231, 455052511}

The only two ways I can think of to get more substantial speedups are

  • Have factors of 2, 3, and 5 pre-sieved (aka wheel sieve). This will reduce the amount of work by 70%!
  • Low level type optimizations you can't really do in Mathematica would provide substantial speed ups :(.

The fastest prime generator I know of is primesieve, and it employs the above ideas. On my machine it's ~60 times faster to generate all primes < 10^9 and ~71.5 times faster to generate all primes < 10^10... Impressive!

enter image description here

$\endgroup$
  • $\begingroup$ Chip, how does your function compare to the built in? $\endgroup$ – MarcoB Jun 12 '15 at 1:15
  • 1
    $\begingroup$ @MarcoB To find all primes less than $n$, this sieve runs in $O(n \log\log n)$ time. Prime[Range[PrimePi[n]]] (the built in version) runs in essentially $O(n^{5/3})$ time for larger $n$ because Prime[x] uses a $O(n^{2/3})$ time algorithm (see reference.wolfram.com/language/tutorial/…). For n == 10^8, this sieve is about 4 times faster, and the larger $n$ grows, the faster the sieve will be compared to the built in way. $\endgroup$ – Chip Hurst Jun 12 '15 at 2:12
  • $\begingroup$ Chip, thank you for the details, and for the link: I had not actually come across that page yet, and it made for interesting reading. $\endgroup$ – MarcoB Jun 12 '15 at 4:31

Not the answer you're looking for? Browse other questions tagged or ask your own question.