Quelle forms.gi Sprache: unbekannt

#############################################################################
##
##  forms.gi              'Forms' package
##                                                              John Bamberg
##                                                              Jan De Beule
##
##  Copyright 2024, Vrije Universiteit Brussel
##  Copyright 2024, The University of Western Austalia
##
##  Implementation of quadratic and sesquilinear forms
##
#############################################################################

#############################################################################
# Fundamental methods:
#############################################################################

InstallMethod( \=, "for two trivial forms",
  [IsTrivialForm and IsFormRep, IsTrivialForm and IsFormRep],
  function( a, b )
    return a!.basefield = b!.basefield;
  end );

InstallMethod( \=, "for two forms",
  [IsForm and IsFormRep, IsForm and IsFormRep],
  function( a, b )
    return (a!.basefield = b!.basefield) and
           (a!.type = b!.type) and
           (a!.matrix = b!.matrix);
  end );

#############################################################################
# Constructor methods
# Form-by-matrix methods
#############################################################################
#############################################################################
#O  FormByMatrix( <mat>, <field>, <string> )
#   A general constructor not for users, not documented.
##
InstallMethod( FormByMatrix, "for a ffe matrix, a field and a string",
  [IsMatrix and IsFFECollColl, IsField and IsFinite, IsString],
  function( m, f, string )
    local el;
    m := ImmutableMatrix(f, m);
    el := rec( matrix := m, basefield := f, type := string );

   ## We follow a certain convention, which is outlined in the manual,
   ## in order to determine from a Gram matrix the type of the constructed
   ## form.

    if string = "hermitian" then
      if not IsInt(Sqrt(Size(f))) then
        Error("No hermitian form exist when the order of <f> is not a square" );
      fi;
      if FORMS_IsHermitianMatrix(m,f) then
        Objectify(NewType( HermitianFormFamily ,  IsFormRep),  el);
        return el;
      else
        Error("Given matrix does not define a hermitian form" );
      fi;
    elif string = "symplectic" then
      if FORMS_IsSymplecticMatrix(m,f) then
        Objectify(NewType( BilinearFormFamily ,  IsFormRep),  el);
        return el;
      else
        Error("Given matrix does not define a symplectic form" );
      fi;
    elif string = "orthogonal" then
      if Characteristic(f) = 2 then
        Error("No orthogonal forms exist in even characteristic" );
      fi;
      if FORMS_IsSymmetricMatrix(m) then
        Objectify(NewType( BilinearFormFamily ,  IsFormRep),  el);
        return el;
      else
        Error("Given matrix does not define an orthogonal form" );
      fi;
    elif string = "pseudo" then
      if IsOddInt(Size(f)) then
        Error("No pseudo forms exist in even characteristic" );
      fi;
      if (FORMS_IsSymmetricMatrix(m) and (not FORMS_IsSymplecticMatrix(m))) then
        Objectify(NewType( BilinearFormFamily ,  IsFormRep),  el);
        return el;
      else
        Error("Given matrix does not define a pseudo-form" );
      fi;
    elif string = "quadratic" then
      el.matrix := ImmutableMatrix(f, Forms_RESET(m));
      Objectify(NewType( QuadraticFormFamily ,  IsFormRep),  el);
      return el;
    else
      Error("Please specify a form properly");
    fi;
  end );

#############################################################################
#O  BilinearFormByMatrixOp( <m>, <f> )
#   A less general constructor not for users, not documented.
# <f> finite field, <m> orthogonal or symplectic matrix over <f>.
# zero matrix is allowed, then the trivial form is returned.
# if <m> is not zero, it is checked that <m> determines a bilinear form which is
# symplectic, orthogonal, pseudo or hermitian
##
InstallMethod( BilinearFormByMatrixOp, "for a ffe matrix and a field",
  [IsMatrix and IsFFECollColl, IsField and IsFinite],
  function( m, f )
    local el, n;
    n := NrRows(m);
    m := ImmutableMatrix(f, m);
    if IsZero(m) then
       el := rec( matrix := m, basefield := f, type := "trivial", vectorspace := FullRowSpace(f,n) );
       Objectify(NewType( TrivialFormFamily ,  IsFormRep),  el);
       return el;
    elif FORMS_IsSymplecticMatrix(m,f) then
       el := rec( matrix := m, basefield := f, type := "symplectic", vectorspace := FullRowSpace(f,n) );
       Objectify(NewType( BilinearFormFamily ,  IsFormRep),  el);
       return el;
    elif FORMS_IsSymmetricMatrix(m) and Characteristic(f) <> 2 then
       el := rec( matrix := m, basefield := f, type := "orthogonal", vectorspace := FullRowSpace(f,n) );
       Objectify(NewType( BilinearFormFamily ,  IsFormRep),  el);
       return el;
    elif FORMS_IsSymmetricMatrix(m) and Characteristic(f) = 2 then
       el := rec( matrix := m, basefield := f, type := "pseudo", vectorspace := FullRowSpace(f,n) );
       Objectify(NewType( BilinearFormFamily ,  IsFormRep),  el);
       return el;
    else
       Error("Invalid Gram matrix");
    fi;
  end );

# 20/01/2016: small change here: added MutableCopyMat, since it is not logical
# that changing the matrix that was used for construction afterwards changes
# the form. See test_forms12.g
#############################################################################
#O BilinearFormByMatrix( <m>, <f>): constructor for users.
# <f> finite field, <m> orthogonal or symplectic matrix over <f>.
# checks whether <m> is over <f>.
# zero matrix is allowed, then the trivial form is returned.
##
InstallMethod( BilinearFormByMatrix, "for a ffe matrix and a field",
  [IsMatrix and IsFFECollColl, IsField and IsFinite],
  function( m, f )
    local gf;
    gf := DefaultFieldOfMatrix(m);
    if not PrimitiveElement(gf) in f then
      Error("<m> is not a matrix over <f>");
    fi;
    return BilinearFormByMatrixOp( m, f);
end );

#############################################################################
#O BilinearFormByMatrix( <m> ): constructor for users. uses above constructor
# finite field <f> is determined from <m>
##
InstallMethod( BilinearFormByMatrix, "for a ffe matrix ",
  [IsMatrix and IsFFECollColl ],
  function( m )
  local f;
  f := DefaultFieldOfMatrix(m);
  return BilinearFormByMatrixOp( m, f);
end );

#############################################################################
#O QuadraticFormByMatrixOp( <m>, <f> ): constructor not for users.
#  Analogous to BilinearFormByMatrixOp
# <f> finite field, <m> matrix over <f>.
# <m> is always "Forms_RESET" to an upper triangle matrix.
# zero matrix is allowed, then the trivial form is returned.
##
InstallMethod( QuadraticFormByMatrixOp, "for a ffe matrix and a field",
  [IsMatrix and IsFFECollColl, IsField and IsFinite],
  function( m, f )
    local el, n;
    n := NrRows(m);
    m := ImmutableMatrix(f, m);
    el := rec( matrix := m, basefield := f, type := "quadratic", vectorspace := FullRowSpace(f,n) );
    if IsZero(m) then
       el.type := "trivial";
       Objectify(NewType( TrivialFormFamily ,  IsFormRep),  el);
       return el;
    else
       el.matrix := ImmutableMatrix(f, Forms_RESET(m));
       Objectify(NewType( QuadraticFormFamily ,  IsFormRep),  el);
       return el;
    fi;
  end );

# 20/01/2016: small change here: added MutableCopyMat, since it is not logical
# that changing the matrix that was used for construction afterwards changes
# the form. See test_forms12.g
#############################################################################
#O QuadraticFormByMatrix( <m>, <f> ): constructor for users.
# <f> finite field, <m> matrix over <f>.
# checks whether <m> is over <f>.
# zero matrix is allowed.
##
InstallMethod( QuadraticFormByMatrix, "for a ffe matrix and a field",
  [IsMatrix and IsFFECollColl, IsField and IsFinite],
  function( m, f )
    local gf;
    gf := DefaultFieldOfMatrix(m);
    if not PrimitiveElement(gf) in f then
      Error("<m> is not a matrix over <f>");
    fi;
    return QuadraticFormByMatrixOp( m, f);
end );

#############################################################################
#O QuadraticFormByMatrix( <m> ): constructor for users. uses above constructor
# finite field <f> is determined from <m>
##
InstallMethod( QuadraticFormByMatrix, "for a ffe matrix ",
  [IsMatrix and IsFFECollColl],
  function( m )
  local f;
  f := DefaultFieldOfMatrix(m);
  return QuadraticFormByMatrixOp( m, f);
end );

