\ 4tH library - FSQRT ZEN - Copyright 2009 J.L. Bezemer
\ You can redistribute this file and/or modify it under
\ the terms of the GNU General Public License
\ include lib/zenfloat.4th
[UNDEFINED] fsqrt [IF]
[UNDEFINED] f+ [IF] [ABORT] [THEN]
\ same number but crank up mantissa
: (scale) ( f1 -- f1 )
max-n 10 / 1+ >r begin over r@ < while 1- swap 10 * swap repeat r> drop
;
\ estimate courtesy of David Johnson
: estimate ( f1 -- f1 f2)
2dup (scale) dup 1 and \ scale for sqrt; crank up mantissa
if 1+ swap 10 / swap then \ if not even, scale down one step
2/ swap max-n over = if 1- then dup \ adjust MAX-N, bugfix
begin over over / over + 2/ tuck - abs 2 < until
nip swap \ take sqrt of mantissa and divide
; \ exponent by two
: fsqrt ( f1 -- f2 )
2dup f0= if exit then \ check for zero
2dup f0< abort" Negative float" \ check for negative float
estimate begin \ get initial estimation
2over 2over f/ 2over f+ (scale) swap 2/ swap \ get next estimate
2swap 2over rot - 0= >r - abs 6 < r> and \ assume exponents are equal
until \ difference half last digit mantissa
2swap 2drop \ babylonian method w/ estimate
; \ clean up after final iteration
[DEFINED] 4TH [IF]
hide estimate
hide (scale)
[THEN]
[THEN]
\ 0 s>f fsqrt f. cr
\ 5 -1 fsqrt f. cr
\ 5 s>f fsqrt f. cr
\ 1 s>f fsqrt f. cr
\ 2 s>f fsqrt f. cr
\ 64 s>f fsqrt f. cr
\ 4 9 fsqrt f. cr
\ 4 -2 fsqrt f. cr
\ 41 -8 fsqrt f. cr
\ 1 10 fsqrt f. cr
\ 1 99 fsqrt f. cr
\ 1 -99 fsqrt f. cr
\ 9224 -4 fsqrt f. cr
\ 120 -4 fsqrt f. cr
\ 9643 -3 fsqrt f. cr
\ max-n 0 fsqrt f. cr
\ max-n 1 fsqrt f. cr
\ -1 s>f fsqrt f. cr
\ 10001 1 do i -4 fsqrt 2drop loop
|
