--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 package CYCLES CycleIndicators
++ Polya-Redfield enumeration by cycle indices.
++ Author: William H. Burge
++ Date Created: 1986
++ Date Last Updated: 11 Feb 1992
++ Keywords:Polya, Redfield, enumeration
++ Examples:
++ References: J.H.Redfield, 'The Theory of Group-Reduced Distributions',
++             American J. Math., 49 (1927) 433-455.
++             G.Polya, 'Kombinatorische Anzahlbestimmungen fur Gruppen,
++               Graphen und chemische Verbindungen', Acta Math. 68
++                (1937) 145-254.
++ Description: Enumeration by cycle indices. 
CycleIndicators: Exports == Implementation where
  I    ==> Integer
  L    ==> List
  B    ==> Boolean
  SPOL ==> SymmetricPolynomial
  PTN  ==> Partition
  RN   ==> Fraction Integer
  FR   ==> Factored Integer
  macro NNI == NonNegativeInteger
  macro PI == PositiveInteger
  Exports ==> with
 
    complete: PI -> SPOL RN
      ++\spad{complete n} is the \spad{n} th complete homogeneous
      ++ symmetric function expressed in terms of power sums.
      ++ Alternatively it is the cycle index of the symmetric
      ++ group of degree n.
 
    powerSum: PI -> SPOL RN
      ++\spad{powerSum n} is the \spad{n} th power sum symmetric
      ++ function.
 
    elementary: PI -> SPOL RN
      ++\spad{elementary n} is the \spad{n} th elementary symmetric
      ++ function expressed in terms of power sums.
 
    alternating: PI -> SPOL RN
      ++\spad{alternating n} is the cycle index of the
      ++ alternating group of degree n.
 
    cyclic: PI -> SPOL RN    --cyclic group
      ++\spad{cyclic n} is the cycle index of the
      ++ cyclic group of degree n.
 
    dihedral: PI -> SPOL RN    --dihedral group
      ++\spad{dihedral n} is the cycle index of the
      ++ dihedral group of degree n.
 
    graphs: PI -> SPOL RN
      ++\spad{graphs n} is the cycle index of the group induced on
      ++ the edges of a graph by applying the symmetric function to the
      ++ n nodes.
 
    cap: (SPOL RN,SPOL RN) -> RN
      ++\spad{cap(s1,s2)}, introduced by Redfield,
      ++ is the scalar product of two cycle indices.
 
    cup: (SPOL RN,SPOL RN) -> SPOL RN
      ++\spad{cup(s1,s2)}, introduced by Redfield,
      ++ is the scalar product of two cycle indices, in which the
      ++ power sums are retained to produce a cycle index.
 
    eval: SPOL RN -> RN
      ++\spad{eval s} is the sum of the coefficients of a cycle index.
 
    wreath: (SPOL RN,SPOL RN) -> SPOL RN
      ++\spad{wreath(s1,s2)} is the cycle index of the wreath product
      ++ of the two groups whose cycle indices are \spad{s1} and
      ++ \spad{s2}.
 
    SFunction:L PI -> SPOL RN
      ++\spad{SFunction(li)} is the S-function of the partition \spad{li}
      ++ expressed in terms of power sum symmetric functions.
 
    skewSFunction:(L I,L I) -> SPOL RN
      ++\spad{skewSFunction(li1,li2)} is the S-function
      ++ of the partition difference \spad{li1 - li2}
      ++ expressed in terms of power sum symmetric functions.
 
  Implementation ==> add
    import IntegerNumberTheoryFunctions
    import Partition
 
    trm: PTN -> SPOL RN
    trm pt == monomial(inv(pdct(pt) :: RN),pt)
 
    list: Stream PTN -> L PTN
    list st == entries complete st
 
    complete i ==
      i=0 => 1
      +/[trm pt for pt in list partitions i]
 
 
    even?: PTN -> B
    even? p == even?( #([i for i in parts p | even? i]))
 
    alternating i ==
      2 * _+/[trm p for p in list partitions i | even? p]

    elementary i ==
       i=0 => 1
       +/[(spol := trm pt; even? pt => spol; -spol)
                          for pt in list partitions i]
 
    divisors: I -> L I
    divisors n ==
      b := factors(n :: FR)
      c := concat(1,"append"/
                 [[a.factor**j for j in 1..a.exponent] for a in b]);
      if #(b) = 1 then c else concat(n,c)
 
    ss: (PI,I) -> SPOL RN
    ss(n,m) ==
      li : L PI := [n for j in 1..m]
      monomial(1,partition li)
 
    powerSum n == ss(n,1)
 
    cyclic n ==
      n = 1 => powerSum 1
      +/[(eulerPhi(i) / n) * ss(i::PI,numer(n/i)) for i in divisors n]
 
    dihedral n ==
      k := n quo 2
      odd? n => (1/2) * cyclic n + (1/2) * ss(2,k) * powerSum 1
      (1/2) * cyclic n + (1/4) * ss(2,k) + (1/4) * ss(2,k-1) * ss(1,2)
 
    trm2: PTN -> SPOL RN
    trm2 li ==
      lli := powers( li)$PTN
      xx := 1/(pdct li)
      prod : SPOL RN := 1
      for ll in lli repeat
        ll0 := first ll; ll1 := second ll
        k := ll0 quo 2
        c :=
          odd? ll0 => ss(ll0,ll1 * k)
          ss(k::PI,ll1) * ss(ll0,ll1 * (k - 1))
        c := c * ss(ll0,ll0 * ((ll1*(ll1 - 1)) quo 2))
        prod2 : SPOL RN := 1
        for r in lli | first(r) < ll0 repeat
          r0 := first r; r1 := second r
          prod2 := ss(lcm(r0,ll0)::PI,gcd(r0,ll0) * r1 * ll1) * prod2
        prod := c * prod2 * prod
      xx * prod
 
    graphs n == +/[trm2 p for p in list(partitions n)]
 
    cupp: (PTN,SPOL RN) -> SPOL RN
    cupp(pt,spol) ==
      zero? spol => 0
      (dg := degree spol) < pt => 0
      dg = pt => (pdct pt) * monomial(leadingCoefficient spol,dg)
      cupp(pt,reductum spol)
 
    cup(spol1,spol2) ==
      zero? spol1 => 0
      p := leadingCoefficient(spol1) * cupp(degree spol1,spol2)
      p + cup(reductum spol1,spol2)
 
    eval spol ==
      zero? spol => 0
      leadingCoefficient(spol) + eval(reductum spol)
 
    cap(spol1,spol2) == eval cup(spol1,spol2)
 
    mtpol: (PI,SPOL RN) -> SPOL RN
    mtpol(n,spol)==
      zero? spol => 0
      deg := partition [n*k for k in (degree spol)::L(PI)]
      monomial(leadingCoefficient spol,deg) + mtpol(n,reductum spol)
 
    evspol: ((PI -> SPOL RN),SPOL RN) -> SPOL RN
    evspol(fn2,spol) ==
      zero? spol => 0
      lc := leadingCoefficient spol
      prod := */[fn2 i for i in (degree spol)::L(PI)]
      lc * prod + evspol(fn2,reductum spol)
 
    wreath(spol1,spol2) == evspol(mtpol(#1,spol2),spol1)
 
    SFunction li==
      a:Matrix SPOL RN :=
        matrix [[complete((k -j+i)::PI) for k in li for j in 1..#li]
                    for i in 1..#li]
      determinant a
 
    roundup:(L I,L I)-> L I
    roundup(li1,li2)==
              #li1 > #li2 => roundup(li1,concat(li2,0))
              li2
 
    skewSFunction(li1,li2)==
      #li1 < #li2 =>
        error "skewSFunction: partition1 does not include partition2"
      li2:=roundup (li1,li2)
      a:Matrix SPOL RN:=matrix [[complete((k-li2.i-j+i)::PI)
               for k in li1 for j in 1..#li1]  for i in 1..#li1]
      determinant a

)abbrev package EVALCYC EvaluateCycleIndicators
++ Author: William H. Burge
++ Date Created: 1986
++ Date Last Updated: Feb 1992
++ Basic Operations:
++ Related Domains:
++ Also See:
++ AMS Classifications:
++ Keywords:
++ Examples:
++ References:
++ Description: This package is to be used in conjuction with
++             the CycleIndicators package. It provides an evaluation
++             function for SymmetricPolynomials.
EvaluateCycleIndicators(F):T==C where
    F:Algebra Fraction Integer
    I==>Integer
    L==>List
    SPOL==SymmetricPolynomial
    RN==>Fraction Integer
    PR==>Polynomial(RN)
    PTN==>Partition()
    lc ==> leadingCoefficient
    red ==> reductum
    T== with
       eval:((I->F),SPOL RN)->F
         ++\spad{eval(f,s)} evaluates the cycle index s by applying
         ++ the function f to each integer in a monomial partition,
         ++ forms their product and sums the results over all monomials.
    C== add
       evp:((I->F),PTN)->F
       fn:I->F
       pt:PTN
       spol:SPOL RN
       evp(fn, pt)== */[fn i for i in pt::L(PositiveInteger)]
 
       eval(fn,spol)==
        if spol=0
        then 0
        else ((lc spol)* evp(fn,degree spol)) + eval(fn,red spol)

