Path: csiph.com!v102.xanadu-bbs.net!xanadu-bbs.net!feeder.erje.net!eu.feeder.erje.net!newsfeed.xs4all.nl!newsfeed3.news.xs4all.nl!xs4all!newsgate.cistron.nl!newsgate.news.xs4all.nl!post.news.xs4all.nl!not-for-mail
To: python-list@python.org
From: Dave Angel <davea@davea.name>
Subject: Re: Working with the set of real numbers
Date: Wed, 5 Mar 2014 07:57:49 -0500 (EST)
Organization: news.gmane.org
References: <8e4c1ab1-e65d-483f-ad9d-6933ae2052c3@googlegroups.com> <CAPTjJmo4qD2B7isr8yAK2t+22KqwaoFrT5b5U66mYe+VgwDkWA@mail.gmail.com> <85r478bv99.fsf_-_@benfinney.id.au> <mailman.6735.1392194885.18130.python-list@python.org> <53153e66$0$24931$e4fe514c@dreader36.news.xs4all.nl> <mailman.7687.1393902132.18130.python-list@python.org> <59dd57ad-39b0-4c71-a58e-b4ae6517b385@googlegroups.com> <mailman.7689.1393904790.18130.python-list@python.org> <53156a42$0$2923$c3e8da3$76491128@news.astraweb.com> <mailman.7696.1393914922.18130.python-list@python.org> <87iortoic0.fsf@elektro.pacujo.net> <CAPTjJmpKASpOEscAjVDC5bpqK1kQt3rLTmCmUp-s4Ov2go6d2g@mail.gmail.com> <CAHVvXxR2TQzoiedi9Dns9RTMr_HVE_A5_tVWPGq_PmQ7zeV9Qw@mail.gmail.com> <CAPTjJmocK7Wpi765JKiNxEtquabm4No9GWomHv-vX7auFkqZ=w@mail.gmail.com> <CAHVvXxTSFnzZbRH2=xRr+yxidmTrUyCYsFPABOQBp3U5_2_-tQ@mail.gmail.com> <lf5mr8$ttr$1@ger.gmane.org> <CAHVvXxQo2-B8t+xmG04Awr1Z9V0o+w9az5UtpSEY=AmaW9pe0A@mail.gmail.com>
Precedence: list
Newsgroups: comp.lang.python
Message-ID: <mailman.7820.1394024033.18130.python-list@python.org>
Lines: 92
NNTP-Posting-Host: 2001:888:2000:d::a6
Xref: csiph.com comp.lang.python:67844

 Oscar Benjamin <oscar.j.benjamin@gmail.com> Wrote in message:
> On 4 March 2014 23:20, Dave Angel <davea@davea.name> wrote:
>>
>> One problem with complexity claims is that it's easy to miss some
>>  contributing time eaters. I haven't done any measuring on modern
>>  machines nor in python, but I'd assume that multiplies take
>>  *much* longer for large integers, and that divides are much
>>  worse. So counting iterations isn't the whole story.
> 
> Agreed but there's a big difference between log(N) iterations and N iterations!
> 
>> On the assumption that division by 2 is very fast,  and that a
>>  general multiply isn't too bad, you could improve on Newton by
>>  observing that the slope is 2.
>>
>>    err = n - guess * guess
>>     guess += err/2
> 
> I gues you mean like this:
> 
> def sqrt(n):
>     err = guess = 1
>     while err > 1e-10:
>         err = n - guess * guess
>         guess += err/2
>     return guess
> 
> This requires log2(10)*N iterations to get N digits. 

No idea how you came up with that,  but I see an error in my
 stated algorithm,  which does surely penalize it. The slope isn't
 2, but 2x.  So the line should have been
        guess += err/(2*guess)

Now if you stop the loop after 3 iterations (or use some other
 approach to get a low-precision estimate,  then you can calculate
 scale = 1/(2*estimate)

and then for remaining iterations, 
        guess += err *scale

> So the penalty
> for using division would have to be extreme in order for this to
> better. Using Decimal to get many digits we can write that as:
> 
> def sqrt2(n, prec=1000):
>     '''Solve x**2 = y'''
>     eps = D(10) ** -(prec + 5)
>     err = guess = D(1)
>     with localcontext() as ctx:
>         ctx.prec = prec + 10
>         while abs(err) > eps:
>             err = n - guess*guess
>             guess += err/2
>     return guess
> 
> This method works out much slower than Newton with division at 10000
> digits: 40s (based on a single trial) vs 80ms (timeit result).

Well considering you did not special-case the divide by 2, I'm not
 surprised it's slower.

> 
>> Some 37 years ago I microcoded  a math package which included
>>  square root.  All the math was decimal, and there was no hardware
>>  multiply or divide. The algorithm I came up with generated the
>>  answer one digit at a time, with no subsequent rounding needed.
>>  And it took just a little less time than two divides. For that
>>  architecture,  Newton's method would've been too
>>  slow.
> 
> If you're working with a fixed small precision then it might be.
> 
>> Incidentally, the algorithm did no divides, not even by 2. No
>>  multiplies either.  Just repeated subtraction,  sorta like divide
>>  was done.
>>
>> If anyone is curious,  I'll be glad to describe the algorithm;
>>  I've never seen it published, before or since.  I got my
>>  inspiration from a method used in mechanical,  non-motorized,
>>  adding machines.  My father had shown me that approach in the
>>  50's.
> 
> I'm curious.
> 

A later message,  I guess.  I can't write that much on the tablet.


-- 
DaveA