Path: csiph.com!v102.xanadu-bbs.net!xanadu-bbs.net!feeder.erje.net!eu.feeder.erje.net!feeds.phibee-telecom.net!newsfeed.xs4all.nl!newsfeed4.news.xs4all.nl!xs4all!post.news.xs4all.nl!not-for-mail
MIME-Version: 1.0
In-Reply-To: <lf5mr8$ttr$1@ger.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>
From: Oscar Benjamin <oscar.j.benjamin@gmail.com>
Date: Wed, 5 Mar 2014 11:59:47 +0000
Subject: Re: Working with the set of real numbers
To: Dave Angel <davea@davea.name>
Content-Type: text/plain; charset=ISO-8859-1
Cc: Python List <python-list@python.org>
Precedence: list
Newsgroups: comp.lang.python
Message-ID: <mailman.7810.1394020816.18130.python-list@python.org>
Lines: 68
NNTP-Posting-Host: 2001:888:2000:d::a6
Xref: csiph.com comp.lang.python:67833

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. 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).

> 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.


Oscar