Sparse Monomial Exponents

This namespace is used by emmy.polynomial.impl and emmy.polynomial to implement a full polynomial data structure.

Polynomials are sums of monomial terms; a monomial is a pair of some non-zero coefficient and a product of some number of variables, each raised to some power.

We represent the exponents of a monomial with an ordered mapping of variable index => the exponent of that variable. Loading... is represented as {0 2, 2 3}, for example. Polynomials are linear combinations of the exponents.

Constructors

(def ^{:arglists '([& i-expt-pairs])}

  make

  "Accepts alternating pairs of integers representing <indeterminate index>,

  <exponent value> and returns a `sorted-map` representing the exponent portion

  of a polynomial term."

  #'sorted-map)

#'clojure.core/sorted-map

(def empty

  "Singleton instance of an empty exponent map."

  (make))

{}

(defn dense->exponents

  "Accepts a sequence of pairs of indeterminate index => power, and returns a

  sparse representation the exponents portion of a monomial.

  For example:

  ```clojure

  (dense->exponents [1 4 0 0 2])

  ;;=> {0 1, 1 4, 4 2}

  ```"

  [idx->pow]

  (reduce-kv (fn [acc i x]

               (if (zero? x)

acc

                 (core/assoc acc i x)))

             empty

             idx->pow))

#object[emmy.polynomial.exponent$dense__GT_exponents 0x756048b9 "

emmy.polynomial.exponent$dense__GT_exponents@756048b9"

]

(defn mul

  "Returns a new exponent vector generated by multiplying all terms in exponent

  vectors `l` and `r`."

  [l r]

  (merge-with + l r))

#object[emmy.polynomial.exponent$mul 0x32facc2d "

emmy.polynomial.exponent$mul@32facc2d"

]

(defn div

  "Returns a new exponent vector generated by dividing terms in exponent vector

  `l` by terms in `r`."

  [l r]

  (dense->exponents

   (merge-with + l (u/map-vals - r))))

#object[emmy.polynomial.exponent$div 0x2d950306 "

emmy.polynomial.exponent$div@2d950306"

]

(defn gcd

  "Returns the exponent vector that is the greatest common divisor of the exponent

  vectors `l` and `r`.

  Calling [[gcd]] with a single argument acts as identity."

  ([l] l)

  ([l r]

   (let [l' (select-keys l (keys r))

         r' (select-keys r (keys l))]

     (merge-with min l' r'))))

#object[emmy.polynomial.exponent$gcd 0x6cc66ea8 "

emmy.polynomial.exponent$gcd@6cc66ea8"

]

(defn lcm

  "Returns the exponent vector that is the least common multiple of the exponent

  vectors `l` and `r`.

  Calling [[lcm]] with no arguments returns the empty exponent vector; calling

  with a single argument acts as identity."

  ([] empty)

  ([l] l)

  ([l r]

   (merge-with max l r)))

#object[emmy.polynomial.exponent$lcm 0x1b14c609 "

emmy.polynomial.exponent$lcm@1b14c609"

]

(defn every-power?

  "Returns true if `f` returns true for every positive non-zero exponent in the

  supplied exponent vector `m`, false otherwise.

  Defaults to `true` if `m` is empty."

  [f m]

  (every? f (vals m)))

#object[emmy.polynomial.exponent$every_power_QMARK_ 0x7f6beae2 "

emmy.polynomial.exponent$every_power_QMARK_@7f6beae2"

]

(defn assoc

  "Replaces the entry for variable `x` in the exponents vector `m` with power `n`.

  If `n` is zero, removes `x`'s entry from the returned exponent vector."

  [m x n]

  (if (zero? n)

    (dissoc m x)

    (core/assoc m x n)))

#object[emmy.polynomial.exponent$assoc 0x7ae35b3a "

emmy.polynomial.exponent$assoc@7ae35b3a"

]

(defn lower

  "Given some exponent vector `expts`, and an optional variable index

  `i` (defaults to `0`), returns a new exponent vector with:

  - the `i`th variable's entry removed

  - decremented index for each variable with index > `i`.

  For example:

  ```clojure

  (lower {0 3, 1 2, 4 4})

  ;;=> {0 2, 3 4}

  (lower {0 3, 1 2, 4 4} 1)

  ;;=> {0 3, 3 4}

  ```"

  ([expts]

   (lower expts 0))

  ([expts i]

   (reduce-kv (fn [acc k v]

                (if (> k i)

                  (core/assoc acc (dec k) v)

                  (core/assoc acc k v)))

              empty

              (dissoc expts i))))

#object[emmy.polynomial.exponent$lower 0x250fc9e0 "

emmy.polynomial.exponent$lower@250fc9e0"

]

(defn raise

  "Given some exponent vector `expts`, an optional variable index `i` (defaults to

  `0`) and an optional exponent power `n` (defaults to `0`), returns a new

  exponent vector with:

  - incremented indices for each variable with index >= `i`

  - a new `i`th variable created with power `n`

  For example:

  ```clojure

  (raise {0 3, 1 2, 4 4})

  ;;=> {1 3, 2 2, 5 4}

  (raise {0 3, 1 2, 4 4} 1)

  ;;=> {0 3, 2 2, 5 4}

  (raise {0 3, 1 2, 4 4} 1)

  ;;=> {0 3, 1 10, 2 2, 5 4}

  ```"

  ([expts]

   (raise expts 0 0))

  ([expts i]

   (raise expts i 0))

  ([expts i n]

   (let [m (reduce-kv (fn [acc k v]

                        (if (>= k i)

                          (core/assoc acc (inc k) v)

                          (core/assoc acc k v)))

                      empty

                      expts)]

     (if (zero? n)

m

       (core/assoc m i n)))))

