Anforderungen  |   Konzepte  |   Entwurf  |   Entwicklung  |   Qualitätssicherung  |   Lebenszyklus  |   Steuerung
 
 
 
 


Quelle  algebra.gi   Sprache: unbekannt

 
Spracherkennung für: .gi vermutete Sprache: Unknown {[0] [0] [0]} [Methode: Schwerpunktbildung, einfache Gewichte, sechs Dimensionen]

#############################################################################
##
##  This file is part of GAP, a system for computational discrete algebra.
##  This file's authors include Thomas Breuer, and Willem de Graaf.
##
##  Copyright of GAP belongs to its developers, whose names are too numerous
##  to list here. Please refer to the COPYRIGHT file for details.
##
##  SPDX-License-Identifier: GPL-2.0-or-later
##
##  This file contains generic methods for algebras and algebras-with-one.
##


#############################################################################
##
#M  Representative( <A> ) . . . . . . . . one element of a left operator ring
##
InstallMethod( Representative,
    "for left operator ring with known generators",
    [ IsLeftOperatorRing and HasGeneratorsOfLeftOperatorRing ],
    RepresentativeFromGenerators( GeneratorsOfLeftOperatorRing ) );


#############################################################################
##
#M  Representative( <A> ) . . .  one element of a left operator ring-with-one
##
InstallMethod( Representative,
    "for left operator ring-with-one with known generators",
    [ IsLeftOperatorRingWithOne and HasGeneratorsOfLeftOperatorRingWithOne ],
    RepresentativeFromGenerators( GeneratorsOfLeftOperatorRingWithOne ) );


#############################################################################
##
#M  FLMLORByGenerators( <R>, <gens> ) . . . .  <R>-FLMLOR generated by <gens>
#M  FLMLORByGenerators( <R>, <gens>, <zero> )
##
InstallMethod( FLMLORByGenerators,
    "for ring and collection",
    [ IsRing, IsCollection ],
    function( R, gens )
    local A;
    A:= Objectify( NewType( FamilyObj( gens ),
                            IsFLMLOR and IsAttributeStoringRep ),
                   rec() );
    SetLeftActingDomain( A, R );
    SetGeneratorsOfLeftOperatorRing( A, AsList( gens ) );

    CheckForHandlingByNiceBasis( R, gens, A, false );
    return A;
    end );

InstallOtherMethod( FLMLORByGenerators,
    "for ring, homogeneous list, and ring element",
    [ IsRing, IsHomogeneousList, IsRingElement ],
    function( R, gens, zero )
    local A;
    A:= Objectify( NewType( CollectionsFamily( FamilyObj( zero ) ),
                            IsFLMLOR and IsAttributeStoringRep ),
                   rec() );
    SetLeftActingDomain( A, R );
    SetGeneratorsOfLeftOperatorRing( A, gens );
    SetZero( A, zero );

    if IsEmpty( gens ) then
      SetDimension( A, 0 );
      SetGeneratorsOfLeftModule( A, gens );
    fi;

    CheckForHandlingByNiceBasis( R, gens, A, zero );
    return A;
    end );


#############################################################################
##
#M  FLMLORWithOneByGenerators( <R>, <gens> )  unit. <R>-FLMLOR gen. by <gens>
#M  FLMLORWithOneByGenerators( <R>, <gens>, <zero> )
##
InstallMethod( FLMLORWithOneByGenerators,
    "for ring and collection",
    [ IsRing, IsCollection ],
    function( R, gens )
    local A;
    A:= Objectify( NewType( FamilyObj( gens ),
                            IsFLMLORWithOne and IsAttributeStoringRep ),
                   rec() );
    SetLeftActingDomain( A, R );
    SetGeneratorsOfLeftOperatorRingWithOne( A, AsList( gens ) );

    CheckForHandlingByNiceBasis( R, gens, A, false );
    return A;
    end );

InstallOtherMethod( FLMLORWithOneByGenerators,
    "for ring, homogeneous list, and ring element",
    [ IsRing, IsHomogeneousList, IsRingElement ],
    function( R, gens, zero )
    local A;
    A:= Objectify( NewType( CollectionsFamily( FamilyObj( zero ) ),
                            IsFLMLORWithOne and IsAttributeStoringRep ),
                   rec() );
    SetLeftActingDomain( A, R );
    SetGeneratorsOfLeftOperatorRingWithOne( A, AsList( gens ) );
    SetZero( A, zero );

    CheckForHandlingByNiceBasis( R, gens, A, zero );
    return A;
    end );


#############################################################################
##
#F  Algebra( <F>, <gens> )
#F  Algebra( <F>, <gens>, <zero> )
#F  Algebra( <F>, <gens>, "basis" )
#F  Algebra( <F>, <gens>, <zero>, "basis" )
##
InstallGlobalFunction( FLMLOR, function( arg )
    local A;

    # ring and list of generators
    if Length( arg ) = 2 and IsRing( arg[1] )
                         and IsList( arg[2] ) and 0 < Length( arg[2] ) then
      A:= FLMLORByGenerators( arg[1], arg[2] );

    # ring, list of generators plus zero
    elif Length( arg ) = 3 and IsRing( arg[1] )
                           and IsList( arg[2] ) then
      if arg[3] = "basis" then
        A:= FLMLORByGenerators( arg[1], arg[2] );
        UseBasis( A, arg[2] );
      else
        A:= FLMLORByGenerators( arg[1], arg[2], arg[3] );
      fi;

    # ring, list of generators plus zero
    elif Length( arg ) = 4 and IsRing( arg[1] )
                           and IsList( arg[2] )
                           and arg[4] = "basis" then
      A:= FLMLORByGenerators( arg[1], arg[2], arg[3] );
      UseBasis( A, arg[2] );

    # no argument given, error
    else
      Error( "usage: FLMLOR( <F>, <gens> ), ",
             "FLMLOR( <F>, <gens>, <zero> )" );
    fi;

    # Return the result.
    return A;
end );


#############################################################################
##
#F  Subalgebra( <A>, <gens> ) . . . . . subalgebra of <A> generated by <gens>
#F  Subalgebra( <A>, <gens>, "basis" )
##
InstallGlobalFunction( SubFLMLOR, function( arg )
    local S;
    if    Length( arg ) <= 1
       or not IsFLMLOR( arg[1] )
       or not IsHomogeneousList( arg[2] ) then
      Error( "first argument must be a FLMLOR,\n",
             "second argument must be a list of generators" );

    elif IsEmpty( arg[2] ) then

      return SubFLMLORNC( arg[1], arg[2] );

    elif     IsIdenticalObj( FamilyObj( arg[1] ),
                          FamilyObj( arg[2] ) )
         and ForAll( arg[2], v -> v in arg[1] ) then

      S:= FLMLORByGenerators( LeftActingDomain( arg[1] ), arg[2] );
      SetParent( S, arg[1] );
      if Length( arg ) = 3 and arg[3] = "basis" then
        UseBasis( S, arg[2] );
      fi;
      return S;

    fi;
    Error( "usage: SubFLMLOR( <V>, <gens> [, \"basis\"] )" );
end );


#############################################################################
##
#F  SubalgebraNC( <A>, <gens>, "basis" )
#F  SubalgebraNC( <A>, <gens> )
##
InstallGlobalFunction( SubFLMLORNC, function( arg )
    local S;
    if IsEmpty( arg[2] ) then
      S:= Objectify( NewType( FamilyObj( arg[1] ),
                                  IsFLMLOR
                              and IsTrivial
                              and IsTwoSidedIdealInParent
                              and IsAttributeStoringRep ),
                     rec() );
      SetDimension(S, 0);
      SetLeftActingDomain( S, LeftActingDomain( arg[1] ) );
      SetGeneratorsOfLeftModule( S, AsList( arg[2] ) );
    else
      S:= FLMLORByGenerators( LeftActingDomain( arg[1] ), arg[2] );
    fi;
    if Length( arg ) = 3 and arg[3] = "basis" then
      UseBasis( S, arg[2] );
    fi;
    SetParent( S, arg[1] );
    return S;
end );


#############################################################################
##
#F  AlgebraWithOne( <F>, <gens> )
#F  AlgebraWithOne( <F>, <gens>, <zero> )
#F  AlgebraWithOne( <F>, <gens>, "basis" )
#F  AlgebraWithOne( <F>, <gens>, <zero>, "basis" )
##
InstallGlobalFunction( FLMLORWithOne, function( arg )
    local A;

    # ring and list of generators
    if Length( arg ) = 2 and IsRing( arg[1] )
                         and IsList( arg[2] ) and 0 < Length( arg[2] ) then
      A:= FLMLORWithOneByGenerators( arg[1], arg[2] );

    # ring, list of generators plus zero
    elif Length( arg ) = 3 and IsRing( arg[1] )
                           and IsList( arg[2] ) then
      if arg[3] = "basis" then
        A:= FLMLORWithOneByGenerators( arg[1], arg[2] );
        UseBasis( A, arg[2] );
      else
        A:= FLMLORWithOneByGenerators( arg[1], arg[2], arg[3] );
      fi;

    # ring, list of generators plus zero
    elif Length( arg ) = 4 and IsRing( arg[1] )
                           and IsList( arg[2] )
                           and arg[4] = "basis" then
      A:= FLMLORWithOneByGenerators( arg[1], arg[2], arg[3] );
      UseBasis( A, arg[2] );

    # no argument given, error
    else
      Error( "usage: FLMLORWithOne( <F>, <gens> ), ",
             "FLMLORWithOne( <F>, <gens>, <zero> )" );
    fi;

    # Return the result.
    return A;
end );


#############################################################################
##
#F  SubalgebraWithOne( <A>, <gens> )   subalg.-with-one of <A> gen. by <gens>
##
InstallGlobalFunction( SubFLMLORWithOne, function( arg )
    local S;
    if Length( arg ) <= 1 or not IsFLMLOR( arg[1] )
                          or not IsHomogeneousList( arg[2] ) then

      Error( "first argument must be a FLMLOR,\n",
             "second argument must be a list of generators" );

    elif IsEmpty( arg[2] ) then

      return SubFLMLORWithOneNC( arg[2], arg[2] );

    elif     IsIdenticalObj( FamilyObj( arg[1] ),
                          FamilyObj( arg[2] ) )
         and ForAll( arg[2], v -> v in arg[1] )
         and ( IsFLMLORWithOne( arg[1] ) or One( arg[1] ) <> fail ) then
      S:= FLMLORWithOneByGenerators( LeftActingDomain( arg[1] ), arg[2] );
      SetParent( S, arg[1] );
      if Length( arg ) = 3 and arg[3] = "basis" then
        UseBasis( S, arg[2] );
      fi;
      return S;
    fi;
    Error( "usage: SubFLMLORWithOne( <V>, <gens> [, \"basis\"] )" );
end );


#############################################################################
##
#F  SubalgebraWithOneNC( <A>, <gens> )
##
InstallGlobalFunction( SubFLMLORWithOneNC, function( arg )
    local S, gens;
    if IsEmpty( arg[2] ) then

      # Note that `S' is in general not trivial,
      # and if we call `Objectify' here then `S' does not get a special
      # representation (e.g., as a matrix algebra).
      # So the argument that special methods would catch this case
      # does not hold!
      gens:= [ One( arg[1] ) ];
      S:= FLMLORWithOneByGenerators( LeftActingDomain( arg[1] ), gens );
      UseBasis( S, gens );

    else

      S:= FLMLORWithOneByGenerators( LeftActingDomain( arg[1] ), arg[2] );
      if Length( arg ) = 3 and arg[3] = "basis" then
        UseBasis( S, arg[2] );
      fi;

    fi;

    SetParent( S, arg[1] );
    return S;
end );


#############################################################################
##
#M  LieAlgebraByDomain( <A> )
##
##  The Lie algebra of the associative algebra <A>
##
InstallMethod( LieAlgebraByDomain,
    "for an algebra",
    [ IsAlgebra ],
    function( A )
       local T, n, zero, nullvec, S, i, j, k, m, cfs, cij, cji;

       if not IsAssociative( A ) then TryNextMethod(); fi;

# We construct a structure constants table for the  Lie algebra
# corresponding to <A>. If the structure constants of <A> are given by
# d_{ij}^k, then the structure constants of the Lie algebra will be given
# by d_{ij}^k - d_{ji}^k.


       T:= StructureConstantsTable( Basis( A ) );
       n:= Dimension( A );
       zero:= Zero( LeftActingDomain( A ) );
       nullvec:= List( [1..n], x -> zero );
       S:= EmptySCTable( n, zero, "antisymmetric" );
       for i in [1..n] do
         for j in [i+1..n] do
           cfs:= ShallowCopy( nullvec );
           cij:= T[i][j]; cji:= T[j][i];
           for m in [1..Length(cij[1])] do
             k:= cij[1][m];
             cfs[k]:= cfs[k] + cij[2][m];
           od;
           for m in [1..Length(cji[1])] do
             k:= cji[1][m];
             cfs[k]:= cfs[k] - cji[2][m];
           od;

           cij:= [ ];

           for m in [1..n] do
             if cfs[m] <> zero then
                 Add( cij, cfs[m] );
                 Add( cij, m );
             fi;
           od;
           SetEntrySCTable( S, i, j, cij );
         od;
       od;

       return LieAlgebraByStructureConstants( LeftActingDomain( A ), S );
  end );


#############################################################################
##
#F  LieAlgebra( <A> )
#F  LieAlgebra( <F>, <gens> )
#F  LieAlgebra( <F>, <gens>, <zero> )
#F  LieAlgebra( <F>, <gens>, "basis" )
#F  LieAlgebra( <F>, <gens>, <zero>, "basis" )
##
InstallGlobalFunction( LieAlgebra, function( arg )
#T check that the families have the same characteristic?
#T `CharacteristicFamily' ?
    local A,gens;

    # In the case of one domain argument,
    # construct the isomorphic Lie algebra.
    if Length( arg ) = 1 and IsDomain( arg[1] ) then
      A:= LieAlgebraByDomain( arg[1] );

    # division ring and list of generators
    elif Length( arg ) >= 2 and IsList( arg[2] ) then

       gens:= List( arg[2], x -> LieObject( x ) );

       if Length( arg ) = 2 and IsDivisionRing( arg[1] )
                            and 0 < Length( arg[2] ) then
         A:= AlgebraByGenerators( arg[1], gens );

       # division ring, list of generators plus zero
       elif Length( arg ) = 3 and IsDivisionRing( arg[1] ) then
          if arg[3] = "basis" then
            A:= AlgebraByGenerators( arg[1], gens );
            UseBasis( A, gens );
          else
            A:= AlgebraByGenerators( arg[1], gens, arg[3] );
          fi;

       # division ring, list of generators plus zero
       elif Length( arg ) = 4 and IsDivisionRing( arg[1] )
                              and arg[4] = "basis" then
          A:= AlgebraByGenerators( arg[1], gens, arg[3] );
          UseBasis( A, gens );

       else
         Error( "usage: LieAlgebra( <F>, <gens> ), ",
              "LieAlgebra( <F>, <gens>, <zero> ), LieAlgebra( <D> )");

       fi;
    # no argument given, error
    else
      Error( "usage: LieAlgebra( <F>, <gens> ), ",
             "LieAlgebra( <F>, <gens>, <zero> ), LieAlgebra( <D> )");
    fi;

    # Return the result.
    return A;
end );


#############################################################################
##
#F  EmptySCTable( <dim>, <zero> )
#F  EmptySCTable( <dim>, <zero>, \"symmetric\" )
#F  EmptySCTable( <dim>, <zero>, \"antisymmetric\" )
##
InstallGlobalFunction( EmptySCTable, function( arg )
    local dim, T, entry, i;

    if     2 <= Length( arg )
       and IsInt( arg[1] ) and IsZero( arg[2] )
       and ( Length( arg ) = 2 or IsString( arg[3] ) ) then

      dim:= arg[1];
      T:= [];
      entry:= Immutable( [ [], [] ] );
      for i in [ 1 .. dim ] do
        T[i]:= List( [ 1 .. dim ], x -> entry );
      od;

      # Store the symmetry flag.
      if Length( arg ) = 3 then
        if   arg[3] = "symmetric" then
          Add( T, 1 );
        elif arg[3] = "antisymmetric" then
          Add( T, -1 );
        else
          Error("third argument must be \"symmetric\" or \"antisymmetric\"");
        fi;
      else
        Add( T, 0 );
      fi;

      # Store the zero coefficient.
      Add( T, arg[2] );

    else
      Error( "usage: EmptySCTable( <dim>, <zero> [,\"symmetric\"] )" );
    fi;

    return T;
end );


