From mboxrd@z Thu Jan 1 00:00:00 1970 X-Spam-Checker-Version: SpamAssassin 3.4.4 (2020-01-24) on polar.synack.me X-Spam-Level: X-Spam-Status: No, score=-1.3 required=5.0 tests=BAYES_00,INVALID_MSGID autolearn=no autolearn_force=no version=3.4.4 X-Google-Language: ENGLISH,ASCII-7-bit X-Google-Thread: 103376,fedc2d05e82c9174 X-Google-Attributes: gid103376,public From: "Howard W. LUDWIG" Subject: Re: Calculating SQRT in ADA Date: 1999/03/26 Message-ID: <36FBAFF6.6C6648CC@lmco.com>#1/1 X-Deja-AN: 459391960 Content-Transfer-Encoding: 7bit References: <7dbv6t$4u5$1@nnrp1.dejanews.com> <19990324201959.00800.00000708@ngol04.aol.com> <7dei9a$dvo$1@nnrp1.dejanews.com> To: "David C. Hoos, Sr." Content-Type: text/plain; charset=us-ascii Organization: Lockheed Martin -- Information Systems Center Mime-Version: 1.0 Reply-To: howard.w.ludwig@lmco.com Newsgroups: comp.lang.ada Date: 1999-03-26T00:00:00+00:00 List-Id: David C. Hoos, Sr. wrote: > Mark A Biggar wrote in message <36FAA3DF.42C31CF@lmco.com>... > >Hans Marqvardsen wrote: > >> > >> John Herro wrote: > > > >No,no,no,no! > > > >This algorithm has quadratic convergence. This means that each time > through the > >loop doubles the number of correct bits in the answer. So with IEEE double > >precession float you need loop at most 6 times, and even with extended > float > >7 times, so don't test against a torerance at all, just loop 6 or 7 times > >and be done with it. > > > No,no,no,no! (to coin a phrase) > > The number of correct bits doubles with every iteration, only if you are > counting place holding bits as well as mantissa bits. So, > if you take the square root of 2.0e100, for example, 173 iterations are > required if the argument is used as the initial estimate. > > similar number is required for 0.5e-100, but only six iterations are > required for 2.0e00 or 0.5e00. No, no, no, no! (to plagiarize a coined phrase) It is a common oversimplified assumption that Newton[-Raphson] has quadratic convergence. There are many well-known cases where the method does not converge at all. Even when the method does converge, the convergence can be linear. Quadratic convergence applies when the quadratic term is negligible compared to the linear term in the Taylor series expansion of the function in question. When the initial guess is far off, "negligible" does not describe the situation. To find the actual a = sqrt(c), where c is known, involves calculating x(n+1) = 0.5 * (x(n) + c/x(n)). Let x(n) = a + d(n), where d(n) measures the deviation of the estimate x(n) from the actual value a. Then d(n+1) = x(n+1) - a = 0.5 * (a + d(n) + a**2/(a + d(n))) - a = 0.5 * d(n)**2/(a + d(n)). When d(n) is small compared to a, the denominator behaves more and more like a, so the deviation is indeed proportional to d(n)**2. However, when d(n) is large compared to a, the denominator for some iterations looks more like d(n), canceling out one of the d(n) factors in the numerator, so the deviation is cut each time by a factor of 2, thus linear convergence, until the iterations get to the place of having the deviation be of the order of a, at which time the convergence transitions to quadratic. In the case of the example given [c = 2E100, a =1.4142...E50], starting with an initial guess of c requires 166 iterations of halving the deviation each time to get to the point of the deviation being smaller than the actual root a. Then quadratic convergence kicks in, apparently requiring only another 7 iterations to achieve the ultimate desired result. If we simply counted place-holder bits, roughly 167, plus the desired 53 bits of precision, then we have about 220 bits to deal with. _Assuming_ the first bit to be fine and doubling the number of good bits each iteration, only 8 iterations would be needed (2^8 = 256), nowhere close to the 173 iterations quoted above. This shows that you cannot assume that the first bit is good. If you want to do a "back of the envelope" calculation based on quadratic convergence, the 2E100 case using 2E100 as the initial guess would say your first guess has less than 1E-50 bit good, in order to double the number of good bits each iteration and to obtain 53 good bits after 173 iterations. A better way of looking at it is that the place-holder bits are lopped off at the rate of one per iteration for the first 166 iterations. HWL