POL_ME_mat2vec := function( mat )
    local a;
    a := mat[1]{[2..8]};
    return a;
end;

POL_ME_vec2mat := function( vec )
    local mat,i;
    mat := IdentityMat(8);
    for i in [1..7] do
        mat[1][i+1] := vec[i];
    od;
    return mat;
end;

POL_ME_S2Lvector := function( vec )
    local n,lvec,i;
    n := Length( vec );
    lvec := [];
    for i in [1..n] do 
        Append( lvec, Coefficients( Basis( CF(7) ), vec[i] ));
    od; 
    return lvec;
end;

POL_ME_S2Lvectors := function( vecs )
    local lvecs,v;
    lvecs := [];
    for v in vecs do 
        Add( lvecs, POL_ME_S2Lvector( v ));
    od; 
    return lvecs;
end;

POL_ME_L2Svector := function( vec )
    local v,i,l;
    v := [];
    for i in [0..6] do
        l := vec{[(i*6)+1..(i*6)+6]}; 
        Add( v, LinearCombination( Basis( CF(7) ), l ));
    od;
    return v;
end;

POL_ME_CPCS_K_p := function( gens )
    local vecs,g,vecsL,latbas;
    vecs := [];
    for g in gens do
        Add( vecs, POL_ME_mat2vec( g ) );
    od;
    vecsL := POL_ME_S2Lvectors( vecs );
    latbas := LatticeBasis( vecsL );
    return latbas;
end;


POL_ME_CPCS_K_pClosedConj := function( gens, gens_K_p)
    local latbas,vecL, gens2,gens_K_p2,g,h,vec,x;
    latbas := POL_ME_CPCS_K_p( gens_K_p );
    gens2 :=  StructuralCopy( gens );
    Append( gens2, List(gens, x-> x^-1) );
    gens_K_p2 := StructuralCopy( gens_K_p );
    for g in gens_K_p2 do
        for h in gens2 do 
            vec := POL_ME_mat2vec( g^h);
            vecL :=  POL_ME_S2Lvector( vec);
            x := SolutionIntMat( latbas, vecL );
            if x = fail then 
                Add( gens_K_p2, g^h );
                latbas := POL_ME_CPCS_K_p( gens_K_p2 );
                #Print( "Extended\n" );
            fi;
        od;
    od;
    return latbas;
end;

POL_ME_2finite := function( mat )
    local n,matf,a,l,i,j,vec;
    # maps 1 to Z(2)^0 and E(7) to Z(2^3)
    n := Length( mat );
    matf := IdentityMat( n , GF(2^3) );
    a := Z(2^3 );
    l := List( [1..6], i-> a^i );
    for i in [1..n] do
        for j in [1..n] do 
            vec :=  Coefficients( Basis( CF(7) ), mat[i][j] );
            matf[i][j] := LinearCombination( l, vec );
        od;
    od;
    return matf;
end;


POL_ME_CPCS := function( gens )
    local pcgs_I_p,gens_K_p,pcgs_K_p,preImgsI_p,r,s,r2,s2,v,vS,mat,
          pcs_K_p, pcs;
    r := gens[1];
    s := gens[2];
    r2 := POL_ME_2finite( r );
    s2 := POL_ME_2finite( s );
    pcgs_I_p := CPCS_finite_word( [r2,s2], 10 );
    gens_K_p := POL_NormalSubgroupGeneratorsOfK_p( pcgs_I_p, [r,s] );
    gens_K_p := Filtered( gens_K_p, x -> not x = IdentityMat(8) );
    pcgs_K_p := POL_ME_CPCS_K_pClosedConj( [r,s], gens_K_p);

    preImgsI_p := POL_PreImagesPcsI_p_G( pcgs_I_p, [r,s]);
    preImgsI_p := Reversed( preImgsI_p );
    pcs_K_p := [];
    for v in pcgs_K_p do
         vS := POL_ME_L2Svector( v );
         mat := POL_ME_vec2mat( vS );
         Add( pcs_K_p, mat ); 
    od;
    pcs := Concatenation( preImgsI_p, pcs_K_p );
    return rec( pcgs_I_p := pcgs_I_p, preImgsI_p := preImgsI_p,
                pcgs_K_p := pcgs_K_p , pcs := pcs);
end;


POL_ME_ExpVector :=  function( pcgs, g)
    local h,exp_h,k,exp_k, vec, vecL; 

    # first we have to compute the part related to preImgs_I_p
        # compute the image in I_p
        h := POL_ME_2finite( g );
        exp_h :=  ExponentvectorPcgs_finite( pcgs.pcgs_I_p, h );
        if exp_h = fail then return fail; fi;

        # divide off to get the part of g which is in K_p
        k := POL_GetPartinK_P( g, exp_h, pcgs.preImgsI_p );
    # now compute the  exponent of the part in K_p
        vec := POL_ME_mat2vec( k );
        vecL :=  POL_ME_S2Lvector( vec);
        exp_k := SolutionIntMat( pcgs.pcgs_K_p, vecL );
   
    exp:=Concatenation(exp_h,exp_k);
    Assert( 2, Exp2Groupelement( pcgs.pcs, exp) = g,
            "error in ExponentVector_CPCS_PRMGroup");
    return exp;