#############################################################################
##
#F  SetEntrySCTable( <T>, <i>, <j>, <list> )
##
InstallGlobalFunction( SetEntrySCTable, function( T, i, j, list )
    local range, zero, Fam, entry, k, val, pos;

    # Check that `i' and `j' are admissible.
    range:= [ 1 .. Length( T ) - 2 ];
    if   not i in range then
      Error( "<i> must lie in ", range );
    elif not j in range then
      Error( "<j> must lie in ", range );
    fi;

    # Check `list', and construct the table entry.
    zero:= Last(T);
    Fam:= FamilyObj( zero );
    entry:= [ [], [] ];
    for k in [ 1, 3 .. Length( list ) -1 ] do

      val:= list[k];
      pos:= list[k+1];

      # Check that `pos' is inside the table,
      # and that its entry is assigned only once.
      if not pos in range then
        Error( "list entry ", list[k+1], " must lie in ", range );
      elif pos in entry[1] then
        Error( "position ", pos, " must occur at most once in <list>" );
      fi;

      # Check that the coefficients either fit to the zero element
      # or are rationals (with suitable denominators).
      if FamilyObj( val ) = Fam then
        if val <> zero then
          Add( entry[1], pos );
          Add( entry[2], val );
        fi;
      elif IsRat( val ) then
        if val <> 0 then
          Add( entry[1], pos );
          Add( entry[2], val * One( zero ) );
        fi;
      else
        Error( "list entry ", list[k], " does not fit to zero element" );
      fi;

    od;

    # Set the table entry.
    SortParallel( entry[1], entry[2] );
    T[i][j]:= Immutable( entry );

    # Add the value `T[j][i]' in the case of (anti-)symmetric tables.
    if   T[ Length(T) - 1 ] =  1 then
      T[j][i]:= T[i][j];
    elif T[ Length(T) - 1 ] = -1 then
      T[j][i]:= Immutable( [ entry[1], -entry[2] ] );
    fi;
end );


#############################################################################
##
#F  ReducedSCTable( <T>, <one> )
##
InstallGlobalFunction( ReducedSCTable, function( T, one )
    local new, n, i, j, entry;

    new:= [];
    n:= Length( T ) - 2;

    # Reduce the entries.
    for i in [ 1 .. n ] do
      new[i]:= [];
      for j in [ 1 .. n ] do
        entry:= T[i][j];
        entry:= [ Immutable( entry[1] ), entry[2] * one ];
        MakeImmutable( entry );
        new[i][j]:= entry;
      od;
    od;

    # Store zero coefficient and symmetry flag.
    new[ n+1 ]:= T[ n+1 ];
    new[ n+2 ]:= T[ n+2 ] * one;

    # Return the immutable new table.
    MakeImmutable( new );
    return new;
end );


#############################################################################
##
#F  GapInputSCTable( <T>, <varnam> )
##
InstallGlobalFunction( GapInputSCTable, function( T, varnam )
    local dim, str, lower, i, j, entry, k;

    # Initialize, and set the ranges for the loops.
    dim:= Length( T ) - 2;
    str:= Concatenation( varnam, ":= EmptySCTable( ",
                         String( dim ), ", ", String( Last(T) ) );
    lower:= [ 1 .. dim ];
    if   T[ dim+1 ] =  1 then
      Append( str, ", \"symmetric\"" );
    elif T[ dim+1 ] = -1 then
      Append( str, ", \"antisymmetric\"" );
    else
      lower:= ListWithIdenticalEntries( dim, 1 );
    fi;
    Append( str, " );\n" );

    # Fill up the table.
    for i in [ 1 .. dim ] do
      for j in [ lower[i] .. dim ] do
        entry:= T[i][j];
        if not IsEmpty( entry[1] ) then
          Append( str, "SetEntrySCTable( " );
          Append( str, varnam );
          Append( str, ", " );
          Append( str, String(i) );
          Append( str, ", " );
          Append( str, String(j) );
          Append( str, ", [" );
          for k in [ 1 .. Length( entry[1] )-1 ] do
            Append( str, String( entry[2][k] ) );
            Add( str, ',' );
            Append( str, String( entry[1][k] ) );
            Add( str, ',' );
          od;
          k:= Length( entry[1] );
          Append( str, String( entry[2][k] ) );
          Add( str, ',' );
          Append( str, String( entry[1][k] ) );
          Append( str, "] );\n" );
        fi;
      od;
    od;

    ConvertToStringRep( str );
    return str;
end );


#############################################################################
##
#F  IdentityFromSCTable( <T> )
##
InstallGlobalFunction( IdentityFromSCTable, function( T )
    local n,       # dimension of the underlying algebra
          equ,     # equation system to solve
          zero,    # zero of the field
          zerovec, # zero vector
          vec,     # right hand side of the equation system
          one,     # identity of the field
          i, j, k, # loop over rows of `equ'
          row,     # one row of the equation system
          Tpos,    #
          Tval,    #
          p,       #
          sol,
          sum;

    n:= Length( T ) - 2;
    zero:= T[ n+2 ];

    # If the table belongs to a trivial algebra,
    # the identity is equal to the zero.
    if n = 0 then
      return EmptyRowVector( FamilyObj( zero ) );
    fi;

    # Set up the equation system,
    # in row $i$ and column $(k-1)*n + j$ we have $c_{ijk}$.
    equ:= [];
    zerovec:= ListWithIdenticalEntries( n^2, zero );
    vec:= ShallowCopy( zerovec );
    one:= One( zero );

    for i in [ 1 .. n ] do
      row:= ShallowCopy( zerovec );
      for j in [ 1 .. n ] do
        Tpos:= T[i][j][1];
        Tval:= T[i][j][2];
        p:= (j-1)*n;
        for k in [ 1 .. Length( Tpos ) ] do
          row[ p + Tpos[k] ]:= Tval[k];
        od;
      od;
      Add( equ, row );
      vec[ (i-1)*n + i ]:= one;
    od;

    sol:= SolutionMat( equ, vec );

    # If we have a candidate and if the algebra is not known
    # to be commutative then check whether the candidate
    # acts trivially also from the right.
    if sol <> fail and T[ n+1 ] <> 1 then

      for j in [ 1 .. n ] do
        for k in [ 1 .. n ] do
          sum:= zero;
          for i in [ 1 .. n ] do
            Tpos:= T[j][i];
            p:= Position( Tpos[1], k );
#T cheaper !!!
            if p <> fail then
              sum:= sum + sol[i] * Tpos[2][p];
            fi;
          od;
          if ( j = k and sum <> one ) or ( j <> k and sum <> zero ) then
            return fail;
          fi;
        od;
      od;

    fi;

    # Return the result.
    return sol;
end );


#############################################################################
##
#F  QuotientFromSCTable( <T>, <num>, <den> )
##
##  We solve the equation system $<num> = x <den>$.
##  If no solution exists, `fail' is returned.
##
##  In terms of the basis $B$ with vectors $b_1, \ldots, b_n$ this means
##  for $<num> = \sum_{i=1}^n a_i b_i$,
##      $<den> = \sum_{i=1}^n c_i b_i$,
##      $x     = \sum_{i=1}^n x_i b_i$ that
##  $a_k = \sum_{i,j} c_i x_j c_{ijk}$ for all $k$.
##  Here $c_{ijk}$ denotes the structure constants w.r.t. $B$.
##  This means $a = x M$ with $M_{ik} = \sum_{j=1}^n c_{ijk} c_j$.
##
InstallGlobalFunction( QuotientFromSCTable, function( T, x, c )
    local M,        # matrix of the equation system
          n,        # dimension of the algebra
          zero,     # zero vector
          i, j,     # loop variables
          row,      # one row of `M'
          entry,    #
          val;      #

    M:= [];
    n:= Length( c );

    # If the algebra is zero dimensional,
    # the zero is also the identity and thus also its inverse.
    if n = 0 then
      return c;
    fi;

    zero:= ListWithIdenticalEntries( n, Last(T) );
    for i in [ 1 .. n ] do
      row:= ShallowCopy( zero );
      for j in [ 1 .. n ] do
        entry:= T[i][j];
        val:= c[j];
        row{ entry[1] }:= row{ entry[1] } + val * entry[2];
#T better!
      od;
      Add( M, row );
    od;

    # Return the quotient, or `fail'.
    return SolutionMat( M, x );
end );


#############################################################################
##
#F  TestJacobi( <T> )
##
##  We check whether for all $1 \leq m \leq n$ the equality
##  $\sum_{l=1}^n c_{jkl} c_{ilm} + c_{kil} c_{jlm} + c_{ijl} c_{klm} = 0$
##  holds.
##
InstallGlobalFunction( TestJacobi, function( T )
    local zero,           # the zero of the field
          n,              # dimension of the algebra
          i, j, k, m,     # loop variables
          cij, cki, cjk,  # structure constant vectors
          sum,
          t;

    zero:= Last(T);
    n:= Length( T ) - 2;

    for i in [ 1 .. n ] do
      for j in [ i+1 .. n ] do
        cij:= T[i][j];
        for k in [ j+1 .. n ] do
          cki:= T[k][i];
          cjk:= T[j][k];
          for m in [ 1 .. n ] do
            sum:= zero;
            for t in [ 1 .. Length( cjk[1] ) ] do
              sum:= sum + cjk[2][t] * SCTableEntry( T, i, cjk[1][t], m );
            od;
            for t in [ 1 .. Length( cki[1] ) ] do
              sum:= sum + cki[2][t] * SCTableEntry( T, j, cki[1][t], m );
            od;
            for t in [ 1 .. Length( cij[1] ) ] do
              sum:= sum + cij[2][t] * SCTableEntry( T, k, cij[1][t], m );
            od;
            if sum <> zero then
              return [ i, j, k ];
            fi;
          od;
        od;
      od;
    od;

    return true;
end );


#############################################################################
##
#M  MultiplicativeNeutralElement( <A> )
##
##  is the multiplicative neutral element of <A> if this exists,
##  otherwise is `fail'.
##
##  Let $(b_1, b_2, \ldots, b_n)$ be a basis of $A$, and $e$ the result of
##  `MultiplicativeNeutralElement( <A> )'.
##  Then $e = \sum_{i=1}^n a_i b_i$, and for $1 \leq k \leq n$ we have
##  $e \cdot b_j = b_j$, or equivalently
##  $\sum_{i=1}^n a_i b_i \cdot b_j = b_j$.
##  Define the structure constants by
##  $b_i \cdot b_j = \sum_{k=1}^n c_{ijk} b_k$.
##
##  Then $\sum_{i=1}^n a_i c_{ijk} = \delta_{jk}$ for $1 \leq k \leq n$.
##
##  This yields $n^2$ linear equations for the $n$ indeterminates $a_i$,
##  and a solution is a left identity.
##  For this we have to test whether it is also a right identity.
##
InstallMethod( MultiplicativeNeutralElement,
    [ IsFLMLOR and IsFiniteDimensional ],
    function( A )
    local B,       # basis of `A'
          one;     # result

    B:= Basis( A );
    one:= IdentityFromSCTable( StructureConstantsTable( B ) );
    if one <> fail then
      one:= LinearCombination( B, one );
    fi;
    return one;
    end );


#############################################################################
##
#M  IsAssociative( <A> )
##
##  We check whether the vectors of a basis satisfy the associativity law.
##  (Bilinearity of the multiplication is of course assumed.)
##
##  If $b_i \cdot b_j = \sum_{l=1}^n c_{ijl} b_l$ then we have
##  $b_i \cdot ( b_j \cdot b_k ) = ( b_i \cdot b_j ) \cdot b_k$
##  if and only if
##  $\sum_{l=1}^n c_{jkl} c_{ilm} = \sum_{l=1}^n c_{ijl} c_{lkm}$ for all
##  $1 \leq m \leq n$.
##
##  We check this equality for all $1 \leq i, j, k \leq n$.
##
InstallMethod( IsAssociative,
    "generic method for a (finite dimensional) FLMLOR",
    [ IsFLMLOR ],
    function( A )
    local T,            # structure constants table w.r.t. a basis of `A'
          zero,
          range,
          i, j, k, l, m,
          Ti,
          Tj,
          cijpos,
          cijval,
          cjkpos,
          cjkval,
          sum,
          x,
          pos;

    if not IsFiniteDimensional( A ) then
      TryNextMethod();
    fi;

    T:= StructureConstantsTable( Basis( A ) );
    zero:= Zero( LeftActingDomain( A ) );
    range:= [ 1 .. Length( T[1] ) ];
    for i in range do
      Ti:= T[i];
      for j in range do
        cijpos:= Ti[j][1];
        cijval:= Ti[j][2];
        Tj:= T[j];
        for k in range do
          cjkpos:= Tj[k][1];
          cjkval:= Tj[k][2];
          for m in range do
            sum:= zero;
            for l in [ 1 .. Length( cjkpos ) ] do
              x:= Ti[ cjkpos[l] ];
              pos:= Position( x[1], m );
              if pos <> fail then
                sum:= sum + cjkval[l] * x[2][ pos ];
              fi;
            od;
            for l in [ 1 .. Length( cijpos ) ] do

              x:= T[ cijpos[l] ][k];
              pos:= Position( x[1], m );
              if pos <> fail then
                sum:= sum - cijval[l] * x[2][ pos ];
              fi;
            od;
            if sum <> zero then
              # $i, j, k$ fail
              Info( InfoAlgebra, 2,
                    "IsAssociative fails for i = ", i, ", j = ", j,
                    ", k = ", k );
              return false;
            fi;
          od;
        od;
      od;
    od;
    return true;
    end );


#############################################################################
##
#M  IsAnticommutative( <A> )  . . . . . . . . . . . . .for a fin.-dim. FLMLOR
##
##  is `true' if the multiplication in <A> is anticommutative,
##  and `false' otherwise.
##
InstallMethod( IsAnticommutative,
    "generic method for a (finite dimensional) FLMLOR",
    [ IsFLMLOR ],
    function( A )
    local n,      # dimension of `A'
          T,      # table of structure constants for `A'
          i, j;   # loop over rows and columns ot `T'

    if not IsFiniteDimensional( A ) then
      TryNextMethod();
    fi;

    n:= Dimension( A );
    T:= StructureConstantsTable( Basis( A ) );
    for i in [ 2 .. n ] do
      for j in [ 1 .. i-1 ] do
        if    T[i][j][1] <> T[j][i][1]
           or ( not IsEmpty( T[i][j][1] )
                and PositionNonZero( T[i][j][2] + T[j][i][2] )
                        <= Length( T[i][j][2] ) ) then
          return false;
        fi;
      od;
    od;

    if Characteristic( A ) <> 2 then

      # The values on the diagonal must be zero.
      for i in [ 1 .. n ] do
        if not IsEmpty( T[i][i][1] ) then
          return false;
        fi;
      od;

    fi;

    return true;
    end );


#############################################################################
##
#M  IsCommutative( <A> )  . . . . . . . . . . . for finite dimensional FLMLOR
##
##  Check whether every basis vector commutes with every basis vector.
##
InstallMethod( IsCommutative,
    "generic method for a finite dimensional FLMLOR",
    [ IsFLMLOR ],
    IsCommutativeFromGenerators( GeneratorsOfVectorSpace ) );
#T use structure constants!


#############################################################################
##
#M  IsCommutative( <A> )  . . . . . . . . . . . . . for an associative FLMLOR
##
##  If <A> is associative then we can restrict the check to a smaller
##  equation system than that for arbitrary algebras, since we have to check
##  $x a = a x$ only for algebra generators $a$ and $x$, not for all vectors
##  of a basis.
##
InstallMethod( IsCommutative,
    "for an associative FLMLOR",
    [ IsFLMLOR and IsAssociative ],
    IsCommutativeFromGenerators( GeneratorsOfAlgebra ) );

InstallMethod( IsCommutative,
    "for an associative FLMLOR-with-one",
    [ IsFLMLORWithOne and IsAssociative ],
    IsCommutativeFromGenerators( GeneratorsOfAlgebraWithOne ) );


