Re: [Edu-sig] Brute force solutions
The fun part here is we can use numerator/denominator syntax with open-ended precision integers, to like express sqrt of 19 as some humongous fraction (as many digits as memory will allow). This lets us far surpass the floating point barrier.
For example, the sqrt of 19 is rougly: 745969084203762918506900171 6630770611758990188906918072 4013315460866431392624885605 5112125808098871177881508095 4010864598095 ---------------------------- 1711370448973571545640915466 3001505416992487326376156603 4301985589644920244770721090 4993017822174818974832607428 966608330904 (that's a numerator over a denominator). Here's the code behind it (still need to replace int(pow(x,0.5)) with a numerical method that doesn't do all the work to find an actual floating point sqrt -- overkill): """ gonzo.py """ from mathteach import cf2 # http://www.4dsolutions.net/ocn/python/mathteach.py def sqrtcf(n): orig = n denom = 1 addterm = 0 cflist = [] while True: m = int(pow(n,0.5)) # <-- replace with function call cflist.append(m) if len(cflist)>1 and denom == 1: break addterm = -(addterm - m*denom) denom = (orig - addterm**2)/denom n = ((pow(orig,0.5) + addterm)/denom)**2 return cflist def getsqrt(n,default=30): thelist = sqrtcf(n) while len(thelist)<default: thelist += thelist[1:] thelist = thelist[:default] return cf2(thelist * 10)
from gonzo import getsqrt getsqrt(19) ...
Kirby _______________________ Edu-sig mailing list Edu-sig@python.org http://mail.python.org/mailman/listinfo/edu-sig -------------------------------------------------------------------- mail2web - Check your email from the web at http://mail2web.com/ .
urnerk@qwest.net wrote:
The fun part here is we can use numerator/denominator syntax with
open-ended
precision integers, to like express sqrt of 19 as some humongous fraction (as many digits as memory will allow). This lets us far surpass the floating point barrier.
OK, here we go: def test_sqrt(numerator, denominator, trial): '''True iff trial (a num,den pair) over-estimates the sqrt(n/d)''' root_n, root_d = trial # return numerator / denominator >= root_n ** 2 / root_d ** 2 return root_d ** 2 * numerator >= root_n ** 2 * denominator # since we don't have 2.5 yet, here's a version of partial: def partial(function, *args): '''Simple no-kwargs version of partial''' def andthen(*final_args): return function(*(args + final_args)) return andthen def farey_trials(tester): '''Binary search for fraction. value must be between 0 and +Inf tester((num, den)) returns True if fract is high, False if Low ''' low_n, low_d = low = 0, 1 # 0/1 = 0 .. high_n, high_d = high = 1, 0 # 1/0 = Infinity while True: mediant_n = low_n + high_n mediant_d = low_d + high_d mediant = mediant_n, mediant_d yield mediant if tester(mediant): low_n, low_d = low = mediant else: high_n, high_d = high = mediant # Now here is a cute reporter that relies on another trick: # n & -n == n only if n is either 0 or a power of 2. for n, fraction in enumerate(farey_trials(partial(test_sqrt, 19, 1))): if n & -n == n: # report increasingly rarely print n, fraction if n > 4096: break -- Scott David Daniels Scott.Daniels@Acm.Org
participants (2)
-
Scott David Daniels -
urnerk＠qwest.net