--Copyright (c) 1991-2002, The Numerical ALgorithms Group Ltd.
--All rights reserved.
--
--Redistribution and use in source and binary forms, with or without
--modification, are permitted provided that the following conditions are
--met:
--
--    - Redistributions of source code must retain the above copyright
--      notice, this list of conditions and the following disclaimer.
--
--    - Redistributions in binary form must reproduce the above copyright
--      notice, this list of conditions and the following disclaimer in
--      the documentation and/or other materials provided with the
--      distribution.
--
--    - Neither the name of The Numerical ALgorithms Group Ltd. nor the
--      names of its contributors may be used to endorse or promote products
--      derived from this software without specific prior written permission.
--
--THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
--IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
--TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
--PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
--OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
--EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
--PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
--PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
--LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
--NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
--SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

)abbrev domain MRING MonoidRing
++ Authors: Stephan M. Watt; revised by Johannes Grabmeier
++ Date Created: January 1986
++ Date Last Updated: 14 December 1995, Mike Dewar
++ Basic Operations: *, +, monomials, coefficients
++ Related Constructors: Polynomial
++ Also See:
++ AMS Classifications:
++ Keywords: monoid ring, group ring, polynomials in non-commuting
++  indeterminates
++ References:
++ Description:
++  \spadtype{MonoidRing}(R,M), implements the algebra
++  of all maps from the monoid M to the commutative ring R with
++  finite support.
++  Multiplication of two maps f and g is defined
++  to map an element c of M to the (convolution) sum over {\em f(a)g(b)}
++  such that {\em ab = c}. Thus M can be identified with a canonical
++  basis and the maps can also be considered as formal linear combinations
++  of the elements in M. Scalar multiples of a basis element are called
++  monomials. A prominent example is the class of polynomials
++  where the monoid is a direct product of the natural numbers
++  with pointwise addition. When M is
++  \spadtype{FreeMonoid Symbol}, one gets polynomials
++  in infinitely many non-commuting variables. Another application
++  area is representation theory of finite groups G, where modules
++  over \spadtype{MonoidRing}(R,G) are studied.