# 20/01/2016: small change here: added MutableCopyMat, since it is not logical
# that changing the matrix that was used for construction afterwards changes
# the form. See test_forms12.g
#############################################################################
#O HermitianFormByMatrix( <m>, <f>): constructor  for users.
# <f> finite field of square order, <m> hermitian matrix over <f>.
# zero matrix is allowed, then the trivial form is returned.
##
InstallMethod( HermitianFormByMatrix, "for a ffe matrix and a field",
  [IsMatrix and IsFFECollColl, IsField and IsFinite],
  function( m, f )
    local el,gf,n;
    n := NrRows(m);
    gf := DefaultFieldOfMatrix(m);
    if not PrimitiveElement(gf) in f then
      Error("<m> is not a matrix over <f>");
    fi;
    if not IsInt(Sqrt(Size(f))) then
        Error("No hermitian form exists when the order of <f> is not a square" );
    fi;
    if FORMS_IsHermitianMatrix(m,f) then
       m := ImmutableMatrix(f, m);
       el := rec( matrix := m, basefield := f, type := "hermitian", vectorspace := FullRowSpace(f,n) );
       Objectify(NewType( HermitianFormFamily ,  IsFormRep),  el);
       return el;
    else
       Error("Given matrix does not define a hermitian form" );
    fi;
  end );

#############################################################################
#O UpperTriangleMatrixByPolynomialForForm( <pol>,<ff>,<int>,<list>),
# <ff>: fin field, <pol>: polynomial over <ff>, <n>: size of the matrix to
# construct; <list>: used variables in <pol>.
# not for users, creates a Gram matrix from a given polynomial to construct a
# quadratic form, and also an orthogonal bilinear form in odd characteristic.
##
InstallMethod( UpperTriangleMatrixByPolynomialForForm,
                      [IsPolynomial, IsField, IsInt, IsList],
  function(poly, gf, n, varlist)
    local vars, mat, i, j, vals, der;

   ## UpperTriangleMatrixByPolynomialForForm returns an upper triangular
   ## matrix in all cases. It can be used for quadratic forms (all chars)
   ## and symmetric bilinear forms (odd char). This is not a restriction,
   ## since in even char, we have no orthogonal bilinear forms, as by
   ## convention, symmetric bilinear forms in odd char, hermitian forms (all char)
   ## and quadratic forms (all char) can be constructed using polynomials.
   ##  Extermely important: the polynomial of a symmetric bilinear form is
   ##  NOT the polynomial of the associated quadratic form!

    mat := NullMat(n, n, gf);
    vars := ShallowCopy(varlist);
    if IsZero(poly) then
       return mat;
    fi;
    n := Length(varlist);
    vals := List([1..n], i -> 0);
    for i in [1..n] do
      vals[i] := 1;
      mat[i,i] := Value(poly, vars, vals);
      vals[i] := 0;
    od;
    for i in [1..n-1] do
      der := Derivative(poly, vars[i]);
      for j in [i+1..n] do
        vals[j] := 1;
        mat[i,j] := Value(der,vars,vals);
        vals[j] := 0;
      od;
    od;
    vars := ShallowCopy(varlist);
    if vars*mat*vars <> poly then
      Error( "<poly> should be a homogeneous polynomial over GF(q)");
    fi;
    return mat;
end);

#############################################################################
#O GramMatrixByPolynomialForHermitianForm( <pol>,<ff>,<int>,<list>),
# <ff>: fin field, <pol>: polynomial over <ff>, <n>: size of the matrix to
# construct; <list>: used variables in <pol>.
# not for users, creates a Gram matrix from a given polynomial to construct a
# hermitian form.
##
InstallMethod( GramMatrixByPolynomialForHermitianForm,
                     [IsPolynomial, IsField, IsInt, IsList],
  function(poly, gf, n, varlist)
    local vars, varst, q, t, a, polarity, i, j, vals, der;
    q := Size(gf);
    t := Sqrt(q);
    vars := ShallowCopy(varlist);
    varst := List(vars, x -> x^t);
    polarity := NullMat(n,n,gf);
    if IsZero(poly) then
      return polarity;
    fi;
    n := Length(varlist);
    vals := [];
    vals := List([1..n], i -> 0);
    for i in [1..n] do
      vals[i] := 1;
      a := Value(poly,vars,vals);
      if a <> a^t then
        Error( "<poly> does not generate a Hermitian matrix" );
      else
        polarity[i,i] := a;
      fi;
      vals[i] := 0;
    od;
    for i in [1..n-1] do
      der := Derivative(poly,vars[i]);
      for j in [i+1..n] do
          vals[j] := 1;
          polarity[i,j] := Value(der,vars,vals);
          polarity[j,i] := polarity[i,j]^t;
          vals[j] := 0;
      od;
    od;
    vars := ShallowCopy(varlist);
      #Check whether the polynomial has the right form.
    if vars*polarity*varst <> poly then
      Error( "<poly> does not generate a Hermitian matrix" );
    fi;
    return polarity;
  end );

#############################################################################
#O BilinearFormByPolynomial( <pol>, <ring>, <int>): constructor for users.
#  <ring>: polynomial ring over finite field, <po>: polyniomial in <ring>, <int>:
#  desired dimension of the vectorspace on which the form acts.
#  an error is returned by UpperTriangleMatrixByPolynomialForForm if the given
#  <pol> is not suitable to define a bilinear form
#  zero <pol > is allowed, then the trivial form is returned.
##
InstallMethod( BilinearFormByPolynomial, "for a polynomial over a field, and a dimension",
  [IsPolynomial, IsFiniteFieldPolynomialRing, IsInt],
  function( pol, pring, n )
    local mat, form, gf, vars, polarity,el;
    gf := CoefficientsRing( pring );
    vars := IndeterminatesOfPolynomialRing( pring );
    if IsZero(pol) then
      mat := NullMat(n, n, gf);
      mat := ImmutableMatrix(gf, mat);
      el := rec( matrix := mat, basefield := gf, type := "trivial" );
      Objectify(NewType( TrivialFormFamily ,  IsFormRep),  el);
      SetPolynomialOfForm(el, pol);
      return el;
    fi;
    if Characteristic(gf) = 2 then
       Error("No orthogonal form can be associated with a quadratic polynomial in even characteristic");
    else
       mat := UpperTriangleMatrixByPolynomialForForm(pol,gf,n,vars);
       polarity := (mat + TransposedMat(mat))/(One(gf)*2);
       form := FormByMatrix(polarity,gf,"orthogonal");
       SetPolynomialOfForm(form, pol);
       return form;
    fi;
  end );

#############################################################################
#O BilinearFormByPolynomial( <pol>, <ring> ): constructor for users.
#  <ring>: polynomial ring over finite field, <po>: polyniomial in <ring>
#  dimension is the length of variables in <ring>.
##
InstallMethod( BilinearFormByPolynomial,  "no dimension",
  [IsPolynomial, IsFiniteFieldPolynomialRing],
  function( pol, pring )
    local n;
    n := Length(IndeterminatesOfPolynomialRing( pring ));
    return BilinearFormByPolynomial( pol, pring, n);
  end );

#############################################################################
#O HermitianFormByPolynomial( <pol>, <ring>, <int>): constructor for users.
#  <ring>: polynomial ring over finite field, <po>: polyniomial in <ring>, <int>:
#  desired dimension of the vectorspace on which the form acts.
#  an error is returned by GramMatrixByPolynomialForHermitianForm if the given
#  <pol> is not suitable to define a hermitian form
#  zero <pol > is allowed, then the trivial form is returned.
##
InstallMethod( HermitianFormByPolynomial, "for a polynomial over a field, and a dimension",
  [IsPolynomial, IsFiniteFieldPolynomialRing, IsInt],
  function(pol, pring, n)
    local mat, form, gf, vars, el;
    gf := CoefficientsRing( pring );
    if not IsInt(Sqrt(Size(gf))) then
        Error("the order of the underlying field of <r> must be a square" );
    fi;
    if IsZero(pol) then
      mat := NullMat(n, n, gf);
      mat := ImmutableMatrix(gf, mat);
      el := rec( matrix := mat, basefield := gf, type := "trivial" );
      Objectify(NewType( TrivialFormFamily ,  IsFormRep),  el);
      SetPolynomialOfForm(el, pol);
      return el;
    fi;
    vars := IndeterminatesOfPolynomialRing( pring );
    mat := GramMatrixByPolynomialForHermitianForm(pol,gf,n,vars);
    form := FormByMatrix(mat,gf,"hermitian");
    SetPolynomialOfForm(form, pol);
    return form;
  end );

#############################################################################
#O HermitianFormByPolynomial( <pol>, <ring> ): constructor for users.
#  <ring>: polynomial ring over finite field, <pol>: polyniomial in <ring>
#  dimension is the length of variables in <ring>.
##
InstallMethod( HermitianFormByPolynomial,  "no dimension",
  [IsPolynomial, IsFiniteFieldPolynomialRing],
  function( pol, pring )
    local n;
    n := Length(IndeterminatesOfPolynomialRing( pring ));
    return HermitianFormByPolynomial( pol, pring, n);
  end );

