/* In this case we let H denote a copy of PSL(2,25), L denote the subgroup 5^2.12, and let G denote F4(13^4).
   We give matrices defining L, an element of order 5 and an element of order 12. L lies in the normalizer of
   a torus of F4(13^4), as constructed using the GroupOfLieType command. We give H as L and a single matrix,
   which were chosen at random. A posteriori, one could produce H as a subgroup of F4 and show that no other
   overgroup of L preserves the F4 trilinear form, but we prefer a less synthetic proof.
*/

F<ww>:=GF(13^4);
l1:=GL(26,F)!GL(26,F)![[ww^22848,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,ww^17136,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,ww^17136,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,ww^5712,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,ww^17136,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,ww^5712,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,ww^17136,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,ww^11424,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,ww^5712,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,ww^5712,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^22848,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^17136,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^22848,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^11424,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^22848,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^11424,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^22848,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^11424,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^11424,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^5712]];

l2:=GL(26,F)![[0,0,12,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,12,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,12,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,12,0,0,0,0,0,0,0,0],
[1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,12,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,12,0,0,0],
[0,0,0,0,0,12,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,12,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,1,12,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,12,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,12,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,12,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,12],
[0,0,0,0,0,0,0,0,0,12,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,12,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,12,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0]];

l3:=GL(26,F)![[0,0,0,0,0,0,0,0,0,ww^8568,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^5712,0,0,0,0,0,0,0],
[0,0,0,0,ww^8568,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^22848,0,0],
[0,0,0,0,0,0,0,0,0,0,0,ww^25704,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,ww^25704,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^22848,0,0,0,0,0],
[1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,ww^11424,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,12,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^25704,0],
[0,0,0,ww^17136,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,1,12,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,12,0,0,0,0,0,0,0,0,0,0,0,0],
[0,ww^2856,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^11424,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^2856,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,12,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,ww^19992,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^17136,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^17136,0,0,0,0,0,0,0,0,0,0],
[0,0,ww^19992,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^19992,0,0,0,0],
[0,0,0,0,0,0,ww^22848,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,ww^5712,0,0,0,0,0,0,0,0]];

L:=sub<GL(26,F)|[l1,l2]>;