#############################################################################
##
#M  IsZeroSquaredRing( <A> )  . . . . . . . . for a finite dimensional FLMLOR
##
InstallMethod( IsZeroSquaredRing,
    "for a finite dimensional FLMLOR",
    [ IsFLMLOR ],
    function( A )

    if not IsAnticommutative( A ) then

      # Every zero squared ring is anticommutative.
      return false;

    elif ForAny( BasisVectors( Basis( A ) ),
                 x -> not IsZero( x*x ) ) then

      # If not all basis vectors are zero squared then we return `false'.
      return false;

    elif IsCommutative( LeftActingDomain( A ) ) then

      # If otherwise the left acting domain is commutative then we return
      # `true' because we know that <A> is anticommutative and the basis
      # vectors are zero squared.
      return true;

    else

      # Otherwise we give up.
      TryNextMethod();

    fi;
    end );


#############################################################################
##
#M  IsJacobianRing( <A> )
##
InstallMethod( IsJacobianRing,
    "for a (finite dimensional) FLMLOR",
    [ IsFLMLOR ],
    function( A )
    local n,   # dimension of `A'
          T,   # table of structure constants for `A'
          i;   # loop over the diagonal of `T'

    if not IsFiniteDimensional( A ) then
      TryNextMethod();
    fi;

    # In characteristic 2 we have to make sure that $a \* a = 0$.
#T really?
    T:= StructureConstantsTable( Basis( A ) );
    if Characteristic( A ) = 2 then
      n:= Dimension( A );
      for i in [ 1 .. n ] do
        if not IsEmpty( T[i][i][1] ) then
          return false;
        fi;
      od;
    fi;

    # Check the Jacobi identity $[a,[b,c]] + [b,[c,a]] + [c,[a,b]] = 0$.
    return TestJacobi( T ) = true;
    end );


#############################################################################
##
#M  Intersection2( <A1>, <A2> ) . . . . . . . . . intersection of two FLMLORs
##
InstallMethod( Intersection2,
    "generic method for two FLMLORs",
    IsIdenticalObj,
    [ IsFLMLOR, IsFLMLOR ],
    Intersection2Spaces( AsFLMLOR, SubFLMLORNC, FLMLOR ) );


#############################################################################
##
#M  Intersection2( <A1>, <A2> ) . . . .  intersection of two FLMLORs-with-one
##
InstallMethod( Intersection2,
    "generic method for two FLMLORs-with-one",
    IsIdenticalObj,
    [ IsFLMLORWithOne, IsFLMLORWithOne ],
    Intersection2Spaces( AsFLMLORWithOne, SubFLMLORWithOneNC,
                         FLMLORWithOne ) );


#############################################################################
##
#M  \/( <A>, <I> )  . . . . . . . . . . . .  factor of an algebra by an ideal
#M  \/( <A>, <relators> ) . . . . . . . . .  factor of an algebra by an ideal
##
##  is the factor algebra of the finite dimensional algebra <A> modulo
##  the ideal <I> or the ideal spanned by the collection <relators>.
##
InstallOtherMethod( \/,
    "for FLMLOR and collection",
    IsIdenticalObj,
    [ IsFLMLOR, IsCollection ],
    function( A, relators )
    if IsFLMLOR( relators ) then
      TryNextMethod();
    else
      return A / TwoSidedIdealByGenerators( A, relators );
    fi;
    end );

InstallOtherMethod( \/,
    "for FLMLOR and empty list",
    [ IsFLMLOR, IsList and IsEmpty ],
    function( A, empty )
    # `NaturalHomomorphismByIdeal( A, TrivialSubFLMLOR( A ) )' is the
    # identity mapping on `A', and `ImagesSource' of it yields `A'.
    return A;
    end );

InstallOtherMethod( \/,
    "generic method for two FLMLORs",
    IsIdenticalObj,
    [ IsFLMLOR, IsFLMLOR ],
    function( A, I )
    return ImagesSource( NaturalHomomorphismByIdeal( A, I ) );
    end );


#############################################################################
##
#M  TrivialSubadditiveMagmaWithZero( <A> )  . . . . . . . . . .  for a FLMLOR
##
InstallMethod( TrivialSubadditiveMagmaWithZero,
    "for a FLMLOR",
    [ IsFLMLOR ],
    A -> SubFLMLORNC( A, [] ) );


#############################################################################
##
#M  AsFLMLOR( <R>, <D> )  . . view a collection as a FLMLOR over the ring <R>
##
InstallMethod( AsFLMLOR,
    "for a ring and a collection",
    [ IsRing, IsCollection ],
    function( F, D )
    local A, L;

    D:= AsSSortedList( D );
    L:= ShallowCopy( D );
    A:= TrivialSubFLMLOR( AsFLMLOR( F, D ) );
    SubtractSet( L, AsSSortedList( A ) );
    while 0 < Length(L)  do
      A:= ClosureLeftOperatorRing( A, L[1] );
#T call explicit function that maintains an elements list?
      SubtractSet( L, AsSSortedList( A ) );
    od;
    if Length( AsList( A ) ) <> Length( D )  then
      return fail;
    fi;
    A:= FLMLOR( F, GeneratorsOfLeftOperatorRing( A ), Zero( D[1] ) );
    SetAsSSortedList( A, D );
    SetSize( A, Length( D ) );
    SetIsFinite( A, true );
#T ?

    # Return the FLMLOR.
    return A;
    end );


#############################################################################
##
#M  AsFLMLOR( <F>, <V> )  . . view a left module as FLMLOR over the field <F>
##
##  is an algebra over <F> that is equal (as set) to <V>.
##  For that, perhaps the field of <A> has to be changed before
##  getting the correct list of generators.
##
InstallMethod( AsFLMLOR,
    "for a division ring and a free left module",
    [ IsDivisionRing, IsFreeLeftModule ],
    function( F, V )
    local L, A;

    if   LeftActingDomain( V ) = F then

      A:= FLMLOR( F, GeneratorsOfLeftModule( V ) );
      if A <> V then
        return fail;
      fi;
      if HasBasis( V ) then
        SetBasis( A, Basis( V ) );
      fi;

    elif IsTrivial( V ) then

      # We need the zero.
      A:= FLMLOR( F, [], Zero( V ) );

    elif IsSubset( LeftActingDomain( V ), F ) then

      # Make sure that the field change does not change the elements.
      L:= BasisVectors( Basis( AsField( F, LeftActingDomain(V) ) ) );
      L:= Concatenation( List( L, x -> List( GeneratorsOfLeftModule( V ),
                                             y -> x * y ) ) );
      A:= FLMLOR( F, L );
      if A <> V then
        return fail;
      fi;

    elif IsSubset( F, LeftActingDomain( V ) ) then

      # Make sure that the field change does not change the elements.
      L:= BasisVectors( Basis( AsField( LeftActingDomain(V), F ) ) );
      if ForAny( L, x -> ForAny( GeneratorsOfLeftModule( V ),
                                 y -> not x * y in V ) ) then
        return fail;
      fi;
      A:= FLMLOR( F, GeneratorsOfLeftModule( V ) );
      if A <> V then
        return fail;
      fi;

    else

      V:= AsFLMLOR( Intersection( F, LeftActingDomain( V ) ), V );
      return AsFLMLOR( F, V );

    fi;

    UseIsomorphismRelation( V, A );
    UseSubsetRelation( V, A );

    return A;
    end );


#############################################################################
##
#M  AsFLMLOR( <F>, <A> )  . . . view an algebra as algebra over the field <F>
##
##  is an algebra over <F> that is equal (as set) to <D>.
##  For that, perhaps the field of <A> has to be changed before
##  getting the correct list of generators.
##
InstallMethod( AsFLMLOR,
    "for a division ring and an algebra",
    [ IsDivisionRing, IsFLMLOR ],
    function( F, D )
    local L, A;

    if   LeftActingDomain( D ) = F then

      return D;

    elif IsTrivial( D ) then

      # We need the zero.
      A:= FLMLOR( F, [], Zero( D ) );

    elif IsSubset( LeftActingDomain( D ), F ) then

      # Make sure that the field change does not change the elements.
      L:= BasisVectors( Basis( AsField( F, LeftActingDomain(D) ) ) );
      L:= Concatenation( List( L, x -> List( GeneratorsOfAlgebra( D ),
                                             y -> x * y ) ) );
      A:= FLMLOR( F, L );

    elif IsSubset( F, LeftActingDomain( D ) ) then

      # Make sure that the field change does not change the elements.
      L:= BasisVectors( Basis( AsField( LeftActingDomain(D), F ) ) );
      if ForAny( L, x -> ForAny( GeneratorsOfAlgebra( D ),
                                 y -> not x * y in D ) ) then
        return fail;
      fi;
      A:= FLMLOR( F, GeneratorsOfAlgebra( D ) );

    else

      D:= AsFLMLOR( Intersection( F, LeftActingDomain( D ) ), D );
      return AsFLMLOR( F, D );

    fi;

    UseIsomorphismRelation( D, A );
    UseSubsetRelation( D, A );
    return A;
    end );


#############################################################################
##
#M  AsFLMLORWithOne( <R>, <D> ) . . . . . . view a coll. as a FLMLOR-with-one
##
InstallMethod( AsFLMLORWithOne,
    "for a ring and a collection",
    [ IsRing, IsCollection ],
    function( F, D )
    return AsFLMLORWithOne( AsFLMLOR( F, D ) );
    end );


#############################################################################
##
#M  AsFLMLORWithOne( <F>, <V> ) . . view a left module as an algebra-with-one
##
InstallMethod( AsFLMLORWithOne,
    "for a division ring and a free left module",
    [ IsDivisionRing, IsFreeLeftModule ],
    function( F, V )
    local L, A;

    # Check that `V' contains the identity.
    if One( V ) = fail then

      return fail;

    elif LeftActingDomain( V ) = F then

      A:= FLMLORWithOne( F, GeneratorsOfLeftModule( V ) );
      if A <> V then
        return fail;
      fi;

      # Left module generators and basis are maintained.
      if HasGeneratorsOfLeftModule( V ) then
        SetGeneratorsOfLeftModule( A, GeneratorsOfLeftModule( V ) );
      fi;
      if HasBasis( V ) then
        SetBasis( A, Basis( V ) );
      fi;

    elif IsSubset( LeftActingDomain( V ), F ) then

      # Make sure that the field change does not change the elements.
      L:= BasisVectors( Basis( AsField( F, LeftActingDomain(V) ) ) );
      L:= Concatenation( List( L, x -> List( GeneratorsOfLeftModule( V ),
                                             y -> x * y ) ) );
      A:= FLMLORWithOne( F, L );
      if A <> V then
        return fail;
      fi;

    elif IsSubset( F, LeftActingDomain( V ) ) then

      # Make sure that the field change does not change the elements.
      L:= BasisVectors( Basis( AsField( LeftActingDomain(V), F ) ) );
      if ForAny( L, x -> ForAny( GeneratorsOfLeftModule( V ),
                                 y -> not x * y in V ) ) then
        return fail;
      fi;
      A:= FLMLORWithOne( F, GeneratorsOfLeftModule( V ) );
      if A <> V then
        return fail;
      fi;

    else

      # Note that we need not use the isomorphism and subset relations
      # (see below) because this is the task of the calls to
      # `AsAlgebraWithOne'.
      A:= AsAlgebraWithOne( Intersection( F, LeftActingDomain( V ) ), V );
      return AsAlgebraWithOne( F, A );

    fi;

    UseIsomorphismRelation( V, A );
    UseSubsetRelation( V, A );
    return A;
    end );


#############################################################################
##
#M  AsFLMLORWithOne( <F>, <D> ) . . .  view an algebra as an algebra-with-one
##
InstallMethod( AsFLMLORWithOne,
    "for a division ring and an algebra",
    [ IsDivisionRing, IsFLMLOR ],
    function( F, D )
    local L, A;

    # Check that `D' contains the identity.
    if One( D ) = fail then

      return fail;

    elif LeftActingDomain( D ) = F then

      A:= FLMLORWithOne( F, GeneratorsOfLeftOperatorRing( D ) );

      # Left module generators and basis are maintained.
      if HasGeneratorsOfLeftModule( D ) then
        SetGeneratorsOfLeftModule( A, GeneratorsOfLeftModule( D ) );
      fi;
      if HasBasis( D ) then
        SetBasis( A, Basis( D ) );
      fi;

    elif IsSubset( LeftActingDomain( D ), F ) then

      # Make sure that the field change does not change the elements.
      L:= BasisVectors( Basis( AsField( F, LeftActingDomain(D) ) ) );
      L:= Concatenation( List( L, x -> List( GeneratorsOfAlgebra( D ),
                                             y -> x * y ) ) );
      A:= FLMLORWithOne( F, L );

    elif IsSubset( F, LeftActingDomain( D ) ) then

      # Make sure that the field change does not change the elements.
      L:= BasisVectors( Basis( AsField( LeftActingDomain(D), F ) ) );
      if ForAny( L, x -> ForAny( GeneratorsOfAlgebra( D ),
                                 y -> not x * y in D ) ) then
        return fail;
      fi;
      A:= FLMLORWithOne( F, GeneratorsOfLeftOperatorRing( D ) );

    else

      # Note that we need not use the isomorphism and subset relations
      # (see below) because this is the task of the calls to
      # `AsAlgebraWithOne'.
      D:= AsAlgebraWithOne( Intersection( F, LeftActingDomain( D ) ), D );
      return AsAlgebraWithOne( F, D );

    fi;

    UseIsomorphismRelation( D, A );
    UseSubsetRelation( D, A );
    return A;
    end );


#############################################################################
##
#M  AsFLMLORWithOne( <F>, <D> ) . . view an alg.-with-one as an alg.-with-one
##
InstallMethod( AsFLMLORWithOne,
    "for a division ring and an algebra-with-one",
    [ IsDivisionRing, IsFLMLORWithOne ],
    function( F, D )
    local L, A;

    if   LeftActingDomain( D ) = F then

      return D;

    elif IsSubset( LeftActingDomain( D ), F ) then

      # Make sure that the field change does not change the elements.
      L:= BasisVectors( Basis( AsField( F, LeftActingDomain(D) ) ) );
      L:= Concatenation( List( L, x -> List( GeneratorsOfAlgebra( D ),
                                             y -> x * y ) ) );
      A:= AlgebraWithOne( F, L );

    elif IsSubset( F, LeftActingDomain( D ) ) then

      # Make sure that the field change does not change the elements.
      L:= BasisVectors( Basis( AsField( LeftActingDomain(D), F ) ) );
      if ForAny( L, x -> ForAny( GeneratorsOfAlgebra( D ),
                                 y -> not x * y in D ) ) then
        return fail;
      fi;
      A:= AlgebraWithOne( F, GeneratorsOfAlgebra( D ) );

    else

      # Note that we need not use the isomorphism and subset relations
      # (see below) because this is the task of the calls to
      # `AsAlgebraWithOne'.
      D:= AsAlgebraWithOne( Intersection( F, LeftActingDomain( D ) ), D );
      return AsAlgebraWithOne( F, D );

    fi;

    UseIsomorphismRelation( D, A );
    UseSubsetRelation( D, A );
    return A;
    end );


#############################################################################
##
#M  ClosureLeftOperatorRing( <A>, <a> ) . . . . . . . closure with an element
##
InstallMethod( ClosureLeftOperatorRing,
    "for a FLMLOR and a ring element",
#T why not a general function for any left operator ring?
#T (need `LeftOperatorRingByGenerators' ?)
    IsCollsElms,
    [ IsFLMLOR, IsRingElement ],
    function( A, a )

    # if possible test if the element lies in the ring already,
    if     HasGeneratorsOfLeftOperatorRing( A )
       and a in GeneratorsOfLeftOperatorRing( A ) then
      return A;

    # otherwise make a new left operator ring
    else
      return FLMLOR( LeftActingDomain( A ),
                 Concatenation( GeneratorsOfLeftOperatorRing( A ), [ a ] ) );
    fi;
    end );

InstallMethod( ClosureLeftOperatorRing,
    "for an FLMLOR with basis, and a ring element",
    IsCollsElms,
    [ IsFLMLOR and HasBasis, IsRingElement ],
    function( A, a )

    # test if the element lies in the FLMLOR already,
    if a in A then
      return A;

    # otherwise make a new FLMLOR
    else
      return FLMLOR( LeftActingDomain( A ),
#T FLMLORByGenerators?
                 Concatenation( BasisVectors( Basis( A ) ), [ a ] ),
                 "basis" );
    fi;
    end );