#############################################################################
#O QuadraticFormByPolynomial( <pol>, <ring>, <int>): constructor for users.
#  <ring>: polynomial ring over finite field, <po>: polyniomial in <ring>, <int>:
#  desired dimension of the vectorspace on which the form acts.
#  an error is returned by UpperTriangleMatrixByPolynomialForForm if the given
#  <pol> is not suitable to define a quadratic form
#  zero <pol > is allowed, then the trivial form is returned.
##
InstallMethod( QuadraticFormByPolynomial,
  "for a polynomial over a field, and a dimension",
  [IsPolynomial, IsFiniteFieldPolynomialRing, IsInt],
  function(pol, pring, n)
    local mat, form, gf, vars;
    gf := CoefficientsRing( pring );
    vars := IndeterminatesOfPolynomialRing( pring );
    mat := UpperTriangleMatrixByPolynomialForForm(pol,gf,n,vars);
    form := FormByMatrix(mat,gf,"quadratic");
    SetPolynomialOfForm(form, pol);
    return form;
  end );

#############################################################################
#O QuadraticFormByPolynomial( <pol>, <ring> ): constructor for users.
#  <ring>: polynomial ring over finite field, <pol>: polyniomial in <ring>
#  dimension is the length of variables in <ring>.
##
InstallMethod( QuadraticFormByPolynomial,  "no dimension",
  [IsPolynomial, IsFiniteFieldPolynomialRing],
  function( pol, pring )
    local n;
    n := Length(IndeterminatesOfPolynomialRing( pring ));
    return QuadraticFormByPolynomial( pol, pring, n);
  end );

#############################################################################
# Form-by-form methods (for users):
# see Package documentation for the interpretation of these functions.
#############################################################################

#############################################################################
#O BilinearFormByQuadraticForm( <form> ): constructor for users.
#  <form>: quadratic form.
##
InstallMethod( BilinearFormByQuadraticForm, [IsQuadraticForm],
  function(f)
    ## This method constructs a bilinear form from a quadratic form.
    ## very important: we use the relation 2Q(v) = f(v,v)
    local m, gf;
    m := f!.matrix;
    gf := f!.basefield;
    if Characteristic(gf) = 2 then
      Error( "No orthogonal bilinear form can be associated with a quadratic form in even characteristic");
    else
      return BilinearFormByMatrix((m+TransposedMat(m))/(One(gf)*2),gf);
    fi;
  end );

#############################################################################
#O BilinearFormByQuadraticForm( <form> ): constructor for users.
#  <form>: bilinear form.
##
InstallMethod( QuadraticFormByBilinearForm, [IsBilinearForm],
  function(f)
    ## This method constructs a quadratic form from a bilinear form.
    ## very important: we use the relation 2Q(v) = f(v,v)
    local m, gf;
    m := f!.matrix;
    gf := f!.basefield;
    if Characteristic(gf) = 2 then
       Error( "No quadratic form can be associated to a symmetric form in even characteristic");
    elif IsAlternatingForm(f) then
      Error( "No quadratic form can be associated properly to an alternating form" );
    else
      return QuadraticFormByMatrix(m,gf);
    fi;
  end );

#############################################################################
#  Attributes ... (for users).
#  see package documentation for more information.
##
InstallMethod( PolynomialOfForm, "for a trivial form",
  [IsTrivialForm],
  function(f)
    ## returns zero polynomial of the corresponding polynomialring.
    local gf, d, r;
    gf := f!.basefield;
    d := NrRows(f!.matrix);
    r := PolynomialRing(gf,d);
    return Zero(r);
  end );

InstallMethod( PolynomialOfForm, "for a quadratic form",
  [IsQuadraticForm],
  function(f)
    ## This method finds the polynomial associated to
    ## the Gram matrix of f.
    local m, gf, d, r, indets, poly;
    m := f!.matrix;
    gf := f!.basefield;
    d := NrRows(m);
    r := PolynomialRing(gf,d);
    indets := IndeterminatesOfPolynomialRing(r);
    poly := indets * m * indets;
    return poly;
  end );

InstallMethod( PolynomialOfForm, "for a hermitian form",
  [IsHermitianForm],
  function(f)
    ## This method finds the polynomial associated to
    ## the Gram matrix of f.
    local m, gf, d, r, indets, poly, q;
    gf := f!.basefield;
    q := Sqrt(Size(gf));
    m := f!.matrix;
    d := NrRows(m);
    r := PolynomialRing(gf,d);
    indets := IndeterminatesOfPolynomialRing(r);
    poly := indets * m * List(indets,t->t^q);
    return poly;
  end );

InstallMethod( PolynomialOfForm, "for a bilinear form",
  [IsBilinearForm],
  function(f)
    ## This method finds the polynomial associated to
    ## the Gram matrix of f. For a symplectic form, we
    ## get the 0 polynomial.
    ## Very important: for an orthogonal form, the polynomial is
    ## NOT the polynomial of the associated quadratic form.
    local m, gf, d, r, indets, poly;
    m := f!.matrix;
    gf := f!.basefield;
    d := NrRows(m);
    if Characteristic(gf) = 2 then
       Error( "No polynomial can be (naturally) associated to a bilinear form in even characteristic");
    fi;
    r := PolynomialRing(gf,d);
    indets := IndeterminatesOfPolynomialRing(r);
    poly := indets * m * indets;
    return poly;
  end );

InstallOtherMethod( BaseField, "for a form", [IsForm],
  function( f )
    return f!.basefield;
  end );

InstallMethod( GramMatrix, "for a form", [IsForm],
  function( f )
    return f!.matrix;
  end );

InstallMethod( CompanionAutomorphism, [ IsSesquilinearForm ],
  function( form )
    local field, aut;
    field := BaseField( form );
    if IsHermitianForm( form ) then
       aut := FrobeniusAutomorphism( field );
       aut := aut ^ (Order(aut)/2);
    else
       aut := IdentityMapping( field );
    fi;
    return aut;
  end );

InstallMethod( AssociatedBilinearForm, "for a quadratic form",
  [ IsQuadraticForm ],
  function( form )
    local gram, f;
    gram := form!.matrix;
    f := form!.basefield;
    return BilinearFormByMatrix( gram + TransposedMat(gram), f);
  end );

##next two operations are not documented. Could be in next versions.

InstallMethod( RadicalOfFormBaseMat, [IsSesquilinearForm],
  function( f )
    local m, gf, d;
    if not IsReflexiveForm( f ) then
       Error( "Form must be reflexive");
    fi;
    m := f!.matrix;
    gf := f!.basefield;
    d := NrRows(m);
    return NullspaceMat( m );
  end );

InstallMethod( RadicalOfFormBaseMat, [IsQuadraticForm],
  function( f )
    local m, null, gf, d;
    m := f!.matrix;
    m := m + TransposedMat(m);
    gf := f!.basefield;
    d := NrRows(m);
    null := NullspaceMat( m );
    if Characteristic(gf) = 2 then
      null := Filtered(SubspaceNC(gf^d,null), x -> IsZero(x^f)); #find vectors vanishing under f
    fi;
    null := Filtered(null,x-> not IsZero(x));
    return null;
  end );

InstallMethod( RadicalOfForm, "for a sesquilinear form",
  [IsSesquilinearForm],
  function( f )
    local m, null, gf, d;
    if not IsReflexiveForm( f ) then
       Error( "Form must be reflexive");
    fi;
    m := f!.matrix;
    gf := f!.basefield;
    d := NrRows(m);
    null := NullspaceMat( m );
    return Subspace( gf^d, null, "basis" );
  end );

InstallMethod( RadicalOfForm, "for a quadratic form",
  [IsQuadraticForm],
  function( f )
    local m, null, gf, d;
    m := f!.matrix;
    m := m + TransposedMat(m);
    gf := f!.basefield;
    d := NrRows(m);
    null := NullspaceMat( m );
    if Characteristic(gf) = 2 then
      null := Filtered(SubspaceNC(gf^d,null), x -> IsZero(x^f)); #find vectors vanishing under f
    fi;
    null := Filtered(null,x-> not IsZero(x));
    return Subspace( gf^d, null );
  end );

InstallMethod( RadicalOfForm, "for a trivial form",
  [IsTrivialForm],
  function( f )
    local m, null, gf, d;
    m := f!.matrix;
    gf := f!.basefield;
    d := Size(m);
    null := NullspaceMat( m );
    return Subspace( gf^d, null, "basis" );
  end );

InstallMethod( DiscriminantOfForm, [ IsQuadraticForm ],
function( f )
   local m, gf, d, det, squares, primroot;
   m := f!.matrix;
   gf := f!.basefield;
   d := Size(m);
   if IsOddInt(d) then
      Error( "Quadratic form must be defined by a matrix of even dimension");
   fi;
   det := Determinant( m + TransposedMat(m) );
   if IsZero( det ) then
      Error( "Form must be nondegenerate" );
   fi;
   if IsOddInt(Size(gf)) then
      primroot := PrimitiveRoot( gf );
      squares := AsList( Group( primroot^2 ) );
      if det in squares then
         return "square";
      else
         return "nonsquare";
      fi;
   else
      return "square";
   fi;
end );

