\documentclass{article}
\usepackage{open-axiom}
\begin{document}
\title{\$SPAD/src/algebra naalg.spad}
\author{Johannes Grabmeier, Robert Wisbauer}
\maketitle
\begin{abstract}
\end{abstract}
\eject
\tableofcontents
\eject
\section{domain ALGSC AlgebraGivenByStructuralConstants}
<<domain ALGSC AlgebraGivenByStructuralConstants>>=
)abbrev domain ALGSC AlgebraGivenByStructuralConstants
++ Authors: J. Grabmeier, R. Wisbauer
++ Date Created: 01 March 1991
++ Date Last Updated: 22 January 1992
++ Basic Operations:
++ Related Constructors:
++ Also See:
++ AMS Classifications:
++ Keywords: algebra, structural constants
++ Reference:
++  R.D. Schafer: An Introduction to Nonassociative Algebras
++  Academic Press, New York, 1966
++ Description:
++  AlgebraGivenByStructuralConstants implements finite rank algebras
++  over a commutative ring, given by the structural constants \spad{gamma}
++  with respect to a fixed  basis \spad{[a1,..,an]}, where
++  \spad{gamma} is an \spad{n}-vector of n by n matrices
++  \spad{[(gammaijk) for k in 1..rank()]} defined by
++  \spad{ai * aj = gammaij1 * a1 + ... + gammaijn * an}.
++  The symbols for the fixed basis
++  have to be given as a list of symbols.
AlgebraGivenByStructuralConstants(R:Field, n : PositiveInteger,_
  ls : List Symbol, gamma: Vector Matrix R ): public == private where

  V  ==> Vector
  M  ==> Matrix
  I  ==> Integer
  NNI  ==> NonNegativeInteger
  REC  ==> Record(particular: Union(V R,"failed"),basis: List V R)
  LSMP ==> LinearSystemMatrixPackage(R,V R,V R, M R)

  --public ==> FramedNonAssociativeAlgebra(R) with
  public ==> Join(FramedNonAssociativeAlgebra(R), _
    LeftModule(SquareMatrix(n,R)) ) with

    coerce : Vector R -> %
      ++ coerce(v) converts a vector to a member of the algebra
      ++ by forming a linear combination with the basis element.
      ++ Note: the vector is assumed to have length equal to the
      ++ dimension of the algebra.

  private ==> DirectProduct(n,R) add

    Rep := DirectProduct(n,R)

    x,y : %
    dp : DirectProduct(n,R)
    v : V R


    recip(x) == recip(x)$FiniteRankNonAssociativeAlgebra_&(%,R)

    (m:SquareMatrix(n,R))*(x:%) == apply((m :: Matrix R),x)
    coerce v == directProduct(v) :: %

    structuralConstants() == gamma

    coordinates(x) == vector(entries(x :: Rep)$Rep)$Vector(R)

    coordinates(x,b) ==
      --not (maxIndex b = n) =>
      --  error("coordinates: your 'basis' has not the right length")
      m : NonNegativeInteger := (maxIndex b) :: NonNegativeInteger
      transitionMatrix   : Matrix R := new(n,m,0$R)$Matrix(R)
      for i in 1..m repeat
        setColumn!(transitionMatrix,i,coordinates(b.i))
      res : REC := solve(transitionMatrix,coordinates(x))$LSMP
      if (not every?(zero?$R,first res.basis)) then
        error("coordinates: warning your 'basis' is linearly dependent")
      (res.particular  case "failed") =>
        error("coordinates: first argument is not in linear span of second argument")
      (res.particular) :: (Vector R)

    basis() == [unitVector(i::PositiveInteger)::% for i in 1..n]

    someBasis() == basis()$%

    rank() == n

    elt(x,i) == elt(x:Rep,i)$Rep

    coerce(x:%):OutputForm ==
      zero?(x::Rep)$Rep => (0$R) :: OutputForm
      le : List OutputForm := nil
      for i in 1..n repeat
        coef : R := elt(x::Rep,i)
        not zero?(coef)$R =>
          one?(coef)$R =>
            -- sy : OutputForm := elt(ls,i)$(List Symbol) :: OutputForm
            le := cons(elt(ls,i)$(List Symbol) :: OutputForm, le)
          le := cons(coef :: OutputForm *  elt(ls,i)$(List Symbol)_
              :: OutputForm, le)
      reduce("+",le)

    x * y ==
      v : Vector R :=  new(n,0)
      for k in 1..n repeat
        h : R := 0
        for i in 1..n repeat
          for j in 1..n repeat
            h := h  +$R elt(x,i) *$R elt(y,j) *$R elt(gamma.k,i,j )
        v.k := h
      directProduct v



    alternative?() ==
      for i in 1..n repeat
        -- expression for check of left alternative is symmetric in i and j:
        -- expression for check of right alternative is symmetric in j and k:
        for j in 1..i-1 repeat
          for k in j..n repeat
            -- right check
            for r in 1..n repeat
              res := 0$R
              for l in 1..n repeat
                res := res - _
                  (elt(gamma.l,j,k)+elt(gamma.l,k,j))*elt(gamma.r,i,l)+_
                  (elt(gamma.l,i,j)*elt(gamma.r,l,k) + elt(gamma.l,i,k)*_
                  elt(gamma.r,l,j) )
              not zero? res =>
                messagePrint("algebra is not right alternative")$OutputForm
                return false
        for j in i..n repeat
          for k in 1..j-1 repeat
            -- left check
            for r in 1..n repeat
              res := 0$R
              for l in 1..n repeat
                res := res + _
                  (elt(gamma.l,i,j)+elt(gamma.l,j,i))*elt(gamma.r,l,k)-_
                  (elt(gamma.l,j,k)*elt(gamma.r,i,l) + elt(gamma.l,i,k)*_
                   elt(gamma.r,j,l) )
              not (zero? res) =>
                messagePrint("algebra is not left alternative")$OutputForm
                return false

          for k in j..n repeat
            -- left check
            for r in 1..n repeat
              res := 0$R
              for l in 1..n repeat
                res := res + _
                  (elt(gamma.l,i,j)+elt(gamma.l,j,i))*elt(gamma.r,l,k)-_
                  (elt(gamma.l,j,k)*elt(gamma.r,i,l) + elt(gamma.l,i,k)*_
                   elt(gamma.r,j,l) )
              not (zero? res) =>
                messagePrint("algebra is not left alternative")$OutputForm
                return false
            -- right check
            for r in 1..n repeat
              res := 0$R
              for l in 1..n repeat
                res := res - _
                  (elt(gamma.l,j,k)+elt(gamma.l,k,j))*elt(gamma.r,i,l)+_
                  (elt(gamma.l,i,j)*elt(gamma.r,l,k) + elt(gamma.l,i,k)*_
                  elt(gamma.r,l,j) )
              not (zero? res) =>
                messagePrint("algebra is not right alternative")$OutputForm
                return false

      messagePrint("algebra satisfies 2*associator(a,b,b) = 0 = 2*associator(a,a,b) = 0")$OutputForm
      true

    -- should be in the category, but is not exported