InstallMethod( ClosureLeftOperatorRing,
    "for a FLMLOR-with-one and a ring element",
    IsCollsElms,
    [ IsFLMLORWithOne, IsRingElement ],
    function( A, a )

    # if possible test if the element lies in the ring already,
    if a in GeneratorsOfLeftOperatorRingWithOne( A ) then
      return A;

    # otherwise make a new left operator ring-with-one
    else
      return FLMLORWithOne( LeftActingDomain( A ),
                 Concatenation( GeneratorsOfLeftOperatorRingWithOne( A ),
                                [ a ] ) );
    fi;
    end );

InstallMethod( ClosureLeftOperatorRing,
    "for a FLMLOR-with-one with basis, and a ring element",
    IsCollsElms,
    [ IsFLMLORWithOne and HasBasis, IsRingElement ],
    function( A, a )

    # test if the element lies in the FLMLOR already,
    if a in A then
        return A;

    # otherwise make a new FLMLOR-with-one
    else
      return FLMLORWithOne( LeftActingDomain( A ),
                 Concatenation( BasisVectors( Basis( A ) ), [ a ] ),
                 "basis" );
    fi;
    end );

InstallMethod( ClosureLeftOperatorRing,
    "for a FLMLOR containing the whole family, and a ring element",
    IsCollsElms,
    [ IsFLMLOR and IsWholeFamily, IsRingElement ],
    SUM_FLAGS, # this is better than everything else
    ReturnFirst);

#############################################################################
##
#M  ClosureLeftOperatorRing( <A>, <U> ) .  closure of two left operator rings
##
InstallMethod( ClosureLeftOperatorRing,
    "for two left operator rings",
    IsIdenticalObj,
    [ IsLeftOperatorRing, IsLeftOperatorRing ],
    function( A, S )
    local   g;          # one generator

    for g in GeneratorsOfLeftOperatorRing( S ) do
      A := ClosureLeftOperatorRing( A, g );
    od;
    return A;
    end );

InstallMethod( ClosureLeftOperatorRing,
    "for two left operator rings-with-one",
    IsIdenticalObj,
    [ IsLeftOperatorRingWithOne, IsLeftOperatorRingWithOne ],
    function( A, S )
    local   g;          # one generator

    for g in GeneratorsOfLeftOperatorRingWithOne( S ) do
      A := ClosureLeftOperatorRing( A, g );
    od;
    return A;
    end );

InstallMethod( ClosureLeftOperatorRing,
    "for a left op. ring cont. the whole family, and a collection",
    IsIdenticalObj,
    [ IsLeftOperatorRing and IsWholeFamily, IsCollection ],
    SUM_FLAGS, # this is better than everything else
    ReturnFirst);

#############################################################################
##
#M  ClosureLeftOperatorRing( <A>, <list> )  . . . .  closure of left op. ring
##
InstallMethod( ClosureLeftOperatorRing,
    "for left operator ring and list of elements",
    IsIdenticalObj,
    [ IsLeftOperatorRing, IsCollection ],
    function( A, list )
    local   g;          # one generator
    for g in list do
      A:= ClosureLeftOperatorRing( A, g );
    od;
    return A;
    end );


#############################################################################
##
#F  MutableBasisOfClosureUnderAction( <F>, <Agens>, <from>, <init>, <opr>,
#F                                    <zero>, <maxdim> )
##
##  This function is used to compute bases of finite dimensional ideals $I$
##  in *associative* algebras that are given by ideal generators $J$ and
##  (generators of) the acting algebra $A$.
##  An important special case is that of the algebra $A$ itself, given by
##  algebra generators.
##
##  The algorithm assumes that it is possible to deal with mutable bases of
##  vector spaces generated by elements of $A$.
##  It proceeds as follows.
##
##  Let $A$ be a finite dimensional algebra over the ring $F$,
##  and $I$ a two-sided ideal in $A$ that is generated (as a two-sided ideal)
##  by the set $J$.
##  (For the cases of one-sided ideals, see below.)
##
##  Let $S$ be a set of algebra generators of $A$.
##  The identity of $A$, if exists, need not be contained in $S$.
##
##  Each element $x$ in $I$ can be written as a linear combination of
##  products $j a_1 a_2 \cdots a_n$, with $j \in J$ and $a_i \in S$ for
##  $1 \leq i \leq n$.
##  Let $l(x)$ denote the minimum of the largest $n$ for involved words
##  $a_1 a_2 \cdots a_n$, taken over all possible expressions for $x$.
##
##  Define $I_i = \{ x \in I \mid l(x) \leq i \}$.
##  Then $I_i$ is an $F$-space, $A_0 = \langle J \rangle_F$,
##  and $I_0 \< I_1 \< I_2 \< \cdots I_k \< I_{k+1} \< \ldots$
##  is an ascending chain that eventually reaches $I$.
##  For $i > 0$ we have
##  $I_{i+1} = \langle I_i\cup\bigcup_{s\in S} ( I_i s\cup s I_i )\rangle_F$.
##
##  (*Note* that the computation of the $I_i$ gives us the smallest value $k$
##  such that every element is a linear combination of words in terms of the
##  algebra generators, of maximal length $k$.)
##
InstallGlobalFunction( MutableBasisOfClosureUnderAction,
    function( F, Agens, from, init, opr, zero, maxdim )
    local MB,        # mutable basis, result
          gen,       # loop over generators
          v,         #
          dim,       # dimension of the actual left module
          right,     # `true' if we have to multiply from the right
          left;      # `true' if we have to multiply from the left

    # Get the side(s) from where to multiply.
    left  := true;
    right := true;
    if   from = "left" then
      right:= false;
    elif from = "right" then
      left:= false;
    fi;

    # $I_0$
    MB  := MutableBasis( F, init, zero );
    dim := 0;

    while dim < NrBasisVectors( MB ) and dim < maxdim do

      # `MB' is a mutable basis of $I_i$.
      dim:= NrBasisVectors( MB );

      if right then

        # Compute $I^{\prime}_i = I_i + \sum_{s \in S} I_i s$.
        for gen in Agens do
          for v in BasisVectors( MB ) do
            CloseMutableBasis( MB, opr( v, gen ) );
          od;
        od;

      fi;
      if left then

        # Compute $I_i + \sum_{s \in S} s I_i$
        # resp. $I^{\prime}_i + \sum_{s \in S} s I_i$.
        for gen in Agens do
          for v in BasisVectors( MB ) do
            CloseMutableBasis( MB, opr( gen, v ) );
          od;
        od;

      fi;

    od;

    # Return the mutable basis.
    return MB;
end );


#############################################################################
##
#F  MutableBasisOfNonassociativeAlgebra( <F>, <Agens>, <zero>, <maxdim> )
##
InstallGlobalFunction( MutableBasisOfNonassociativeAlgebra,
    function( F, Agens, zero, maxdim )
    local MB,        # mutable basis, result
          dim,       # dimension of the current left module
          bv,        # current basis vectors
          v, w;      # loop over basis vectors found already

    MB  := MutableBasis( F, Agens, zero );
    dim := 0;

    while dim < NrBasisVectors( MB ) and dim < maxdim do

      dim := NrBasisVectors( MB );
      bv  := BasisVectors( MB );

      for v in bv do
        for w in bv do
          CloseMutableBasis( MB, v * w );
          CloseMutableBasis( MB, w * v );
        od;
      od;

    od;

    # Return the mutable basis.
    return MB;
end );


#############################################################################
##
#F  MutableBasisOfIdealInNonassociativeAlgebra( <F>, <Vgens>, <Igens>,
#F                                              <zero>, <from>, <maxdim> )
##
InstallGlobalFunction( MutableBasisOfIdealInNonassociativeAlgebra,
    function( F, Vgens, Igens, zero, from, maxdim )
    local MB,        # mutable basis, result
          dim,       # dimension of the current left module
          bv,        # current basis vectors
          right,     # `true' if we have to multiply from the right
          left,      # `true' if we have to multiply from the left
          v, gen;    # loop over basis vectors found already

    # Get the side(s) from where to multiply.
    left  := true;
    right := true;
    if   from = "left" then
      right:= false;
    elif from = "right" then
      left:= false;
    fi;

    dim := 0;
    MB  := MutableBasis( F, Igens, zero );

    while dim < NrBasisVectors( MB ) and dim < maxdim do

      dim := NrBasisVectors( MB );
      bv  := BasisVectors( MB );

      for v in bv do
        for gen in Vgens do

          if left then
            CloseMutableBasis( MB, gen * v );
          fi;
          if right then
            CloseMutableBasis( MB, v * gen );
          fi;
        od;
      od;

    od;

    # Return the mutable basis.
    return MB;
end );


#############################################################################
##
#M  IsSubset( <A>, <B> )  . . . . . . . . . . . .  test for subset of FLMLORs
##
##  These methods are preferable to that for free left modules because they
##  use algebra generators.
##
##  We assume that generators of an extension of the left acting domains can
##  be computed if they are fields; note that infinite field extensions do
##  not (yet) occur in {\GAP} as `LeftActingDomain' values.
##
InstallMethod( IsSubset,
    "for two FLMLORs",
    IsIdenticalObj,
    [ IsFLMLOR, IsFLMLOR ],
    function( D1, D2 )
    local F1, F2;

    F1:= LeftActingDomain( D1 );
    F2:= LeftActingDomain( D2 );
    if not ( HasIsDivisionRing( F1 ) and IsDivisionRing( F1 ) and
             HasIsDivisionRing( F2 ) and IsDivisionRing( F2 ) ) then
      TryNextMethod();
    fi;

    # catch trivial case
    if IsSubset(GeneratorsOfLeftOperatorRing(D1),
                  GeneratorsOfLeftOperatorRing(D2)) then
        return true;
    fi;
    return IsSubset( D1, GeneratorsOverIntersection( D2,
                             GeneratorsOfLeftOperatorRing( D2 ),
                             F2, F1 ) );
    end );

InstallMethod( IsSubset,
    "for two FLMLORs-with-one",
    IsIdenticalObj,
    [ IsFLMLORWithOne, IsFLMLORWithOne ],
    function( D1, D2 )
    local F1, F2;

    F1:= LeftActingDomain( D1 );
    F2:= LeftActingDomain( D2 );
    if not ( HasIsDivisionRing( F1 ) and IsDivisionRing( F1 ) and
             HasIsDivisionRing( F2 ) and IsDivisionRing( F2 ) ) then
      TryNextMethod();
    fi;
    return IsSubset( D1, GeneratorsOverIntersection( D2,
                             GeneratorsOfLeftOperatorRingWithOne( D2 ),
                             F2, F1 ) );
    end );


#############################################################################
##
#M  ViewObj( <A> )  . . . . . . . . . . . . . . . . . . . . . . view a FLMLOR
##
##  print left acting domain, if known also dimension or no. of generators
##
InstallMethod( ViewObj,
    "for a FLMLOR",
    [ IsFLMLOR ],
    function( A )
    Print( "<free left module over ", LeftActingDomain( A ), ", and ring>" );
    end );

InstallMethod( ViewObj,
    "for a FLMLOR with known dimension",
    [ IsFLMLOR and HasDimension ], 1,   # override method requiring gens.
    function( A )
    Print( "<free left module of dimension ", Dimension( A ),
           " over ", LeftActingDomain( A ), ", and ring>" );
    end );

InstallMethod( ViewObj,
    "for a FLMLOR with known generators",
    [ IsFLMLOR and HasGeneratorsOfAlgebra ],
    function( A )
    Print( "<free left module over ", LeftActingDomain( A ),
           ", and ring, with ",
           Pluralize( Length( GeneratorsOfFLMLOR( A ) ), "generator" ), ">" );
    end );


#############################################################################
##
#M  PrintObj( <A> ) . . . . . . . . . . . . . . . . . . . . .  print a FLMLOR
##
InstallMethod( PrintObj,
    "for a FLMLOR",
    [ IsFLMLOR ],
    function( A )
    Print( "FLMLOR( ", LeftActingDomain( A ), ", ... )" );
    end );

InstallMethod( PrintObj,
    "for a FLMLOR with known generators",
    [ IsFLMLOR and HasGeneratorsOfFLMLOR ],
    function( A )
    if IsEmpty( GeneratorsOfFLMLOR( A ) ) then
      Print( "FLMLOR( ", LeftActingDomain( A ), ", [], ", Zero( A ), " )" );
    else
      Print( "FLMLOR( ", LeftActingDomain( A ), ", ",
             GeneratorsOfFLMLOR( A ), " )" );
    fi;
    end );


#############################################################################
##
#M  ViewObj( <A> )  . . . . . . . . . . . . . . . . .  view a FLMLOR-with-one
##
##  print left acting domain, if known also dimension or no. of generators
##
InstallMethod( ViewObj,
    "for a FLMLOR-with-one",
    [ IsFLMLORWithOne ],
    function( A )
    Print( "<free left module over ", LeftActingDomain( A ),
           ", and ring-with-one>" );
    end );

InstallMethod( ViewObj,
    "for a FLMLOR-with-one with known dimension",
    [ IsFLMLORWithOne and HasDimension ], 1,   # override method requ. gens.
    function( A )
    Print( "<free left module of dimension ", Dimension( A ),
           " over ", LeftActingDomain( A ), ", and ring-with-one>" );
    end );

InstallMethod( ViewObj,
    "for a FLMLOR-with-one with known generators",
    [ IsFLMLORWithOne and HasGeneratorsOfFLMLORWithOne ],
    function( A )
    Print( "<free left module over ", LeftActingDomain( A ),
           ", and ring-with-one, with ",
           Pluralize( Length( GeneratorsOfAlgebraWithOne( A ) ), "generator" ),
           ">" );
    end );


#############################################################################
##
#M  PrintObj( <A> ) . . . . . . . . . . . . . . . . . print a FLMLOR-with-one
##
InstallMethod( PrintObj,
    "for a FLMLOR-with-one",
    [ IsFLMLORWithOne ],
    function( A )
    Print( "FLMLORWithOne( ", LeftActingDomain( A ), ", ... )" );
    end );

InstallMethod( PrintObj,
    "for a FLMLOR-with-one with known generators",
    [ IsFLMLORWithOne and HasGeneratorsOfFLMLOR ],
    function( A )
    if IsEmpty( GeneratorsOfFLMLORWithOne( A ) ) then
      Print( "FLMLORWithOne( ", LeftActingDomain( A ), ", [], ",
             Zero( A ), " )" );
    else
      Print( "FLMLORWithOne( ", LeftActingDomain( A ), ", ",
             GeneratorsOfFLMLORWithOne( A ), " )" );
    fi;
    end );


#############################################################################
##
#M  ViewObj( <A> )  . . . . . . . . . . . . . . . . . . . . . view an algebra
##
##  print left acting domain, if known also dimension or no. of generators
##
InstallMethod( ViewObj,
    "for an algebra",
    [ IsAlgebra ],
    function( A )
    Print( "<algebra over ", LeftActingDomain( A ), ">" );
    end );

InstallMethod( ViewObj,
    "for an algebra with known dimension",
    [ IsAlgebra and HasDimension ], 1,   # override method requiring gens.
    function( A )
    Print( "<algebra of dimension ", Dimension( A ),
           " over ", LeftActingDomain( A ), ">" );
    end );

InstallMethod( ViewObj,
    "for an algebra with known generators",
    [ IsAlgebra and HasGeneratorsOfAlgebra ],
    function( A )
    Print( "<algebra over ", LeftActingDomain( A ), ", with ",
           Pluralize( Length( GeneratorsOfAlgebra( A ) ), "generator" ), ">" );
    end );


#############################################################################
##
#M  PrintObj( <A> ) . . . . . . . . . . . . . . . . . . . .  print an algebra
##
InstallMethod( PrintObj,
    "for an algebra",
    [ IsAlgebra ],
    function( A )
    Print( "Algebra( ", LeftActingDomain( A ), ", ... )" );
    end );

InstallMethod( PrintObj,
    "for an algebra with known generators",
    [ IsAlgebra and HasGeneratorsOfAlgebra ],
    function( A )
    if IsEmpty( GeneratorsOfAlgebra( A ) ) then
      Print( "Algebra( ", LeftActingDomain( A ), ", [], ", Zero( A ), " )" );
    else
      Print( "Algebra( ", LeftActingDomain( A ), ", ",
             GeneratorsOfAlgebra( A ), " )" );
    fi;
    end );


#############################################################################
##
#M  ViewObj( <A> )  . . . . . . . . . . . . . . . .  view an algebra-with-one
##
##  print left acting domain, if known also dimension or no. of generators
##
InstallMethod( ViewObj, "for an algebra-with-one", [ IsAlgebraWithOne ],
function( A )
  if IsIdenticalObj(A,LeftActingDomain(A)) then
    Print( "<algebra-with-one over itself>" );
  else
    Print( "<algebra-with-one over ", LeftActingDomain( A ), ">" );
  fi;
end );

