4
$\begingroup$

I want to generate a graphic using

DiscretePlot[primePercent[x], {x, 10^10}, 
  ImageSize -> Full, 
  Filling -> None, 
  Joined -> True] 

where

primePercent[x_Integer] := 100 PrimePi[x]/x

When I generate it for 10^6, my computer becomes slow as there are simply too many data points. What I'm wondering is if I can get some load off my RAM by continually saving the evaluated data to a file (or using a different evaluation method/graphic options, of course)?

Does Export[] save the file only once it has been generated?

$\endgroup$

4 Answers 4

7
$\begingroup$

For a nice plot you don't need all the million points in your graphic. 10^10 points is a bit heavy for this curve, so why not using Plot? It has the advantage of adaptive sampling and it seems to do it absolutely correct in this case. The runtime is only about 1 second:

Plot[100 PrimePi[x]/x, {x, 1, 10^10}, 
  PlotRange -> {Automatic, {0, 10}}]

enter image description here

And even if you reduce your range, it catches all features

Plot[100 PrimePi[x]/x, {x, 1, 50}, PlotRange -> Automatic]

enter image description here

Update on Accuracy

In the comment you asked how accurate this is. I assume you mean "are there enough points in the plot to represent important features?" Consider the following Manipulate where you can visually inspect the sampling points in your plot with different settings for PlotPoints and MaxRecursion

Manipulate[
 Plot[100 PrimePi[x]/x, {x, 1, xEnd}, PlotRange -> Automatic, 
  Mesh -> All, MeshStyle -> Directive[PointSize[.005], Red], 
  MaxRecursion -> maxRecursion, PlotPoints -> plotPoints, 
  AxesOrigin -> {0, 0}],
 {xEnd, 10, 1000}, {maxRecursion, 0, 5, 1}, {plotPoints, 5, 30, 1}]

enter image description here

Basically those options work like that: PlotPoints defines how many sampling points are in your region initially. With the MaxRecursion setting Plot has the opportunity to subdivide the region at places where the function is not represented properly. This can be deduced by some curvature measure, but what exactly Mathematica is using for this decision is afaik not documented.

Therefore, MaxRecursion is more important to catch features than many PlotPoints alone. On the other hand when you want a plot from 1 to 10^10 than you are clearly more interested in the overall tendency of the function than in the structure of the sharp edges between 4 and 10, which couln't be displayed anyway.

So when you are asking about accuracy, I'm not really sure what to answer. Every (red) point you see above is as accurate as it would be in DiscretePlot and usually the default settings for Plot work nicely and catch most features you want to see. But to have a feeling for that you can now use the Manipulate and watch what happens.

In return I could ask, in your DiscretePlot do you see the difference between point 10^5 and point 10^5+1? Or is this maybe the wrong scale and you are more interested in the global bend from 0 to 2*10^9, because then you could leave out many points and exactly this is, what Plot tries to do.

$\endgroup$
4
  • $\begingroup$ That is unimaginably faster. Thanks. How accurate is it, though? $\endgroup$
    – CHM
    Feb 14, 2012 at 15:06
  • $\begingroup$ @CHM I updated my post. $\endgroup$
    – halirutan
    Feb 15, 2012 at 3:08
  • $\begingroup$ Great. The Manipulate helped. My concern was if Plot would too much approximation, but of course this can be dealt with by tweaking the options. Thanks. $\endgroup$
    – CHM
    Feb 15, 2012 at 3:29
  • $\begingroup$ @halirutan would you please join me in chat if you are here. $\endgroup$
    – Mr.Wizard
    Feb 15, 2012 at 12:43
5
$\begingroup$

Keep in mind that 10^10 points is an outrageously large amount of computation no matter what hoops you jump through. Even something much more modest, say 10^7 points, is going to eat up your RAM when using any Graphics object to store it all in.

I suggest creating a series of rasterized images and using some basic image processing to join them together. The following code does that, though I warn you, it is slow. The plus side is that it won't run you out of memory!

primePercent[x_Integer] := 100 PrimePi[x]/x

Here step, max, and isize are the number of points to plot in intermediate images, the total number of desired points, and the image size respectively.

step = 10^5;
max = 10^7;
isize = 500;

I store an initial image img0.

img0 = Rasterize[
   DiscretePlot[primePercent[k], {k, step}, 
    PlotRange -> {{0, max}, {0, 60}}, ImageSize -> Full, 
    Filling -> None, Joined -> True], ImageSize -> isize];

Here, I set up a dynamic progress indicator. Its nice to know how long you have to wait if you are going to be waiting a while!

Dynamic[ProgressIndicator[(i/max)]]

This is where I put it all together. I'm using ImageMultiply to combine rasterized images. This will overemphasize repeated components which is why AxesStyle is set to 0 opacity.

Do[img0 = 
  ImageMultiply[img0, 
   Rasterize[
    DiscretePlot[primePercent[k], {k, i - step + 1, i}, 
     PlotRange -> {{0, max}, {0, 60}}, ImageSize -> Full, 
     Filling -> None, Joined -> True, 
     AxesStyle -> Directive[Opacity[0]]], ImageSize -> isize]]
     , {i, 2 step, max, step}
]

The resulting image looks like...

enter image description here

The easiest first step to getting a higher quality image would be to set isize to something larger.

Edit: As David points out in his response, it is better for this particular example to specify a step size in DiscretePlot. Only use the method I've provided if David's approach isn't reasonable for your particular application.

$\endgroup$
4
$\begingroup$

Have you already tried Export["picture.png", DiscretePlot[100 PrimePi[k]/k, {k, 10^10}, ImageSize -> Full, Filling -> None, Joined -> True]]? This directly saves to a file instead of you having to wait for the front end to display the image in addition to the computational effort expended by evaluating your function.

$\endgroup$
1
  • $\begingroup$ I did. It works up to a point (10^6) well. after that, good luck. :/ $\endgroup$
    – CHM
    Feb 14, 2012 at 14:52
4
$\begingroup$

In this specific case, you don't want to have all the $10^{11}$ data points DiscretePlot generated. Luckily, there's an option for that: Instead of specifying {x, 10^11} like you did, there's an optional parameter dealing with the step size the plot uses. The bracket is actually of the form {x, xMin, xMax, dx}, where dx is the step width to be taken, just like for example in the Table function.

Benchmark:

AbsoluteTiming@DiscretePlot[
    100 PrimePi[x]/x,
    {x, 1, 10^11, 10^9},
    ImageSize -> 600,
    Filling -> None,
    Joined -> True
]

enter image description here

Way better, isn't it? 100 points are easier to plot than 100 billion ;-)

A shortcoming here is that DiscretePlot does not provide adaptive sampling, i.e. automatically generates denser point sets where the data fluctuates. This is something you'll have to take care of on your own somehow.

$\endgroup$
1
  • $\begingroup$ This is also very fast. I think the slowest function is PrimePi: both your method and halirutan's are very quick up to 10^11. Once you get to 10^12, evaluation time takes long enough for me not to bother evaluating it. As this is only for a qualitative purpose, I don't mind not having 10^12. $\endgroup$
    – CHM
    Feb 14, 2012 at 15:23

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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