--
--
-- Arithmetic Operation
--
--

toMathExpr arg := mathNormalize (toMathExpr' arg)

(+') $x $y := b.+ x y

(-') $x $y := b.- x y

(*') $x $y := b.* x y

(/') $x $y := b./ x y

(+) $x $y :=
  match (isFloat x, isFloat y) as eq with
    | #(True, True)  -> f.+ x y
    | #(True, False) -> f.+ x (itof y)
    | #(False, True) -> f.+ (itof x) y
    | _              -> mathNormalize (x +' y)

(-) $x $y :=
  match (isFloat x, isFloat y) as eq with
    | #(True, True)  -> f.- x y
    | #(True, False) -> f.- x (itof y)
    | #(False, True) -> f.- (itof x) y
    | _              -> mathNormalize (x -' y)

(*) $x $y :=
  match (isFloat x, isFloat y) as eq with
    | #(True, True)  -> f.* x y
    | #(True, False) -> f.* x (itof y)
    | #(False, True) -> f.* (itof x) y
    | _              -> mathNormalize (x *' y)

(/) $x $y :=
  match (isFloat x, isFloat y) as eq with
    | #(True, True)  -> f./ x y
    | #(True, False) -> f./ x (itof y)
    | #(False, True) -> f./ (itof x) y
    | _              -> x /' y

reduceFraction := id

sum xs := foldl (+) 0 xs

sum' xs := foldl (+') 0 xs

product xs := foldl (*) 1 xs

product' xs := foldl (*') 1 xs

power $x $n := mathNormalize (power' x n)

power' $x $n := foldl (*') 1 (take n (repeat1 x))

(^) $x $n :=
  if x = e
    then exp n
    else if isRational n
      then if n >= 0
        then if isInteger n then power x n else `(^) x n
        else 1 / x ^ neg n
      else `(^) x n

(^') $x $n :=
  if x = e
    then exp n
    else if isRational n
      then if n >= 0
        then if isInteger n then power' x n else `(^) x n
        else 1 /' x ^' neg n
      else `(^) x n

gcd $x $y :=
  match (x, y) as (termExpr, termExpr) with
    | (_, #0) -> x
    | (#0, _) -> y
    | (term $a $xs, term $b $ys) ->
      gcd' (abs a) (abs b) *' foldl (*') 1 (map (^') (AC.intersect xs ys))

gcd' $x $y :=
  match (x, y) as (integer, integer) with
    | (_, #0) -> x
    | (#0, _) -> y
    | (_, ?(>= x)) -> gcd' (modulo y x) x
    | (_, _) -> gcd' y x

P./ fx $gx $x :=
  let xs := reverse (coefficients fx x)
      ys := reverse (coefficients gx x)
      (zs, rs) := L./ xs ys
   in ( sum' (map2 2#(%1 *' x ^' %2) (reverse zs) nats0)
      , sum' (map2 2#(%1 *' x ^' %2) (reverse rs) nats0) )