InstallMethod( ViewObj,
    "for an algebra-with-one with known dimension",
    [ IsAlgebraWithOne and HasDimension ], 1,   # override method requ. gens.
    function( A )
    Print( "<algebra-with-one of dimension ", Dimension( A ),
           " over ", LeftActingDomain( A ), ">" );
    end );

InstallMethod( ViewObj,
    "for an algebra-with-one with known generators",
    [ IsAlgebraWithOne and HasGeneratorsOfAlgebraWithOne ],
    function( A )
    Print( "<algebra-with-one over ", LeftActingDomain( A ), ", with ",
           Pluralize( Length( GeneratorsOfAlgebraWithOne( A ) ), "generator" ),
           ">" );
    end );


#############################################################################
##
#M  PrintObj( <A> ) . . . . . . . . . . . . . . . . print an algebra-with-one
##
InstallMethod( PrintObj,
    "for an algebra-with-one",
    [ IsAlgebraWithOne ],
    function( A )
    Print( "AlgebraWithOne( ", LeftActingDomain( A ), ", ... )" );
    end );

InstallMethod( PrintObj,
    "for an algebra-with-one with known generators",
    [ IsAlgebraWithOne and HasGeneratorsOfAlgebra ],
    function( A )
    if IsEmpty( GeneratorsOfAlgebraWithOne( A ) ) then
      Print( "AlgebraWithOne( ", LeftActingDomain( A ), ", [], ",
             Zero( A ), " )" );
    else
      Print( "AlgebraWithOne( ", LeftActingDomain( A ), ", ",
             GeneratorsOfAlgebraWithOne( A ), " )" );
    fi;
    end );


#############################################################################
##
#M  ViewObj( <A> )  . . . . . . . . . . . . . . . . . . view a Lie algebra
##
##  print left acting domain, if known also dimension or no. of generators
##
InstallMethod( ViewObj,
    "for a Lie algebra",
    [ IsLieAlgebra ],
    function( A )
    Print( "<Lie algebra over ", LeftActingDomain( A ), ">" );
    end );

InstallMethod( ViewObj,
    "for a Lie algebra with known dimension",
    [ IsLieAlgebra and HasDimension ], 1,       # override method requ. gens.
    function( A )
    Print( "<Lie algebra of dimension ", Dimension( A ),
           " over ", LeftActingDomain( A ), ">" );
    end );

InstallMethod( ViewObj,
    "for a Lie algebra with known generators",
    [ IsLieAlgebra and HasGeneratorsOfAlgebra ],
    function( A )
    Print( "<Lie algebra over ", LeftActingDomain( A ), ", with ",
           Pluralize( Length( GeneratorsOfAlgebra( A ) ), "generator" ), ">" );
end );


#############################################################################
##
#M  AsSubalgebra(<A>, <U>) . view an algebra as subalgebra of another algebra
##
InstallMethod( AsSubalgebra,
    "for two algebras",
    IsIdenticalObj,
    [ IsAlgebra, IsAlgebra ],
    function( A, U )
    local samecoeffs, S;
    if not IsSubset( A, U ) then
      return fail;
    fi;

    # Construct the generators list.
    samecoeffs:= LeftActingDomain( A ) = LeftActingDomain( U );
    if not samecoeffs then
      U:= AsAlgebra( LeftActingDomain( A ), U );
    fi;

    # Construct the subalgebra.
    S:= SubalgebraNC( A, GeneratorsOfAlgebra( U ) );

    # Maintain useful information.
    UseIsomorphismRelation( U, S );
    UseSubsetRelation( U, S );
    if samecoeffs and HasDimension( U ) then
      SetDimension( S, Dimension( U ) );
    fi;

    # Return the subalgebra.
    return S;
    end );

InstallMethod( AsSubalgebra,
    "for an algebra and an algebra-with-one",
    IsIdenticalObj,
    [ IsAlgebra, IsAlgebraWithOne ],
    function( A, U )
    local samecoeffs, S;
    if not IsSubset( A, U ) then
      return fail;
    fi;

    # Construct the generators list.
    samecoeffs:= LeftActingDomain( A ) = LeftActingDomain( U );
    if not samecoeffs then
      U:= AsAlgebraWithOne( LeftActingDomain( A ), U );
    fi;

    # Construct the subalgebra.
    S:= SubalgebraWithOneNC( A, GeneratorsOfAlgebraWithOne( U ) );

    # Maintain useful information.
    UseIsomorphismRelation( U, S );
    UseSubsetRelation( U, S );
    if samecoeffs and HasDimension( U ) then
      SetDimension( S, Dimension( U ) );
    fi;

    # Return the subalgebra.
    return S;
    end );


#############################################################################
##
#M  AsSubalgebraWithOne(<A>, <U>) . . . view algebra as subalgebra of another
##
InstallMethod( AsSubalgebraWithOne,
    "for two algebras",
    IsIdenticalObj,
    [ IsAlgebra, IsAlgebra ],
    function( A, U )
    local S;
    if not IsSubset( A, U ) or One( U ) = fail then
      return fail;
    fi;

    if LeftActingDomain( A ) <> LeftActingDomain( U ) then
      U:= AsAlgebraWithOne( LeftActingDomain( A ), U );
    fi;

    # Construct and return the subalgebra.
    S:= SubalgebraWithOneNC( A, GeneratorsOfAlgebra( U ) );

    UseIsomorphismRelation( U, S );
    UseSubsetRelation( U, S );

    return S;
    end );

InstallMethod( AsSubalgebraWithOne,
    "for an algebra and an algebra-with-one",
    IsIdenticalObj,
    [ IsAlgebra, IsAlgebraWithOne ],
    function( A, U )
    local S;
    if not IsSubset( A, U ) then
      return fail;
    fi;

    if LeftActingDomain( A ) <> LeftActingDomain( U ) then
      U:= AsAlgebraWithOne( LeftActingDomain( A ), U );
    fi;

    # Construct and return the subalgebra.
    S:= SubalgebraWithOneNC( A, GeneratorsOfAlgebraWithOne( U ) );

    UseIsomorphismRelation( U, S );
    UseSubsetRelation( U, S );

    return S;
    end );


#############################################################################
##
#F  CentralizerInFiniteDimensionalAlgebra( <A>, <S>, <issubset> )
##
##  Let $(b_1, \ldots, b_n)$ be a basis of <A>, and $(s_1, \ldots, s_m)$
##  be a basis of <S>, with $s_j = \sum_{l=1}^n v_{jl} b_l$.
##  The structure constants of <A> are $c_{ijk}$ with
##  $b_i b_j = \sum_{k=1}^n c_{ijk} b_k$.
##  Then we compute a basis of the solution space of the system
##  $\sum_{i=1}^n a_i \sum_{l=1}^n v_{jl} ( c_{ilk} - c_{lik} )$ for
##  $1 \leq j \leq m$ and $1 \leq k \leq n$.
##
##  (left null space of an $n \times (nm)$ matrix)
##
##  If the multiplication in <A> is known to be anticommutative this is used.
##  (Note that the case of commutative multiplication is handled in a more
##  general way.)
##
#T We will have problems with this approach if centralizers of algebras over
#T non-commutative coefficients domains are considered.
#T In this case, <S> might stand for generators w.r.t. a larger coefficients
#T domain that also must be centralized ...
##
InstallGlobalFunction( CentralizerInFiniteDimensionalAlgebra,
    function( A, S, issubset )

    local B,           # basis of `A'
          T,           # structure constants table w. r. to `B'
          n,           # dimension of `A'
          m,           # length of `S'
          M,           # matrix of the equation system
          v,           # coefficients of basis vectors of `S' w.r. to `B'
          zerovector,  # initialize one row of `M'
          row,         # one row of `M'
          i, j, k, l,  # loop variables
          cil, cli,    #
          offset,
          pos;

    # Handle the case that `S' may be not contained in `A'.
    # (If `A' knows a basis and `S' is a subset of `A' then
    # there are methods that return `A' itself as closure.)
    if not issubset then
      M:= ClosureAlgebra( A, S );
      return Intersection2( A,
                 CentralizerInFiniteDimensionalAlgebra( M, S, true ) );
    fi;

    # Now `S' is known to be contained in `A'.
    B:= Basis( A );
    T:= StructureConstantsTable( B );
    n:= Dimension( A );
    m:= Length( S );
    M:= [];
    v:= List( S, x -> Coefficients( B, x ) );

    zerovector:= [ 1 .. n*m ] * Zero( LeftActingDomain( A ) );

    if HasIsAnticommutative( A ) and IsAnticommutative( A ) then

      # Column $(j-1)*n + k$ contains in row $i$ the value
      # $\sum_{l=1}^n v_{jl} c_{ilk}$

      for i in [ 1 .. n ] do
        row:= ShallowCopy( zerovector );
        for l in [ 1 .. n ] do
          cil:= T[i][l];
          for j in [ 1 .. m ] do
            offset:= (j-1)*n;
            for k in [ 1 .. Length( cil[1] ) ] do
              pos:= cil[1][k] + offset;
              row[ pos ]:= row[ pos ] + v[j][l] * cil[2][k];
            od;
          od;
        od;
        Add( M, row );
      od;

    else

      # Column $(j-1)*n + k$ contains in row $i$ the value
      # $\sum_{l=1}^n v_{jl} ( c_{ilk} - c_{lik} )$

      for i in [ 1 .. n ] do
        row:= ShallowCopy( zerovector );
        for l in [ 1 .. n ] do
          cil:= T[i][l];
          cli:= T[l][i];
          for j in [ 1 .. m ] do
            offset:= (j-1)*n;
            for k in [ 1 .. Length( cil[1] ) ] do
              pos:= cil[1][k] + offset;
              row[ pos ]:= row[ pos ] + v[j][l] * cil[2][k];
            od;
            for k in [ 1 .. Length( cli[1] ) ] do
              pos:= cli[1][k] + offset;
              row[ pos ]:= row[ pos ] - v[j][l] * cli[2][k];
            od;
          od;
        od;
        Add( M, row );
      od;

    fi;

    # Solve the equation system.
    M:= NullspaceMat( M );

    # Construct the generators from the coefficient vectors.
    M:= List( M, x -> LinearCombination( B, x ) );

    # Return the subalgebra.
    if IsFLMLORWithOne( A ) then
      return SubalgebraWithOneNC( A, M, "basis" );
    else
      return SubalgebraNC( A, M, "basis" );
    fi;
end );


#############################################################################
##
#M  CentralizerOp( <A>, <S> ) . . . . . cent. of a vector space in an algebra
##
InstallMethod( CentralizerOp,
    "for a finite dimensional algebra and a vector space with parent",
    IsIdenticalObj,
    [ IsAlgebra, IsVectorSpace and HasParent ],
    function( A, S )
    if    not IsIdenticalObj( A, Parent( S ) )
       or not IsFiniteDimensional( A ) then
      TryNextMethod();
    fi;
    return CentralizerInFiniteDimensionalAlgebra( A,
               BasisVectors( Basis( S ) ),
               true );
    end );


#############################################################################
##
#M  CentralizerOp( <A>, <S> ) . . .  cent. of an algebra in an assoc. algebra
##
InstallMethod( CentralizerOp,
    "for a fin. dim. assoc. algebra and an algebra with parent",
    IsIdenticalObj,
    [ IsAlgebra and IsAssociative, IsAlgebra and HasParent ],
    function( A, S )
    if    not IsIdenticalObj( A, Parent( S ) )
       or not IsFiniteDimensional( A ) then
      TryNextMethod();
    fi;
    return CentralizerInFiniteDimensionalAlgebra( A,
               GeneratorsOfAlgebra( S ),
               true );
    end );


#############################################################################
##
#M  CentralizerOp( <A>, <S> ) . . . . . cent. of a vector space in an algebra
##
InstallMethod( CentralizerOp,
    "for a finite dimensional algebra and a vector space",
    IsIdenticalObj,
    [ IsAlgebra, IsVectorSpace ],
    function( A, S )
    if not IsFiniteDimensional( A ) then
      TryNextMethod();
    fi;
    return CentralizerInFiniteDimensionalAlgebra( A,
               BasisVectors( Basis( S ) ),
               false );
    end );


#############################################################################
##
#M  CentralizerOp( <A>, <S> ) . . .  cent. of an algebra in an assoc. algebra
##
InstallMethod( CentralizerOp,
    "for a fin. dim. assoc. algebra and an algebra",
    IsIdenticalObj,
    [ IsAlgebra and IsAssociative, IsAlgebra ],
    function( A, S )
    if not IsFiniteDimensional( A ) then
      TryNextMethod();
    fi;
    return CentralizerInFiniteDimensionalAlgebra( A,
               GeneratorsOfAlgebra( S ),
               false );
    end );


#############################################################################
##
#M  CentralizerOp( <A>, <elm> ) . . . . . . cent. of an element in an algebra
##
InstallMethod( CentralizerOp,
    "for an algebra and an element",
    IsCollsElms,
    [ IsAlgebra, IsObject ],
    function( A, elm )
    return Centralizer( A, Algebra( LeftActingDomain( A ), [ elm ] ) );
    end );


#############################################################################
##
#F  CentreFromSCTable( <T> )
##
##  Given a structure constants table <T> w.r.t. a basis $B$,
##  `CentreFromSCTable' returns a list of $B$-coefficients vectors
##  of a basis for the centre of an $F$-algebra with s.c. table <T>,
##  where $F$ is assumed to be commutative.
##
##  We have to solve the system
##  $\sum_{i=1}^n a_i ( c_{ijk} - c_{jik} ) = 0$ for $1 \leq j, k \leq n$.
##
BindGlobal( "CentreFromSCTable", function( T )
    local n,       # the dimension
          M,       # matrix of the equation system
          i, j, k, # loop variables
          row,     # one row in `M'
          offset,  # offset between entry in `T' and column in `M'
          pos,     # nonzero positions in $c_{ij}$
          val;     # loop over structure constants in $c_{ij}$

    n:= Length( T ) - 2;
    M:= NullMat( n, n*n, Last(T) );
    for i in [ 1 .. n ] do
      row:= M[i];
      for j in [ 1 .. n ] do

        offset:= (j-1)*n;
        row{ offset + T[i][j][1] }:= T[i][j][2];
        pos:= T[j][i][1];
        val:= T[j][i][2];
        for k in [ 1 .. Length( pos ) ] do
          row[ offset + pos[k] ]:= row[ offset + pos[k] ] - val[k];
#T cheaper!
        od;

      od;
    od;

    # Solve the equation system.
    return NullspaceMat( M );
    end );


#############################################################################
##
#M  Centre( <A> )
##
InstallMethod( Centre,
    "for a finite dimensional FLMLOR",
    [ IsFLMLOR ],
    function( A )
    local   C,        # centre of `A', result
            B,        # a basis of `A'
            M;        # matrix of the equation system

    if not IsFiniteDimensional( A ) then
      TryNextMethod();
    fi;

    # If necessary convert `A' to a FLMLOR over a commutative ring.
    if not IsCommutative( LeftActingDomain( A ) ) then
      A:= AsFLMLOR( Centre( LeftActingDomain( A ) ), A );
    fi;

    # Solve the equation system.
    # If a s.c. table is already known, we use it since this allows us to
    # avoid multiplications.
    # Only if no s.c. table is known and if the algebra is known to be
    # associative and if the number of algebra generators is less than the
    # number of basis vectors, we do not force the computation of a
    # s.c. table.
    # (Note that for associative algebras,
    # we have to check $x a = a x$ only for algebra generators $a$,
    # not for all vectors of a basis.)
    B:= Basis( A );
    if    HasStructureConstantsTable( B )
       or not ( HasIsAssociative( A ) and IsAssociative( A ) ) then
      M:= CentreFromSCTable( StructureConstantsTable( B ) );
    else
      if HasGeneratorsOfAlgebraWithOne( A ) then
        M:= GeneratorsOfAlgebraWithOne( A );
      else
        M:= GeneratorsOfAlgebra( A );
      fi;
      if Dimension( A ) <= 2 * Length( M ) then
        M:= CentreFromSCTable( StructureConstantsTable( B ) );
      else
        M:= List( BasisVectors( B ), bi -> Concatenation( List( M,
                             a -> Coefficients( B, bi * a - a * bi ) ) ) );
        M:= NullspaceMat( M );
      fi;
    fi;

    # Construct the generators from the coefficient vectors.
    M:= List( M, x -> LinearCombination( B, x ) );

    # Construct the centre.
    if IsFLMLORWithOne( A ) then
      C:= SubalgebraWithOneNC( A, M, "basis" );
    else
      C:= SubalgebraNC( A, M, "basis" );
    fi;
    Assert( 1, IsAbelian( C ) );
    SetIsAbelian( C, true );

    # Return the centre.
    return C;
    end );


