newton_sqrt.ml (1554B)
1 (* 15 2 3 NEW IDEA: recursion that stops when it is close enough, because floats cannot 4 be compared for equality. 5 6 newton_sqrt x approximates the square root of x >= 0 by Newton's method. Start 7 from a guess and repeatedly improve it with 8 9 next = (guess +. x /. guess) /. 2. 10 11 Stop when the guess is good enough, when abs_float (guess *. guess -. x) is 12 below a small tolerance such as 1e-10. Note that x itself works as a starting 13 guess, except for x = 0. 14 15 Do not write guess *. guess = x. Floats are approximations, that test would 16 loop forever, and it is the same reason 0.1 +. 0.2 = 0.3 is false in every 17 language with binary floats. 18 19 newton_sqrt 2.0 is about 1.4142135, newton_sqrt 81.0 = 9.0, 20 newton_sqrt 0.0 = 0.0. 21 22 Afterwards try newton_sqrt 1e12 and work out why an absolute tolerance of 23 1e-10 cannot terminate there. A double has about 16 significant digits, so near 24 1e12 the gap between neighbouring representable values is already larger than 25 1e-10. Fixing it means comparing relatively, for instance stopping when 26 successive guesses barely change. 27 28 *) 29 30 #use "topfind";; 31 #require "qcheck";; 32 33 let abs_diff x y = abs_float (x -. y);; 34 35 let newton_sqrt x = 36 let rec go guess = 37 let diff = abs_diff x (guess *. guess) in 38 if diff < 0.0000000001 39 then guess 40 else go ((guess +. x /. guess) /. 2.) 41 in go x;; 42 43 44 let reverse = QCheck.Test.make 45 (QCheck.float_range 0. 100000.) 46 (fun n -> (abs_diff n ((newton_sqrt n) *. (newton_sqrt n))) < 0.00001);; 47 48 QCheck_base_runner.run_tests [ reverse ];;