InstallMethod( DiscriminantOfForm, [ IsSesquilinearForm ],
function( f )
   local m, gf, d, det, squares, primroot;
   m := f!.matrix;
   gf := f!.basefield;
   d := Size(m);
   if IsOddInt(d) then
      Error( "Sesquilinear form must be defined by a matrix of even dimension");
   fi;
   det := Determinant( m );
   if IsZero( det ) then
      Error( "Form must be nondegenerate" );
   fi;
   if IsOddInt(Size(gf)) then
      primroot := PrimitiveRoot( gf );
      squares := AsList( Group( primroot^2 ) );
      if det in squares then
         return "square";
      else
         return "nonsquare";
      fi;
   else
      return "square";
   fi;
end );

InstallMethod( DiscriminantOfForm, [ IsHermitianForm ],
function( f )
   Error( "The discriminant of a hermitian form is not defined");
end );

InstallMethod( DiscriminantOfForm, [ IsTrivialForm ],
function( f )
   Error( "<form> must be sesquilinear or quadratic");
end );

####
#  ... and Properties (for users).
#  see package documentation for more information.
####

InstallMethod( IsDegenerateForm, [IsSesquilinearForm],
  function( f )
    return not IsTrivial( RadicalOfForm( f ) );
  end );

InstallMethod( IsSingularForm, [IsQuadraticForm],
  function( f )
    return not IsTrivial( RadicalOfForm( f ) );
  end );

InstallMethod( IsSingularForm, [IsTrivialForm], #new in 1.2.1
  function( f )
  return true;
  end );

InstallMethod( IsDegenerateForm, [IsQuadraticForm],
  function( f )
    Info(InfoWarning,1,"Testing degeneracy of the *associated bilinear form*");
    return not IsTrivial( RadicalOfForm( AssociatedBilinearForm( f ) ) );
  end );

InstallMethod( IsDegenerateForm, [IsTrivialForm],
  function( f )
    return true;
  end );

InstallMethod( IsReflexiveForm, [IsBilinearForm],
  function( f )
    ## simply check if form is symmetric or alternating
    local m;
    m := f!.matrix;
    return m = TransposedMat(m) or m = -TransposedMat(m);
  end );

InstallMethod( IsReflexiveForm, [IsHermitianForm],
  function( f )
    return FORMS_IsHermitianMatrix( f!.matrix, f!.basefield);
  end );

InstallMethod( IsReflexiveForm, [IsTrivialForm], #new in 1.2.1
function( f )
   return true;
end );

InstallMethod( IsReflexiveForm, [IsQuadraticForm],
  function( f )
    Error( "<form> must be sesquilinear" );
end );

InstallMethod( IsAlternatingForm, [IsBilinearForm],
  function( f )
    return FORMS_IsSymplecticMatrix( f!.matrix, f!.basefield );
end );

InstallMethod( IsAlternatingForm, [IsHermitianForm],
  function( f )
    return false;
end );

InstallMethod( IsAlternatingForm, [IsTrivialForm], #new in 1.2.1
function( f )
   return true;
end );

InstallMethod( IsAlternatingForm, [IsQuadraticForm],
  function( f )
    Error( "<form> must be sesquilinear" );
end );

InstallMethod( IsSymmetricForm, [IsBilinearForm],
  function( f )
    return FORMS_IsSymmetricMatrix( f!.matrix );
end );

InstallMethod( IsSymmetricForm, [IsHermitianForm],
  function( f )
    return false;
end );

InstallMethod( IsSymmetricForm, [IsTrivialForm], #new in 1.2.1
function( f )
   return true;
end );

InstallMethod( IsSymmetricForm, [IsQuadraticForm],
  function( f )
    Error( "<form> must be sesquilinear" );
end );

InstallMethod( IsSymplecticForm, "for sesquilinear forms",
  [IsSesquilinearForm and IsFormRep],
  function( f )
    local string;
    string := f!.type;
    if string = "symplectic" then
       if IsAlternatingForm( f ) then
          return true;
       else
          Error( "<form> was incorrectly specified" );
       fi;
    else
       return false;
    fi;
  end );

InstallMethod( IsSymplecticForm, [IsTrivialForm], #new in 1.2.1
function( f )
   return true;
end );

InstallMethod( IsOrthogonalForm, "for sesquilinear forms",
  [IsSesquilinearForm and IsFormRep],
  function( f )
    local string;
    string := f!.type;
    if string = "orthogonal" then
       if IsSymmetricForm( f ) then
          return true;
       else
          Error( "Form was incorrectly specified" );
       fi;
    else
       return false;
    fi;
  end );

InstallMethod( IsOrthogonalForm, [IsTrivialForm], #new in 1.2.1
function( f )
   return false;
end );

InstallMethod( IsOrthogonalForm, [IsQuadraticForm],
  function( f )
    Error( "<form> must be sesquilinear" );
end );

InstallMethod( IsPseudoForm, "for sesquilinear forms",
  [IsSesquilinearForm and IsFormRep],
  function( f )
    local string;
    string := f!.type;
    if string = "pseudo" then
       if (IsSymmetricForm( f ) and (not IsAlternatingForm(f))) then
          return true;
       else
          Error( "Form was incorrectly specified" );
       fi;
    else
       return false;
    fi;
  end );

InstallMethod( IsPseudoForm, [IsTrivialForm], #new in 1.2.1
function( f )
   return false;
end );

InstallMethod( IsPseudoForm, [IsQuadraticForm],
  function( f )
    Error( "<form> must be sesquilinear" );
end );

##
#############################################################################

#############################################################################
# Base change methods (for users):
#############################################################################

#############################################################################
#O BaseChangeToCanonical( <form> ) <form>: trivial form.
#  returns warning and identity matrix.
##
InstallMethod( BaseChangeToCanonical, "for a trivial form",
  [IsTrivialForm],
  function(f)
    local b,m,n;
    m := f!.matrix;
    n := NrRows(m);
    b := IdentityMat(n,m);
    Info(InfoWarning,1,"<form> is trivial, trivial base change is returned");
    return b;
end );

#############################################################################
#O BaseChangeToCanonical( <form> ) <form>: sesquilinear form.
#  this function checks the type of the given form and calls the appropriate
#  main base change operation.
##
InstallMethod( BaseChangeToCanonical, "for a sesquilinear form",
  [IsSesquilinearForm and IsFormRep],
  function(f)
    local string,m,gf,b;
    string := f!.type;
    m := f!.matrix;
    gf := f!.basefield;
    if IsOrthogonalForm(f) then
      b := BaseChangeOrthogonalBilinear(m, gf);
      SetWittIndex(f, (b[2]+b[3]-1) / 2);
      SetIsEllipticForm(f,b[3]=0);
      SetIsParabolicForm(f,b[3]=1);
      SetIsHyperbolicForm(f,b[3]=2);
      return b[1];
    elif IsSymplecticForm(f) then
      b := BaseChangeSymplectic(m, gf);
      SetWittIndex(f, b[2]/2 );
      return b[1];
    elif string = "hermitian" then
      b := BaseChangeHermitian(m, gf);
      SetWittIndex(f, Int( (b[2]+1)/2 ));
      return b[1];
    elif string = "pseudo" then
      Error("BaseChangeToCanonical not yet implemented for pseudo forms");
    fi;
end );

#############################################################################
#O BaseChangeToCanonical( <form> ) <form>: quadratic form.
#  this function checks the type of the given form and calls the appropriate
#  main base change operation.
##
InstallMethod( BaseChangeToCanonical, "for a quadratic form",
  [IsQuadraticForm and IsFormRep],
  function(f)
    local m, gf, b, polarity;
    m := f!.matrix;
    gf := f!.basefield;
    if IsOddInt(Size(gf)) then
      polarity := (m+TransposedMat(m))/2;
      b := BaseChangeOrthogonalBilinear(polarity, gf);
    else
      b := BaseChangeOrthogonalQuadratic(m, gf);
    fi;
    SetWittIndex(f, (b[2]+b[3]-1) / 2);
    SetIsEllipticForm(f,b[3]=0);
    SetIsParabolicForm(f,b[3]=1);
    SetIsHyperbolicForm(f,b[3]=2);
    return b[1];
  end );

#############################################################################
# Overloading: Frobenius Automorphisms
#############################################################################

InstallOtherMethod( \^, "for a FFE vector and a Frobenius automorphism",
  [ IsVector and IsFFECollection, IsFrobeniusAutomorphism ],
  function( v, f )
    return List(v,x->x^f);
  end );

InstallOtherMethod( \^, "for a FFE vector and a trivial Frobenius automorphism",
  [ IsVector and IsFFECollection, IsMapping and IsOne ],
  function( v, f )
    return v;
  end );

InstallOtherMethod( \^,
  "for a compressed GF2 vector and a Frobenius automorphism",
  [ IsVector and IsFFECollection and IsGF2VectorRep, IsFrobeniusAutomorphism ],
  function( v, f )
    local w;
    w := List(v,x->x^f);
    ConvertToVectorRepNC(w,2);
    return w;
  end );