#############################################################################
##
#F  MutableBasisOfProductSpace( <U>, <V> )
##
##  Once we have the basis vectors for the product space,
##  we only have to decide whether the result of `ProductSpace' is an ideal,
##  an algebra, or just a vector space.
##  This decision is left to the methods of `ProductSpace',
##  the computation of basis vectors is done by `MutableBasisOfProductSpace'.
##
BindGlobal( "MutableBasisOfProductSpace", function( U, V )
    local inter, # intersection of left acting domains
          u, v,  # loop over the bases
          MB;    # mutable basis of the commutator subspace, result

    if LeftActingDomain( U ) = LeftActingDomain( V ) then
      inter:= LeftActingDomain( U );
    else
      inter:= Intersection2( LeftActingDomain( U ), LeftActingDomain( V ) );
      U:= AsVectorSpace( inter, U );
      V:= AsVectorSpace( inter, V );
    fi;

    MB:= MutableBasis( inter, [], Zero( U ) );
    V:= BasisVectors( Basis( V ) );
    for u in BasisVectors( Basis( U ) ) do
      for v in V do
        CloseMutableBasis( MB, u * v );
      od;
    od;

    # Return the result.
    return [ MB, inter ];
end );


#############################################################################
##
#M  ProductSpace( <U>, <V> )  . . . . . . . . . . . for two free left modules
##
InstallMethod( ProductSpace,
    "for two free left modules",
    IsIdenticalObj,
    [ IsFreeLeftModule, IsFreeLeftModule ],
    function( U, V )
    local MB, vectors, C;

    # Compute basis vectors.
    MB:= MutableBasisOfProductSpace( U, V );

    vectors:= BasisVectors( MB[1] );
    if IsEmpty( vectors ) then
      return TrivialSubspace( U );
    fi;

    # Create the appropriate domain.
    if HasParent( U ) and HasParent( V )
                      and IsIdenticalObj( Parent( U ), Parent( V ) ) then
      C:= SubmoduleNC( Parent( U ), vectors, "basis" );
    else
      C:= FreeLeftModule( MB[2], vectors, "basis" );
    fi;

    # Insert the basis.
    SetBasis( C, ImmutableBasis( MB[1] ) );

    # Return the result.
    return C;
    end );


#############################################################################
##
#M  ProductSpace( <U>, <V> )  . . . . . . . . . . . . . . .  for two algebras
##
##  If $<U> = <V>$ is known to be an algebra then the product space is also
##  an algebra, moreover it is an ideal in <U>.
##  If <U> and <V> are known to be ideals in an algebra $A$
##  then the product space is known to be an algebra and an ideal in $A$.
##
InstallMethod( ProductSpace,
    "for two algebras",
    IsIdenticalObj,
    [ IsAlgebra, IsAlgebra ],
    function( U, V )
    local P, MB, C;

    # Look for the ideal relation that allows one to construct an ideal.
    if IsIdenticalObj( U, V ) then
      P:= U;
    elif HasParent( V ) and IsIdenticalObj( Parent( V ), U )
                        and HasIsTwoSidedIdealInParent( V )
                        and IsTwoSidedIdealInParent( V ) then
      P:= U;
    elif HasParent( U ) and IsIdenticalObj( Parent( U ), V )
                        and HasIsTwoSidedIdealInParent( U )
                        and IsTwoSidedIdealInParent( U ) then
      P:= V;
    else
      TryNextMethod();
    fi;

    # Compute basis vectors for the result.
    MB:= MutableBasisOfProductSpace( U, V )[1];

    # The result is an ideal in `U'.
    C:= SubalgebraNC( P, BasisVectors( MB ), "basis" );

    SetIsTwoSidedIdealInParent( C, true );
    SetBasis( C, ImmutableBasis( MB, C ) );

    # Return the result.
    return C;
    end );


#############################################################################
##
#M  ProductSpace( <U>, <V> )  . . . . . . . . for two ideals with same parent
##
InstallMethod( ProductSpace,
    "for two ideals with same parent",
    IsIdenticalObj,
    [ IsAlgebra and HasParent and IsTwoSidedIdealInParent,
      IsAlgebra and HasParent and IsTwoSidedIdealInParent ],
    function( U, V )
    local MB, C;

    if not IsIdenticalObj( Parent( U ), Parent( V ) ) then
      TryNextMethod();
    fi;

    # The result is an ideal in the parent of `U'.
    MB:= MutableBasisOfProductSpace( U, V )[1];
    C:= SubalgebraNC( Parent( U ), BasisVectors( MB ), "basis" );

    SetIsTwoSidedIdealInParent( C, true );
    SetBasis( C, ImmutableBasis( MB, C ) );

    # Return the result.
    return C;
    end );


#############################################################################
##
#M  RadicalOfAlgebra( <A> ) . . . . . . . . radical of an associative algebra
##
##  `RadicalOfAlgebra' computes the radical (maximal nilpotent ideal)
##  of an associative algebra <A> by first constructing a faithful
##  matrix representation.
##  (Note that there is a special implementation for associative matrix
##  algebras.)
##
##  If <A> contains an identity element then the adjoint representation is
##  already faithful.
##
##  Otherwise we add an identity element (Dorroh extension).
##  More precisely we consider the space $B = \{ (x,t) | x\in A, t\in F }$.
##  We let <A> act on this space via $a (x,t) = (ax+ta,0)$.
##  Then it is easily seen that this representation is faithful.
##
InstallMethod( RadicalOfAlgebra,
    "for an associative algebra",
    [ IsAlgebra ],
    function( A )
    local bb,    # list of matrices representing the basis elements of <A>
          n,     # dimension of <A>
          BA,    # basis of `A'
          bv,    # basis vectors of `BA'
          F,     # left acting domain of `A'
          M,     # (n+1) x (n+1) matrix
          i,j,   # loop variables
          col,   # a column of `M'
          B,     # the representation of <A>
          R,     # the radical of `B'
          bas,   # a basis of `B' (corresponding to the basis of <A>)
          rad;   # a basis of the radical of <A>

    # Make sure that the algebra is associative and not a matrix algebra.
    if not IsAssociative( A ) then
      TryNextMethod();
    fi;

    n:= Dimension( A );

    if n = 0 then
      return A;
    fi;

    BA:= Basis( A );
    bv:= BasisVectors( BA );
    F:= LeftActingDomain( A );

    if One( A ) <> fail then

      bb:= List( bv, x -> AdjointMatrix( BA, x ) );

    else

      bb:= [];

      for i in [1..n] do
        M:=[];
        for j in [1..n] do
          col:= Coefficients( BA, bv[i] * bv[j] );
          if not IsMutable( col ) then
            col:= ShallowCopy( col );
          fi;
          col[n+1]:= Zero( F );
          Add( M, col );
        od;
        col:= [ 1 .. n+1 ] * Zero( F );
        col[i]:= One( F );
        Add( M, col );
        Add( bb, M );
      od;

    fi;

    # Compute the radical of the isomorphic matrix algebra.
    B:= Algebra( F, bb, "basis" );
    R:= RadicalOfAlgebra( B );

    # Transfer the radical back to the original algebra.
    bas:= BasisNC( B, bb );
    rad:= List( BasisVectors( Basis( R ) ),
                x -> LinearCombination( BA, Coefficients( bas, x ) ) );

    return SubalgebraNC( A, rad, "basis" );
    end );


#############################################################################
##
#M  IsTrivial( <A> )  . . . . . . . . . . . . . . . . . . . . .  for a FLMLOR
##
InstallMethod( IsTrivial,
    "for a FLMLOR",
    [ IsFLMLOR ],
    A -> ForAll( GeneratorsOfLeftOperatorRing( A ), IsZero ) );


#############################################################################
##
#M  IsTrivial( <A> )  . . . . . . . . . . . . . . . . . for a FLMLOR-with-one
##
##  A FLMLOR-with-one is trivial if and only if its identity is equal to its
##  zero.
##
InstallMethod( IsTrivial,
    "for a FLMLOR-with-one",
    [ IsFLMLORWithOne ],
    A -> IsZero( One( A ) ) );


#############################################################################
##
#M  GeneratorsOfLeftModule( <A> )
##
##  We assume that it is possible to construct a basis for the algebra <A>.
##  If <A> is finite dimensional and if we know algebra generators
##  then the process of successive closure under the action of <A> on itself
##  yields this.
##
InstallMethod( GeneratorsOfLeftModule,
    "for a FLMLOR",
    [ IsFLMLOR ],
    A -> BasisVectors( Basis( A ) ) );


#############################################################################
##
#M  Basis( <A> )  . . . . . . . . . . . .  basis from FLMLOR gens. for FLMLOR
##
InstallMethod( Basis,
    "for a FLMLOR",
    [ IsFLMLOR ],
    function( A )

    # If generators as left module are known
    # we do not need to multiply at all.
    if HasGeneratorsOfLeftModule( A ) then
      TryNextMethod();
    fi;

    return ImmutableBasis( MutableBasisOfNonassociativeAlgebra(
             LeftActingDomain( A ),
             GeneratorsOfLeftOperatorRing( A ),
             Zero( A ),
             infinity ), A );
    end );


#############################################################################
##
#M  Basis( <A> )  . . . . . .  basis from FLMLOR gens. for associative FLMLOR
##
BindGlobal("BasisLieAlgebraAndAssociativeFLMOR",
    function( A )
    local  mb;

    # If generators as left module are known
    # we do not need to multiply at all.
    if HasGeneratorsOfLeftModule( A ) then
      TryNextMethod();
    fi;

    mb:= MutableBasisOfClosureUnderAction(
             LeftActingDomain( A ),
             GeneratorsOfLeftOperatorRing( A ),
             "left",
             GeneratorsOfLeftOperatorRing( A ),
             \*,
             Zero( A ),
             infinity );
    return ImmutableBasis( mb, A );
    end );
InstallMethod( Basis,
    "for an associative FLMLOR",
    [ IsFLMLOR and IsAssociative ],
    BasisLieAlgebraAndAssociativeFLMOR );


#############################################################################
##
#M  Basis( <A> )   .  basis from FLMLOR gens. for associative FLMLOR-with-one
##
InstallMethod( Basis,
    "for an associative FLMLOR-with-one",
    [ IsFLMLORWithOne and IsAssociative ],
    function( A )

    # If generators as left module are known
    # we do not need to multiply at all.
    if HasGeneratorsOfLeftModule( A ) then
      TryNextMethod();
    fi;

    return ImmutableBasis( MutableBasisOfClosureUnderAction(
             LeftActingDomain( A ),
             GeneratorsOfLeftOperatorRing( A ),
             "left",
             [ One( A ) ],
             \*,
             Zero( A ),
             infinity ), A );
    end );


#############################################################################
##
#M  Basis( <A> )  . . . . . . . . . . basis from FLMLOR gens. for Lie algebra
##
##  In a Lie algebra, every word (with brackets) in terms of algebra
##  generators is a linear combination of left-normed words;
##  this means that it is sufficient to multiply with generators from one
##  side.
##
InstallMethod( Basis,
    "for a Lie algebra",
    [ IsLieAlgebra ],
    BasisLieAlgebraAndAssociativeFLMOR );


#############################################################################
##
#M  PowerSubalgebraSeries( <A> )
##
InstallOtherMethod( PowerSubalgebraSeries,
    "for an algebra",
    [ IsAlgebra ],
    function ( A )
    local   S,          # power subalgebra series of <A>, result
            D;          # power subalgebras

    # Compute the series by repeated calling of `ProductSpace'.
    S := [ A ];
    D := ProductSpace( A, A );
    while D <> Last(S)  do
      Add( S, D );
      D:= ProductSpace( D, D );
    od;

    # Return the series when it becomes stable.
    return S;
    end );


#############################################################################
##
#M  IsNilpotentElement( <L>, <x> )  . . . . . . for an algebra and an element
##
##  <x> is nilpotent if its adjoint matrix $A$ (i.e. the linear map coming
##  from left multiplication by <x>) is nilpotent.
##  To check this, we only need to check whether $A^n$ (or a smaller power)
##  is zero, where $n$ denotes the dimension of <L>.
##
InstallMethod( IsNilpotentElement,
    "for an algebra, and an element",
    IsCollsElms,
    [ IsAlgebra, IsRingElement ],
    function( L, x )
    local B,     # a basis of `L'
          A,     # adjoint matrix of `x w.r. to `B'
          n,     # dimension of `L'
          i;     # loop variable

    B := Basis( L );
    A := AdjointMatrix( B, x );
    n := Dimension( L );
    i := 1;

    if ForAll( A, x -> n < PositionNonZero( x ) ) then
#T better ask IsZero?
      return true;
    fi;

    while i < n do
      i:= 2 * i;
      A:= A * A;
      if ForAll( A, x -> n < PositionNonZero( x ) ) then
#T better ask IsZero?
        return true;
      fi;
    od;

    return false;
    end );


#############################################################################
##
#M  GeneratorsOfLeftOperatorRing( <A> ) . . . . . . . . for a FLMLOR-with-one
##
InstallMethod( GeneratorsOfLeftOperatorRing,
    "for a FLMLOR-with-one",
    [ IsFLMLORWithOne ],
    A -> Concatenation( [ One( A ) ],
                        GeneratorsOfLeftOperatorRingWithOne( A ) ) );


#############################################################################
##
#M  GeneratorsOfLeftOperatorRing( <A> ) . . . .  for FLMLOR with module gens.
##
InstallMethod( GeneratorsOfLeftOperatorRing,
    "for a FLMLOR with known left module generators",
    [ IsFLMLOR and HasGeneratorsOfLeftModule ],
    GeneratorsOfLeftModule );


#############################################################################
##
#M  GeneratorsOfLeftOperatorRingWithOne( <A> ) . for FLMLOR with module gens.
##
InstallMethod( GeneratorsOfLeftOperatorRingWithOne,
    "for a FLMLOR-with-one with known left module generators",
    [ IsFLMLORWithOne and HasGeneratorsOfLeftModule ],
    GeneratorsOfLeftModule );