h1:=GL(26,F)![[ww^17952,ww^13464,ww^13464,ww^17377,ww^15133,0,ww^8809,0,0,ww^3097,ww^17989,ww^15708,ww^22784,ww^22978,ww^28189,ww^12852,0,ww^14521,ww^23089,0,0
,ww^16765,ww^14521,ww^15096,ww^816,ww^24888],
[ww^7752,ww^18564,0,ww^3097,ww^853,ww^24888,ww^8809,ww^853,0,ww^3097,0,ww^6528,ww^8698,ww^20834,0,ww^7752,0,ww^14521,ww^23089,ww^16765,ww^3672,ww^2485
,ww^241,0,ww^9996,ww^6528],
[ww^7752,0,ww^27744,0,ww^20233,ww^20808,ww^17989,ww^5953,ww^9996,ww^22477,ww^8809,0,ww^20834,ww^8504,ww^23089,0,ww^4284,ww^9421,ww^28189,ww^25945,ww^7752
,ww^25945,0,ww^816,0,ww^6528],
[ww^25463,ww^16895,0,0,ww^13464,ww^4859,0,ww^4284,ww^3227,ww^10608,ww^12240,ww^14039,ww^6835,ww^7029,ww^2040,ww^11183,ww^21995,ww^17952,0,ww^24276
,ww^20363,ww^15096,0,0,ww^8327,ww^28319],
[ww^6083,ww^26075,ww^16895,ww^24888,0,ww^28319,ww^12240,ww^23664,ww^22607,ww^15708,0,0,ww^21115,ww^21309,0,0,ww^16895,ww^27132,ww^2040,ww^19176,ww^11183
,0,ww^17952,ww^22607,ww^13427,ww^4859],
[0,ww^13464,ww^9384,ww^8197,ww^20233,ww^1428,0,ww^20233,ww^4896,ww^8197,ww^8809,0,ww^22978,ww^6554,ww^23089,0,ww^23664,ww^23701,0,ww^11665,ww^12852
,ww^25945,ww^9421,ww^19176,ww^816,0],
[ww^11183,ww^16895,ww^26075,0,ww^23664,0,ww^12240,ww^27744,ww^17507,0,5,ww^14039,ww^21309,ww^4885,8,ww^11183,ww^7715,0,ww^2040,ww^15096,0,ww^4896,0
,ww^27707,ww^8327,ww^14039],
[0,ww^26075,ww^2615,ww^15708,ww^23664,ww^28319,ww^16320,0,ww^8327,ww^24888,0,ww^4859,ww^4885,ww^21115,0,ww^6083,ww^2615,ww^3672,ww^26520,0,ww^11183
,ww^4896,ww^12852,ww^8327,ww^13427,0],
[0,0,ww^4284,ww^12277,ww^20233,ww^10608,ww^3709,ww^5953,ww^4896,0,ww^23089,ww^24888,ww^22784,ww^22978,ww^8809,ww^17952,ww^23664,0,ww^13909,ww^25945
,ww^17952,ww^25945,ww^5341,ww^9996,0,0],
[ww^11183,ww^16895,ww^7715,ww^10608,ww^4284,ww^4859,0,ww^13464,0,0,ww^16320,ww^28319,ww^19165,ww^6835,ww^26520,ww^25463,0,0,0,ww^816,ww^20363,ww^9996
,ww^3672,ww^17507,ww^8327,ww^14039],
[ww^20363,0,ww^16895,ww^6528,0,ww^28319,5,0,ww^8327,ww^10608,ww^16320,ww^23219,ww^7029,ww^19165,ww^26520,ww^16283,ww^2615,ww^17952,8,0,ww^11183,0,ww^7752
,ww^22607,0,ww^19139],
[ww^27132,ww^23664,0,ww^17377,0,0,ww^23089,ww^25333,ww^19176,ww^3097,ww^3709,ww^6528,ww^20834,ww^8504,ww^13909,ww^7752,ww^9384,ww^14521,ww^8809,ww^20845
,0,0,ww^14521,0,ww^19176,ww^1428],
[ww^7552,ww^9154,ww^27266,ww^12077,ww^14933,ww^6298,ww^13679,ww^14655,ww^18976,ww^11799,ww^27959,ww^10130,ww^17206,ww^11424,ww^13679,ww^7274,ww^13264
,ww^23223,ww^27959,ww^6087,ww^17722,ww^20645,ww^9221,ww^4418,ww^586,ww^10408],
[ww^21554,ww^13264,ww^23434,ww^26079,ww^375,ww^10408,ww^17789,ww^10823,ww^4418,ww^7967,ww^3509,ww^6298,ww^3398,ww^2926,ww^17789,ww^3442,ww^27266,ww^19391
,ww^3509,ww^2255,ww^21832,ww^6087,ww^23223,ww^586,ww^4696,ww^24410],
[ww^2003,0,ww^2615,ww^24888,0,ww^14039,8,0,ww^22607,ww^20808,ww^26520,ww^4859,ww^21309,ww^4885,ww^16320,ww^6083,ww^16895,ww^7752,5,0,ww^25463,0,ww^17952
,ww^8327,0,ww^8939],
[ww^12852,ww^13464,0,ww^3097,0,0,ww^8809,ww^15133,ww^816,ww^17377,ww^13909,ww^24888,ww^6554,ww^22784,ww^3709,ww^17952,ww^27744,ww^241,ww^23089,ww^2485,0
,0,ww^241,0,ww^816,ww^15708],
[0,0,ww^4284,ww^8197,ww^20233,ww^6528,ww^28189,ww^5953,ww^816,0,ww^23089,ww^20808,ww^22784,ww^22978,ww^8809,ww^22032,ww^27744,0,ww^17989,ww^25945
,ww^22032,ww^25945,ww^9421,ww^9996,0,0],
[ww^11183,ww^16895,ww^11795,ww^6528,ww^4284,ww^8939,0,ww^9384,0,0,ww^12240,ww^28319,ww^19165,ww^6835,ww^2040,ww^25463,0,0,0,ww^4896,ww^16283,ww^9996
,ww^7752,ww^13427,ww^8327,ww^14039],
[ww^25463,ww^2615,ww^7715,0,ww^13464,0,ww^2040,ww^9384,ww^27707,0,8,ww^28319,ww^7029,ww^19165,5,ww^25463,ww^26075,0,ww^12240,ww^4896,0,ww^15096,0
,ww^17507,ww^22607,ww^28319],
[0,ww^7715,ww^16895,ww^1428,ww^13464,ww^14039,ww^26520,0,ww^22607,ww^6528,0,ww^23219,ww^19165,ww^6835,0,ww^16283,ww^16895,ww^22032,ww^16320,0,ww^25463
,ww^15096,ww^27132,ww^22607,ww^3227,0],
[0,ww^9384,ww^13464,ww^12277,ww^20233,ww^1428,0,ww^20233,ww^816,ww^12277,ww^8809,0,ww^22978,ww^6554,ww^23089,0,ww^27744,ww^19621,0,ww^11665,ww^12852
,ww^25945,ww^5341,ww^15096,ww^4896,0],
[ww^2003,ww^21995,ww^16895,ww^20808,0,ww^28319,ww^16320,ww^27744,ww^22607,ww^15708,0,0,ww^21115,ww^21309,0,0,ww^16895,ww^27132,ww^26520,ww^15096,ww^11183
,0,ww^22032,ww^22607,ww^17507,ww^8939],
[ww^11183,ww^2615,0,0,ww^23664,ww^23219,0,ww^18564,ww^13427,ww^20808,ww^2040,ww^28319,ww^21115,ww^21309,ww^12240,ww^25463,ww^11795,ww^7752,0,ww^9996
,ww^2003,ww^4896,0,0,ww^22607,ww^14039],
[ww^3672,0,ww^23664,0,ww^20233,ww^24888,ww^13909,ww^5953,ww^9996,ww^26557,ww^8809,0,ww^20834,ww^8504,ww^23089,0,ww^4284,ww^5341,ww^3709,ww^25945,ww^3672
,ww^25945,0,ww^4896,0,ww^10608],
[ww^17952,ww^4284,0,ww^17377,ww^11053,ww^6528,ww^23089,ww^11053,0,ww^17377,0,ww^24888,ww^22978,ww^6554,0,ww^17952,0,ww^241,ww^8809,ww^6565,ww^22032
,ww^20845,ww^14521,0,ww^24276,ww^24888],
[ww^7752,ww^23664,ww^23664,ww^3097,ww^25333,0,ww^23089,0,0,ww^17377,ww^28189,ww^1428,ww^8504,ww^8698,ww^17989,ww^27132,0,ww^241,ww^8809,0,0,ww^6565
,ww^241,ww^4896,ww^19176,ww^6528]];