InstallOtherMethod( \^,
  "for a compressed GF2 vector and a trivial Frobenius automorphism",
  [ IsVector and IsFFECollection and IsGF2VectorRep, IsMapping and IsOne ],
  function( v, f )
    return v;
  end );

InstallOtherMethod( \^,
  "for a compressed 8bit vector and a Frobenius automorphism",
  [ IsVector and IsFFECollection and Is8BitVectorRep, IsFrobeniusAutomorphism ],
  function( v, f )
    local w;
    w := List(v,x->x^f);
    ConvertToVectorRepNC(w,Q_VEC8BIT(v));
    return w;
  end );

InstallOtherMethod( \^,
  "for a compressed 8bit vector and a trivial Frobenius automorphism",
  [ IsVector and IsFFECollection and Is8BitVectorRep, IsMapping and IsOne ],
  function( v, f )
    return v;
  end );

InstallOtherMethod( \^, "for a FFE matrix and a Frobenius automorphism",
  [ IsMatrix and IsFFECollColl, IsFrobeniusAutomorphism ],
  function( m, f )
    return List(m,v->List(v,x->x^f));
  end );

InstallOtherMethod( \^, "for a FFE matrix and a trivial Frobenius automorphism",
  [ IsMatrix and IsFFECollColl, IsMapping and IsOne ],
  function( m, f )
    return m;
  end );

InstallOtherMethod( \^,
  "for a compressed GF2 matrix and a Frobenius automorphism",
  [ IsMatrix and IsFFECollColl and IsGF2MatrixRep, IsFrobeniusAutomorphism ],
  function( m, f )
    local w,l,i;
    l := [];
    for i in [1..NrRows(m)] do
        w := List(m[i],x->x^f);
        ConvertToVectorRepNC(w,2);
        Add(l,w);
    od;
    ConvertToMatrixRepNC(l,2);
    return l;
  end );

InstallOtherMethod( \^,
  "for a compressed GF2 matrix and a trivial Frobenius automorphism",
  [ IsMatrix and IsFFECollColl and IsGF2MatrixRep, IsMapping and IsOne ],
  function( m, f )
    return m;
  end );

InstallOtherMethod( \^,
  "for a compressed 8bit matrix and a Frobenius automorphism",
  [ IsMatrix and IsFFECollColl and Is8BitMatrixRep, IsFrobeniusAutomorphism ],
  function( m, f )
    local w,l,i,q;
    l := [];
    q := Q_VEC8BIT(m[1]);
    for i in [1..NrRows(m)] do
        w := List(m[i],x->x^f);
        ConvertToVectorRepNC(w,q);
        Add(l,w);
    od;
    ConvertToMatrixRepNC(l,q);
    return l;
  end );

InstallOtherMethod( \^,
  "for a compressed 8bit matrix and a trivial Frobenius automorphism",
  [ IsMatrix and IsFFECollColl and Is8BitMatrixRep, IsMapping and IsOne ],
  function( m, f )
    return m;
  end );

#############################################################################
# Overloading: Forms
#############################################################################

InstallOtherMethod( \^, "for a pair of FFE vectors and a sesquilinear form",
  [ IsVectorList and IsFFECollColl, IsBilinearForm ],
  function( pair, f )
    if Size(pair) <> 2 then
       Error("The first argument must be a pair of vectors");
    fi;
    return pair[1] * f!.matrix * pair[2];
  end );

InstallOtherMethod( \^, "for a pair of FFE matrices and a sesquilinear form",
  [ IsFFECollCollColl, IsBilinearForm ],
  function( pair, f )
    if Size(pair) <> 2 then
       Error("The first argument must be a pair of vectors");
    fi;
    return pair[1] * f!.matrix * TransposedMat(pair[2]);
  end );

InstallOtherMethod( \^, "for a pair of FFE vectors and an hermitian form",
  [ IsVectorList and IsFFECollColl, IsHermitianForm ],
  function( pair, f )
    local frob,hh,bf;
    if Size(pair) <> 2 then
       Error("The first argument must be a pair of vectors");
    fi;
    bf := f!.basefield;
    hh := DegreeOverPrimeField(bf) / 2;
    frob := FrobeniusAutomorphism(bf)^hh;
    return pair[1] * f!.matrix * (pair[2]^frob);
  end );

InstallOtherMethod( \^, "for a pair of FFE matrices and an hermitian form",
  [ IsFFECollCollColl, IsHermitianForm ],
  function( pair, f )
    local frob,hh,bf;
    if Size(pair) <> 2 then
       Error("The first argument must be a pair of vectors");
    fi;
    bf := f!.basefield;
    hh := DegreeOverPrimeField(bf) / 2;
    frob := FrobeniusAutomorphism(bf)^hh;
    return pair[1] * f!.matrix * (TransposedMat(pair[2])^frob);
  end );

InstallOtherMethod( \^, "for a pair of FFE matrices and a trivial form", #new in 1.2.1
  [ IsVectorList and IsFFECollColl, IsTrivialForm ],
  function( pair, f )
    if Size(pair) <> 2 then
       Error("The first argument must be a pair of vectors or a vector");
    fi;
    return Zero(BaseField(f));
  end );

InstallOtherMethod( \^, "for a FFE vector and a quadratic form",
  [ IsVector and IsFFECollection, IsQuadraticForm ],
  function( v, f )
    return v * f!.matrix * v;
  end );

InstallOtherMethod( \^, "for a FFE matrix and a quadratic form",
  [ IsMatrix and IsFFECollColl, IsQuadraticForm ],
  function( m, f )
    return m * f!.matrix * TransposedMat(m);
  end );

InstallOtherMethod( \^, "for a FFE vector and a quadratic form", #new in 1.2.1
  [ IsVector and IsFFECollection, IsTrivialForm ],
  function( m, f )
    return Zero(BaseField(f));
  end );

#############################################################################
# Viewing methods:
#############################################################################
#
InstallMethod( ViewObj, [ IsTrivialForm ],
  function( f )
    Print("< trivial form >");
  end );

InstallMethod( PrintObj, [ IsTrivialForm ],
  function( f )
    Print("Trivial form\n");
    Print("Gram Matrix:\n",f!.matrix,"\n");
  end );

InstallMethod( Display, [ IsTrivialForm ],
  function( f )
    Print("Trivial form\n");
    Print("Gram Matrix:\n");
    Display(f!.matrix);
  end);

InstallMethod( ViewObj, [ IsHermitianForm ],
  function( f )
    if HasIsDegenerateForm(f) then
      if IsDegenerateForm(f) then
        Print(" < degenerate hermitian form >");
      else
        Print(" < non-degenerate hermitian form >");
      fi;
    else
      Print("< hermitian form >");
    fi;
  end );

InstallMethod( PrintObj, [ IsHermitianForm ],
  function( f )
    Print("Hermitian form\n");
    Print("Gram Matrix:\n",f!.matrix,"\n");
    if HasPolynomialOfForm( f ) then
       Print("Polynomial: ", PolynomialOfForm, "\n");
    fi;
    if HasWittIndex( f ) then
       Print("Witt Index: ", WittIndex(f), "\n");
    fi;
  end );

InstallMethod( Display, [ IsHermitianForm ],
  function( f )
    Print("Hermitian form\n");
    Print("Gram Matrix:\n");
    Display(f!.matrix);
    if HasPolynomialOfForm( f ) then
       Print("Polynomial: ");
       Display(PolynomialOfForm(f));
       Print("\n");
    fi;
    if HasWittIndex( f ) then
       Print("Witt Index: ", WittIndex(f), "\n");
    fi;
  end);

InstallMethod( ViewObj, [ IsQuadraticForm ],
  function( f )
    local string;
    string := ["< "];
    if HasIsSingularForm(f) then
      if IsSingularForm(f) then
        Add(string,"singular ");
      else
        Add(string,"non-singular ");
      fi;
    fi;
    if HasIsEllipticForm( f ) or HasIsHyperbolicForm( f ) or
      HasIsParabolicForm( f ) then
      if IsEllipticForm( f ) then
          Add(string,"elliptic ");
      elif IsHyperbolicForm( f ) then
          Add(string,"hyperbolic ");
      elif IsParabolicForm( f)  then
          Add(string,"parabolic ");
      fi;
    fi;
    Add(string,"quadratic form >");
    string := Concatenation(string);
    Print(string);
  end );

InstallMethod( PrintObj, [ IsQuadraticForm ],
  function( f )
    local string;
    string := [];
    if HasIsSingularForm(f) then
       if IsSingularForm(f) then
          Add(string,"Singular ");
       else
          Add(string,"Non-singular ");
       fi;
    fi;
    if HasIsEllipticForm( f ) or
       HasIsHyperbolicForm( f ) or
       HasIsParabolicForm( f ) then

       if IsEllipticForm( f ) then
          Add(string,"Elliptic ");
       elif IsHyperbolicForm( f ) then
          Add(string,"Hyperbolic ");
       elif IsParabolicForm( f)  then
          Add(string,"Parabolic ");
       fi;
     fi;
     Add(string,"Quadratic form\n");
     string := Concatenation(string[1],LowercaseString(Concatenation(string{[2..Length(string)]})));
     Print(string);
     Print("Gram Matrix:\n",f!.matrix,"\n");
     if HasPolynomialOfForm( f ) then
        Print("Polynomial: ", PolynomialOfForm(f), "\n");
     fi;
     if HasWittIndex( f ) then
        Print("Witt Index: ", WittIndex(f), "\n");
     fi;
  end );