--    conditionsForIdempotents b  ==
--      n := rank()
--      --gamma : Vector Matrix R := structuralConstants b
--      listOfNumbers : List String :=  [string(q)$String for q in 1..n]
--      symbolsForCoef : Vector Symbol :=
--        [concat("%", concat("x", i))::Symbol  for i in listOfNumbers]
--      conditions : List Polynomial R := []
 --     for k in 1..n repeat
 --       xk := symbolsForCoef.k
 --       p : Polynomial R :=  monomial( - 1$Polynomial(R), [xk], [1] )
 --       for i in 1..n repeat
 --         for j in 1..n repeat
 --           xi := symbolsForCoef.i
 --           xj := symbolsForCoef.j
 --           p := p + monomial(_
 --             elt((gamma.k),i,j) :: Polynomial(R), [xi,xj], [1,1])
 --       conditions := cons(p,conditions)
 --     conditions

    associative?() ==
      for i in 1..n repeat
       for j in 1..n repeat
        for k in 1..n repeat
         for r in 1..n repeat
           res := 0$R
           for l in 1..n repeat
             res := res + elt(gamma.l,i,j)*elt(gamma.r,l,k)-_
                          elt(gamma.l,j,k)*elt(gamma.r,i,l)
           not (zero? res) =>
            messagePrint("algebra is not associative")$OutputForm
            return false
      messagePrint("algebra is associative")$OutputForm
      true


    antiAssociative?() ==
      for i in 1..n repeat
       for j in 1..n repeat
        for k in 1..n repeat
         for r in 1..n repeat
           res := 0$R
           for l in 1..n repeat
             res := res + elt(gamma.l,i,j)*elt(gamma.r,l,k)+_
                          elt(gamma.l,j,k)*elt(gamma.r,i,l)
           not (zero? res) =>
            messagePrint("algebra is not anti-associative")$OutputForm
            return false
      messagePrint("algebra is anti-associative")$OutputForm
      true

    commutative?() ==
      for i in 1..n repeat
       for j in (i+1)..n repeat
        for k in 1..n repeat
           not ( elt(gamma.k,i,j)=elt(gamma.k,j,i) ) =>
            messagePrint("algebra is not commutative")$OutputForm
            return false
      messagePrint("algebra is commutative")$OutputForm
      true

    antiCommutative?() ==
      for i in 1..n repeat
       for j in i..n repeat
        for k in 1..n repeat
          not zero? (i=j => elt(gamma.k,i,i); elt(gamma.k,i,j)+elt(gamma.k,j,i) ) =>
            messagePrint("algebra is not anti-commutative")$OutputForm
            return false
      messagePrint("algebra is anti-commutative")$OutputForm
      true

    leftAlternative?() ==
      for i in 1..n repeat
       -- expression is symmetric in i and j:
       for j in i..n repeat
        for k in 1..n repeat
         for r in 1..n repeat
           res := 0$R
           for l in 1..n repeat
             res := res + (elt(gamma.l,i,j)+elt(gamma.l,j,i))*elt(gamma.r,l,k)-_
               (elt(gamma.l,j,k)*elt(gamma.r,i,l) + elt(gamma.l,i,k)*elt(gamma.r,j,l) )
           not (zero? res) =>
            messagePrint("algebra is not left alternative")$OutputForm
            return false
      messagePrint("algebra is left alternative")$OutputForm
      true


    rightAlternative?() ==
      for i in 1..n repeat
       for j in 1..n repeat
       -- expression is symmetric in j and k:
        for k in j..n repeat
         for r in 1..n repeat
           res := 0$R
           for l in 1..n repeat
             res := res - (elt(gamma.l,j,k)+elt(gamma.l,k,j))*elt(gamma.r,i,l)+_
               (elt(gamma.l,i,j)*elt(gamma.r,l,k) + elt(gamma.l,i,k)*elt(gamma.r,l,j) )
           not (zero? res) =>
            messagePrint("algebra is not right alternative")$OutputForm
            return false
      messagePrint("algebra is right alternative")$OutputForm
      true


    flexible?() ==
      for i in 1..n repeat
       for j in 1..n repeat
       -- expression is symmetric in i and k:
        for k in i..n repeat
         for r in 1..n repeat
           res := 0$R
           for l in 1..n repeat
             res := res + elt(gamma.l,i,j)*elt(gamma.r,l,k)-_
                          elt(gamma.l,j,k)*elt(gamma.r,i,l)+_
                          elt(gamma.l,k,j)*elt(gamma.r,l,i)-_
                          elt(gamma.l,j,i)*elt(gamma.r,k,l)
           not (zero? res) =>
            messagePrint("algebra is not flexible")$OutputForm
            return false
      messagePrint("algebra is flexible")$OutputForm
      true

    lieAdmissible?() ==
      for i in 1..n repeat
       for j in 1..n repeat
        for k in 1..n repeat
         for r in 1..n repeat
           res := 0$R
           for l in 1..n repeat
             res := res_
              + (elt(gamma.l,i,j)-elt(gamma.l,j,i))*(elt(gamma.r,l,k)-elt(gamma.r,k,l)) _
              + (elt(gamma.l,j,k)-elt(gamma.l,k,j))*(elt(gamma.r,l,i)-elt(gamma.r,i,l)) _
              + (elt(gamma.l,k,i)-elt(gamma.l,i,k))*(elt(gamma.r,l,j)-elt(gamma.r,j,l))
           not (zero? res) =>
            messagePrint("algebra is not Lie admissible")$OutputForm
            return false
      messagePrint("algebra is Lie admissible")$OutputForm
      true

    jordanAdmissible?()  ==
      recip(2 * 1$R) case "failed" =>
        messagePrint("this algebra is not Jordan admissible, as 2 is not invertible in the ground ring")$OutputForm
        false
      for i in 1..n repeat
       for j in 1..n repeat
        for k in 1..n repeat
         for w in 1..n repeat
          for t in 1..n repeat
           res := 0$R
           for l in 1..n repeat
            for r in 1..n repeat
             res := res_
              + (elt(gamma.l,i,j)+elt(gamma.l,j,i))_
                * (elt(gamma.r,w,k)+elt(gamma.r,k,w))_
                * (elt(gamma.t,l,r)+elt(gamma.t,r,l))_
              - (elt(gamma.r,w,k)+elt(gamma.r,k,w))_
                * (elt(gamma.l,j,r)+elt(gamma.l,r,j))_
                * (elt(gamma.t,i,l)+elt(gamma.t,l,i))_
              + (elt(gamma.l,w,j)+elt(gamma.l,j,w))_
                * (elt(gamma.r,k,i)+elt(gamma.r,i,k))_
                * (elt(gamma.t,l,r)+elt(gamma.t,r,l))_
              - (elt(gamma.r,k,i)+elt(gamma.r,k,i))_
                * (elt(gamma.l,j,r)+elt(gamma.l,r,j))_
                * (elt(gamma.t,w,l)+elt(gamma.t,l,w))_
              + (elt(gamma.l,k,j)+elt(gamma.l,j,k))_
                * (elt(gamma.r,i,w)+elt(gamma.r,w,i))_
                * (elt(gamma.t,l,r)+elt(gamma.t,r,l))_
              - (elt(gamma.r,i,w)+elt(gamma.r,w,i))_
                * (elt(gamma.l,j,r)+elt(gamma.l,r,j))_
                * (elt(gamma.t,k,l)+elt(gamma.t,l,k))
           not (zero? res) =>
             messagePrint("algebra is not Jordan admissible")$OutputForm
             return false
      messagePrint("algebra is Jordan admissible")$OutputForm
      true

    jordanAlgebra?()  ==
      recip(2 * 1$R) case "failed" =>
        messagePrint("this is not a Jordan algebra, as 2 is not invertible in the ground ring")$OutputForm
        false
      not commutative?() =>
        messagePrint("this is not a Jordan algebra")$OutputForm
        false
      for i in 1..n repeat
       for j in 1..n repeat
        for k in 1..n repeat
         for l in 1..n repeat
           for t in 1..n repeat
             res := 0$R
             for r in 1..n repeat
               for s in 1..n repeat
                 res := res +  _
                   elt(gamma.r,i,j)*elt(gamma.s,l,k)*elt(gamma.t,r,s) - _
                   elt(gamma.r,l,k)*elt(gamma.s,j,r)*elt(gamma.t,i,s) + _
                   elt(gamma.r,l,j)*elt(gamma.s,k,k)*elt(gamma.t,r,s) - _
                   elt(gamma.r,k,i)*elt(gamma.s,j,r)*elt(gamma.t,l,s) + _
                   elt(gamma.r,k,j)*elt(gamma.s,i,k)*elt(gamma.t,r,s) - _
                   elt(gamma.r,i,l)*elt(gamma.s,j,r)*elt(gamma.t,k,s)
                 not zero? res =>
                   messagePrint("this is not a Jordan algebra")$OutputForm
                   return false
      messagePrint("this is a Jordan algebra")$OutputForm
      true


    jacobiIdentity?()  ==
      for i in 1..n repeat
       for j in 1..n repeat
        for k in 1..n repeat
         for r in 1..n repeat
           res := 0$R
           for s in 1..n repeat
                 res := res +  elt(gamma.r,i,j)*elt(gamma.s,j,k) +_
                               elt(gamma.r,j,k)*elt(gamma.s,k,i) +_
                               elt(gamma.r,k,i)*elt(gamma.s,i,j)
           not zero? res =>
                 messagePrint("Jacobi identity does not hold")$OutputForm
                 return false
      messagePrint("Jacobi identity holds")$OutputForm
      true