#object[emmy.polynomial.exponent$raise 0x181667c6 "

emmy.polynomial.exponent$raise@181667c6"

]

(defn monomial-degree

  "Returns the [monomial degree](https://en.wikipedia.org/wiki/Monomial#Degree) of

  the exponent vector `m`, i.e., the sum of the powers of all variables in `m`.

  If the optional `i` is supplied, returns the degree of the `i`th variable, i.e.,

  the entry for `i` in `m`, defaulting to `0`."

  ([m]

   (apply + (vals m)))

  ([m i]

   (m i 0)))

#object[emmy.polynomial.exponent$monomial_degree 0x380aa2dc "

emmy.polynomial.exponent$monomial_degree@380aa2dc"

]

(defn ->sort+unsort

  "Given a power product `m`, returns a pair of `sort` and `unsort` functions of a

  single power product argument.

  `sort` rearranges the indices of its argument to match the order of increasing

  variable degree in `m`. `unsort` undoes this transformation."

[m]

  (let [indices (range (count m))

        order (into [] (sort-by m (keys m)))]

    (letfn [(sort [m']

              (into empty

                    (mapcat (fn [i]

                              (when-let [v (m' (order i))]

                                [[i v]])))

                    indices))

            (unsort [m']

              (into empty

                    (mapcat (fn [i]

                              (when-let [v (m' i)]

                                [[(order i) v]])))

                    indices))]

      [sort unsort])))

#object[emmy.polynomial.exponent$__GT_sort_PLUS_unsort 0x5850146c "

emmy.polynomial.exponent$__GT_sort_PLUS_unsort@5850146c"

]

Monomial Orderings

This section implements a number of Monomial orderings that are useful for making various polynomial algorithms efficient.

Each of the following functions matches Java's "comparator" interface:

If the left argument is greater than the right argument, the function returns a negative number
equal values return 0
if the right value is greater, returns a positive number

(defn lex-order

  "Comparator that responds based on the [lexicographic

  order](https://en.wikipedia.org/wiki/Monomial_order#Lexicographic_order), or

  'lex order', of the exponent vectors `xs` and `ys`. Accepts any sequence of

  pairs of the form `[variable, power]` for `xs` and `ys`.

  Lex order first compares the power of variable `0`, then, in case of equality,

  variable `1` and so on. If all powers match, returns `0`.

  If `:reverse? true` is passed, the inputs are compared in reverse lexicographic order.

  Reverse here means two things:

  1. the variables are considered in reverse order; the variable with the

     largest index is compared first, then the next-largest, and on down the line.

  2. The _smaller_ exponent is grevlex-greater in this comparison."

  ([xs ys & {:keys [reverse?]}]

   (let [xs (vec (if reverse? (rseq xs) xs))

         ys (vec (if reverse? (rseq ys) ys))]

     (loop [i (long 0)]

       (let [x (nth xs i nil)

             y (nth ys i nil)]

         (cond (and (not x) (not y)) 0

               (not x) -1

               (not y)  1

               :else (let [bit (compare (nth x 0) (nth y 0))]

                       (cond (zero? bit)

                             (let [xv (nth x 1)

                                   yv (nth y 1)]

                               (if (= xv yv)

                                 (recur (inc i))

                                 (if reverse?

                                   (- yv xv)

                                   (- xv yv))))

                             (neg? bit) 1

                             :else -1))))))))

#object[emmy.polynomial.exponent$lex_order 0x34b964f "

emmy.polynomial.exponent$lex_order@34b964f"

]

(defn graded-lex-order

  "Comparator that responds based on the [graded lexicographic

  order](https://en.wikipedia.org/wiki/Monomial_order#Graded_lexicographic_order),

  or 'grlex order', of the exponent vectors `xs` and `ys`. Accepts any sequence

  of pairs of the form `[variable, power]` for `xs` and `ys`.

  grlex order first compares the total degree of `xs` and `ys`, and falls back

  to [[lex-order]] in case of a tie. See [[lex-order]] for details on this

  case."

  [xs ys]

  (let [xd (monomial-degree xs)

        yd (monomial-degree ys)]

    (if (= xd yd)

      (lex-order xs ys)

      (- xd yd))))

#object[emmy.polynomial.exponent$graded_lex_order 0x28c8e9b7 "

emmy.polynomial.exponent$graded_lex_order@28c8e9b7"

]

(defn graded-reverse-lex-order

  "Comparator that responds based on the [graded reverse lexicographic

  order](https://en.wikipedia.org/wiki/Monomial_order#Graded_reverse_lexicographic_order),

  or 'grevlex order', of the exponent vectors `xs` and `ys`. Accepts any

  sequence of pairs of the form `[variable, power]` for `xs` and `ys`.

  grevlex order first compares the total degree of `xs` and `ys`, and falls back

  to reverse [[lex-order]] in case of a tie. Reverse here means two things:

  1. the variables are considered in reverse order; the variable with the

     largest index is compared first, then the next-largest, and on down the line.

  2. The _smaller_ exponent is grevlex-greater in this comparison."

  [xs ys]

  (let [xd (monomial-degree xs)

        yd (monomial-degree ys)]

    (if (= xd yd)

      (lex-order xs ys :reverse? true)

      (- xd yd))))

#object[emmy.polynomial.exponent$graded_reverse_lex_order 0x45b9a230 "

emmy.polynomial.exponent$graded_reverse_lex_order@45b9a230"

]