end ;

POL_ME_RelativeOrders :=  function( pcgs )
    local relOrders_I_p,relOrders_K_p,n,relOrders,l;
    relOrders_I_p := RelativeOrdersPcgs_finite( pcgs.pcgs_I_p );
    l := Length( pcgs.pcgs_K_p );
    relOrders_K_p := List( [1..l], x-> 0 );
    # merge the relOrders
    relOrders := Concatenation(relOrders_I_p,relOrders_K_p);
    return relOrders;
end ;

POL_ME_TestExp := function( pcgs, gens )
    local numberOfTests, i, g, exp;
    Print( "Start testing\n" );
    SetAssertionLevel( 2 );
    numberOfTests := 20;
    for i in [1..numberOfTests] do
        Print(i);
        g := POL_RandomGroupElement( gens );
        #Print( "g is equal to ", g ,"\n");
        exp := POL_ME_ExpVector(  pcgs, g);
     od;
    Print( "\n" );
end;


POL_ME_SetPcPresentation:= function(pcgs)
    local genList,ftl,n,ro,i,j,exp,conj,f_i,f_j,r_i,pcsInv;
    # Setup
    n:=Length(pcgs.pcs);
    ftl:=FromTheLeftCollector(n);
    pcsInv:=[];
    for i in [1..n] do
        pcsInv[i]:=pcgs.pcs[i]^-1;
    od;
    # the relative orders
    ro:= POL_ME_RelativeOrders(pcgs);
    for i in [1..n] do
        if ro[i]<>0 then
            SetRelativeOrder(ftl,i,ro[i]);
        fi;
    od;
    # Set power relations
    for i in [1..n] do
        if ro[i]<>0 then
            f_i:=pcgs.pcs[i];
            r_i:=ro[i];
            exp:=POL_ME_ExpVector( pcgs, f_i^r_i );
            genList := POL_Exp2GenList(exp);
            SetPower(ftl,i,genList);
        fi;
    od;
    # Set the conjugation relations
    
    for i in [1..n] do
        for j in [1..(i-1)] do
            f_i:=pcgs.pcs[i];
            f_j:=pcgs.pcs[j];
            conj:=(pcsInv[j])*f_i*f_j;
            if not conj = f_i then
                exp:=POL_ME_ExpVector(pcgs, conj);
                genList:=POL_Exp2GenList(exp);
                SetConjugate(ftl,i,j,genList);
            fi;
        od;
    od;
    UpdatePolycyclicCollector(ftl);
    return ftl;
end;

PcpGroupOfMarcosExam := function()
    local CPCS, pcp, K,r,s;
    r :=
    [ [ 1, 0, 0, 0, 0, 0, 0, 0 ], [ 0, 0, 1, 0, 0, 0, 0, 0 ],
      [ 0, 0, 0, 1, 0, 0, 0, 0 ], [ 0, 0, 0, 0, 1, 0, 0, 0 ],
      [ 0, 0, 0, 0, 0, 1, 0, 0 ], [ 0, 0, 0, 0, 0, 0, 1, 0 ],
      [ 0, 0, 0, 0, 0, 0, 0, 1 ], [ 0, 1, 0, 0, 0, 0, 0, 0 ] ];

    s :=
    [ [ 1, 0, 0, 0, 0, 0, 0, 1 ], [ 0, E(7)^6, 0, 0, 0, 0, 0, 0 ],
      [ 0, 0, E(7), 0, 0, 0, 0, 0 ], [ 0, 0, 0, 1, 0, 0, 0, 0 ],
      [ 0, 0, 0, 0, 1, 0, 0, 0 ], [ 0, 0, 0, 0, 0, 1, 0, 0 ],
      [ 0, 0, 0, 0, 0, 0, 1, 0 ], [ 0, 0, 0, 0, 0, 0, 0, 1 ] ];
    CPCS := POL_ME_CPCS( [r,s] ); 
    pcp := POL_ME_SetPcPresentation( CPCS );
    K := PcpGroupByCollector( pcp );
    return K;
end;

PCentralSeriesSteps := function(G,p,n)
    local ser, U, C, pcp, new, igs,i;
    ser := [G];
    U := G;
    G!.isNormal := true;
    for i in [1..n] do
        C := CommutatorSubgroup( U, G );
        pcp := Pcp( U, C );
        new := List( pcp, x -> x^p );
        igs := AddToIgs(Igs(C), new);
        U := SubgroupByIgs(G, Cgs(igs));
        Add( ser, U );
    od;
    Unbind( G!.isNormal );
    return ser;
end;