InstallMethod( Display,  [ IsQuadraticForm ],
  function( f )
    local string;
    string := [];
    if HasIsSingularForm(f) then
       if IsSingularForm(f) then
          Add(string,"Singular ");
       else
          Add(string,"Non-singular ");
       fi;
    fi;
    if HasIsEllipticForm( f ) or
       HasIsHyperbolicForm( f ) or
       HasIsParabolicForm( f ) then

       if IsEllipticForm( f ) then
          Add(string,"Elliptic ");
       elif IsHyperbolicForm( f ) then
          Add(string,"Hyperbolic ");
       elif IsParabolicForm( f)  then
          Add(string,"Parabolic ");
       fi;
    fi;
    Add(string,"Quadratic form\n");
    string := Concatenation(string[1],LowercaseString(Concatenation(string{[2..Length(string)]})));
    Print(string);
    Print("Gram Matrix:\n");
    Display(f!.matrix);
    if HasPolynomialOfForm( f ) then
       Print("Polynomial: ");
       Display(PolynomialOfForm(f));
       Print("\n");
    fi;
    if HasWittIndex( f ) then
       Print("Witt Index: ", WittIndex(f), "\n");
    fi;
  end );

InstallMethod( ViewObj, [ IsBilinearForm ],
  function( f )
    local string;
    string := ["< "];
    if HasIsDegenerateForm(f) then
      if IsDegenerateForm(f) then
        Add(string,"degenerate ");
      else
        Add(string,"non-degenerate ");
      fi;
    fi;
    if HasIsOrthogonalForm(f) then
      if HasIsEllipticForm( f ) or
         HasIsHyperbolicForm( f ) or
         HasIsParabolicForm( f ) then

         if IsEllipticForm( f ) then
            Add(string,"elliptic bilinear ");
         elif IsHyperbolicForm( f ) then
            Add(string,"hyperbolic bilinear ");
         elif IsParabolicForm( f)  then
            Add(string,"parabolic bilinear ");
         fi;
      elif IsOrthogonalForm(f) then
         Add(string,"orthogonal ");
      fi;
    elif HasIsSymplecticForm(f) then
      if IsSymplecticForm(f) then
         Add(string,"symplectic ");
      fi;
    elif HasIsPseudoForm(f) then
      if IsPseudoForm(f) then
         Add(string,"pseudo ");
      fi;
    else
      Add(string,"bilinear ");
    fi;
    Add(string,"form >");
    string := Concatenation(string);
    Print(string);
  end );

InstallMethod( PrintObj, [ IsBilinearForm ],
  function( f )
    local string;
    string := [];
    if HasIsDegenerateForm(f) then
      if IsDegenerateForm(f) then
         Add(string,"Degenerate ");
      else
         Add(string,"Non-degenerate ");
      fi;
    fi;
    if HasIsOrthogonalForm(f) then
       if HasIsEllipticForm( f ) or
          HasIsHyperbolicForm( f ) or
          HasIsParabolicForm( f ) then
        if IsEllipticForm( f ) then
           Add(string,"Elliptic bilinear ");
        elif IsHyperbolicForm( f ) then
           Add(string,"Hyperbolic bilinear ");
        elif IsParabolicForm( f)  then
           Add(string,"Parabolic bilinear ");
        fi;
      elif IsOrthogonalForm(f) then
           Add(string,"Orthogonal ");
      fi;
    elif HasIsSymplecticForm(f) then
      if IsSymplecticForm(f) then
         Add(string,"Symplectic ");
      fi;
    elif HasIsPseudoForm(f) then
      if IsPseudoForm(f) then
         Add(string,"Pseudo ");
      fi;
    else
         Add(string,"Bilinear ");
    fi;
    Add(string,"form\n");
    string := Concatenation(string[1],LowercaseString(Concatenation(string{[2..Length(string)]})));
    Print(string);
    Print("Gram Matrix:\n",f!.matrix,"\n");
    if HasPolynomialOfForm( f ) then
       Print("Polynomial: ", PolynomialOfForm(f), "\n");
    fi;
    if HasWittIndex( f ) then
       Print("Witt Index: ", WittIndex(f), "\n");
    fi;
  end );

InstallMethod( Display, [ IsBilinearForm ],
  function( f )
    local string;
    string := [];
    if HasIsDegenerateForm(f) then
       if IsDegenerateForm(f) then
          Add(string,"Degenerate ");
       else
          Add(string,"Non-degenerate ");
       fi;
    fi;
    if HasIsOrthogonalForm(f) then
       if HasIsEllipticForm( f ) or
          HasIsHyperbolicForm( f ) or
          HasIsParabolicForm( f ) then

          if IsEllipticForm( f ) then
             Add(string,"Elliptic bilinear ");
          elif IsHyperbolicForm( f ) then
             Add(string,"Hyperbolic bilinear ");
          elif IsParabolicForm( f)  then
             Add(string,"Parabolic bilinear ");
          fi;
       elif IsOrthogonalForm(f) then
            Add(string,"Orthogonal ");
       fi;
    elif HasIsSymplecticForm(f) then
       if IsSymplecticForm(f) then
          Add(string,"Symplectic ");
        fi;
    elif HasIsPseudoForm(f) then
       if IsPseudoForm(f) then
          Add(string,"Pseudo ");
       fi;
    else
       Add(string,"Bilinear ");
    fi;
    Add(string,"form\n");
    string := Concatenation(string[1],LowercaseString(Concatenation(string{[2..Length(string)]})));
    Print(string);
    Print("Gram Matrix:\n");
    Display(f!.matrix);
    if HasPolynomialOfForm( f ) then
       Print("Polynomial: ");
       Display(PolynomialOfForm(f));
       Print("\n");
    fi;
    if HasWittIndex( f ) then
       Print("Witt Index: ", WittIndex(f), "\n");
    fi;
  end );

# end of Viewing methods.
#############################################################################

#############################################################################
# Functions to support Base Change operations (not for the user):
##
if IsBound(SwapMatrixColumns) and IsBound(SwapMatrixRows) then
  # For GAP >= 4.12
  BindGlobal("Forms_SwapRows", SwapMatrixRows);
  BindGlobal("Forms_SwapCols", SwapMatrixColumns);
  BindGlobal("Forms_AddRows", AddMatrixRows);
  BindGlobal("Forms_AddCols", AddMatrixColumns);
  BindGlobal("Forms_MultRow", MultMatrixRow);
  BindGlobal("Forms_MultCol", MultMatrixColumn);
else
  # For GAP <= 4.11
  BindGlobal("Forms_SwapRows", function(mat, i, j)
    mat{[i,j]} := mat{[j,i]};
  end);

  BindGlobal("Forms_SwapCols", function(mat, i, j)
    local row;
    for row in mat do
      row{[i,j]} := row{[j,i]};
    od;
  end);

  BindGlobal("Forms_AddRows", function(mat, i, j, scalar)
    mat[i] := mat[i] + mat[j] * scalar;
  end);

  BindGlobal("Forms_AddCols", function(mat, i, j, scalar)
    local row;
    for row in mat do
      row[i] := row[i] + row[j] * scalar;
    od;
  end);

  BindGlobal("Forms_MultRow", function(mat, i, scalar)
    mat[i] := mat[i] * scalar;
  end);

  BindGlobal("Forms_MultCol", function(mat, i, scalar)
    local row;
    for row in mat do
      row[i] := row[i] * scalar;
    od;
  end);

fi;

# express v as sum of two squares of elements in gf
InstallGlobalFunction(Forms_SUM_OF_SQUARES,
  function(v,gf)
    local dummy,i,v1,v2, primroot;
    primroot := PrimitiveRoot(gf);
    i := 0;
    repeat
      dummy := LogFFE(v - primroot^(2*i), primroot);
      if dummy mod 2 = 0 then
        v1 := primroot^i;
        v2 := primroot^(dummy/2);
        break;
      else
        i := i + 1;
      fi;
    until false;
    return [v1,v2];
  end );

# Apply a 2x2 transformation matrix to rows 'p1' and 'p2' of the matrix 'D'.
# The matrix 'D' may be changed in-place by this. For efficiency the matrix
# entries are passed in as arguments 'a11', 'a12', 'a21', 'a22'.
BindGlobal("Forms_TRANSFORM_2_BY_2",
  function(D,p1,p2,a11,a12,a21,a22)
    local r1,r2;
    r1 := D[p1];
    r2 := D[p2];
    D[p1] := a11 * r1 + a12 * r2;
    D[p2] := a21 * r1 + a22 * r2;
  end );

