gap-packages · AliaBonnet · Nov 30, 2025 · Oct 30, 2025 · Nov 12, 2025 · Nov 12, 2025
diff --git a/PackageInfo.g b/PackageInfo.g
@@ -77,7 +77,7 @@ Dependencies := rec(
 
 AvailabilityTest := ReturnTrue,
 
-TestFile := "tst/testall.g",
+TestFile := "tst/testall.tst",
 
 #Keywords := [ "TODO" ],
 

diff --git a/doc/nofoma.xml b/doc/nofoma.xml
@@ -0,0 +1,17 @@
+<?xml version="1.0" encoding="UTF-8"?>
+
+<!-- This is an automatically generated file. -->
+<!DOCTYPE Book SYSTEM "gapdoc.dtd"
+[
+    [<#Include SYSTEM "_entities.xml">
+]
+>
+<Book Name="nofoma">
+<#Include SYSTEM "title.xml">
+<TableOfContents/>
+<Body>
+<#Include SYSTEM "_AutoDocMainFile.xml">
+</Body>
+<Bibliography Databases="nofoma.bib"/>
+<TheIndex/>
+</Book>
diff --git a/doc/title.xml b/doc/title.xml
@@ -0,0 +1,26 @@
+<?xml version="1.0" encoding="UTF-8"?>
+
+<!-- This is an automatically generated file. -->
+<TitlePage>
+  <Title>
+    nofoma
+  </Title>
+  <Subtitle>
+    Normal forms of matrices
+  </Subtitle>
+  <Version>
+    1.0
+  </Version>
+  <Author>
+    Meinolf Geck<Alt Only="LaTeX"><Br/></Alt>
+<Address>
+Fachbereich Mathematik, Pfaffenwaldring 57, 70569 Stuttgart, Germany<Br/>
+</Address>
+<Email>meinolf.geck@mathematik.uni-stuttgart.de</Email>
+<Homepage>https://pnp.mathematik.uni-stuttgart.de/idsr/idsr1/geckmf/</Homepage>
+
+  </Author>
+  <Date>
+    1 June 2022
+  </Date>
+  </TitlePage>
diff --git a/gap/nofoma.gd b/gap/nofoma.gd
@@ -80,10 +80,61 @@ DeclareGlobalFunction("nfmCoeffsPol");
 DeclareGlobalFunction("nfmPolCoeffs");
 DeclareGlobalFunction("nfmGcd");
 DeclareGlobalFunction("nfmLcm");
+
+#! @Arguments a,b
+#! @Description 
+#! 'GcdCoprimeSplit' computes a divisor <M>a_1</M> of the polynomial <M>a</M> and a
-#! 'GcdCoprimeSplit' computes a divisor <M>a_1</M> of the polynomial <M>a</M> and a
+#! Computes a divisor <M>a_1</M> of the polynomial <M>a</M> and a
-#! 'GcdCoprimeSplit' computes a divisor <M>a_1</M> of the polynomial <M>a</M> and a
+#! Computes a divisor <M>a_1</M> of the polynomial <M>a</M> and a
+#! divisor <M>b_1</M> of the polynomial <M>b</M> such that <M>a_1</M> and <M>b_1</M> are coprime
+#! and the lcm of <M>a</M>, <M>b</M> is <M>a_1</M>*<M>b_1</M>.  This is based on Lemma 5 in <Cite Key = "Bon14"/>.
-#! and the lcm of <M>a</M>, <M>b</M> is <M>a_1</M>*<M>b_1</M>.  This is based on Lemma 5 in <Cite Key = "Bon14"/>.
+#! and the lcm of <M>a</M>, <M>b</M> is <M>a_1</M>*<M>b_1</M>.  This is based on Lemma 5 in <Cite Key="Bon14"/>.
-#! and the lcm of <M>a</M>, <M>b</M> is <M>a_1</M>*<M>b_1</M>.  This is based on Lemma 5 in <Cite Key = "Bon14"/>.
+#! and the lcm of <M>a</M>, <M>b</M> is <M>a_1</M>*<M>b_1</M>.  This is based on Lemma 5 in <Cite Key="Bon14"/>.
+#! (see also Lemma 4.3 in <Cite Key = "Gec20"/>).
+#!
+#! (Note that it does not use the prime factorisation of polynomials but 
+#! only gcd computations.)
+#! 
+#! @BeginExampleSession
+#! gap> a:=x^2*(x-1)^3*(x-2)*(x-3);
+#! x^7-8*x^6+24*x^5-34*x^4+23*x^3-6*x^2
+#! gap> b:=x^2*(x-1)^2*(x-2)^4*(x-4);
+#! x^9-14*x^8+81*x^7-252*x^6+456*x^5-480*x^4+272*x^3-64*x^2
+#! gap> GcdCoprimeSplit(a,b);
+#! [ x^5-4*x^4+5*x^3-2*x^2,               # the (monic) gcd
+#! x^4-6*x^3+12*x^2-10*x+3,               # a1
+#! x^7-12*x^6+56*x^5-128*x^4+144*x^3-64*x^2 ]  # b1
-#! [ x^5-4*x^4+5*x^3-2*x^2,               # the (monic) gcd
-#! x^4-6*x^3+12*x^2-10*x+3,               # a1
-#! x^7-12*x^6+56*x^5-128*x^4+144*x^3-64*x^2 ]  # b1
+#! [ x^5-4*x^4+5*x^3-2*x^2,               # the (monic) gcd
+#!   x^4-6*x^3+12*x^2-10*x+3,               # a1
+#!   x^7-12*x^6+56*x^5-128*x^4+144*x^3-64*x^2 ]  # b1
-#! [ x^5-4*x^4+5*x^3-2*x^2,               # the (monic) gcd
-#! x^4-6*x^3+12*x^2-10*x+3,               # a1
-#! x^7-12*x^6+56*x^5-128*x^4+144*x^3-64*x^2 ]  # b1
+#! [ x^5-4*x^4+5*x^3-2*x^2,               # the (monic) gcd
+#!   x^4-6*x^3+12*x^2-10*x+3,               # a1
+#!   x^7-12*x^6+56*x^5-128*x^4+144*x^3-64*x^2 ]  # b1
+#! @EndExampleSession
 DeclareGlobalFunction("GcdCoprimeSplit");
+
+#! @Arguments A,pol,v
+#! @Description
+#! 'PolynomialToMatVec' returns  the row vector  obtained  by multiplying 
+#! the row vector <M>v</M> with the matrix <M>pol</M>(<M>A</M>), where <M>pol</M> is the list 
-#! the row vector <M>v</M> with the matrix <M>pol</M>(<M>A</M>), where <M>pol</M> is the list 
+#! the row vector <A>v</A> with the matrix <M><A>pol</A>(<A>A</A>)</M>, where <A>pol</A> is the list 
-#! the row vector <M>v</M> with the matrix <M>pol</M>(<M>A</M>), where <M>pol</M> is the list 
+#! the row vector <A>v</A> with the matrix <M><A>pol</A>(<A>A</A>)</M>, where <A>pol</A> is the list 
+#! of coefficients of a polynomial.
+#!
+#! @BeginExampleSession
+#! gap> A:=([ [ 0, 1, 0, 1 ],
+#! gap>       [ 0, 0, 0, 0 ],
+#! gap>       [ 0, 1, 0, 1 ],
+#! gap>       [ 1, 1, 1, 1 ] ]);;
-#! gap> A:=([ [ 0, 1, 0, 1 ],
-#! gap>       [ 0, 0, 0, 0 ],
-#! gap>       [ 0, 1, 0, 1 ],
-#! gap>       [ 1, 1, 1, 1 ] ]);;
+#! gap> A:= [ [ 0, 1, 0, 1 ],
+#! gap>       [ 0, 0, 0, 0 ],
+#! gap>       [ 0, 1, 0, 1 ],
+#! gap>       [ 1, 1, 1, 1 ] ];;
-#! gap> A:=([ [ 0, 1, 0, 1 ],
-#! gap>       [ 0, 0, 0, 0 ],
-#! gap>       [ 0, 1, 0, 1 ],
-#! gap>       [ 1, 1, 1, 1 ] ]);;
+#! gap> A:= [ [ 0, 1, 0, 1 ],
+#! gap>       [ 0, 0, 0, 0 ],
+#! gap>       [ 0, 1, 0, 1 ],
+#! gap>       [ 1, 1, 1, 1 ] ];;
+#! gap> f:=x^6-6*x^5+12*x^4-10*x^3+3*x^2;;
+#! gap> v:=[ 1, 1, 1, 1];;
+#! gap> l:=nfmCoeffsPol(f);
+#! gap> [ 0, 0, 3, -10, 12, -6, 1 ]
+#! gap> PolynomialToMat(A,last,v);
+#! [ 8, -16, 8, -16 ]
+#! @EndExampleSession
 DeclareGlobalFunction("PolynomialToMatVec");
+
 DeclareGlobalFunction("PolynomialToMat");
+
+#! @Arguments A,v1,v2,pol1,pol2
+#! @Description
+#!  'LcmPolynomialToMatVec' returns,  given  a matrix  <M>A</M>,  vectors <M>v1</M>,
+#!  <M>v2</M> with minimal polynomials <M>pol1</M>, <M>pol2</M>,  a new pair [<M>v</M>,<M>pol</M>],  
+#!  where <M>v</M> has minimal polynomial <M>pol</M>, and <M>pol</M> is the least common
+#!  multiple of <M>pol1</M> and <M>pol2</M>.  
+#!  This crucially relies on  'GcdCoprimeSplit' to avoid  factorisation of 
+#!  polynomials.
 DeclareGlobalFunction("LcmMaximalVectorMat");
+
 DeclareGlobalFunction("SpinMatVector1");
 
 #! @Arguments A,v
@@ -122,7 +173,32 @@ DeclareGlobalFunction("SpinMatVector1");
 #! @EndExampleSession
 DeclareGlobalFunction("SpinMatVector");
 
+#! @Arguments A
+#! @Description
+#!  'CyclicChainMat' repeatedly applies 'SpinMatVector1' (relative version
+#!  of 'SpinMatVector') to compute a chain of cyclic subspaces. The output
+#!  is a triple [B,C,svec]  where  C is such that  C*<M>A</M>*C^-1  has a block
-#!  is a triple [B,C,svec]  where  C is such that  C*<M>A</M>*C^-1  has a block
+#!  is a triple <C>[B,C,svec]</C> where <M>C</M> is such that  <M>C <A>A</A> C^-1</M>  has a block
-#!  is a triple [B,C,svec]  where  C is such that  C*<M>A</M>*C^-1  has a block
+#!  is a triple <C>[B,C,svec]</C> where <M>C</M> is such that  <M>C <A>A</A> C^-1</M>  has a block
+#!  triangular shape with companian matrices along the diagonal), B is the
-#!  triangular shape with companian matrices along the diagonal), B is the
+#!  triangular shape with companion matrices along the diagonal), B is the
-#!  triangular shape with companian matrices along the diagonal), B is the
+#!  triangular shape with companion matrices along the diagonal), B is the
+#!  row echelon form of C and svec is the list of indices where the blocks
+#!  begin.
+#!
+#! @BeginExampleSession
+#! gap> A:=[ [ 0, 1, 0, 1 ],
+#! gap>      [ 0, 0, 1, 0 ],
+#! gap>      [ 0, 1, 0, 1 ],
+#! gap>      [ 1, 1, 1, 1 ] ]);;
-#! gap>      [ 1, 1, 1, 1 ] ]);;
+#! gap>      [ 1, 1, 1, 1 ] ];;
-#! gap>      [ 1, 1, 1, 1 ] ]);;
+#! gap>      [ 1, 1, 1, 1 ] ];;
+#! gap> sp:=CyclicChainMat(A);
+#! [ [ [ 1, 0, 0, 0 ], [ 0, 1, 0, 1 ], [ 0, 0, 1, 0 ], [ 0, 0, 0, 1 ] ],
+#!   [ [ 1, 0, 0, 0 ], [ 0, 1, 0, 1 ], [ 1, 1, 2, 1 ], [ 0, 0, 0, 1 ] ],
+#!   [ 1, 4, 5 ] ]
+#! gap> PrintArray(sp[2]*A*sp[2]^-1);
+#! [ [    0,    1,    0,    0 ],  #There are 2 diagonal blocks,
+#! [    0,    0,    1,    0 ],    #one (size 3x3) starting at index 1,
+#! [    0,    3,    1,    0 ],    #one (size 1x1) starting at index 4.
+#! [  1/2,  1/2,  1/2,    0 ] ]
+#! @EndExampleSession
 DeclareGlobalFunction("CyclicChainMat");
+
 DeclareGlobalFunction("nfmRelMinPols");
 DeclareGlobalFunction("nfmOrderPolM");
 DeclareGlobalFunction("MinPolyMat");
@@ -133,15 +209,13 @@ DeclareGlobalFunction("MinPolyMat");
 #!  of the matrix <A>A</A>, that is, a vector whose local minimal polynomial 
 #!  is that of <A>A</A>. This is done by repeatedly spinning up vectors until 
 #!  a maximal one is found. The exact algorithm is a combination of 
-#!  * the minimal polynomial algorithm by Neunhoeffer-Praeger; see 
-#!      J. Comput. Math. 11, 2008 
-#!      (<URL>http://doi.org/10.1112/S1461157000000590</URL>); and  
+#!  * the minimal polynomial algorithm by Neunhoeffer-Praeger; see <Cite Key = "Neu08"/>; and  
 #!  * the method described by Bongartz 
-#!      (see <URL>https://arxiv.org/abs/1410.1683</URL>) for computing 
+#!      (see <Cite Key = "Bon14"/>) for computing 
 #!      maximal vectors.
 #!
 #!  See also the article by Geck at
-#!  <URL>https://doi.org/10.13001/ela.2020.5055</URL>. 
+#!  <Cite Key = "Gec20"/>. 
 #!
 #! @BeginExampleSession
 #! gap> A:=[ [  2,  2,  0,  1,  0,  2,  1 ],
@@ -177,11 +251,51 @@ DeclareGlobalFunction("MinPolyMat");
 #! @EndExampleSession
 DeclareGlobalFunction("MaximalVectorMat");
 
+#! @Arguments T,d
+#! @Description
+#! 'JacobMatComplement' modifies an already given  complementary subspace 
+#! to the  complementary subspace defined by  Jacob;  concretely, this is 
+#! realized by assuming that  <M>T</M>  is a matrix in block triangular shape, 
+#! where the upper left diagonal block is a companion matrix (as returned 
+#! by 'RatFormStep1'; the variable <M>d</M> gives the size of that block.  
+#! (If <M>T</M> gives a  maximal cyclic subspace,  then  Jacob's complement is  
+#! also  <M>T</M>-invariant;  but even if not,  it appears  to be  very useful 
+#! because it produces many zeroes.)
 DeclareGlobalFunction("JacobMatComplement");
+
 DeclareGlobalFunction("BuildBlockDiagonalMat");
 DeclareGlobalFunction("BuildBlockDiagonalMat1");
+
+#! @Arguments A,v
+#! @Description
+#! 'RatFormStep1J' spins up a vector  <M>v</M> under a  matrix  <M>A</M>,  computes
+#! a complementary subspace  (using  Jacob's construction),  and performs
+#! the base change. The output is a quadruple  [A1,P,pol,str] where A1 is
+#! the new matrix,  P is the base change,  pol is  the minimal polynomial
+#! and str is either 'split' or 'not', according to whether the extension
+#! is split or not. The second form repeatedly applies 'RatFormStep1J' in
+#! order to obtain an invariant complement.
+#!
+#! @BeginExampleSession
+#! gap> v:=[ 1, 1, 1, 1 ];;
+#! gap> A:=[ [ 0, 1, 0, 1 ],
+#! gap>      [ 0, 0, 1, 0 ],
+#! gap>      [ 0, 1, 0, 1 ],
+#! gap>      [ 1, 1, 1, 1 ] ];;
+#! gap> PrintArray(RatFormStep1J(A,v)[1])
+#! [ [  0,  1,  0,  0 ],    #There are 2 diagonal blocks but
+#!   [  0,  0,  1,  0 ],    #(because of the (4,1)-entry 1)
+#!   [  0,  3,  1,  0 ],    #the extension is not split.
+#!   [  1,  0,  0,  0 ] ]   
+#! gap> PrintArray(RatFormStep1Js(A,v)[1])";
+#! [ [  0,  1,  0,  0 ],    #Now we actually see that
+#!   [  0,  0,  1,  0 ],    #the matrix is cyclic.
+#!   [  0,  0,  0,  1 ],     
+#!   [  0,  0,  3,  1 ] ]   
+#! @EndExampleSession
 DeclareGlobalFunction("RatFormStep1");
 DeclareGlobalFunction("RatFormStep1J");
+
 DeclareGlobalFunction("nfmCompanionMat");
 DeclareGlobalFunction("nfmCompanionMat1");
 
@@ -191,14 +305,17 @@ DeclareGlobalFunction("nfmCompanionMat1");
 #!  and an invertible matrix <M>P</M> such that <M>P.A.P^{{-1}}</M> is the 
 #!  Frobenius normal form of <A>A</A>. The algorithm first computes a maximal 
 #!  vector and an <A>A</A>-invariant complement following Jacob's construction
-#!  (as described in matrix language in Geck, Electron. J. Linear Algebra 36, 
-#!  2020, <URL>https://doi.org/10.13001/ela.2020.5055</URL>); then the 
+#!  (as described in matrix language in <Cite Key = "Gec20"/>); then the 
 #!  algorithm continues recursively. It works for matrices over any field 
 #!  that is available in GAP. The output is a triple with 
 #!  * 1st component  = list of invariant factors; 
 #!  * 2nd component = base change matrix <M>P</M>; and 
 #!  * 3rd component = indices where the various blocks in the normal form 
 #!       begin.
+#!  You can also use  'CreateNormalForm( f[1] );' to produce the Frobenius
+#!  normal form. (This function just builds the block diagonal matrix with 
+#!  diagonal blocks given by the companion matrices corresponding to the 
+#!  various invariant factors of <A>A</A>.) 
 #! 
 #! @BeginExampleSession
 #! gap> A:=[ [  2,  2,  0,  1,  0,  2,  1 ],
@@ -231,16 +348,32 @@ DeclareGlobalFunction("nfmCompanionMat1");
 #! (This is the Frobenius normal form; there are 3 diagonal blocks,
 #!  one of size 4, one of size 2 and one of size 1.)
 #! @EndExampleSession
-#! 
-#! You can also use  'CreateNormalForm( f[1] );' to produce the Frobenius
-#! normal form. (This function just builds the block diagonal matrix with 
-#! diagonal blocks given by the companion matrices corresponding to the 
-#! various invariant factors of <A>A</A>.) 
 DeclareGlobalFunction("FrobeniusNormalForm");
 
 DeclareGlobalFunction("CreateNormalForm");
 DeclareGlobalFunction("FrobeniusNormalForm1");
+
+#! @Arguments A
+#! @Description
+#! 'InvariantFactorsMat' returns the invariant factors of the matrix <M>A</M>,
+#! i.e.,  the minimal polynomials of the  diagonal blocks in the rational 
+#! canonical form  of <M>A</M>. Thus, 'InvariantFactorsMat' also specifies the
+#! rational canonical form of <M>A</M>, but without computing the base change.
+#!
+#! @BeginExampleSession
+#! gap> InvariantFactorsMat([ [ 2,  2, 0, 1, 0,  2, 1 ],
+#!                            [ 0,  4, 0, 0, 0,  1, 0 ],
+#!                            [ 0,  1, 1, 0, 0,  1, 1 ],
+#!                            [ 0, -1, 0, 1, 0, -1, 0 ],
+#!                            [ 0, -7, 0, 0, 1, -5, 0 ],
+#!                            [ 0, -2, 0, 0, 0,  1, 0 ],
+#!                            [ 0, -1, 0, 0, 0, -1, 1 ] ]);
+#!   #I Degree of minimal polynomial is 4
+#!   #I Degree of minimal polynomial is 2
+#!   [ x^4-7*x^3+17*x^2-17*x+6, x^2-3*x+2, x-1 ]
+#! @EndExampleSession
 DeclareGlobalFunction("InvariantFactorsMat");
+
 DeclareGlobalFunction("nfmFrobInv");
 DeclareGlobalFunction("SquareFreePol");
 
@@ -253,8 +386,7 @@ DeclareGlobalFunction("SquareFreePol");
 #!  <A>A</A>), such that <M>D.N=N.D</M>; the argument <A>f</A> is a 
 #!  polynomial such that <M>f(A)=0</M> (e.g., the minimal polynomial of 
 #!  <A>A</A>). This is called the Jordan-Chevalley decomposition of <A>A</A>; 
-#!  the algorithm is based on the preprint at 
-#!  <URL>https://arxiv.org/abs/2205.05432</URL>. Note that this 
+#!  the algorithm is based on <Cite Key = "Gec22"/>. Note that this 
 #!  algorithm does not require the knowledge of the eigenvalues of <A>A</A>; 
 #!  it works over any perfect field that is available in GAP.
 #!
@@ -289,19 +421,23 @@ DeclareGlobalFunction("SquareFreePol");
 #!  'JordanChevalleyDecMat(<A>A</A>);)'
 DeclareGlobalFunction("JordanChevalleyDecMat");
 
+#! @Arguments A
+#! @Description
+#!  'JordanChevalleyDecMatF'  first computes the Frobenius normal form and
+#!  then applies 'JordanChevalleyDecMat' to each diagonal block.
 DeclareGlobalFunction("JordanChevalleyDecMatF");
+
 DeclareGlobalFunction("CheckFrobForm");
 DeclareGlobalFunction("CheckJordanChev");
 
-#! @Arguments lev
-#! @Description
-#!  'Testnofoma' runs a number of tests on the functions in this package;
-#!  the argument <A>lev</A> is a positive integer specifying the InfoLevel.
-DeclareGlobalFunction("Testnofoma");
+DeclareGlobalFunction("nfmmat1");
 
 #! @Section Further documentation
 #! The above functions, as well as a number of further auxiliary functions, 
 #! are all contained and defined in the file 'pgk/nofoma-1.0/gap/nofoma.gi'; 
 #! in that file, you can also find further inline documentation for the 
 #! auxiliary functions.
 
+<Bibliography Databases="../doc/nofoma.bib"/>
+
+
-<Bibliography Databases="../doc/nofoma.bib"/>
-<Bibliography Databases="../doc/nofoma.bib"/>