MonoidRing(R: Ring, M: Monoid): MRcategory == MRdefinition where
    Term ==> Record(coef: R, monom: M)

    MRcategory ==> Join(Ring, RetractableTo M, RetractableTo R) with
        monomial         : (R, M) -> %
          ++ monomial(r,m) creates a scalar multiple of the basis element m.
        coefficient : (%, M) -> R
          ++ coefficient(f,m) extracts the coefficient of m in f with respect
          ++ to the canonical basis M.
        coerce:   List Term -> %
          ++ coerce(lt) converts a list of terms and coefficients to a member of the domain.
        terms       : % -> List Term
          ++ terms(f) gives the list of non-zero coefficients combined
          ++ with their corresponding basis element as records.
          ++ This is the internal representation.
        map         : (R -> R, %) -> %
          ++ map(fn,u) maps function fn onto the coefficients
          ++ of the non-zero monomials of u.
        monomial?   : % -> Boolean
          ++ monomial?(f) tests if f is a single monomial.
        coefficients: % -> List R
          ++ coefficients(f) lists all non-zero coefficients.
        monomials: % -> List %
           ++ monomials(f) gives the list of all monomials whose
           ++ sum is f.
        numberOfMonomials: % -> NonNegativeInteger
           ++ numberOfMonomials(f) is the number of non-zero coefficients
           ++ with respect to the canonical basis.
        if R has CharacteristicZero then CharacteristicZero
        if R has CharacteristicNonZero then CharacteristicNonZero
        if R has CommutativeRing then Algebra(R)
        if (R has Finite and M has Finite) then Finite
        if M has OrderedSet then
          leadingMonomial   : % -> M
            ++ leadingMonomial(f) gives the monomial of f whose
            ++ corresponding monoid element is the greatest
            ++ among all those with non-zero coefficients.
          leadingCoefficient: % -> R
            ++ leadingCoefficient(f) gives the coefficient of f, whose
            ++ corresponding monoid element is the greatest
            ++ among all those with non-zero coefficients.
          reductum          : % -> %
            ++ reductum(f) is f minus its leading monomial.

    MRdefinition ==> add
        Ex ==> OutputForm
        Cf ==> coef
        Mn ==> monom

        Rep  := List Term

        coerce(x: List Term): % == x :: %

        monomial(r:R, m:M)  ==
          r = 0 => empty()
          [[r, m]]

        if (R has Finite and M has Finite) then
          size() == size()$R ** size()$M

          index k ==
            -- use p-adic decomposition of k
            -- coefficient of p**j determines coefficient of index(i+p)$M
            i:Integer := k rem size()
            p:Integer := size()$R
            n:Integer := size()$M
            ans:% := 0
            for j in 0.. while positive? i repeat
              h := i rem p
              -- we use index(p) = 0$R
              if h ~= 0 then
                c : R := index(h :: PositiveInteger)$R
                m : M := index((j+n) :: PositiveInteger)$M
                --ans := ans + c *$% m
                ans := ans + monomial(c, m)$%
              i := i quo p
            ans

          lookup(z : %) : PositiveInteger ==
            -- could be improved, if M has OrderedSet
            -- z = index lookup z, n = lookup index n
            -- use p-adic decomposition of k
            -- coefficient of p**j determines coefficient of index(i+p)$M
            zero?(z) => size()$% pretend PositiveInteger
            liTe : List Term := terms z  -- all non-zero coefficients
            p  : Integer := size()$R
            n  : Integer := size()$M
            res : Integer := 0
            for te in liTe repeat
              -- assume that lookup(p)$R = 0
              l:NonNegativeInteger:=lookup(te.Mn)$M
              ex : NonNegativeInteger := (n=l => 0;l)
              co : Integer := lookup(te.Cf)$R
              res := res + co * p ** ex
            res pretend PositiveInteger

          random() == index( (1+(random()$Integer rem size()$%) )_
            pretend PositiveInteger)$%

        0                   == empty()
        1                   == [[1, 1]]
        terms a             == (copy a) pretend List(Term)
        monomials a         == [[t] for t in a]
        coefficients a      == [t.Cf for t in a]
        coerce(m:M):%       == [[1, m]]
        coerce(r:R): % ==
        -- coerce of ring
          r = 0 => 0
          [[r,    1]]
        coerce(n:Integer): % ==
        -- coerce of integers
          n = 0 => 0
          [[n::R, 1]]
        - a                 == [[ -t.Cf, t.Mn] for t in a]
        if R has noZeroDivisors
           then
            (r:R) * (a:%) ==
              r = 0 => 0
              [[r*t.Cf, t.Mn] for t in a]
           else
            (r:R) * (a:%) ==
              r = 0 => 0
              [[rt, t.Mn] for t in a | (rt:=r*t.Cf) ~= 0]
        if R has noZeroDivisors
           then
            (n:Integer) * (a:%) ==
              n = 0 => 0
              [[n*t.Cf, t.Mn] for t in a]
           else
            (n:Integer) * (a:%) ==
              n = 0 => 0
              [[nt, t.Mn] for t in a | (nt:=n*t.Cf) ~= 0]
        map(f, a)           == [[ft, t.Mn] for t in a | (ft:=f(t.Cf)) ~= 0]
        numberOfMonomials a == #a

        retractIfCan(a:%):Union(M, "failed") ==
          one?(#a) and one?(a.first.Cf) => a.first.Mn
          "failed"

        retractIfCan(a:%):Union(R, "failed") ==
          one?(#a) and one?(a.first.Mn) => a.first.Cf
          "failed"

        if R has noZeroDivisors then
          if M has Group then
            recip a ==
              lt := terms a
              not one?(#lt) => "failed"
              (u := recip lt.first.Cf) case "failed" => "failed"
              --(u::R) * inv lt.first.Mn
              monomial((u::R), inv lt.first.Mn)$%
          else
            recip a ==
              not one?(#a) or not one?(a.first.Mn) => "failed"
              (u := recip a.first.Cf) case "failed" => "failed"
              u::R::%

        mkTerm(r:R, m:M):Ex ==
            r=1 => m::Ex
            r=0 or m=1 => r::Ex
            r::Ex * m::Ex

        coerce(a:%):Ex ==
            empty? a => (0$Integer)::Ex
            empty? rest a => mkTerm(a.first.Cf, a.first.Mn)
            reduce(_+, [mkTerm(t.Cf, t.Mn) for t in a])$List(Ex)

        if M has OrderedSet then -- we mean totally ordered
            -- Terms are stored in decending order.
            leadingCoefficient a == (empty? a => 0; a.first.Cf)
            leadingMonomial a    == (empty? a => 1; a.first.Mn)
            reductum a           == (empty? a => a; rest a)

            a = b ==
                #a ~= #b => false
                for ta in a for tb in b repeat
                    ta.Cf ~= tb.Cf or ta.Mn ~= tb.Mn => return false
                true

            a + b ==
                c:% := empty()
                while not empty? a and not empty? b repeat
                  ta := first a; tb := first b
                  ra := rest a;  rb := rest b
                  c :=
                    ta.Mn > tb.Mn => (a := ra; concat!(c, ta))
                    ta.Mn < tb.Mn => (b := rb; concat!(c, tb))
                    a := ra; b := rb
                    not zero?(r := ta.Cf+tb.Cf) =>
                                        concat!(c, [r, ta.Mn])
                    c
                concat!(c, concat(a, b))

            coefficient(a, m) ==
                for t in a repeat
                    if t.Mn = m then return t.Cf
                    if t.Mn < m then return 0
                0


            if M has OrderedMonoid then

            -- we use that multiplying an ordered list of monoid elements
            -- by a single element respects the ordering

              if R has noZeroDivisors then
                a:% * b:% ==
                  +/[[[ta.Cf*tb.Cf, ta.Mn*tb.Mn]$Term
                    for tb in b ] for ta in reverse a]
              else
                a:% * b:% ==
                  +/[[[r, ta.Mn*tb.Mn]$Term
                    for tb in b | not zero?(r := ta.Cf*tb.Cf)]
                      for ta in reverse a]
            else -- M hasn't OrderedMonoid

            -- we cannot assume that mutiplying an ordered list of
            -- monoid elements by a single element respects the ordering:
            -- we have to order and to collect equal terms
              ge : (Term,Term) -> Boolean
              ge(s,t) == t.Mn <= s.Mn

              sortAndAdd : List Term -> List Term
              sortAndAdd(liTe) ==  -- assume liTe not empty
                liTe := sort(ge,liTe)
                m : M :=  (first liTe).Mn
                cf : R := (first liTe).Cf
                res : List Term := []
                for te in rest liTe repeat
                  if m = te.Mn then
                    cf := cf + te.Cf
                  else
                    if not zero? cf then res := cons([cf,m]$Term, res)
                    m := te.Mn
                    cf := te.Cf
                if not zero? cf then res := cons([cf,m]$Term, res)
                reverse res


              if R has noZeroDivisors then
                a:% * b:% ==
                  zero? a => a
                  zero? b => b  -- avoid calling sortAndAdd with []
                  +/[sortAndAdd [[ta.Cf*tb.Cf, ta.Mn*tb.Mn]$Term
                    for tb in b ] for ta in reverse a]
              else
                a:% * b:% ==
                  zero? a => a
                  zero? b => b  -- avoid calling sortAndAdd with []
                  +/[sortAndAdd [[r, ta.Mn*tb.Mn]$Term
                    for tb in b | not zero?(r := ta.Cf*tb.Cf)]
                      for ta in reverse a]


        else -- M hasn't OrderedSet
            -- Terms are stored in random order.
          a = b ==
            #a ~= #b => false
            brace(a pretend List(Term)) =$Set(Term) brace(b pretend List(Term))

          coefficient(a, m) ==
            for t in a repeat
              t.Mn = m => return t.Cf
            0

          addterm(Tabl: AssociationList(M,R), r:R, m:M):R ==
              (u := search(m, Tabl)) case "failed" => Tabl.m := r
              zero?(r := r + u::R) => (remove!(m, Tabl); 0)
              Tabl.m := r

          a + b ==
              Tabl := table()$AssociationList(M,R)
              for t in a repeat
                  Tabl t.Mn := t.Cf
              for t in b repeat
                  addterm(Tabl, t.Cf, t.Mn)
              [[Tabl m, m]$Term for m in keys Tabl]

          a:% * b:% ==
              Tabl := table()$AssociationList(M,R)
              for ta in a repeat
                  for tb in (b pretend List(Term)) repeat
                      addterm(Tabl, ta.Cf*tb.Cf, ta.Mn*tb.Mn)
              [[Tabl.m, m]$Term for m in keys Tabl]

)abbrev package MRF2 MonoidRingFunctions2
++ Author: Johannes Grabmeier
++ Date Created: 14 May 1991
++ Date Last Updated: 14 May 1991
++ Basic Operations: map
++ Related Constructors:
++ Also See:
++ AMS Classifications:
++ Keywords: monoid ring, group ring, change of coefficient domain
++ References:
++ Description:
++  MonoidRingFunctions2 implements functions between
++  two monoid rings defined with the same monoid over different rings.
MonoidRingFunctions2(R,S,M) : Exports == Implementation where
    R  : Ring
    S  : Ring
    M  : Monoid
    Exports ==> with
      map: (R -> S, MonoidRing(R,M)) -> MonoidRing(S,M)
        ++ map(f,u) maps f onto the coefficients f the element
        ++ u of the monoid ring to create an element of a monoid
        ++ ring with the same monoid b.
    Implementation ==> add
      map(fn, u) ==
        res : MonoidRing(S,M) := 0
        for te in terms u repeat
          res := res + monomial(fn(te.coef), te.monom)
        res