InstallGlobalFunction(Forms_REDUCE2,
  function(D,start,stop,gf)
    local n,t,i,half,primroot;
    n := NrRows(D);
    primroot := PrimitiveRoot(gf);
    half := One(gf) / 2;
    t := primroot^(LogFFE(-One(gf),primroot)/2) / 2;
    i := start;
    while i < stop do
      Forms_TRANSFORM_2_BY_2(D,i,i+1,half,t,half,-t);
      i := i + 2;
    od;
  end );

InstallGlobalFunction(Forms_REDUCE4,
  function(D,start,stop,gf)
    local n,c,d,i,dummy;
    n := NrRows(D);
    i := start;
    dummy := Forms_SUM_OF_SQUARES(-One(gf),gf);
    c := dummy[1];
    d := dummy[2];
    while i < stop do
      Forms_TRANSFORM_2_BY_2(D,i+1,i+3,c,d,d,-c);
      i := i + 4;
    od;
  end );

InstallGlobalFunction(Forms_DIFF_2_S,
  function(D,start,stop)
    local n,i,half;
    n := NrRows(D);
    i := start;
    half := One(D[1,1]) / 2;
    while i < stop do
      Forms_TRANSFORM_2_BY_2(D,i,i+1,half,half,half,-half);
      i := i + 2;
    od;
  end );

InstallGlobalFunction(Forms_HERM_CONJ,
  function(mat,t)
    local n,i,j,dummy;
    n := NrRows(mat);
    dummy := MutableTransposedMat(mat);
    for i in  [1..n] do
      for j in [1..n] do
        dummy[i,j] := dummy[i,j]^t;
      od;
    od;
    return dummy;
  end );

#Forms_RESET: given a matrix mat, computes an upper triangular matrix A such that mat and A determine
#the same quadratic form. The convention about quadratic forms in this package
#is that their Gram matrix is always an upper triangular matrix, although the user
#is free to use any matrix to construct the form.

BindGlobal("Forms_RESET_inplace",
  function(A)
    local i,j,t;
    t := Zero(A[1,1]);
    for i in [2..NrRows(A)] do
      for j in [1..i-1] do
        A[j,i] := A[j,i] + A[i,j];
        A[i,j] := t;
      od;
    od;
  end );

# HACK: ignore extra arguments to Forms_RESET, which used to
# take three arguments; the other two were redundant, so we dropped them.
# But the FinInG function also calls Forms_RESET, so for its sake, we
# ignore the other arguments
InstallGlobalFunction(Forms_RESET,
  function(A,extra...)
    A := MutableCopyMat(A);
    Forms_RESET_inplace(A);
    return A;
  end );

InstallGlobalFunction(Forms_SQRT2,
  function(a, gf)
    local z, q;
    if IsZero( a ) then
      return a;
    fi;
    q:=Size(gf);
    if q mod 2 = 0 then
      return a^(q/2);
    fi;
    z:=PrimitiveRoot(gf);
    return z^(LogFFE(a,z)/2);
  end );

InstallGlobalFunction(Forms_PERM_VAR,
  function(D,r)
    local i;
    i := Remove(D, r);
    Add(D, i, 1);
  end );

InstallGlobalFunction(Forms_C1,
  function(gf, h)
    local i, primroot;
    if h mod 2 = 0 then
      primroot := PrimitiveRoot(gf);
      i := 1;
      while i <= h - 1 do
        if not IsZero( Trace(gf, primroot^i) ) then
          return primroot^i;
        else
          i := i + 1;
        fi;
      od;
    else
      return One(gf);
    fi;
  end );

InstallGlobalFunction(Forms_QUAD_EQ,
  function(delta, gf, h)
    local i,k,dummy,result;
    k := Forms_C1(gf,h);
    dummy := Zero(gf);
    result := Zero(gf);
    for i in [1..h-1] do
      dummy := dummy + k^(2^(i-1));
      result := result + dummy*(delta^(2^i));
    od;
    return result;
  end );
##
#end support functions base change
#############################################################################

#############################################################################
# Operations to check input (most likely not for the user):
#############################################################################

InstallMethod( FORMS_IsSymplecticMatrix, [IsFFECollColl, IsField],
  function(m,f)
    if Characteristic(f) = 2 then
       if not ForAll([1..NrRows(m)], i -> IsZero(m[i,i]) ) then
           return false;
       fi;
    fi;
    return m = -TransposedMat(m);
  end );

InstallMethod( FORMS_IsSymmetricMatrix, [IsFFECollColl],
  function(m)
    return m=TransposedMat(m);
  end );

InstallMethod( FORMS_IsHermitianMatrix, [IsFFECollColl, IsField],
  function(m,f)
    return m=Forms_HERM_CONJ(m,Sqrt(Size(f)));
  end );