@
\section{package SCPKG StructuralConstantsPackage}
<<package SCPKG StructuralConstantsPackage>>=
)abbrev package SCPKG StructuralConstantsPackage
++ Authors: J. Grabmeier
++ Date Created: 02 April 1992
++ Date Last Updated: 14 April 1992
++ Basic Operations:
++ Related Constructors: AlgebraPackage, AlgebraGivenByStructuralConstants
++ Also See:
++ AMS Classifications:
++ Keywords: structural constants
++ Reference:
++ Description:
++  StructuralConstantsPackage provides functions creating
++  structural constants from a multiplication tables or a basis
++  of a matrix algebra and other useful functions in this context.
StructuralConstantsPackage(R:Field): public == private where

  L  ==> List
  S  ==> Symbol
  FRAC ==> Fraction
  POLY ==> Polynomial
  V  ==> Vector
  M  ==> Matrix
  REC  ==> Record(particular: Union(V R,"failed"),basis: List V R)
  LSMP ==> LinearSystemMatrixPackage(R,V R,V R, M R)

  public ==>  with
      -- what we really want to have here is a matrix over
      -- linear polynomials in the list of symbols, having arbitrary
      -- coefficients from a ring extension of R, e.g. FRAC POLY R.
      structuralConstants : (L S, M FRAC POLY R) -> V M FRAC POLY R
        ++ structuralConstants(ls,mt) determines the structural constants
        ++ of an algebra with generators ls and multiplication table mt, the
        ++ entries of which must be given as linear polynomials in the
        ++ indeterminates given by ls. The result is in particular useful
        ++  as fourth argument for \spadtype{AlgebraGivenByStructuralConstants}
        ++  and \spadtype{GenericNonAssociativeAlgebra}.
      structuralConstants : (L S, M POLY R) -> V M POLY R
        ++ structuralConstants(ls,mt) determines the structural constants
        ++ of an algebra with generators ls and multiplication table mt, the
        ++ entries of which must be given as linear polynomials in the
        ++ indeterminates given by ls. The result is in particular useful
        ++  as fourth argument for \spadtype{AlgebraGivenByStructuralConstants}
        ++  and \spadtype{GenericNonAssociativeAlgebra}.
      structuralConstants: L M R -> V M R
        ++ structuralConstants(basis)  takes the basis of a matrix
        ++ algebra, e.g. the result of \spadfun{basisOfCentroid} and calculates
        ++ the structural constants.
        ++ Note, that the it is not checked, whether basis really is a
        ++ basis of a matrix algebra.
      coordinates: (M R, L M R) -> V R
        ++ coordinates(a,[v1,...,vn]) returns the coordinates of \spad{a}
        ++ with respect to the \spad{R}-module basis \spad{v1},...,\spad{vn}.

  private ==> add

      matrix2Vector: M R -> V R
      matrix2Vector m ==
        lili : L L R := listOfLists m
        --li : L R  := reduce(concat, listOfLists m)
        li : L R  := reduce(concat, lili)
        construct(li)$(V R)

      coordinates(x,b) ==
        m : NonNegativeInteger := (maxIndex b) :: NonNegativeInteger
        n : NonNegativeInteger := nrows(b.1) * ncols(b.1)
        transitionMatrix   : Matrix R := new(n,m,0$R)$Matrix(R)
        for i in 1..m repeat
          setColumn!(transitionMatrix,i,matrix2Vector(b.i))
        res : REC := solve(transitionMatrix,matrix2Vector(x))$LSMP
        if (not every?(zero?$R,first res.basis)) then
          error("coordinates: the second argument is linearly dependent")
        (res.particular  case "failed") =>
          error("coordinates: first argument is not in linear span of _
second argument")
        (res.particular) :: (Vector R)

      structuralConstants b ==
        --n := rank()
        -- be careful with the possibility that b is not a basis
        m : NonNegativeInteger := (maxIndex b) :: NonNegativeInteger
        sC : Vector Matrix R := [new(m,m,0$R) for k in 1..m]
        for i in 1..m repeat
          for j in 1..m repeat
            covec : Vector R := coordinates(b.i * b.j, b)$%
            for k in 1..m repeat
               setelt( sC.k, i, j, covec.k )
        sC

      structuralConstants(ls:L S, mt: M POLY R)  ==
        nn := #(ls)
        nrows(mt) ~= nn or ncols(mt) ~= nn =>
          error "structuralConstants: size of second argument does not _
agree with number of generators"
        gamma : L M POLY R := []
        lscopy : L S := copy ls
        while not null lscopy repeat
          mat : M POLY R := new(nn,nn,0)
          s : S := first lscopy
          for i in 1..nn repeat
            for j in 1..nn repeat
              p := qelt(mt,i,j)
              totalDegree(p,ls) > 1 =>
                error "structuralConstants: entries of second argument _
must be linear polynomials in the generators"
              if (c := coefficient(p, s, 1) ) ~= 0 then qsetelt!(mat,i,j,c)
          gamma := cons(mat, gamma)
          lscopy := rest lscopy
        vector reverse gamma

      structuralConstants(ls:L S, mt: M FRAC POLY R)  ==
        nn := #(ls)
        nrows(mt) ~= nn or ncols(mt) ~= nn =>
          error "structuralConstants: size of second argument does not _
agree with number of generators"
        gamma : L M FRAC(POLY R) := []
        lscopy : L S := copy ls
        while not null lscopy repeat
          mat : M FRAC(POLY R) := new(nn,nn,0)
          s : S := first lscopy
          for i in 1..nn repeat
            for j in 1..nn repeat
              r := qelt(mt,i,j)
              q := denom(r)
              totalDegree(q,ls) ~= 0 =>
                error "structuralConstants: entries of second argument _
must be (linear) polynomials in the generators"
              p := numer(r)
              totalDegree(p,ls) > 1 =>
                error "structuralConstants: entries of second argument _
must be linear polynomials in the generators"
              if (c := coefficient(p, s, 1) ) ~= 0 then qsetelt!(mat,i,j,c/q)
          gamma := cons(mat, gamma)
          lscopy := rest lscopy
        vector reverse gamma

@
\section{package ALGPKG AlgebraPackage}
<<package ALGPKG AlgebraPackage>>=
)abbrev package ALGPKG AlgebraPackage
++ Authors: J. Grabmeier, R. Wisbauer
++ Date Created: 04 March 1991
++ Date Last Updated: 04 April 1992
++ Basic Operations:
++ Related Constructors:
++ Also See:
++ AMS Classifications:
++ Keywords: rank, nucleus, nucloid, structural constants
++ Reference:
++  R.S. Pierce: Associative Algebras
++  Graduate Texts in Mathematics 88
++  Springer-Verlag,  Heidelberg, 1982, ISBN 0-387-90693-2
++
++  R.D. Schafer: An Introduction to Nonassociative Algebras
++  Academic Press, New York, 1966
++
++  A. Woerz-Busekros: Algebra in Genetics
++  Lectures Notes in Biomathematics 36,
++  Springer-Verlag,  Heidelberg, 1980
++ Description:
++  AlgebraPackage assembles a variety of useful functions for
++  general algebras.
AlgebraPackage(R:IntegralDomain, A: FramedNonAssociativeAlgebra(R)): _
   public == private where

  V  ==> Vector
  M  ==> Matrix
  I  ==> Integer
  NNI  ==> NonNegativeInteger
  REC  ==> Record(particular: Union(V R,"failed"),basis: List V R)
  LSMP ==> LinearSystemMatrixPackage(R,V R,V R, M R)

  public ==>  with

      leftRank: A -> NonNegativeInteger
        ++ leftRank(x) determines the number of linearly independent elements
        ++ in \spad{x*b1},...,\spad{x*bn},
        ++ where \spad{b=[b1,...,bn]} is a basis.
      rightRank: A -> NonNegativeInteger
        ++ rightRank(x) determines the number of linearly independent elements
        ++ in \spad{b1*x},...,\spad{bn*x},
        ++ where \spad{b=[b1,...,bn]} is a basis.
      doubleRank: A -> NonNegativeInteger
        ++ doubleRank(x) determines the number of linearly
        ++ independent elements
        ++ in \spad{b1*x},...,\spad{x*bn},
        ++ where \spad{b=[b1,...,bn]} is a basis.
      weakBiRank: A -> NonNegativeInteger
        ++ weakBiRank(x) determines the number of
        ++ linearly independent elements
        ++ in the \spad{bi*x*bj}, \spad{i,j=1,...,n},
        ++ where \spad{b=[b1,...,bn]} is a basis.
      biRank: A -> NonNegativeInteger
        ++ biRank(x) determines the number of linearly independent elements
        ++ in \spad{x}, \spad{x*bi}, \spad{bi*x}, \spad{bi*x*bj},
        ++ \spad{i,j=1,...,n},
        ++ where \spad{b=[b1,...,bn]} is a basis.
        ++ Note: if \spad{A} has a unit,
        ++ then \spadfunFrom{doubleRank}{AlgebraPackage},
        ++ \spadfunFrom{weakBiRank}{AlgebraPackage}
        ++ and \spadfunFrom{biRank}{AlgebraPackage} coincide.
      basisOfCommutingElements: () -> List A
        ++ basisOfCommutingElements() returns a basis of the space of
        ++ all x of \spad{A} satisfying \spad{0 = commutator(x,a)} for all
        ++ \spad{a} in \spad{A}.
      basisOfLeftAnnihilator: A -> List A
        ++ basisOfLeftAnnihilator(a) returns a basis of the space of
        ++ all x of \spad{A} satisfying \spad{0 = x*a}.
      basisOfRightAnnihilator: A -> List A
        ++ basisOfRightAnnihilator(a) returns a basis of the space of
        ++ all x of \spad{A} satisfying \spad{0 = a*x}.
      basisOfLeftNucleus: () -> List A
        ++ basisOfLeftNucleus() returns a basis of the space of
        ++ all x of \spad{A} satisfying \spad{0 = associator(x,a,b)}
        ++ for all \spad{a},b in \spad{A}.
      basisOfRightNucleus: () -> List A
        ++ basisOfRightNucleus() returns a basis of the space of
        ++ all x of \spad{A} satisfying \spad{0 = associator(a,b,x)}
        ++ for all \spad{a},b in \spad{A}.
      basisOfMiddleNucleus: () -> List A
        ++ basisOfMiddleNucleus() returns a basis of the space of
        ++ all x of \spad{A} satisfying \spad{0 = associator(a,x,b)}
        ++ for all \spad{a},b in \spad{A}.
      basisOfNucleus: () -> List A
        ++ basisOfNucleus() returns a basis of the space of all x of \spad{A} satisfying
        ++ \spad{associator(x,a,b) = associator(a,x,b) = associator(a,b,x) = 0}
        ++ for all \spad{a},b in \spad{A}.
      basisOfCenter: () -> List A
        ++ basisOfCenter() returns a basis of the space of
        ++ all x of \spad{A} satisfying \spad{commutator(x,a) = 0} and
        ++ \spad{associator(x,a,b) = associator(a,x,b) = associator(a,b,x) = 0}
        ++ for all \spad{a},b in \spad{A}.
      basisOfLeftNucloid:()-> List Matrix R
        ++ basisOfLeftNucloid() returns a basis of the space of
        ++ endomorphisms of \spad{A} as right module.
        ++ Note: left nucloid coincides with left nucleus if \spad{A} has a unit.
      basisOfRightNucloid:()-> List Matrix R
        ++ basisOfRightNucloid() returns a basis of the space of
        ++ endomorphisms of \spad{A} as left module.
        ++ Note: right nucloid coincides with right nucleus if \spad{A} has a unit.
      basisOfCentroid:()-> List Matrix R
        ++ basisOfCentroid() returns a basis of the centroid, i.e. the
        ++ endomorphism ring of \spad{A} considered as \spad{(A,A)}-bimodule.
      radicalOfLeftTraceForm: () -> List A
        ++ radicalOfLeftTraceForm() returns basis for null space of
        ++ \spad{leftTraceMatrix()}, if the algebra is
        ++ associative, alternative or a Jordan algebra, then this
        ++ space equals the radical (maximal nil ideal) of the algebra.
      if R has EuclideanDomain then
        basis : V A ->  V A
          ++ basis(va) selects a basis from the elements of va.


  private ==>  add

      -- constants

      n  : PositiveInteger := rank()$A
      n2 : PositiveInteger := n*n
      n3 : PositiveInteger := n*n2
      gamma : Vector Matrix R  := structuralConstants()$A


      -- local functions

      convVM : Vector R -> Matrix R
        -- converts n2-vector to (n,n)-matrix row by row
      convMV : Matrix R -> Vector R
        -- converts n-square matrix to  n2-vector row by row
      convVM v  ==
        cond : Matrix(R) := new(n,n,0$R)$M(R)
        z : Integer := 0
        for i in 1..n repeat
          for j in 1..n  repeat
            z := z+1
            setelt(cond,i,j,v.z)
        cond


      -- convMV m ==
      --     vec : Vector(R) := new(n*n,0$R)
      --     z : Integer := 0
      --     for i in 1..n repeat
      --       for j in 1..n  repeat
      --         z := z+1
      --         setelt(vec,z,elt(m,i,j))
      --     vec


      radicalOfLeftTraceForm() ==
        ma : M R := leftTraceMatrix()$A
        map(represents, nullSpace ma)$ListFunctions2(Vector R, A)


      basisOfLeftAnnihilator a ==
        ca : M R := transpose (coordinates(a) :: M R)
        cond : M R := reduce(vertConcat$(M R),
          [ca*transpose(gamma.i) for i in 1..#gamma])
        map(represents, nullSpace cond)$ListFunctions2(Vector R, A)

      basisOfRightAnnihilator a ==
        ca : M R := transpose (coordinates(a) :: M R)
        cond : M R := reduce(vertConcat$(M R),
          [ca*(gamma.i) for i in 1..#gamma])
        map(represents, nullSpace cond)$ListFunctions2(Vector R, A)

      basisOfLeftNucloid() ==
        cond : Matrix(R) := new(n3,n2,0$R)$M(R)
        condo: Matrix(R) := new(n3,n2,0$R)$M(R)
        z : Integer := 0
        for i in 1..n repeat
          for j in 1..n repeat
            r1  : Integer := 0
            for k in 1..n repeat
              z := z + 1
              -- z equals (i-1)*n*n+(j-1)*n+k (loop-invariant)
              r2 : Integer := i
              for r in 1..n repeat
                r1 := r1 + 1
                -- here r1 equals (k-1)*n+r (loop-invariant)
                setelt(cond,z,r1,elt(gamma.r,i,j))
                -- here r2 equals (r-1)*n+i (loop-invariant)
                setelt(condo,z,r2,-elt(gamma.k,r,j))
                r2 := r2 + n
        [convVM(sol) for sol in nullSpace(cond+condo)]

      basisOfCommutingElements() ==
        --gamma1 := first gamma
        --gamma1 := gamma1 - transpose gamma1
        --cond : Matrix(R) := gamma1 :: Matrix(R)
        --for  i in  2..n repeat
        --  gammak := gamma.i
        --  gammak := gammak - transpose gammak
        --  cond :=  vertConcat(cond, gammak :: Matrix(R))$Matrix(R)
        --map(represents, nullSpace cond)$ListFunctions2(Vector R, A)

        cond : M R := reduce(vertConcat$(M R),
          [(gam := gamma.i) - transpose gam for i in 1..#gamma])
        map(represents, nullSpace cond)$ListFunctions2(Vector R, A)

      basisOfLeftNucleus() ==
        condi: Matrix(R) := new(n3,n,0$R)$Matrix(R)
        z : Integer := 0
        for k in 1..n repeat
         for j in 1..n repeat
          for s in 1..n repeat
            z := z+1
            for i in 1..n repeat
              entry : R := 0
              for l in 1..n repeat
                entry :=  entry+elt(gamma.l,j,k)*elt(gamma.s,i,l)_
                               -elt(gamma.l,i,j)*elt(gamma.s,l,k)
              setelt(condi,z,i,entry)$Matrix(R)
        map(represents, nullSpace condi)$ListFunctions2(Vector R,A)

      basisOfRightNucleus() ==
        condo : Matrix(R) := new(n3,n,0$R)$Matrix(R)
        z : Integer := 0
        for k in 1..n repeat
         for j in 1..n repeat
          for s in 1..n repeat
            z := z+1
            for i in 1..n repeat
              entry : R := 0
              for l in 1..n repeat
                entry :=  entry+elt(gamma.l,k,i)*elt(gamma.s,j,l) _
                               -elt(gamma.l,j,k)*elt(gamma.s,l,i)
              setelt(condo,z,i,entry)$Matrix(R)
        map(represents, nullSpace condo)$ListFunctions2(Vector R,A)

      basisOfMiddleNucleus() ==
        conda : Matrix(R) := new(n3,n,0$R)$Matrix(R)
        z : Integer := 0
        for k in 1..n repeat
         for j in 1..n repeat
          for s in 1..n repeat
            z := z+1
            for i in 1..n repeat
              entry : R := 0
              for l in 1..n repeat
                entry :=  entry+elt(gamma.l,j,i)*elt(gamma.s,l,k)
                               -elt(gamma.l,i,k)*elt(gamma.s,j,l)
              setelt(conda,z,i,entry)$Matrix(R)
        map(represents, nullSpace conda)$ListFunctions2(Vector R,A)


      basisOfNucleus() ==
        condi: Matrix(R) := new(3*n3,n,0$R)$Matrix(R)
        z : Integer := 0
        u : Integer := n3
        w : Integer := 2*n3
        for k in 1..n repeat
         for j in 1..n repeat
          for s in 1..n repeat
            z := z+1
            u := u+1
            w := w+1
            for i in 1..n repeat
              entry : R := 0
              enter : R := 0
              ent   : R := 0
              for l in 1..n repeat
                entry :=  entry + elt(gamma.l,j,k)*elt(gamma.s,i,l) _
                                - elt(gamma.l,i,j)*elt(gamma.s,l,k)
                enter :=  enter + elt(gamma.l,k,i)*elt(gamma.s,j,l) _
                                - elt(gamma.l,j,k)*elt(gamma.s,l,i)
                ent :=  ent  +  elt(gamma.l,j,k)*elt(gamma.s,i,l) _
                             -  elt(gamma.l,j,i)*elt(gamma.s,l,k)
              setelt(condi,z,i,entry)$Matrix(R)
              setelt(condi,u,i,enter)$Matrix(R)
              setelt(condi,w,i,ent)$Matrix(R)
        map(represents, nullSpace condi)$ListFunctions2(Vector R,A)

      basisOfCenter() ==
        gamma1 := first gamma
        gamma1 := gamma1 - transpose gamma1
        cond : Matrix(R) := gamma1 :: Matrix(R)
        for  i in  2..n repeat
          gammak := gamma.i
          gammak := gammak - transpose gammak
          cond :=  vertConcat(cond, gammak :: Matrix(R))$Matrix(R)
        B := cond :: Matrix(R)
        condi: Matrix(R) := new(2*n3,n,0$R)$Matrix(R)
        z : Integer := 0
        u : Integer := n3
        for k in 1..n repeat
         for j in 1..n repeat
          for s in 1..n repeat
            z := z+1
            u := u+1
            for i in 1..n repeat
              entry : R := 0
              enter : R := 0
              for l in 1..n repeat
                entry :=  entry + elt(gamma.l,j,k)*elt(gamma.s,i,l) _
                                - elt(gamma.l,i,j)*elt(gamma.s,l,k)
                enter :=  enter + elt(gamma.l,k,i)*elt(gamma.s,j,l) _
                                - elt(gamma.l,j,k)*elt(gamma.s,l,i)
              setelt(condi,z,i,entry)$Matrix(R)
              setelt(condi,u,i,enter)$Matrix(R)
        D := vertConcat(condi,B)$Matrix(R)
        map(represents, nullSpace D)$ListFunctions2(Vector R, A)

      basisOfRightNucloid() ==
        cond : Matrix(R) := new(n3,n2,0$R)$M(R)
        condo: Matrix(R) := new(n3,n2,0$R)$M(R)
        z : Integer := 0
        for i in 1..n repeat
          for j in 1..n repeat
            r1  : Integer := 0
            for k in 1..n repeat
              z := z + 1
              -- z equals (i-1)*n*n+(j-1)*n+k (loop-invariant)
              r2 : Integer := i
              for r in 1..n repeat
                r1 := r1 + 1
                -- here r1 equals (k-1)*n+r (loop-invariant)
                setelt(cond,z,r1,elt(gamma.r,j,i))
                -- here r2 equals (r-1)*n+i (loop-invariant)
                setelt(condo,z,r2,-elt(gamma.k,j,r))
                r2 := r2 + n
        [convVM(sol) for sol in nullSpace(cond+condo)]

      basisOfCentroid() ==
        cond : Matrix(R) := new(2*n3,n2,0$R)$M(R)
        condo: Matrix(R) := new(2*n3,n2,0$R)$M(R)
        z : Integer := 0
        u : Integer := n3
        for i in 1..n repeat
          for j in 1..n repeat
            r1  : Integer := 0
            for k in 1..n repeat
              z := z + 1
              u := u + 1
              -- z equals (i-1)*n*n+(j-1)*n+k (loop-invariant)
              -- u equals n**3 + (i-1)*n*n+(j-1)*n+k (loop-invariant)
              r2 : Integer := i
              for r in 1..n repeat
                r1 := r1 + 1
                -- here r1 equals (k-1)*n+r (loop-invariant)
                setelt(cond,z,r1,elt(gamma.r,i,j))
                setelt(cond,u,r1,elt(gamma.r,j,i))
                -- here r2 equals (r-1)*n+i (loop-invariant)
                setelt(condo,z,r2,-elt(gamma.k,r,j))
                setelt(condo,u,r2,-elt(gamma.k,j,r))
                r2 := r2 + n
        [convVM(sol) for sol in nullSpace(cond+condo)]


      doubleRank x ==
        cond : Matrix(R) := new(2*n,n,0$R)
        for k in 1..n repeat
         z : Integer := 0
         u : Integer := n
         for j in 1..n repeat
           z := z+1
           u := u+1
           entry : R := 0
           enter : R := 0
           for i in 1..n repeat
             entry := entry + elt(x,i)*elt(gamma.k,j,i)
             enter := enter + elt(x,i)*elt(gamma.k,i,j)
           setelt(cond,z,k,entry)$Matrix(R)
           setelt(cond,u,k,enter)$Matrix(R)
        rank(cond)$(M R)

      weakBiRank(x) ==
        cond : Matrix(R) := new(n2,n,0$R)$Matrix(R)
        z : Integer := 0
        for i in 1..n repeat
          for j in 1..n repeat
            z := z+1
            for k in 1..n repeat
              entry : R := 0
              for l in 1..n repeat
               for s in 1..n repeat
                entry:=entry+elt(x,l)*elt(gamma.s,i,l)*elt(gamma.k,s,j)
              setelt(cond,z,k,entry)$Matrix(R)
        rank(cond)$(M R)

      biRank(x) ==
        cond : Matrix(R) := new(n2+2*n+1,n,0$R)$Matrix(R)
        z : Integer := 0
        for j in 1..n repeat
          for i in 1..n repeat
            z := z+1
            for k in 1..n repeat
              entry : R := 0
              for l in 1..n repeat
               for s in 1..n repeat
                entry:=entry+elt(x,l)*elt(gamma.s,i,l)*elt(gamma.k,s,j)
              setelt(cond,z,k,entry)$Matrix(R)
        u : Integer := n*n
        w : Integer := n*(n+1)
        c := n2 + 2*n + 1
        for j in 1..n repeat
           u := u+1
           w := w+1
           for k in 1..n repeat
             entry : R := 0
             enter : R := 0
             for i in 1..n repeat
               entry := entry + elt(x,i)*elt(gamma.k,j,i)
               enter := enter + elt(x,i)*elt(gamma.k,i,j)
             setelt(cond,u,k,entry)$Matrix(R)
             setelt(cond,w,k,enter)$Matrix(R)
           setelt(cond,c,j, elt(x,j))
        rank(cond)$(M R)

      leftRank x ==
        cond : Matrix(R) := new(n,n,0$R)
        for k in 1..n repeat
         for j in 1..n repeat
           entry : R := 0
           for i in 1..n repeat
             entry := entry + elt(x,i)*elt(gamma.k,i,j)
           setelt(cond,j,k,entry)$Matrix(R)
        rank(cond)$(M R)

      rightRank x ==
        cond : Matrix(R) := new(n,n,0$R)
        for k in 1..n repeat
         for j in 1..n repeat
           entry : R := 0
           for i in 1..n repeat
             entry := entry + elt(x,i)*elt(gamma.k,j,i)
           setelt(cond,j,k,entry)$Matrix(R)
        rank(cond)$(M R)


      if R has EuclideanDomain then
        basis va ==
          v : V A := remove(zero?, va)$(V A)
          v : V A := removeDuplicates v
          empty? v =>  [0$A]
          m : Matrix R := coerce(coordinates(v.1))$(Matrix R)
          for i in 2..maxIndex v repeat
            m := horizConcat(m,coerce(coordinates(v.i))$(Matrix R) )
          m := rowEchelon m
          lj : List Integer := []
          h : Integer := 1
          mRI : Integer := maxRowIndex m
          mCI : Integer := maxColIndex m
          finished? : Boolean := false
          j : Integer := 1
          while not finished? repeat
            not zero? m(h,j) =>  -- corner found
              lj := cons(j,lj)
              h := mRI
              while zero? m(h,j) repeat h := h-1
              finished? := (h = mRI)
              if not finished? then h := h+1
            if j < mCI then
              j := j + 1
            else
              finished? := true
          [v.j for j in reverse lj]

@
\section{package FRNAAF2 FramedNonAssociativeAlgebraFunctions2}
<<package FRNAAF2 FramedNonAssociativeAlgebraFunctions2>>=
)abbrev package FRNAAF2 FramedNonAssociativeAlgebraFunctions2
++ Author: Johannes Grabmeier
++ Date Created: 28 February 1992
++ Date Last Updated: 28 February 1992
++ Basic Operations: map
++ Related Constructors:
++ Also See:
++ AMS Classifications:
++ Keywords: non-associative algebra
++ References:
++ Description:
++  FramedNonAssociativeAlgebraFunctions2 implements functions between
++  two framed non associative algebra domains defined over different rings.
++  The function map is used to coerce between algebras over different
++  domains having the same structural constants.

FramedNonAssociativeAlgebraFunctions2(AR,R,AS,S) : Exports ==
  Implementation where
    R  : CommutativeRing
    S  : CommutativeRing
    AR : FramedNonAssociativeAlgebra R
    AS : FramedNonAssociativeAlgebra S
    V ==> Vector
    Exports ==> with
      map:     (R -> S, AR) -> AS
        ++ map(f,u) maps f onto the coordinates of u to get an element
        ++ in \spad{AS} via identification of the basis of \spad{AR}
        ++ as beginning part of the basis of \spad{AS}.
    Implementation ==> add
      map(fn : R -> S, u : AR): AS ==
        rank()$AR > rank()$AS => error("map: ranks of algebras do not fit")
        vr : V R := coordinates u
        vs : V S := map(fn,vr)$VectorFunctions2(R,S)
@
This line used to read:
\begin{verbatim}
        rank()$AR = rank()$AR => represents(vs)$AS
\end{verbatim}
but the test is clearly always true and cannot be what was intended.
Gregory Vanuxem supplied the fix below.
<<package FRNAAF2 FramedNonAssociativeAlgebraFunctions2>>=
        rank()$AR = rank()$AS => represents(vs)$AS
        ba := basis()$AS
        represents(vs,[ba.i for i in 1..rank()$AR])

@
\section{License}
<<license>>=
--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.
@
<<*>>=
<<license>>

<<domain ALGSC AlgebraGivenByStructuralConstants>>
<<package ALGPKG AlgebraPackage>>
<<package SCPKG StructuralConstantsPackage>>
<<package FRNAAF2 FramedNonAssociativeAlgebraFunctions2>>
@
\eject
\begin{thebibliography}{99}
\bibitem{1} nothing
\end{thebibliography}
\end{document}