H:=sub<GL(26,F)|[l1,l2,h1]>;

// We need the indices of the trilinear form as well.
seqs:=[[i,j,k]:i,j,k in [1..NumberOfRows(l1)]|i le j and j le k];

ents:=[273,295,313,327,332,571,584,597,608,628,644,840,855,883,926,942,1081,1098,1158,1181,1209,1295,1314,1410,1433,1452,1581,1613,
1640,1653,1674,1773,1807,1836,1861,2012,2039,2052,2072,2167,2196,2220,2361,2374,2384,2483,2507,2615,2718,2731];
vals:=[1,12,1,12,1,1,12,12,1,12,1,12,1,12,1,12,1,12,1,12,1,12,1,12,1,12,1,12,1,1,12,12,1,12,1,1,12,12,1,12,1,12,1,1,12,12,1,1,11,11];

f26:=[F!0:i in [1..#seqs]];
for i in [1..#ents] do f26[ents[i]]:=vals[i]; end for;

// Now set up the centralizer of L in GL(26,F), which is given by all possible matrices scal below.

R<a,b,c,d>:=PolynomialRing(F,4);
scal:=DiagonalMatrix([a,a,a,b,b,a,b,b,a,b,b,a,c+3*d,d+3*c,b,a,a,b,b,b,a,b,b,a,a,a]);
mats:=MatrixRing(R,26);

scal:=DiagonalMatrix([a,a,a,b,b,a,b,b,a,b,b,a,8*c+6*d,6*c+8*d,b,a,a,b,b,b,a,b,b,a,a,a]);
scal[13,14]:=2*d-2*c;
scal[14,13]:=2*c-2*d;

mats!l1*scal eq scal*mats!l1 and mats!l2*scal eq scal*mats!l2;

function Matricise(v)
return Matrix(1,NumberOfColumns(v),[R!v[i]:i in [1..NumberOfColumns(v)]]);
end function;

// We now show that g, which is the second generator of H (the first lies in L, so is fine) satisfies the equations mentioned
V:=GModule(L);
g:=mats!h1;

intseqs1:=[seqs[i]:i in ents];
intseqs:=&join{{[i[1],i[2],i[3]],[i[1],i[3],i[2]],[i[2],i[1],i[3]],[i[2],i[3],i[1]],[i[3],i[1],i[2]],[i[3],i[2],i[1]]}:i in intseqs1};

function ProdRel(seq)
u:=V.seq[1]; v:=V.seq[2]; w:=V.seq[3];
u1:=Matricise(u)*scal; v1:=Matricise(v)*scal; w1:=Matricise(w)*scal;
u2:=Matricise(u)*g*scal; v2:=Matricise(v)*g*scal; w2:=Matricise(w)*g*scal;
e1:=&+[u1[1,i[1]]*v1[1,i[2]]*w1[1,i[3]]*f26[Position(seqs,Sort([i[1],i[2],i[3]]))]:i in intseqs];
e2:=&+[u2[1,i[1]]*v2[1,i[2]]*w2[1,i[3]]*f26[Position(seqs,Sort([i[1],i[2],i[3]]))]:i in intseqs];

return e1-e2;
end function;

function CheckRel(seq,gg)
u:=V.seq[1]; v:=V.seq[2]; w:=V.seq[3];
u2:=u^gg; v2:=v^gg; w2:=w^gg;
e1:=&+[u[i[1]]*v[i[2]]*w[i[3]]*f26[Position(seqs,Sort([i[1],i[2],i[3]]))]:i in intseqs];
e2:=&+[u2[i[1]]*v2[i[2]]*w2[i[3]]*f26[Position(seqs,Sort([i[1],i[2],i[3]]))]:i in intseqs];
return e1-e2;
end function;

function CheckDetermination()

ProdRel([1,10,21])-ProdRel([1,4,9]) eq ww^22131*a*(a+ww^1669*b)*(a+ww^8809*b);
ProdRel([1,5,23])-ProdRel([1,7,8]) eq ww^711*b*(a+ww^1669*b)*(a+ww^23089*b);
// Since a,b ne 0, we must have a+ww^1669*b=0, so a-ww^15949*b

coeffs1:=[ww^15949*b,b,c,d];
Evaluate(ProdRel([1,1,23])/ww^19998-ProdRel([1,4,23])/ww^4049,coeffs1) eq ww^12664*b^2*(c+ww^9834*d);
// Since b ne 0 we have c + ww^9834*d=0, so c=ww^24114*d

coeffs2:=[ww^15949*b,b,ww^24114*d,d];
Evaluate(ProdRel([1,1,22]),coeffs2) eq ww^13936*b^2*(d-ww^6062*b);

coeffs3:=[ww^15949*b,b,ww^24114*ww^6062*b,ww^6062*b];
scal3:=Evaluate(scal,coeffs3);
sca:=GL(26,F)!Evaluate(scal3,[1:i in [1..4]]);
gg:=h1^sca;

// Then that gg does lie in F4.
{CheckRel(i,gg):i in seqs} eq {0};

// Now check that the normalizer of L normalizes H1.
H1:=sub<GL(26,F)|l1,l2,gg>;
NH:=sub<GL(26,F)|l1,l2,gg,l3>;
Index(NH,H1) eq 2 and IsIsomorphic(H1,PSL(2,25));
// We will check that <L,l3> is the normalizer of L in NGT, which is the normalizer of L in G, in CheckF4Form.

return "";
end function;

// We now give code that can prove that the trilinear form we gave is the correct one. The normalizer of a torus
// has a small space of trilinear forms for F4. We use two random generators to quickly reduce to that space, and
// and then choose a random element of F4 to reduce the space of trilinear forms down to 1 dimension.

function CheckF4Form()

q:=13^4;
F:=GF(q);
z:=PrimitiveElement(F)^5712;
GG:=GroupOfLieType("F4",q);
W:=VectorSpace(GF(q),4);
rho:=StandardRepresentation(GG);
Over:=GL(26,q);
g1:=elt<GG|W![z,1,1,1]>;
g2:=elt<GG|W![1,z,1,1]>;
g3:=elt<GG|W![1,1,z,1]>;
g4:=elt<GG|W![1,1,1,z]>;
W:=Reflections(GG);
Mats1:=[rho(i):i in [g1,g2,g3,g4]];
MatsE:=[rho(i):i in W];
NGT:=sub<Over|Mats1 cat MatsE>;
T:=sub<NGT|Mats1>;
E:=sub<NGT|MatsE>;

// Confirm that L is a subgroup of F4:

l1 in NGT and l2 in NGT and l3 in NGT;
Normalizer(NGT,L) eq sub<NGT|L,l3>;
// Check our assertions from the paper:

// Note that the character of M(F4) on elements of order 5 is 1.
subs:=[i`subgroup:i in Subgroups(NGT:OrderEqual:=25)|{BrauerCharacter(GModule(i`subgroup))[j]:j in {2..25}} eq {1}];
#subs eq 1;
Bool:=IsConjugate(NGT,DerivedSubgroup(L),subs[1]); Bool;
Ns:=Normalizer(NGT,subs[1]);
subs2:=Subgroups(Ns:OrderEqual:=12,IsCyclic:=true);
Ltest:=sub<NGT|subs[1],subs2[1]`subgroup>;
Bool:=IsConjugate(NGT,L,Ltest); Bool;
Index(NGT,Centralizer(NGT,subs[1])) eq 1152; //the order of W(F4). This shows that C_G(L0)=T.

// Now build the matrix of linear relations that the trilinear form must satisfy.
mat:=[];
repeat y1:=Random(NGT); y2:=Random(NGT); until NGT eq sub<NGT|y1,y2>;
for h in [y1,y2] do
  for nn in [1..#seqs] do a:=seqs[nn,1]; b:=seqs[nn,2]; c:=seqs[nn,3];
    val:=[F!0:i in [1..#seqs]];
    for i in [1..NumberOfRows(h)] do if(h[a,i] ne 0) then for j in [1..NumberOfRows(h)] do if(h[b,j] ne 0) then for k in [1..NumberOfRows(h)] do
      val[Position(seqs,Sort([i,j,k]))]+:=h[a,i]*h[b,j]*h[c,k];
    end for; end if; end for; end if; end for;
    val[nn]-:=1; Append(~mat,val); delete val;
  end for;
  delete h;
end for;

repeat
  h:=rho(Random(GG));
  for iii in [1..5] do nn:=Random([1..#seqs]); a:=seqs[nn,1]; b:=seqs[nn,2]; c:=seqs[nn,3];
    val:=[F!0:i in [1..#seqs]];
    for i in [1..NumberOfRows(h)] do if(h[a,i] ne 0) then for j in [1..NumberOfRows(h)] do if(h[b,j] ne 0) then for k in [1..NumberOfRows(h)] do
      val[Position(seqs,Sort([i,j,k]))]+:=h[a,i]*h[b,j]*h[c,k];
    end for; end if; end for; end if; end for;
    val[nn]-:=1; Append(~mat,val); delete val;
  end for;
  delete h;
until Dimension(Nullspace(Transpose(Matrix(F,mat)))) eq 1;

ftest:=Nullspace(Transpose(Matrix(F,mat))).1;
"Is the trilinear form given in this document correct?";
return &and[ftest[i] eq f26[i]:i in [1..#f26]];

end function;

"The function CheckF4Form() checks that the trilinear form f26 from F4 is correct";
"The function CheckDetermination() checks that the determination of the possible conjugates into F4 is correct.";