3
$\begingroup$

I'd like to speed up this calculation:

 First[Sort[Table[{k,Abs[Round[k/Pi]Pi - k]1.{k,1,100000}],#1[[2]]<#2[[2]] &]]

The code above generates a list of approximations of multiples of $\pi$ out of the first 100000 integers, and sorts through it to find the best approximation. I'd like to calculate this for values larger than 100000, but it takes too long using the above code. So I'm wondering: Since sorting seems to be the bottleneck, is there a faster way to do this? I'm guessing there is some way to run through the list once, without sorting, and find the best fit.

Can anyone help? Any help would be greatly appreciated!

$\endgroup$
  • 2
    $\begingroup$ I don't understand what you're doing, but this does the same quite fast Ordering[Abs[Round[#/Pi] Pi - #] 1. & /@ Range[1000000], 1] $\endgroup$ – Dr. belisarius Sep 28 '14 at 12:48
  • $\begingroup$ @Belisarius: Wow! It is very fast! I'm technically finding the best approximation for a multiple of $\pi$ out of the integers from 1 to 100,000, or any large number. Thanks! $\endgroup$ – Matt Groff Sep 28 '14 at 12:52
  • $\begingroup$ Can Rationalize[Pi, tolerance] helps you? $\endgroup$ – ybeltukov Sep 28 '14 at 14:40
  • 2
    $\begingroup$ The best rational approximations to a real number are the continued fractions, a classic result. In fact, you can only get a better approximation than the convergent if you increase the denominator. There are numerous simple proofs/theorems about this, the earliest, to my knowledge, from Lagrange (unless you count Euclid in 300 BC, because the Euclidean algorithm gives you the c. f. approximation). So the comment from ybeltukov above and the answer from ubpqdn below are dead-on (as Rationalize gives you the c. f. approximation). $\endgroup$ – Andreas Lauschke Sep 28 '14 at 14:49
  • $\begingroup$ @AndreasLauschke Note that Rationalize[Pi, 1/1000] does not give a c.f. convergent. (From the docs on Convergents.) $\endgroup$ – Michael E2 Sep 28 '14 at 15:11
3
$\begingroup$

This is not meant as an "answer" in the traditional sense, but as a series of comments to previous comments.

@Michael E2: I guess I should have been more precise in my phrasing. Most of the fractions in the output list returned by Rationalize are indeed convergents, see also Michael Trott's Numerics book, page 74 -- although that was written for 5.1/5.2, but it seems Rationalize hasn't changed. And the approximation criterion for Rationalize (see documentation) is pretty much the lattice criterion for diophantine approximation, that means several of the list elements Rationalize will output are members of the convergents list. Indeed, if you want convergents, use Convergents. Both are useful, Convergents limits with n and Rationalize limits with dx, and these are both useful. For a future version of M it would be great if one could specify which convergent type to use, because there are dozens of convergent types (even, odd, etc.) for c. f.s, and Convergents only outputs the canonical convergents (but those are the most useful for general purposes).

@Matt Groff: I may misunderstand you, but if you want the closest integer to a multiple of Pi, couldn't you use Round?

Two more notes:

a) shameless self-promotion: I have written a demonstration on fast c. f. convergence towards Pi that may be instructive in this context: http://demonstrations.wolfram.com/ApproximatingPiWithContinuedFractions/

It includes the comparison to a contraction of the fastest-converging c. f. expansion known to me which I have developed for very fast convergence acceleration from the even branch of convergents.

b) With that said, I want to point out that shortly after I had uploaded that demonstration in 2011, Jon McLoone wrote a good blog post that makes the point that all rational approximations to Pi are useless. http://blog.wolfram.com/2011/06/30/all-rational-approximations-of-pi-are-useless/

Paraphrasing what he wrote, all rational best approximations (and those are the c. f. convergents) to Pi need to carry so much "information" or "data" with them, that this "cost" does not justify the increase in precision of the approximation. I agree with Jon. Rational approximations to Pi are "expensive" in terms of digit count (in fact, you can clearly see how fast the numerators and denominators of the Convergents output grow), so in order to get more ("data") efficient approximations to Pi, you would have to resort to other methods (exp, log, sqrt, trig, inverse trig, hypergeom, sums, prods, etc.). For practical applications rational approximations to Pi are indeed a bad choice.

$\endgroup$
  • $\begingroup$ +1 I wouldn't call your answer "shameless self-promotion" but very instructive $\endgroup$ – eldo Sep 28 '14 at 18:04
  • $\begingroup$ I'm looking for the integer (in a given range) that best approximates any multiple of pi. In other words, Out of the integers from 1 to 1000, 355 is the closest to a multiple of pi, with an error of 0.0000301444. No other integer in this range comes closer to a multiple of pi. So that's really what I'm after. I want to go through all integers in a range, and find the integer that best approximates a multiple of pi. $\endgroup$ – Matt Groff Sep 28 '14 at 18:04
  • $\begingroup$ All I meant to point out was that Rationalize is not guaranteed to give a convergent. But I have another question you may know the answer to. Given a real x and a convergent p/q, then is Abs[q x - p] <= Abs[b x - a] for all a, b with b <= q? An answer to this tells us whether Convergents can be used to answer this question. I think this is what Matt was wondering in a comment. (I know it is so for Abs[x - p/q] < Abs[x - a/b], but that is different.) $\endgroup$ – Michael E2 Sep 28 '14 at 18:06
2
$\begingroup$

Update: Integers closest to a multiple of π in various ranges in [1, 100 000 000]:

Quiet[FixedPointList[imPi[# - 1] &, 100000000][[-3 ;; 2 ;;-1]]] // AbsoluteTiming
(* {8.175765, {1, 3, 22, 333, 355, 103993, 104348, 208341, 312689, 833719, 
               1146408, 4272943, 5419351, 80143857}} *)

... and between 100 000 000 and 300 000 000 a single number (245 850 922) is added to this list.


Two modifications of the method suggested by belisarius in the comments make it even faster:

  1. Since both Abs and Round are Listable we can avoid using Map.

  2. Using the second argument of Round we can avoid the division and multiplication in the original.

Ordering[Abs[# - Round[#, (1.) π]] &@(Range[1000000]), 1] // AbsoluteTiming
(* {0.044031, {833719}} *)
Ordering[Abs[Round[#/ π] π - #] (1.) & /@ Range[1000000], 1] // AbsoluteTiming
(* {0.108076, {833719}} *)

Ordering[Abs[# - Round[#, (1.) π]] &@(Range[100000]), 1] // AbsoluteTiming
(* {0.004005, {355}} *)
Ordering[Abs[Round[#/π] π - #] (1.) & /@ Range[100000], 1] // AbsoluteTiming
(* {0.009005, {355}} *)

ClearAll[imPi];
imPi = Ordering[Abs[# - Round[#, (1.) π]] &@(Range[#]),1][[1]] &; 


ListPlot[Table[{n, First@imPi[n]}, {n, 1, 100}], Joined -> True]

enter image description here

ListLogLinearPlot[Table[{n, First@imPi[n]}, {n, 1, 1000, 10}], Joined -> True]

enter image description here

ListLogLogPlot[Table[{n, First@imPi[n]}, {n, 10000, 1000000, 10000}], Joined -> True]

enter image description here

$\endgroup$
1
$\begingroup$

Perhaps, the convergents of continued fractions is helpful to you:

con = Convergents[Pi, 20];
ListPlot[con[[1 ;; 10]], GridLines -> {None, {Pi}}]

enter image description here

TableForm[{N[#, 30], N[(# - Pi), 30]} & /@ con, 
 TableHeadings -> {Range[20], {"Convergent", "Con-Pi"}}]

enter image description here

$\endgroup$
  • $\begingroup$ @MichaelE2 I did not (and do not) sufficiently understand the question and only posted Convergents in event may have been helpful to aim. Fully accept may be a misdirected answer. $\endgroup$ – ubpdqn Sep 28 '14 at 13:13
  • $\begingroup$ I think he wants the fraction, so Convergents is definitely helpful. My deleted suggestion gave the fraction, for max. integer 100000, but didn't work for other bounds. $\endgroup$ – Michael E2 Sep 28 '14 at 13:15
  • $\begingroup$ @MichaelE2: I've updated the question, in the hopes of avoiding confusion. This answer is still interesting, but I'm trying to get the integer in the range (0,10000) that best approximates a multiple of pi. I'm hoping that I can find a faster way to calculate this so that I can extend the range. I hope this helps! $\endgroup$ – Matt Groff Sep 28 '14 at 13:47

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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