#############################################################################
##
#M  DirectSumOfAlgebras( <A1>, <A2> )
##
##  Construct a s.c. algebra.
##  (There are special methods for the sum of appropriate matrix algebras.)
##
#T embeddings/projections should be provided!
##
InstallOtherMethod( DirectSumOfAlgebras,
    "for two algebras",
    [ IsAlgebra, IsAlgebra ],
    function( A1, A2 )
    local n,     # The dimension of the resulting algebra.
          i,j,   # Loop variables.
          T,     # The table of structure constants of the direct sum.
          scT,   #
          n1,    # The dimension of A1.
          n2,    # The dimension of A2.
          ll,    # A list of structure constants.
          L,     # result.
          sym,   # if both products are (anti)symmetric, then the result
                 # will have the same property.
          R1,R2, # Root systems of A1,A2.
          f1,f2, # Embeddings of A1,A2 in L.
          R,     # Root system of L.
          RV,    # List of various things.
          r,
          pos;   # List of positions.

    if LeftActingDomain( A1 ) <> LeftActingDomain( A2 ) then
      Error( "<A1> and <A2> must be written over the same field" );
    fi;

    n1:= Dimension( A1 );
    n2:= Dimension( A2 );
    n:= n1+n2;
    T:= [];

    # Initialize the s.c. table.
    T:= EmptySCTable( n, Zero( LeftActingDomain( A1 ) ) );

    # Enter the structure constants for the first algebra.
    scT:= StructureConstantsTable( Basis( A1 ) );
    sym:= scT[n1+1];
    T{ [ 1 .. n1 ] }{ [ 1 .. n1 ] }:= scT{ [ 1 .. n1 ] };

    scT:= StructureConstantsTable( Basis( A2 ) );

    for i in [1..n2] do
      for j in [1..n2] do
        ll:= ShallowCopy( scT[i][j] );
        ll[1]:= ll[1] + n1;
        T[n1+i][n1+j]:= ll;
      od;
    od;


    # Set the (anti)symmetric flag
    if scT[n2 + 1] = sym  then
        T[n + 1] := sym;
    fi;
    if Characteristic( LeftActingDomain( A1 ) ) = 2 and sym in [ 1, -1 ]
        and scT[n2 + 1] in [ 1, -1 ]  then
        T[n + 1] := 1;
    fi;


    L:= AlgebraByStructureConstants( LeftActingDomain( A1 ), T );

    # Maintain useful information.
    if     HasIsLieAlgebra( A1 ) and HasIsLieAlgebra( A2 )
       and IsLieAlgebra( A1 ) and IsLieAlgebra( A2 ) then
       SetIsLieAlgebra( L, true );
       if HasRootSystem( A1 ) and HasRootSystem( A2 ) then
           # We can easily compute the root system of `L'.
           R1:= RootSystem( A1 ); R2:= RootSystem( A2 );
           f1:= LeftModuleGeneralMappingByImages( A1, L, CanonicalBasis(A1),
                        CanonicalBasis(L){[1..Dimension(A1)]} );
           f2:= LeftModuleGeneralMappingByImages( A2, L, CanonicalBasis(A2),
                        CanonicalBasis(L){[Dimension(A1)+1..Dimension(L)]} );
           R:= Objectify( NewType( NewFamily( "RootSystemFam", IsObject ),
                       IsAttributeStoringRep and IsRootSystemFromLieAlgebra ),
                       rec() );
           RV:= List( PositiveRootVectors( R1 ), x -> Image( f1, x ) );
           Append( RV, List( PositiveRootVectors( R2 ),
                                               x -> Image( f2, x ) ) );
           SetPositiveRootVectors( R, RV );
           RV:= List( NegativeRootVectors( R1 ), x -> Image( f1, x ) );
           Append( RV, List( NegativeRootVectors( R2 ),
                                               x -> Image( f2, x ) ) );
           SetNegativeRootVectors( R, RV );
           RV:= List( PositiveRoots( R1 ), ShallowCopy );
           for i in [1..Length(RV)] do
               Append( RV[i], ListWithIdenticalEntries(
                       Length( CartanMatrix( R2 ) ),
                       Zero( LeftActingDomain( A2 ) ) ) );
           od;
           for i in PositiveRoots( R2 ) do
               r:= ListWithIdenticalEntries( Length( CartanMatrix( R1 ) ),
                           Zero( LeftActingDomain( A1 ) ) );
               Append( r, i );
               Add( RV, r );
           od;
           SetPositiveRoots( R, RV );
           RV:= List( NegativeRoots( R1 ), ShallowCopy );
           for i in [1..Length(RV)] do
               Append( RV[i], ListWithIdenticalEntries(
                       Length( CartanMatrix( R2 ) ),
                       Zero( LeftActingDomain( A2 ) ) ) );
           od;
           for i in NegativeRoots( R2 ) do
               r:= ListWithIdenticalEntries( Length( CartanMatrix( R1 ) ),
                           Zero( LeftActingDomain( A1 ) ) );
               Append( r, i );
               Add( RV, r );
           od;
           SetNegativeRoots( R, RV );
           pos:= List( SimpleSystem( R1 ), x -> Position(
                         PositiveRoots(R1), x ) );
           RV:= PositiveRoots( R ){pos};
           pos:= List( SimpleSystem( R2 ), x -> Position(
                       PositiveRoots(R2), x ) + Length( PositiveRoots(R1) ));
           Append( RV, PositiveRoots( R ){pos} );
           SetSimpleSystem( R, RV );
           RV:= [ ];
           for i in [1,2,3] do
               RV[i]:= List( CanonicalGenerators( R1 )[i],
                             x -> Image( f1, x ) );
               Append( RV[i], List( CanonicalGenerators( R2 )[i],
                              x -> Image( f2, x ) ) );
           od;
           SetCanonicalGenerators( R, RV );
           SetCartanMatrix( R, DirectSumMat( CartanMatrix( R1 ),
                   CartanMatrix( R2 ) ) );
           SetUnderlyingLieAlgebra( R, L );
           SetRootSystem( L, R );
           if HasChevalleyBasis( A1 ) and HasChevalleyBasis( A2 ) then
               RV:= [ ];
               for i in [1,2,3] do
                   RV[i]:= List( ChevalleyBasis( A1 )[i],
                                 x -> Image( f1, x ) );
                   Append( RV[i], List( ChevalleyBasis( A2 )[i],
                           x -> Image( f2, x ) ) );
               od;
               SetChevalleyBasis( L, RV );
           fi;
       fi;

       if HasSemiSimpleType( A1 ) and HasSemiSimpleType( A2 ) then
           SetSemiSimpleType( L, Concatenation( SemiSimpleType(A1)," ",
                   SemiSimpleType( A2 ) ) );
       fi;

       if HasIsRestrictedLieAlgebra ( A1 ) and HasIsRestrictedLieAlgebra ( A2 )
         and IsRestrictedLieAlgebra ( A1 ) and IsRestrictedLieAlgebra ( A2 )
       then
          SetIsRestrictedLieAlgebra( L, true );
          if HasPthPowerImages( Basis ( A1 ) )
             and HasPthPowerImages( Basis ( A2 ) ) then
             if not IsBound (f1) then
                f1:= LeftModuleGeneralMappingByImages( A1, L,
                     CanonicalBasis(A1),
                     CanonicalBasis(L){[1..Dimension(A1)]} );
                f2:= LeftModuleGeneralMappingByImages( A2, L,
                     CanonicalBasis(A2),
                     CanonicalBasis(L){[Dimension(A1)+1..Dimension(L)]});
             fi;

             SetPthPowerImages( Basis ( L ),
               Concatenation ( List (PthPowerImages( Basis ( A1 ) ), x->x^f1),
                          List (PthPowerImages( Basis ( A2 ) ), x->x^f2) ) );

          fi;
       fi;

    fi;
    if     HasIsAssociative( A1 ) and HasIsAssociative( A2 )
       and IsAssociative( A1 ) and IsAssociative( A2 ) then
      SetIsAssociative( L, true );
    fi;

    # Return the result.
    return L;
    end );


#############################################################################
##
#M  DirectSumOfAlgebras( <list> ) . . . . . . .  for a dense list of algebras
##
InstallMethod( DirectSumOfAlgebras,
    "for list of algebras",
    [ IsDenseList ],
    function( list )
    local R, A, i;

    if IsEmpty( list ) then
      Error( "<list> must be nonempty" );
    fi;

    R:= LeftActingDomain( list[1] );
    for A in list do
      if not IsFLMLOR( A ) or LeftActingDomain( A ) <> R then
        Error( "all entries in <list> must be FLMLORs over <R>" );
      fi;
    od;

    A:= list[1];
    for i in [ 2 .. Length( list ) ] do
      A:= DirectSumOfAlgebras( A, list[i] );
    od;

    return A;
    end );


#############################################################################
##
#O  IsCentral( <A>, <U> )  . . . . . . . .  test if <U> is centralized by <A>
##
##  Check whether every basis vector of <A> commutes with every basis vector
##  of the subset <U>.
##
##  For associative algebras, we have to check $u a = a u$ only for algebra
##  generators $a$ and $u$ of $A$ and $U$, respectively,
##  not for all vectors of a basis.
##
InstallMethod( IsCentral,
    "for two FLMLORs",
    IsIdenticalObj,
    [ IsFLMLOR, IsFLMLOR ],
    IsCentralFromGenerators( GeneratorsOfLeftModule,
                             GeneratorsOfLeftModule ) );

InstallMethod( IsCentral,
    "for two associative FLMLORs",
    IsIdenticalObj,
    [ IsFLMLOR and IsAssociative, IsFLMLOR and IsAssociative ],
    IsCentralFromGenerators( GeneratorsOfAlgebra,
                             GeneratorsOfAlgebra ) );

InstallMethod( IsCentral,
    "for two associative FLMLORs-with-one",
    IsIdenticalObj,
    [ IsFLMLORWithOne and IsAssociative,
      IsFLMLORWithOne and IsAssociative ],
    IsCentralFromGenerators( GeneratorsOfAlgebraWithOne,
                             GeneratorsOfAlgebraWithOne ) );

#############################################################################
##
#O  IsCentral( <A>, <x> )  . . . . . . . .  test if <x> is centralized by <A>
##
InstallMethod( IsCentral,
    "for an FLMLOR and an element",
    IsCollsElms,
    [ IsFLMLOR, IsObject ],
    IsCentralElementFromGenerators( GeneratorsOfLeftModule ) );

InstallMethod( IsCentral,
    "for an associative FLMLOR and an element",
    IsCollsElms,
    [ IsFLMLOR and IsAssociative, IsObject ],
    IsCentralElementFromGenerators( GeneratorsOfAlgebra ) );

InstallMethod( IsCentral,
    "for an associative FLMLOs-with-one and an element",
    IsCollsElms,
    [ IsFLMLORWithOne and IsAssociative, IsObject ],
    IsCentralElementFromGenerators( GeneratorsOfAlgebraWithOne ) );


#############################################################################
##
#F  FreeAlgebraConstructor( <name>, <magma> )
##
##  is used for a uniform treatment of free (associative) algebras(-with-one)
##
BindGlobal( "FreeAlgebraConstructor", function( name, magma )
    return function( arg )
    local   R,          # coefficients ring
            names,      # names of the algebra generators
            M,          # magma
            A,          # algebra
            F,          # family
            i;

    # Check the argument list.
    if Length( arg ) = 0 or not IsRing( arg[1] ) then
      Error( "first argument must be a ring" );
    fi;

    R:= arg[1];

    # Construct names of generators.
    if   Length( arg ) = 2 and IsInt( arg[2] ) then
      names:= List( [ 1 .. arg[2] ],
                    i -> Concatenation( "x.", String(i) ) );
      MakeImmutable( names );
    elif     Length( arg ) = 2
         and IsList( arg[2] )
         and ForAll( arg[2], IsString ) then
      names:= arg[2];
    elif Length( arg ) = 3 and IsInt( arg[2] ) and IsString( arg[3] ) then
      names:= List( [ 1 .. arg[2] ],
                    x -> Concatenation( arg[3], ".", String(x) ) );
      MakeImmutable( names );
    elif ForAll( arg{ [ 2 .. Length( arg ) ] }, IsString ) then
      names:= arg{ [ 2 .. Length( arg ) ] };
    else
      Error( "usage: ", name, "( <R>, <rank> )\n",
                 "or ", name, "( <R>, <name1>, ... )" );
    fi;

    M := magma( names );

    # Construct the algebra as free magma algebra of a free magma over `R'.
    A := FreeMagmaRing( R, M );

    # Store the names.
    F := ElementsFamily( FamilyObj( A ) );
    F!.names:= names;

    # Install grading
    if HasOne(M) then i := 0; else i := 1; fi;
    SetGrading( A, rec(min_degree := i,
                                     max_degree := infinity,
                                     source := Integers,
                                     hom_components := function(degree)
        local i, d, B, x, y;
        if HasOne(M) then
            B := [[One(M)],GeneratorsOfMagmaWithOne(M)];
        else
            B := [[],GeneratorsOfMagma(M)];
        fi;
        for d in [2..degree] do
            Add(B,[]);
            if IsAssociative(M) then
                for x in B[2] do for y in B[d] do
                    Add(B[d+1],x*y);
                od; od;
            else
                for i in [2..d] do
                    for x in B[i] do for y in B[d+2-i] do
                        Add(B[d+1],x*y);
                    od; od;
                od;
            fi;
        od;
        x := Zero(R);
        y := [One(R)];
        return FreeLeftModule(R, List(B[degree+1],
                       p->ElementOfMagmaRing( F, x, y, [p] )), Zero(A));
    end));

    # Return the result.
    return A;
    end;
end );


#############################################################################
##
#F  FreeAlgebra( <R>, <rank> ) . . . . . . . . . . free algebra of given rank
#F  FreeAlgebra( <R>, <rank>, <name> )
#F  FreeAlgebra( <R>, <name1>, <name2>, ... )
##
InstallGlobalFunction( FreeAlgebra,
    FreeAlgebraConstructor( "FreeAlgebra", FreeMagma ) );


#############################################################################
##
#F  FreeAlgebraWithOne( <R>, <rank> ) . . free algebra-with-one of given rank
#F  FreeAlgebraWithOne( <R>, <rank>, <name> )
#F  FreeAlgebraWithOne( <R>, <name1>, <name2>, ... )
##
InstallGlobalFunction( FreeAlgebraWithOne,
    FreeAlgebraConstructor( "FreeAlgebraWithOne", FreeMagmaWithOne ) );


#############################################################################
##
#F  FreeAssociativeAlgebra( <R>, <rank> )
#F  FreeAssociativeAlgebra( <R>, <rank>, <name> )
#F  FreeAssociativeAlgebra( <R>, <name1>, <name2>, ... )
##
InstallGlobalFunction( FreeAssociativeAlgebra,
    FreeAlgebraConstructor( "FreeAssociativeAlgebra", FreeSemigroup ) );


#############################################################################
##
#F  FreeAssociativeAlgebraWithOne( <R>, <rank> )
#F  FreeAssociativeAlgebraWithOne( <R>, <rank>, <name> )
#F  FreeAssociativeAlgebraWithOne( <R>, <name1>, <name2>, ... )
##
InstallGlobalFunction( FreeAssociativeAlgebraWithOne,
    FreeAlgebraConstructor( "FreeAssociativeAlgebraWithOne", FreeMonoid ) );


#############################################################################
##
#M  \.( <F>, <n> )  . . . . . . . . .  access to generators of a free algebra
##
InstallAccessToGenerators( IsMagmaRingModuloRelations,
                           "magma ring containing the whole family",
                           GeneratorsOfAlgebra );

InstallAccessToGenerators( IsMagmaRingModuloRelations and IsRingWithOne,
                           "magma ring-with-one containing the whole family",
                           GeneratorsOfAlgebraWithOne );


#############################################################################
##
#M  CentralIdempotentsOfAlgebra( <A> )
##
##   Let A be an associative algebra with one. We construct a maximal
##   system of orthogonal primitive idemoptents in the centre of A.
##   First we let B be the centre of A and Q the
##   the semisimple commutative associative algebra A/Rad(A).
##   We calculate a complete set of orthogonal idempotents in `Q'
##   and then lift them to A.
##   The orthogonal idempotents in `Q' correspond to the decomposition
##   of `Q' as a direct sum of simple ideals. Now `ideals' will contain
##   a list of ideals of `Q' such that the direct sum of these equals
##   `Q'. The variable `ids' will contain the idempotents corresponding
##   to the ideals in `ids'.
##   The algorithms has two parts: one for small fields (of size less than
##   `2*Dimension( Q )', and one for big fields.
##   If the field is big, then using a Las Vegas algorithm we find a
##   splitting element (this is an element that generates `Q'). By
##   factoring the minimal polynomial of such element we can find a
##   complete set of orthogonal idempotents in one step.
##   However, if the field is small splitting elements might not exist.
##   In this case we use decomposable elements (of which the minimum
##   polynomial factors into two (or more) relatively prime factors.
##   Then using the same procedure as for splitting elements we can
##   find some idempotents. But in this case the corresponding ideals
##   might split further. So we have to find decomposable elements in
##   these and so on.
##   Decomposable elements are found as follows: first we calculate
##   the subalgebra of all elements x such that x^q=x
##   (where `q=Size( F )').
##   This subalgebra is a number of copies of the ground field. So any
##   element independent from 1 of this subalgebra will have a minimum
##   polynomial that splits completely. On the other hand, if 1 is the
##   only basis vector of this subalgebra than the original algebra was
##   simple.
##   For a more elaborate description we refer to "W. Eberly and M.
##   Giesbrecht, Efficient Decomposition of Associative Algebras,
##   Proceedings of ISSAC 1996."
##

