arithmetic.egi

--
--
-- Arithmetic Operation
--
--
declare symbol i, w, e, π: MathValue

def (+') : MathValue -> MathValue -> MathValue := i.+
def (-') : MathValue -> MathValue -> MathValue := i.-
def (*') : MathValue -> MathValue -> MathValue := i.*
def (/') : MathValue -> MathValue -> MathValue := i./

def plusForMathValue (x: MathValue) (y: MathValue) : MathValue :=
  mathNormalize (x +' y)

def minusForMathValue (x: MathValue) (y: MathValue) : MathValue :=
  mathNormalize (x -' y)

def multForMathValue (x: MathValue) (y: MathValue) : MathValue :=
  mathNormalize (x *' y)

def divForMathValue (x: MathValue) (y: MathValue) : MathValue :=
  x /' y

def sum {AddMonoid a} (xs: [a]) : a := foldl (+) zero xs
def sum' (xs: [MathValue]) : MathValue := foldl (+') 0 xs

def product {MulMonoid a} (xs: [a]) : a := foldl (*) one xs
def product' (xs: [MathValue]) : MathValue := foldl (*') 1 xs

def power (x: MathValue) (n: MathValue) : MathValue := mathNormalize (power' x n)
-- power' must avoid any operator that dispatches through mathNormalize so it
-- can be safely called from within `mathNormalize` itself (e.g. inside a
-- declare-rule RHS). foldl/take used `n - 1` whose `-` resolves to
-- `minusForMathValue` -> `mathNormalize` and creates a cycle. Use direct
-- recursion with `i.-` (integer subtraction primitive) and `*'` (the
-- un-normalised multiplication).
def power' (x: MathValue) (n: MathValue) : MathValue :=
  if n = 0
    then 1
    else x *' power' x (i.- n 1)

def exp (x: MathValue) : MathValue := 'exp x

def (^) (x: MathValue) (n: MathValue) : MathValue :=
  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 ^ i.neg n
      else '(^) x n

def (^') (x: MathValue) (n: MathValue) : MathValue :=
  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 ^' i.neg n
      else '(^) x n

def gcdForMathValue (x: MathValue) (y: MathValue) : MathValue :=
  match (x, y) as (termExpr, termExpr) with
    | (_, #0) -> x
    | (#0, _) -> y
    | (term $a $xs, term $b $ys) ->
      -- After Term widening (2026-05-06), `$a` and `$b` are statically MathValue.
      -- Only compute integer gcd when both are integers; otherwise fall back to 1.
      (if isInteger a && isInteger b then gcd' (i.abs a) (i.abs b) else 1)
        *' foldl (*') 1 (map (\(s, n) -> s ^' n) (AC.intersect xs ys))
    | _ -> 1  -- fallback when neither is in term form (e.g., level-3/4 raw poly after tower fix)

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

def P./ fx gx x :=
  let xs := reverse (coefficients fx x)
      ys := reverse (coefficients gx x)
      (zs, rs) := L./ xs ys
   in ( sum' (map2 (\c n -> c *' x ^' n) (reverse zs) nats0)
      , sum' (map2 (\c n -> c *' x ^' n) (reverse rs) nats0) )