#############################################################################
# Main Base-change operations
# (helping methods for BaseChangeToCanonical, not for the users):
#############################################################################
#O BaseChangeOrthogonalBilinear( <mat>, <f> )
#  <f>: finite field, <mat>: matrix over <f>
# output: [D,r,w]; D = base change matrix,
#                r = number of non zero rows in D*mat*TransposedMat(D)
#                using r, it follows immediately whether the form is degenerate
#                w = character (0=elliptic, 2=hyperbolic, 1=parabolic).
##
InstallMethod( BaseChangeOrthogonalBilinear,
    [ IsMatrix and IsFFECollColl, IsField and IsFinite ],
  function(mat, gf)
    local row,i,j,k,A,b,c,d,P,D,dummy,r,w,s,v,v1,v2,
          n,q,primroot,one;
    Assert(1, mat = TransposedMat(mat));
    n := NrRows(mat);
    q := Size(gf);
    one := One(gf);

    A := MutableCopyMat(mat);
    ConvertToMatrixRep(A, gf);
    D := IdentityMat(n, gf);
    ConvertToMatrixRep(D, gf);
    row := 0;

    # Diagonalize A

    while row < n - 1 do
      row := row + 1;

      # We look for a nonzero element on the main diagonal, starting
      # from row
      i := row;
      while i <= n and IsZero(A[i,i]) do
        i := i + 1;
      od;

      if i = row then
        # do nothing since A[row,row] <> 0
      elif i <= n then
        # swap things around to ensure A[row,row] <> 0
        Forms_SwapCols(A, row, i);
        Forms_SwapRows(A, row, i);
        Forms_SwapRows(D, row, i);
      else
        # All entries on the main diagonal are zero. We now search for a
        # nonzero element off the main diagonal.
        i := row;
        while i < n do
          k := i + 1;
          while k <= n and IsZero(A[i,k]) do
            k := k + 1;
          od;
          if k = n + 1 then
             i := i + 1;
          else
             break;
          fi;
        od;

        # if i is n, then they are all zero and we can stop.
        if i = n then
          row := row - 1;
          r := row;
          break;
        fi;

        # Otherwise: Go and fetch...
        # Put it on A[row,row+1]
        if i <> row then
          Forms_SwapCols(A, row, i);
          Forms_SwapRows(A, row, i);
          Forms_SwapRows(D, row, i);
        fi;

        Forms_SwapCols(A, row + 1, k);
        Forms_SwapRows(A, row + 1, k);
        Forms_SwapRows(D, row + 1, k);

        b := 1/A[row+1,row];
        Forms_AddCols(A, row, row+1, b);
        Forms_AddRows(A, row, row+1, b);
        Forms_AddRows(D, row, row+1, b);
      fi;   # end if i = row ... elif  i <= n ... else ... fi

      # There is no zero element on the main diagonal, make the rest zero.
      c := -A[row,row]^-1;
      for i in [row+1..n] do
        b := A[i,row] * c;
        if IsZero(b) then continue; fi;
        Forms_AddCols(A, i, row, b);
        Forms_AddRows(A, i, row, b);
        Forms_AddRows(D, i, row, b);
      od;
    od;

    # Count how many variables are used.

    if row = n - 1 then
      if not IsZero(A[n,n]) then
        r := n;
      else
        r := n - 1;
      fi;
    fi;

    # Here we distinguish the quadratic from the non-quadratic

    i := 1;
    s := 0;
    primroot := PrimitiveRoot(gf);
    while i < r do
      if IsOddInt( LogFFE(A[i,i], primroot) ) then
         j := i + 1;
         repeat
            if IsEvenInt( LogFFE(A[j,j], primroot) ) then
               dummy := A[j,j];
               A[j,j] := A[i,i];
               A[i,i] := dummy;
               Forms_SwapRows(D, i, j);
               i := i + 1;
               s := s + 1;
               break;
            else
               j := j + 1;
            fi;
         until j = r + 1;
         if j = r + 1 then
            break;
         fi;
      else
        i := i + 1;
        s := s + 1;
      fi;
    od;
    if IsEvenInt( LogFFE(A[r,r], primroot) ) then
       s := s + 1;
    fi;

    # We do the form x_0^2 + ... + x_s^2 + v(x_s+1^2 + ... + x_r^2)
    # with v not quadratic.

    v := ShallowCopy(primroot);
    for i in [1..s] do
      Forms_MultRow(D,i,(primroot^(LogFFE(A[i,i], primroot)/2))^-1);
    od;
    for i in [s+1..r] do
      Forms_MultRow(D,i,(primroot^(LogFFE(A[i,i]/primroot,primroot)/2))^-1);
    od;

    # We keep as much quadratic part as we can:

    s := s - 1;
    r := r - 1;

    # Case by case:

    if not (s = -1 or r = s )  then
      # We write first v=v1^2 + v2^2
      dummy := Forms_SUM_OF_SQUARES(v,gf);
      v1 := dummy[1];
      v2 := dummy[2];

      if (r - s) mod 2 = 0 then
        v1 := v1/v;
        v2 := v2/v;
        i := s + 2;
        repeat
          Forms_TRANSFORM_2_BY_2(D,i,i+1,v1,-v2,v2,v1);
          i := i + 2;
        until i = r + 2;
        s := r;
      else
        if r mod 2 = 0 then
          i := 1;
          repeat
            Forms_TRANSFORM_2_BY_2(D,i,i+1,v1,v2,-v2,v1);
            i := i + 2;
          until i = s + 2;
          s := -1;
        elif not (s = r - 1) then
          v1 := v1/v;
          v2 := v2/v;
          i := s + 2;
          repeat
            Forms_TRANSFORM_2_BY_2(D,i,i+1,v1,-v2,v2,v1);
            i := i + 2;
          until i = r + 1;
          s := r - 1;
        fi;
      fi;
    fi;

    # Towards standard forms (uses the standard proof of
    # the classification of quadratic forms):

    if r mod 2 <> 0 then
      if s = -1 or s = r then
        if q mod 4 = 1 then
          Forms_REDUCE2(D,1,r+1,gf);
          w := 2;
        else
          if ((r-1)/2) mod 2 <> 0 then
            Forms_REDUCE4(D,1,r+1,gf);
            Forms_DIFF_2_S(D,1,r+1);
            w := 2;
          else
            Forms_REDUCE4(D,3,r+1,gf);
            Forms_DIFF_2_S(D,3,r+1);
            w := 0;
          fi;
        fi;
      else
        if q mod 4 = 1 then
          if 1 < r then
            Forms_SwapRows(D,2,r+1);
            Forms_REDUCE2(D,3,r+1,gf);
          fi;
          w := 0;
        else
          if ((r-1)/2) mod 2 <> 0 then
            Forms_SwapRows(D,4,r+1);
            if 3 < r then
              Forms_REDUCE4(D,5,r+1,gf);
            fi;
            b := primroot^(LogFFE(-v,primroot)/2);
            b := 1/(2*b);
            Forms_TRANSFORM_2_BY_2(D,3,4,one/2,-b,one/2,b);
            if 3 < r then
              Forms_DIFF_2_S(D,5,r+1);
            fi;
            w := 0;
          else
            Forms_SwapRows(D,2,r+1);
            if 1 < r then
              Forms_REDUCE4(D,3,r+1,gf);
            fi;
            b := primroot^(LogFFE(-v,primroot)/2);
            b := 1/(2*b);
            Forms_TRANSFORM_2_BY_2(D,1,2,one/2,-b,one/2,b);
            if 1 < r then
              Forms_DIFF_2_S(D,3,r+1);
            fi;
            w := 2;
          fi;
        fi;
      fi;
    elif r <> 0 then
      w := 1;
      if q mod 4 = 1 then
        Forms_REDUCE2(D,2,r+1,gf);
      else
        if r mod 4 = 0 then
          Forms_REDUCE4(D,2,r+1,gf);
          Forms_DIFF_2_S(D,2,r+1);
        else
          if 3 < r then
            Forms_REDUCE4(D,4,r+1,gf);
          fi;
          dummy := Forms_SUM_OF_SQUARES(-1,gf);
          c := dummy[1];
          d := dummy[2];
          Forms_TRANSFORM_2_BY_2(D,1,3,c,d,d,-c);
          Forms_DIFF_2_S(D,2,r+1);
          i := 3;
          while i <= r + 1 do
            Forms_MultRow(D,i,-one);
            i := i + 2;
          od;
        fi;
      fi;
    else
      w := 1;
    fi;

    return [D,r,w];
end);

#############################################################################
#O BaseChangeOrthogonalQuadratic( <mat>, <f> )
#  <f>: finite field, <mat>: matrix over <f>
# output: [D,r,w]; D = base change matrix,
#                r = number of non zero rows in D*mat*TransposedMat(D)
#                using r, it follows immediately whether the form is degenerate
#                w = character (0=elliptic, 2=hyperbolic, 1=parabolic).
##
InstallMethod(BaseChangeOrthogonalQuadratic, [ IsMatrix and IsFFECollColl, IsField and IsFinite ],
    function(mat, gf)
    local A,r,w,row,dummy,i,j,h,D,P,t,a,b,c,d,e,s,
      zeros,posr,posk,n,zero,one;
    n := NrRows(mat);
    r := n;
    row := 1;
    zero := Zero(gf);
    one := One(gf);
    A := MutableCopyMat(mat);
    D := IdentityMat(n, gf);
    ConvertToMatrixRep(D, gf);
    zeros := [];
    for i in [1..n] do
      zeros[i] := zero;
    od;
    h := DegreeOverPrimeField(gf);
    while row + 2 <= r do
      if not IsZero( A[row,row] ) then
        i := row + 1;
        # check on the main diagonal; we look for a zero
        while i <= r and not IsZero(A[i,i]) do
          i := i + 1;
        od;

        # if there is a zero somewhere, then we go and get it.

        if i <= r then
          Forms_SwapCols(A,row,i);
          Forms_SwapRows(A,row,i);
          Forms_SwapRows(D,row,i);
          Forms_RESET_inplace(A);

        # Otherwise: look in other places.

        else
          dummy := true;
          i := row;
          while i <= r - 1 and dummy do
            j := i + 1;
            while j <= r do
              if not IsZero( A[i,j] ) then
                posr := i;
                posk := j;
                dummy := false;
                break;
              else
                j := j + 1;
              fi;
            od;
            i := i + 1;
          od;

          # If all is zero, STOP
          if dummy then
            t := Forms_SQRT2(A[row,row],gf);
            Forms_MultRow(D,row,1/t);
            for i in [row + 1..r] do
              Forms_AddRows(D,i,row,Forms_SQRT2(A[i,i],gf));
            od;
            # Permutation of the variables, it is a parabolic
            r := row;
            Forms_PERM_VAR(D,r);
            w := 1;
            r := r - 1;
            return [D,r,w];
          # Otherwise: A basischange
          else
            if IsZero( A[row+1,row+2] ) then
               if posr = row + 1 then
                  Forms_SwapCols(A,posk,row+2);
                  Forms_SwapRows(A,posk,row+2);
                  Forms_SwapRows(D,posk,row+2);
               elif posk = row + 2 then
                  Forms_SwapCols(A,posr,row+1);
                  Forms_SwapRows(A,posr,row+1);
                  Forms_SwapRows(D,posr,row+1);
               elif posr = row + 2 then
                  # TODO: Does this case ever occur? I failed to find examples
                  # that trigger it
                  P := TransposedMat(PermutationMat((posk,posr,row+1),n));
                  A := P*A*TransposedMat(P);
                  D := P*D;
               else
                  Forms_SwapCols(A,posk,row+2);
                  Forms_SwapRows(A,posk,row+2);
                  Forms_SwapRows(D,posk,row+2);

                  Forms_SwapCols(A,posr,row+1);
                  Forms_SwapRows(A,posr,row+1);
                  Forms_SwapRows(D,posr,row+1);
               fi;
               Forms_RESET_inplace(A);
            fi;
            #A[row+1,row+2] <> 0
            t := A[row+1,row+2];

            b := 1/t;
            Forms_MultCol(A,row+2,b);
            Forms_MultRow(A,row+2,b);
            Forms_MultRow(D,row+2,b);

            b := A[row,row+1];
            Forms_AddCols(A,row,row+2,b);
            Forms_AddRows(A,row,row+2,b);
            Forms_AddRows(D,row,row+2,b);

            for i in [row+3..n] do
              b := A[row+1,i];
              Forms_AddCols(A,i,row+2,b);
              Forms_AddRows(A,i,row+2,b);
              Forms_AddRows(D,i,row+2,b);
            od;

            Forms_RESET_inplace(A);

            # A has now that special form a_11*X_1^2+X_1*X_2 + G(X_0,X_2,...,X_n);
            b := A[row,row];
            t :=  Forms_SQRT2(b/A[row+1,row+1],gf);
            Forms_AddCols(A,row,row+1,t);
            Forms_AddRows(A,row,row+1,t);
--> --------------------

--> maximum size reached

--> --------------------

Quelle forms.gi Sprache: unbekannt

[ 0.67Quellennavigators Projekt ]