axiom-source 20170501-3 / usr / share / axiom-20170501 / src / algebra / SGCF.spad

)abbrev package SGCF SymmetricGroupCombinatoricFunctions
++ Authors: Johannes Grabmeier, Thorsten Werther
++ Date Created: 03 September 1988
++ Date Last Updated: 07 June 1990
++ References:
++   G. James/ A. Kerber: The Representation Theory of the Symmetric
++    Group. Encycl. of Math. and its Appl., Vol. 16., Cambridge
++    Univ. Press 1981, ISBN 0-521-30236-6.
++   S.G. Williamson: Combinatorics for Computer Science,
++    Computer Science Press, Rockville, Maryland, USA, ISBN 0-88175-020-4.
++   A. Nijenhuis / H.S. Wilf: Combinatoral Algorithms, Academic Press 1978.
++    ISBN 0-12-519260-6.
++   H. Gollan, J. Grabmeier: Algorithms in Representation Theory and
++    their Realization in the Computer Algebra System Scratchpad,
++    Bayreuther Mathematische Schriften, Heft 33, 1990, 1-23.
++ Description:
++ SymmetricGroupCombinatoricFunctions contains combinatoric
++ functions concerning symmetric groups and representation
++ theory: list young tableaus, improper partitions, subsets
++ bijection of Coleman.
 
