(* FILE: numlib.joy *) LIBRA _numlib == true; (* predicates *) positive == 0 >; negative == 0 <; even == 2 rem null; odd == even not; prime == 2 [ [dup * >] nullary [rem 0 >] dip and ] [ succ ] while dup * < ; (* functions *) fact == [1 1] dip [dup [*] dip succ] times pop; fib == [1 0] dip [swap [+] unary] times popd; (* nfib == [1 1] dip [dup [+ succ] dip swap] times pop; *) gcd == [0 >] [dup rollup rem] while pop; fahrenheit == 9 * 5 / 32 + ; celsius == 32 - 5 * 9 / ; pi == 3.14159265; e == 1.0 exp; radians == pi * 180 / ; sindeg == radians sin; cosdeg == radians cos; tandeg == radians tan; (* find roots of the quadratic equation with coefficients a b c : a * X^2 + b * X + c = 0 *) qroots == (* a b c *) [ pop pop null ] (* a = 0 ? *) (* degenerate cases: *) [ [ pop null ] (* b = 0 ? *) [ [ null ] (* c = 0 ? *) [ [_INF] ] (* => [_INF] *) [ [] ] (* => [] *) ifte [ pop pop pop ] dip ] [ 0 swap - swap 1.0 * / (* float divisor *) [] cons popd ] (* => [ -c/b ] *) ifte ] (* standard cases: *) [ [ [ dup * swap ] dip 4 * * - ] unary (* b^2 - 4ac *) [ 0 < ] (* b^2 - 4ac negative ? *) [ pop pop pop [_COMPLEX] ] (* => [_COMPLEX] *) [ [ 0 swap - 1.0 * (* -b (floated) *) swap 2 * ] dip (* 2a *) [ 0 = ] (* b^2 - 4ac zero ? *) [ pop / [] cons ] (* => [-b / 2a] *) [ sqrt swapd dup2 - [+] dip (* -b+s -b-s *) [] cons cons (* [ -b+s -b-s ]*) swap [/] cons map ] (* => [(-b+s)/2a (-b-s)/2a]*) ifte ] ifte ] ifte; (* a simpler version, FEB 05 : *) quadratic-formula == # a b c => [root1 root2 ] [ [ [ pop pop 2 * ] # divisor [ pop 0 swap - ] # minusb [ swap dup * rollup * 4 * - sqrt ] ] # radical [i] map ] ternary i [ [ [ + swap / ] # root1 [ - swap / ] ] # root2 [i] map ] ternary; (* combinators *) deriv == [unary2 swap - 0.001 /] cons [dup 0.001 +] swoncat; (* Newton's method for finding roots of equations in one variable, for example to find the temperature at which Fahrenheit = Celsius: "using 5 as a guess, find the X for which X = fahrenheit(X) or X - fahrenheit(X) = 0 in Joy: 5 [dup fahrenheit -] newton => -40 or -100 [dup celsius -] newton => -40 So -40 is the fixpoint for the two conversion functions. *) newton == (* Usage: guess [F] newton *) dup deriv (* guess [F] [D] *) [ pop i abs 0.0001 > ] (* too big ? *) [ [dupd] dip (* guess guess [F] [D] *) dup2 (* guess guess [F] [D] [F] [D] *) [[cleave / - ] dip] dip ] (* newguess [F] [D] *) while pop pop; use-newton == [[-] cons] dip swoncat 1 swap newton; cube-root == [dup dup * *] use-newton; NUMLIB == "numlib.joy - numerical library\n". (* end LIBRA *) "numlib is loaded\n" putchars. (* END numlib.joy *)