InstallMethod( CentralIdempotentsOfAlgebra,
    "for an associative algebra",
    [ IsAlgebra ],
    function( A )
    local F,B,Rad,Q,bQ,ids,ideals,id,i,j,k,l,set,cf,e,vv,sp,x,f,q,sol,
          eq,facs,hlist,c,p,g,gcd,bb,E,ei,ni,hom,qq;


      if MultiplicativeNeutralElement( A ) = fail then TryNextMethod(); fi;

      F:=LeftActingDomain(A);
      B:=Centre(A);

      Rad:= RadicalOfAlgebra( B );
      hom:= NaturalHomomorphismByIdeal( B, Rad );
      Q:= ImagesSource( hom );
      bQ:= BasisVectors( Basis( Q ) );
      ids:= [ MultiplicativeNeutralElement( Q ) ];
      ideals:= [ Q ];

      # The variable `k' will point to the first element of `ideals' that
      # still has to be decomposed.

      k:=1;

      if Size(F) > 2*Dimension( Q )^2 then
        set:= [ 0 .. 2*Dimension(Q)^2 ]*One( F );
      else
        set:= [ ];
      fi;

      repeat

        if Length( set ) > 1 then

          # We are in the case of a big field.

          repeat

            # We try to find an element of `Q' that generates it.
            # If we take the coefficients of such an element randomly
            # from a set of `2*Dimension(Q)^2' elements,
            # then this element generates `Q' with probability > 1/2

            bQ:= BasisVectors( Basis( ideals[k] ) );
            cf:= List( [ 1 .. Length(bQ) ], x -> Random( set ) );
            e:= LinearCombination( bQ, cf );

            # Now we calculate the minimum polynomial of `e'.

            vv:= [ MultiplicativeNeutralElement( ideals[k] ) ];
            sp:= MutableBasis( F, vv );
            x:= ShallowCopy( e );

            while not IsContainedInSpan( sp, x ) do
              Add( vv, x );
              CloseMutableBasis( sp, x );
              x:= x*e;
            od;
            sp:= UnderlyingLeftModule( ImmutableBasis( sp ) );
            cf:= ShallowCopy(
                   - Coefficients( BasisNC( sp, vv ), x )
                 );
            Add( cf, One( F ) );
            f:= ElementsFamily( FamilyObj( F ) );
            f:= LaurentPolynomialByCoefficients( f, cf, 0 );

          until DegreeOfLaurentPolynomial( f ) = Dimension( Q );

        else

          # Here the field is small.

          q:= Size( F );

        # `sol' will be a basis of the subalgebra of the k-th ideal
        # consisting of all elements x such that x^q=x.
        # If the length of this list is 1,
        # then the ideal is simple and we proceed to the next one. If all
        # ideals are simple then we quit the loop.

          sol:= [ ];
          while Length( sol ) < 2 and k <= Length( ideals ) do
            bQ:= BasisVectors( Basis( ideals[k] ) );
            eq:= [ ];
            for i in [1..Dimension( ideals[k] )] do
              Add( eq, Coefficients( Basis( ideals[k] ), bQ[i]^q-bQ[i] ) );
            od;
            sol:= List( NullspaceMat( eq ),
                        x -> LinearCombination( bQ, x ) );
            if Length(sol) = 1 then k:=k+1; fi;
          od;

          if k>Length(ideals) then break; fi;

          vv:= [ MultiplicativeNeutralElement( ideals[k] ) ];
          sp:= MutableBasis( F, vv );

          e:= sol[1];
          if IsContainedInSpan( sp, e ) then e:=sol[2]; fi;

        # We calculate the minimum polynomial of `e'.

          x:= ShallowCopy( e );
          while not IsContainedInSpan( sp, x ) do
            Add( vv, x );
            CloseMutableBasis( sp, x );
            x:= x*e;
          od;
          sp:= UnderlyingLeftModule( ImmutableBasis( sp ) );
          cf:=  ShallowCopy(
                  - Coefficients( BasisNC( sp, vv ), x )
                );
          Add( cf, One( F ) );

          f:= ElementsFamily( FamilyObj( F ) );
          f:= LaurentPolynomialByCoefficients( f, cf, 0 );

        fi;

        facs:= Factors( PolynomialRing( F ), f );

      # Now we find elements h1,...,hs such that `hi = 1 mod facs[i]' and
      # `hi = 0 mod facs[j]' if `i<>j'.
      # This is possible due to the Chinese remainder theorem.

        hlist:= [ ];
        for i in [1..Length( facs )] do
          cf:= List( [ 1..Length( facs ) ], x -> Zero( F ) );
          cf[i]:= One(F);
          j:= 1;
          c:= cf[1];
          p:= facs[1];
          while j < Length(facs) do
            j:= j + 1;
            g:= GcdRepresentation( p, facs[j] );
            gcd:= g[1]*p+g[2]*facs[j];
            qq:= g[1]*(cf[j]-c)/gcd;
            if qq<>0*qq then
              c:= p*EuclideanRemainder( qq, facs[j] )
                              + c;
            fi;
            p:= p*facs[j] / gcd;
          od;

          Add( hlist, EuclideanRemainder( c*facs[i]^0 , p ) );

        od;

      # Now a set of orthogonal idempotents is given by `hi(e)'.
      # We evaluate `hi(e)' in a rather strange way; this in order to make
      # sure that the one is the one of `ideals[ k ]' ('e^0' will be the
      # one of the big algebra `Q').

        id:= List( hlist, x -> Value( x, e,
                        MultiplicativeNeutralElement( ideals[k] ) ) );

        if Length(set) = 0 then

        # We are in the case of a small field;
        # so we append the new idempotents and ideals
        # (and erase the old ones). (If `E' is an idempotent,
        # then the corresponding ideal is given by `E*Q*E'.)

          Append(ids,id);

          for l in [1..Length(id)] do
            bb:=List(BasisVectors(Basis(ideals[k])),x->id[l]*x*id[l]);
            Add(ideals,Subalgebra(Q,bb));
          od;

          ideals:=Filtered(ideals,x->x<>ideals[k]);
          ids:=Filtered(ids,x->x<>ids[k]);
        else

        # Here the field is big so we found the complete list of idempotents
        # in one step.

          ids:= id;
          k:=Length(ideals)+1;
        fi;

        while k<=Length(ideals) and Dimension( ideals[k] ) = 1 do k:=k+1; od;

      until k>Length(ideals);

      id:= List( ids, e -> PreImagesRepresentative( hom, e ) );

      # Now we lift the idempotents to the big algebra `A'. The
      # first idempotent is lifted as follows:
      # We have that `id[1]^2-id[1]' is an element of `Rad'.
      # We construct the sequences e_{i+1} = e_i + n_i - 2e_in_i,
      # and n_{i+1}=e_{i+1}^2-e_{i+1}, starting with e_0=id[1].
      # It can be proved by induction that e_q is an idempotent in `A'
      # because n_0^{2^q}=0.
      # Now `E' will be the sum of all idempotents lifted so far.
      # Then the next lifted idempotent is obtained by setting
      # `ei:=id[i]-E*id[i]-id[i]*E+E*id[i]*E;'
      # and lifting as above. It can be proved that in this manner we
      # get an orthogonal system of primitive idempotents.

      E:= Zero( F )*id[1];

      for i in [1..Length(id)] do
        ei:= id[i]-E*id[i]-id[i]*E+E*id[i]*E;
        q:= 0;
        while 2^q <= Dimension( Rad ) do
          q:= q+1;
        od;
        ni:= ei*ei-ei;
        for j in [1..q] do
          ei:= ei+ni-2*ei*ni;
          ni:= ei*ei-ei;
        od;
        id[i]:= ei;
        E:= E+ei;
      od;

      return AsSSortedList(id);
end );


##############################################################################
##
#M  IsSimpleAlgebra( <A> )  . . . . . . . . . . . .for an associative algebra
##
##  A test whether <A> is simple.
##
InstallMethod( IsSimpleAlgebra,
    "for an associative algebra",
    [ IsAlgebra ],
    function( A )

    if not IsAssociative( A ) then
      TryNextMethod();
    fi;

    if Dimension( RadicalOfAlgebra( A ) ) <> 0 then
      return false;
    else
      return Length( CentralIdempotentsOfAlgebra( A ) ) = 1;
    fi;
    end );


###############################################################################
##
#M  LeviMalcevDecomposition( <L> )
##
##  A Levi-Malcev subalgebra of `L' is a semisimple subalgebra complementary to
##  the radical `R'. We find a Levi-Malcev subalgebra of `L' by first
##  computing a complementary subspace to `R'. This subspace is a Levi-Malcev
##  subalgebra modulo `R'. Then we change the basis vectors such that they
##  form a basis of a Levi-Malcev subalgebra modulo the second term of the
##  derived series of `R' after that we consider the third term of the
##  derived series, and so on.
##
InstallMethod( LeviMalcevDecomposition,
    "for an associative or a Lie algebra",
    [ IsAlgebra ],
    function( L )
    local R,             # The solvable radical of `L'.
          s,             # The dimension of the Levi subalgebra.
          F,             # coefficients domain of `L'
          bas,bb,        # Lists of basiselements.
          sp,            # A vector space.
          subalg,        # Boolean.
          a,i,j,k,l,m,   # Loop variables.
          x,             # Element of `L'.
          ser,           # The derived series of `R'.
          p,             # The length of the derived series.
          Rbas,          # A special basis of `R'.
          levi,          # Basis of a Levi complement.
          T,             # Structure constants table of `L', w.r.t. a
                         # particular basis.
          cf,cf1,cf2,    # Coefficient vectors.
          klist,         # List of integers.
          comp,          # List of basis vectors of a complement.
          dim,           # The length of `comp'.
          B,             # A basis.
          cij,           # Entry of the table of structure constants.
          eqs,           # Matrix of equation set.
          rl,            # Right hand side of the equation system.
          eqno,          # Number of the equation.
          sol,           # Solution set to the equations.
          r,offset;      # Integers.

    if IsLieAlgebra( L ) then
      R:= LieSolvableRadical( L );
      offset:= 1;
    elif IsAssociative( L ) then
      R:= RadicalOfAlgebra( L );
      offset:=0;
    fi;

    if Dimension( R ) = 0 then
      return [ L, R ];
    elif Dimension( R ) = Dimension( L ) then
      return [ TrivialSubalgebra( L ), R ];
    fi;

    s:= Dimension( L ) - Dimension( R );

    # `bb' will be a basis of a complement to `R' in `L'.

    bas:= ShallowCopy( BasisVectors( Basis( R ) ) );
    F:= LeftActingDomain( L );
    sp:= MutableBasis( F, bas );
    bb:= [ ];
    for k in BasisVectors( Basis( L ) ) do
      if Length( bb ) = s then
        break;
      elif not IsContainedInSpan( sp, k ) then
        Add( bb, k );
        CloseMutableBasis( sp, k );
      fi;
    od;

    sp:= MutableBasis( F, bb );
    subalg:= true;
    for i in [1..Length(bb)] do
      for j in [offset*i+1..Length(bb)] do
        if not IsContainedInSpan( sp, bb[i]*bb[j] ) then
          subalg:= false;
          break;
        fi;
      od;
    od;
    if subalg then
      Info( InfoAlgebra, 1,
            "LeviDecomposition: subalgebra test successful" );
      return [ SubalgebraNC( L, bb, "basis" ), R ];
    fi;

    ser:= PowerSubalgebraSeries( R );

    # We now calculate a basis of `R' such that the first k1 elements
    # form a basis of the last nonzero term of the derived series `ser',
    # the first k2 ( k2>k1 ) elements form a basis of the next to last
    # element of the derived series, and so on.

    p:= Length( ser );
    i:= p-1;
    Rbas:= ShallowCopy( BasisVectors( Basis( ser[p-1] ) ) );
    sp:= MutableBasis( F, Rbas );
    while Length(Rbas) < Dimension(R) do
      if Length(Rbas) = Dimension(ser[i]) then
        i:= i-1;
        k:= 1;
      else
        x:= BasisVectors( Basis( ser[i] ) )[k];
        if not IsContainedInSpan( sp, x ) then
          Add( Rbas, x );
          CloseMutableBasis( sp, x );
        fi;
        k:= k+1;
      fi;
    od;

    levi:= ShallowCopy( bb );
    Append(bb,Rbas);

    # So now `bb' is a list of basis vectors of `L' such that
    # the first elements form a basis of a complement to `R'
    # and the remaining elements are a basis for `R' of the form
    # described above.
    # We now calculate a structure constants table of `L' w.r.t. this basis.

    sp:= VectorSpace( F, bb );
    B:= BasisNC( sp, bb );
    T:= List([1..s],x->[]);
    for i in [1..s] do
      for j in [offset*i+1..s] do
        cf:= Coefficients( B, levi[i]*levi[j] ){[1..s]};
        klist:= Filtered([1..s],i->cf[i]<>0);
        cf:= Filtered(cf,x->x<>0);
        T[i][j]:= [klist,cf];
      od;
    od;

    # Now `levi' is a Levi-Malcev subalgebra modulo `R'.
    # The loop modifies this modulo statement.
    # After the first round `levi' will be a Levi-Malcev subalgebra modulo
    # the second element of the derived series.
    # After the second step `levi' will be a Levi-Malcev subalgebra modulo
    # the third element of the derived series.
    # And so on.

    for a in [1..p-1] do

      # `comp' will be a basis of the complement of the `a+1'-th element
      # of the derived series in the `a'-th element of the derived series.
      # `B' will be a basis of the `a'-th term of the derived series,
      # such that the basis elements of the complement come first.
      # So if we have an element v of the `a'-th term of the derived series,
      # then by taking the coefficients w.r.t. `B', we can easily find
      # the part that belongs to `comp'.
      # The equations we have are vector equations in the space `comp',
      # i.e., in the quotient of two elements of the derived series.
      # But we do not want to work with this quotient directly.

      comp:= Rbas{ [ Dimension(ser[a+1])+1 .. Dimension(ser[a]) ] };
      dim:= Length(comp);
      bb:= ShallowCopy( comp );
      for i in [1..Dimension(ser[a+1])] do
        Add(bb,Rbas[i]);
      od;
      sp:= VectorSpace( F, bb );
      B:= BasisNC( sp, bb );

      cf:= List( comp, x -> Coefficients( B, x ){[1..dim]} );

      eqs:= NullMat( s*dim, dim*s*(s-offset)/(offset+1), F );
      rl:= [];
      for i in [1..s] do
        for j in [offset*i+1..s] do
          cij:= T[i][j];
          for k in [1..dim] do

            cf1:= Coefficients(B,levi[i]*comp[k]){[1..dim]};
            cf2:= Coefficients(B,comp[k]*levi[j]){[1..dim]};

            for l in [1..dim] do
              if IsAssociative( L ) then
                eqno:= (i-1)*s*dim+(j-1)*dim+l;
              else
                eqno:= (i-1)*(2*s-i)*dim/2+(j-i-1)*dim+l;
              fi;

              eqs[(j-1)*dim+k][eqno]:= eqs[(j-1)*dim+k][eqno]+cf1[l];
              eqs[(i-1)*dim+k][eqno]:= eqs[(i-1)*dim+k][eqno]+cf2[l];

              for m in [1..Length(cij[1])] do
                r:=cij[1][m];
                if r <= s then
                  eqs[(r-1)*dim+k][eqno]:= eqs[(r-1)*dim+k][eqno]-
                                           cij[2][m]*cf[k][l];
                fi;
              od;
            od;
          od;

          x:= Zero(L);
          for m in [1..Length(cij[1])] do
            if cij[1][m] <= s then
              x:= x+cij[2][m]*levi[cij[1][m]];
            fi;
          od;
          x:= x-levi[i]*levi[j];
          Append(rl,Coefficients(B,x){[1..dim]});

        od;
      od;

      sol:= SolutionMat( eqs, rl );

      if sol = fail then return sol; fi;

      for i in [1..s] do
        for j in [1..dim] do
          levi[i]:=levi[i]+sol[(i-1)*dim+j]*comp[j];
        od;
      od;
    od;

    return [ SubalgebraNC( L, levi, "basis" ), R ];
    end );


############################################################################
##
#M  DirectSumDecomposition( <A> )   ........direct sum decomposition of <A>
##
InstallMethod( DirectSumDecomposition,
        "for semisimple associative algebras",
        [ IsAlgebra and IsAssociative ],
    function( A )
    local R;

    R:= RadicalOfAlgebra( A );
    if Dimension( R ) > 0 then
        TryNextMethod();
    fi;
    return List( CentralIdempotentsOfAlgebra( A ), x -> Ideal( A, [ x ] ) );
end );

[zur Elbe Produktseite wechseln0.82QuellennavigatorsAnalyse erneut starten2026-05-06]

                                                                                                                                                                                                                                                                                                                                                                                                     


Neuigkeiten

     Aktuelles
     Motto des Tages

Software

     Produkte
     Quellcodebibliothek

Aktivitäten

     Artikel über Sicherheit
     Anleitung zur Aktivierung von SSL

Muße

     Gedichte
     Musik
     Bilder

Jenseits des Üblichen ....
    

Besucherstatistik

Besucherstatistik

Monitoring

Montastic status badge