SymmetricGroupCombinatoricFunctions() : SIG == CODE where
 
  NNI  ==> NonNegativeInteger
  I    ==> Integer
  L    ==> List
  M    ==> Matrix
  V    ==> Vector
  B    ==> Boolean
  ICF  ==> IntegerCombinatoricFunctions Integer
 
  SIG ==> with
 
    coleman : (L I, L I, L I) -> M I
      ++ coleman(alpha,beta,pi):
      ++ there is a bijection from the set of matrices having nonnegative
      ++ entries and row sums alpha, column sums beta
      ++ to the set of Salpha - Sbeta double cosets of the
      ++ symmetric group Sn. (Salpha is the Young subgroup
      ++ corresponding to the improper partition alpha).
      ++ For a representing element pi of such a double coset,
      ++ coleman(alpha,beta,pi) generates the Coleman-matrix
      ++ corresponding to alpha, beta, pi.
      ++ Note that The permutation pi of {1,2,...,n} has to be given
      ++ in list form.
      ++ Note that the inverse of this map is inverseColeman
      ++ (if pi is the lexicographical smallest permutation
      ++ in the coset). For details see James/Kerber.

    inverseColeman : (L I, L I, M I) -> L I
      ++ inverseColeman(alpha,beta,C):
      ++ there is a bijection from the set of matrices having nonnegative
      ++ entries and row sums alpha, column sums beta
      ++ to the set of Salpha - Sbeta double cosets of the
      ++ symmetric group Sn. (Salpha is the Young subgroup
      ++ corresponding to the improper partition alpha).
      ++ For such a matrix C, inverseColeman(alpha,beta,C)
      ++ calculates the lexicographical smallest pi in the
      ++ corresponding double coset.
      ++ Note that the resulting permutation pi of {1,2,...,n} 
      ++ is given in list form.
      ++ Notes: the inverse of this map is coleman.
      ++ For details, see James/Kerber.

    listYoungTableaus : (L I) -> L M I
      ++ listYoungTableaus(lambda) where lambda is a proper partition
      ++ generates the list of all standard tableaus of shape lambda
      ++ by means of lattice permutations. The numbers of the lattice
      ++ permutation are interpreted as column labels. Hence the
      ++ contents of these lattice permutations are the conjugate of
      ++ lambda.
      ++ Notes: the functions nextLatticePermutation and
      ++ makeYoungTableau are used.
      ++ The entries are from 0,...,n-1.

    makeYoungTableau : (L I,L I) -> M I
      ++ makeYoungTableau(lambda,gitter) computes for a given lattice
      ++ permutation gitter and for an improper partition lambda
      ++ the corresponding standard tableau of shape lambda.
      ++ Notes: see listYoungTableaus.
      ++ The entries are from 0,...,n-1.

    nextColeman : (L I, L I, M I) -> M I
      ++ nextColeman(alpha,beta,C) generates the next Coleman matrix
      ++ of column sums alpha and row sums beta according
      ++ to the lexicographical order from bottom-to-top.
      ++ The first Coleman matrix is achieved by C=new(1,1,0).
      ++ Also, new(1,1,0) indicates that C is the last Coleman matrix.

    nextLatticePermutation : (L I, L I, B) -> L I
      ++ nextLatticePermutation(lambda,lattP,constructNotFirst) generates
      ++ the lattice permutation according to the proper partition
      ++ lambda succeeding the lattice permutation lattP in
      ++ lexicographical order as long as constructNotFirst is true.
      ++ If constructNotFirst is false, the first lattice permutation
      ++ is returned.
      ++ The result nil indicates that lattP has no successor.

    nextPartition : (V I, V I, I) -> V I
      ++ nextPartition(gamma,part,number) generates the partition of
      ++ number which follows part according to the right-to-left
      ++ lexicographical order. The partition has the property that
      ++ its components do not exceed the corresponding components of
      ++ gamma. The first partition is achieved by part=[].
      ++ Also, [] indicates that part is the last partition.

    nextPartition : (L I, V I, I) -> V I
      ++ nextPartition(gamma,part,number) generates the partition of
      ++ number which follows part according to the right-to-left
      ++ lexicographical order. The partition has the property that
      ++ its components do not exceed the corresponding components of
      ++ gamma. the first partition is achieved by part=[].
      ++ Also, [] indicates that part is the last partition.

    numberOfImproperPartitions: (I,I) -> I
      ++ numberOfImproperPartitions(n,m) computes the number of partitions
      ++ of the nonnegative integer n in m nonnegative parts with regarding
      ++ the order (improper partitions).
      ++ Example: numberOfImproperPartitions (3,3) is 10,
      ++ since [0,0,3], [0,1,2], [0,2,1], [0,3,0], [1,0,2], [1,1,1],
      ++ [1,2,0], [2,0,1], [2,1,0], [3,0,0] are the possibilities.
      ++ Note that this operation has a recursive implementation.

    subSet : (I,I,I) -> L I
      ++ subSet(n,m,k) calculates the k-th m-subset of the set
      ++ 0,1,...,(n-1) in the lexicographic order considered as
      ++ a decreasing map from 0,...,(m-1) into 0,...,(n-1).
      ++ See S.G. Williamson: Theorem 1.60.
      ++ Error: if not (0 <= m <= n and 0 < = k < (n choose m)).

    unrankImproperPartitions0 : (I,I,I) -> L I
      ++ unrankImproperPartitions0(n,m,k) computes the k-th improper
      ++ partition of nonnegative n in m nonnegative parts in reverse
      ++ lexicographical order.
      ++ Example: [0,0,3] < [0,1,2] < [0,2,1] < [0,3,0] <
      ++ [1,0,2] < [1,1,1] < [1,2,0] < [2,0,1] < [2,1,0] < [3,0,0].
      ++ Error: if k is negative or too big.
      ++ Note that counting of subtrees is done by 
      ++ numberOfImproperPartitions

    unrankImproperPartitions1: (I,I,I) -> L I
      ++ unrankImproperPartitions1(n,m,k) computes the k-th improper
      ++ partition of nonnegative n in at most m nonnegative parts 
      ++ ordered as follows: first, in reverse
      ++ lexicographically according to their non-zero parts, then
      ++ according to their positions (lexicographical order
      ++ using subSet: [3,0,0] < [0,3,0] < [0,0,3] < [2,1,0] <
      ++ [2,0,1] < [0,2,1] < [1,2,0] < [1,0,2] < [0,1,2] < [1,1,1].
      ++ Note that counting of subtrees is done by
      ++ numberOfImproperPartitionsInternal.
 
  CODE ==> add
 
    import Set I
 
    -- declaration of local functions
 
    numberOfImproperPartitionsInternal: (I,I,I) -> I
      -- this is used as subtree counting function in
      -- "unrankImproperPartitions1". For (n,m,cm) it counts
      -- the following set of m-tuples: The  first (from left
      -- to right) m-cm non-zero entries are equal, the remaining
      -- positions sum up to n. Example: (3,3,2) counts
      -- [x,3,0], [x,0,3], [0,x,3], [x,2,1], [x,1,2], x non-zero.
 
    -- definition of local functions
 
    numberOfImproperPartitionsInternal(n,m,cm) ==
      n = 0 => binomial(m,cm)$ICF
      cm = 0 and n > 0 => 0
      s := 0
      for i in 0..n-1 repeat
        s := s + numberOfImproperPartitionsInternal(i,m,cm-1)
      s
 
    -- definition of exported functions
 
    numberOfImproperPartitions(n,m) ==
      if n < 0 or m < 1 then return 0
      if m = 1 or n = 0 then return 1
      s := 0
      for i in 0..n repeat
        s := s + numberOfImproperPartitions(n-i,m-1)
      s
 
    unrankImproperPartitions0(n,m,k) ==
      l : L I  := nil$(L I)
      k < 0 => error"counting of partitions is started at 0"
      k >= numberOfImproperPartitions(n,m) =>
        error"there are not so many partitions"
      for t in 0..(m-2) repeat
        s : I := 0
        for y in 0..n repeat
          sOld := s
          s := s + numberOfImproperPartitions(n-y,m-t-1)
          if s > k then leave
        l := append(l,list(y)$(L I))$(L I)
        k := k - sOld
        n := n - y
      l := append(l,list(n)$(L I))$(L I)
      l
 
    unrankImproperPartitions1(n,m,k) ==
      -- we use the counting procedure of the leaves in a tree
      -- having the following structure: First of all non-zero
      -- labels for the sons. If addition along a path gives n,
      -- then we go on creating the subtree for (n choose cm)
      -- where cm is the length of the path. These subsets determine
      -- the positions for the non-zero labels for the partition
      -- to be formeded. The remaining positions are filled by zeros.
      nonZeros   : L I := nil$(L I)
      partition  : V I :=  new(m::NNI,0$I)$(V I)
      k < 0 => nonZeros
      k >= numberOfImproperPartitions(n,m) => nonZeros
      cm : I := m    --cm gives the depth of the tree
      while n ^= 0 repeat
        s : I := 0
        cm := cm - 1
        for y in n..1 by -1 repeat   --determination of the next son
          sOld := s  -- remember old s
          -- this functions counts the number of elements in a subtree
          s := s + numberOfImproperPartitionsInternal(n-y,m,cm)
          if s > k then leave
        -- y is the next son, so put it into the pathlist "nonZero"
        nonZeros := append(nonZeros,list(y)$(L I))$(L I)
        k := k - sOld    --updating
        n := n - y       --updating
      --having found all m-cm non-zero entries we change the structure
      --of the tree and determine the non-zero positions
      nonZeroPos : L I := reverse subSet(m,m-cm,k)
      --building the partition
      for i in 1..m-cm  repeat partition.(1+nonZeroPos.i) := nonZeros.i
      entries partition
 
    subSet(n,m,k) ==
      k < 0 or n < 0 or m < 0 or m > n =>
        error "improper argument to subSet"
      bin : I := binomial$ICF (n,m)
      k >= bin =>
        error "there are not so many subsets"
      l : L I  := []
      n = 0 => l
      mm : I := k
      s  : I := m
      for t in 0..(m-1) repeat
         for y in (s-1)..(n+1) repeat
            if binomial$ICF (y,s) > mm then leave
         l := append (l,list(y-1)$(L I))
         mm := mm - binomial$ICF (y-1,s)
         s := s-1
      l
 
    nextLatticePermutation(lambda, lattP, constructNotFirst) ==
      lprime  : L I  := conjugate(lambda)$PartitionsAndPermutations
      columns : NNI := (first(lambda)$(L I))::NNI
      rows    : NNI := (first(lprime)$(L I))::NNI
      n       : NNI :=(+/lambda)::NNI
      not constructNotFirst =>   -- first lattice permutation
        lattP := nil$(L I)
        for i in columns..1 by -1 repeat
          for l in 1..lprime(i) repeat
            lattP := cons(i,lattP)
        lattP
      help : V I := new(columns,0) -- entry help(i) stores the number
      -- of occurences of number i on our way from right to left
      rightPosition  : NNI := n
      leftEntry : NNI := lattP(rightPosition)::NNI
      ready  : B  := false
      until (ready or (not constructNotFirst)) repeat
        rightEntry : NNI := leftEntry
        leftEntry := lattP(rightPosition-1)::NNI
        help(rightEntry) := help(rightEntry) + 1
        -- search backward decreasing neighbour elements
        if rightEntry > leftEntry then
          if ((lprime(leftEntry)-help(leftEntry)) >_
            (lprime(rightEntry)-help(rightEntry)+1)) then
            -- the elements may be swapped because the number of occurances
            -- of leftEntry would still be greater than those of rightEntry
            ready := true
            j : NNI := leftEntry + 1
            -- search among the numbers leftEntry+1..rightEntry for the
            -- smallest one which can take the place of leftEntry.
            -- negation of condition above:
            while (help(j)=0) or ((lprime(leftEntry)-lprime(j))
              < (help(leftEntry)-help(j)+2)) repeat j := j + 1
            lattP(rightPosition-1) := j
            help(j) := help(j)-1
            help(leftEntry) := help(leftEntry) + 1
            -- reconstruct the rest of the list in increasing order
            for l in rightPosition..n repeat
              j := 0
              while help(1+j) = 0 repeat j := j + 1
              lattP(l::NNI) := j+1
              help(1+j) := help(1+j) - 1
        -- end of "if rightEntry > leftEntry"
        rightPosition := (rightPosition-1)::NNI
        if rightPosition = 1 then constructNotFirst := false
      -- end of repeat-loop
      not constructNotFirst =>  nil$(L I)
      lattP
 
    makeYoungTableau(lambda,gitter) ==
      lprime  : L I  := conjugate(lambda)$PartitionsAndPermutations
      columns : NNI := (first(lambda)$(L I))::NNI
      rows    : NNI := (first(lprime)$(L I))::NNI
      ytab    : M I  := new(rows,columns,0)
      help    : V I  := new(columns,1)
      i : I := -1     -- this makes the entries ranging from 0,..,n-1
                      -- i := 0 would make it from 1,..,n.
      j : I := 0
      for l in 1..maxIndex gitter repeat
        j := gitter(l)
        i := i + 1
        ytab(help(j),j) := i
        help(j) := help(j) + 1
      ytab
 
    listYoungTableaus(lambda) ==
      lattice   : L I
      ytab      : M I
      younglist : L M I := nil$(L M I)
      lattice   := nextLatticePermutation(lambda,lattice,false)
      until null lattice repeat
        ytab      := makeYoungTableau(lambda,lattice)
        younglist := append(younglist,[ytab]$(L M I))$(L M I)
        lattice   := nextLatticePermutation(lambda,lattice,true)
      younglist
 
    nextColeman(alpha,beta,C) ==
      nrow  : NNI := #beta
      ncol  : NNI := #alpha
      vnull : V I  := vector(nil()$(L I)) -- empty vector
      vzero : V I  := new(ncol,0)
      vrest : V I  := new(ncol,0)
      cnull : M I  := new(1,1,0)
      coleman := copy C
      if coleman ^= cnull then
        -- look for the first row of "coleman" that has a succeeding
        -- partition, this can be atmost row nrow-1
        i : NNI := (nrow-1)::NNI
        vrest := row(coleman,i) + row(coleman,nrow)
        --for k in 1..ncol repeat
        --  vrest(k) := coleman(i,k) + coleman(nrow,k)
        succ := nextPartition(vrest,row(coleman, i),beta(i))
        while (succ = vnull) repeat
          if i = 1 then return cnull -- part is last partition
          i := (i - 1)::NNI
          --for k in 1..ncol repeat
          --  vrest(k) := vrest(k) + coleman(i,k)
          vrest := vrest + row(coleman,i)
          succ := nextPartition(vrest, row(coleman, i), beta(i))
        j : I := i
        coleman := setRow_!(coleman, i, succ)
        --for k in 1..ncol repeat
        --  vrest(k) := vrest(k) - coleman(i,k)
        vrest := vrest - row(coleman,i)
      else
        vrest := vector alpha
        -- for k in 1..ncol repeat
        --  vrest(k) := alpha(k)
        coleman := new(nrow,ncol,0)
        j : I := 0
      for i in (j+1)::NNI..nrow-1 repeat
        succ := nextPartition(vrest,vnull,beta(i))
        coleman := setRow_!(coleman, i, succ)
        vrest := vrest - succ
        --for k in 1..ncol repeat
        --  vrest(k) := vrest(k) - succ(k)
      setRow_!(coleman, nrow, vrest)
 
    nextPartition(gamma:V I, part:V I, number:I) ==
      nextPartition(entries gamma, part, number)
 
    nextPartition(gamma:L I,part:V I,number:I) ==
      n : NNI := #gamma
      vnull : V I := vector(nil()$(L I)) -- empty vector
      if part ^= vnull then
        i : NNI := 2
        sum := part(1)
        while (part(i) = gamma(i)) or (sum = 0) repeat
          sum := sum + part(i)
          i := i + 1
          if i = 1+n then return vnull -- part is last partition
        sum := sum - 1
        part(i) := part(i) + 1
      else
        sum := number
        part := new(n,0)
        i := 1+n
      j : NNI := 1
      while sum > gamma(j) repeat
        part(j) := gamma(j)
        sum := sum - gamma(j)
        j := j + 1
      part(j) := sum
      for k in j+1..i-1 repeat
        part(k) := 0
      part
 
    inverseColeman(alpha,beta,C) ==
      pi   : L I  := nil$(L I)
      nrow : NNI := #beta
      ncol : NNI := #alpha
      help : V I  := new(nrow,0)
      sum  : I   := 1
      for i in 1..nrow repeat
        help(i) := sum
        sum := sum + beta(i)
      for j in 1..ncol repeat
        for i in 1..nrow repeat
          for k in 2..1+C(i,j) repeat
            pi := append(pi,list(help(i))$(L I))
            help(i) := help(i) + 1
      pi
 
    coleman(alpha,beta,pi) ==
      nrow : NNI := #beta
      ncol : NNI := #alpha
      temp : L L I := nil$(L L I)
      help : L I  := nil$(L I)
      colematrix : M I := new(nrow,ncol,0)
      betasum  : NNI := 0
      alphasum : NNI := 0
      for i in 1..ncol repeat
        help := nil$(L I)
        for j in alpha(i)..1 by-1 repeat
          help := cons(pi(j::NNI+alphasum),help)
        alphasum := (alphasum + alpha(i))::NNI
        temp := append(temp,list(help)$(L L I))
      for i in 1..nrow repeat
        help := nil$(L I)
        for j in beta(i)..1 by-1 repeat
          help := cons(j::NNI+betasum, help)
        betasum := (betasum + beta(i))::NNI
        for j in 1..ncol repeat
          colematrix(i,j) := #intersect(brace(help),brace(temp(j)))
      colematrix