/*
 * Calculate the maximal length of a strictly increasing chain in the subgroup
 * lattice of a finite group.
 *
 * SubgroupLength( grp :: Group : Memory :: Assoc ) -> RngIntElt, Assoc
 *
 * Memory is an associative array that caches the values on simple groups.
 *
 * Example:
 *
 * > len, mem := SubgroupLength( Alt( 18 ) ); len;
 * Calculating length of  Alt(6)
 * Calculating length of  Alt(9)
 * Calculating length of  A(1, 8)[= L(2, 8)]
 * Calculating length of  Alt(7)
 * Calculating length of  Alt(11)
 * Calculating length of  Alt(13)
 * Calculating length of  Alt(14)
 * Calculating length of  Alt(15)
 * Calculating length of  Alt(16)
 * Calculating length of  A(1, 16)[= L(2, 16)]
 * 23
 * > len, mem := SubgroupLength( Alt( 18 ) : Memory := mem ); len;
 * 23
 * > len, mem := SubgroupLength( GL(5,3) : Memory := mem ); len;
 * Calculating length of  A(4, 3)[= L(5, 3)]
 * Calculating length of  C(2, 3)[= S(4, 3)]
 * Calculating length of  A(3, 3)[= L(4, 3)]
 * 21
 *
 */
intrinsic SubgroupLength( grp :: Grp : Memory := AssociativeArray() ) -> RngIntElt, Assoc
{The maximal length of a chain of subgroups of grp. Results for simple groups
are cached in the associative array Memory which is the second return value.}
  local sollen, sublen;
  
  sollen := func< grp | grp eq 1 select 0 else &+[ pn[2] : pn in Factorization( grp ) ] >;

  sublen := procedure( grp, ~res, ~mem : Naive := false )
    len := $$;
    if Type( grp ) eq RngIntElt then res := sollen( grp );
    elif Type( grp ) eq Tup then 
      if not IsDefined( mem, grp ) then
	if grp[1] lt 19 then "Calculating length of ", SimpleGroupName( grp ); end if;
	f := grp[1];
	n := grp[2];
	q := grp[3];
	case f: // One should check the values of top and bot
	  when 1: ext := SL( n+1, q ); top := 1; bot := Gcd( n+1, q-1 );
	  when 2: ext := Omega(2*n+1,q); top := 1; bot := Gcd( 2, q-1 );
	  when 3: ext := Sp(2*n,q); top := 1; bot := Gcd( 2, q-1 );
	  when 4: ext := OmegaPlus(2*n,q); top := 1; bot := Gcd( 2, q-1 );
          when 5: ext := ChevalleyGroup( "G", n, q ); top := 1; bot := 1;
          when 6: ext := ChevalleyGroup( "F", n, q ); top := 1; bot := 1;
          when 7: ext := ChevalleyGroup( "E", n, q ); top := 1; bot := Gcd( 3, q-1 );
          when 8: ext := ChevalleyGroup( "E", n, q ); top := 1; bot := Gcd( 2, q-1 );
          when 9: ext := ChevalleyGroup( "E", n, q ); top := 1; bot := 1;
	  when 10: ext := SU(n+1,q); top := 1; bot := Gcd( n+1, q+1 );
	  when 11: ext := Sz(q); top := 1; bot := 1;
	  when 12: ext := OmegaMinus(2*n,q); top := 1; bot := Gcd( 2, q+1 );
          when 13: ext := ChevalleyGroup( "3D", n, q ); top := 1; bot := 1;
	  when 14: ext := Ree(q); top := 1; bot := 1;
	  when 15: ext := LargeRee(q); top := 1; bot := 1;
          when 16: ext := ChevalleyGroup( "E", n, q ); top := 1; bot := Gcd( 3, q+1 );
	  when 17: ext := Alt(n); top:=1; bot:=1;
          when 18: tmp :=  [ [ 11, 6 ], [ 220, 1 ], [ 231, 1 ], [ 253, 4 ], [ 276, 3 ], [ 266, 1 ], [ 100, 35 ], [ 280, 13 ], [ 275, 1 ], [ 1782, 1 ], [], [], [ 2300, 1 ], [ 276, 4 ], [ 2058, 1 ] ]; ext := PrimitiveGroup( tmp[n][1], tmp[n][2] ); top:=1; bot:=1;
          when 19: ext := q; top:=1; bot:=1;
	end case;
        if not assigned ext then "Error!"; end if;
	len( ext, ~res, ~mem : Naive );
	mem[grp] := res - sollen( top ) - sollen( bot );
      end if;
      res := mem[grp];
    elif IsSolvable( grp ) then res := sollen( #grp );
    elif IsSimple(grp) or Naive then
      res := 0;
      for max in MaximalSubgroups( grp ) do
	len( max`subgroup, ~tmp, ~mem );
	if tmp gt res then res := tmp; end if;
      end for;
      res := res + 1;
      grp := CompositionFactors( grp );
      if #grp eq 1 then mem[grp[1]] := res; end if;
    elif true then 
      res := 0;
      for cf in CompositionFactors( grp ) do
	len( cf, ~tmp, ~mem );
	res +:= tmp;
      end for;
    end if;
  end procedure;

  sublen( grp, ~res, ~Memory );
  return res, Memory;

end intrinsic;

/* Implementation notes:
 *
 * SubgroupLength is additive over composition factors.  This is obvious for
 * solvable groups by considering a composition series which is itself a
 * subgroup chain of maximal length.  In general, it should follow by a version
 * of the Schreier refinement lemma.
 *
 * In particular, for solvable groups we just need to factor the order and in
 * general we just need a table for simple groups.
 *
 * However, no table is available for simple groups (I believe the values for
 * the alternating groups are unknown in general, though asymptotics may be
 * known).  hence we need to cache values as we compute them.
 *
 * To compute for an non-solvable group we try to remove any easy solvable
 * parts. This can be done with SolvableQuotient and SolvableRadical or just
 * with CompositionFactors.  The former has the advantage of retaining groups
 * as the datatype, and the latter has the advantage of being a bit quicker and
 * less memory hungry.
 *
 * Once we get to a group with no easy to remove solvable parts (for instance a
 * simple group), then we use len(G) = 1 + max( len(M) : M max in G ) and
 * Magma's table of maximal subgroups of many simple groups.
 *
 * Using CompositionFactors is a bit of a hassle since we don't have a table of
 * lengths to begin with, and so we really do need groups to compute with.
 * Hence the long switch/case statements that construct a group.
 *
 * These long statements are particularly troublesome in that they take a
 * specification of a simple group as input, and create a "nearly simple" group
 * as output. Rather than use PSL we use SL because magma is good with matrix
 * groups and sometimes can work more efficiently in SL than PSL.  However,
 * that means we have to remove the contribution of Z(SL).  Similarly for the
 * other families. However, not all of the families are completely specified by
 * the documentation, and so removing the "nearly" involves a bit of guess
 * work. I hope I have guessed correctly, and have checked many cases.  In many
 * cases simple things like IsPerfect(grp) or AbelianQuotient(grp) fail, so I
 * cannot check that magma's conventions agree with my own.
 */
