Gröbner bases calculations of the paper 
"Conformally Einstein Lorentzian Lie groups: the almost Abelian case"

E. Calviño-Louzao, E. García-Río, I.  Gutiérrez-Rodríguez, R. Vázquez-Lorenzo

ABSTRACT. This notebook is a complementary file of the paper "Conformally Einstein Lorentzian Lie groups: the almost Abelian case". It contains the  expressions of the Bach tensor and the calculations of the different Gröbner bases, aimed to provide evidence of the ideal membership of the polynomials used through the paper. Headings and numbering correspond to those in the paper.

Quit[];

Code needed for the notebook

Options and memory function

SetOptions[MatrixForm,TableAlignments->Left];

$OutputSizeLimit=Infinity;

memory:=
    Module[{i, j, k, l,inv},   
           For[i = 1, i <= dim, i++, j = i + 1;
               While[j <= dim, cor[[j, i]] = -cor[[i, j]];
                     g[[j, i]] = g[[i, j]]; j++]];           
           If[Length[ref] < dim  && Length[ref] == 1,
              ref = Table[Subscript[ref[[1]],i], {i, dim}]];
           dimtemp= dim;
           reftemp= ref;
           cortemp= cor;
           gtemp= g;
           If[(Det[g])===0, meminv=False;Print["Error: inverse of g"],
              meminv=True;                   
              Clear[estmem];
              Clear[chrmem];
              Clear[curcompmem];
              estmem = Table[cte[i, j, k], {i, dim}, {j, dim}, {k, dim}];
              chrmem = Table[Γaux[i, j][[k]], {i, dim}, {j, dim}, {k, dim}];
              curcompmem = Table[R13c[i, j, k, l], {i, dim}, {j, dim}, {k, dim}, {l, dim}];
           ]   
          ];

Geometric functions

cte[i_, j_, k_] := 
    Module[{},
           Coefficient[cor[[i, j]], ref[[k]]]
          ];

(*
Levi-Civita connection
*)
gconex[i_Integer, j_Integer, k_Integer] :=
    Module[{r},
           (1/2) Sum[estmem[[i, j, r]] g[[r, k]] -
                   estmem[[j, k, r]] g[[r, i]] +
                   estmem[[k, i, r]] g[[r, j]], {r,dim}]
          ];

Γaux[i_Integer, j_Integer] :=
    Module[{r, inv, gconexij},
           inv = Inverse[g];
           gconexij = Table[gconex[i, j, r], {r, dim}];
           inv.gconexij
          ];

chr := 
    Module[{i, j, k},
           Table[Γaux[i, j][[k]], {i, dim}, {j, dim}, {k, dim}]
          ];

(*
Curvature tensor
*)
R13[i_Integer, j_Integer, k_Integer] :=
    Module[{a, b},
           Collect[ Sum[ Sum[estmem[[i, j, a]] chrmem[[a, k, b]] -
                             chrmem[[j, k, a]] chrmem[[i, a, b]] +
                             chrmem[[i, k, a]] chrmem[[j, a, b]],{a,dim}] ref[[b]],
                         {b,dim}],
                    ref]
          ];

R13c[i_Integer, j_Integer, k_Integer, l_Integer] := 
    Module[{},
           Coefficient[R13[i, j, k], ref[[l]]]
          ];

(*
Ricci tensor
*)
ρ[i_Integer, j_Integer] :=
    Module[{r},
           Sum[curcompmem[[i, r, j, r]],{r,dim}] 
          ];

ric02 :=
    Module[{i, j},
           Table[ρ[i, j], {i, dim}, {j, dim}]
          ];

(*
Scalar curvature
*)
sc := 
    Module[{inv, i, j}, 
           inv = Inverse[g]; 
           Sum[inv[[i, j]] ρ[i, j],{i,dim},{j,dim}]
          ]  ;

(*
Weyl tensor
*)
W := 
    Module[{i, j, k, l, r, curvaturacomp, ricci, escalar}, 
           ricci = ric02; 
           escalar = sc; 
           Table[Sum[curcompmem[[i, j, k, r]] g[[r,l]],{r,dim}] + 
                 (escalar/((dim - 1) (dim - 2))) *
                    (g[[i, k]] g[[j, l]] - g[[j, k]] g[[i, l]]) - 
                 (1/(dim - 2)) (ricci[[i, k]] g[[j, l]] - 
                    ricci[[j, k]] g[[i, l]] + 
                    ricci[[j, l]] g[[i, k]] - 
                    ricci[[i, l]] g[[j, k]]), 
                  {i, dim}, {j, dim}, {k, dim}, {l, dim}
                 ]
          ];

Bach tensor

SchoutenTensor:=Module[{ric02temp,sctemp},
ric02temp=ric02//Simplify;
sctemp=sc//Simplify;
ric02temp-(sctemp/(2(dim-1)))g
];

CottonTensor03:=Module[{schoutentemp,chrtemp,calctemp,i,j,k,α,β},
schoutentemp=SchoutenTensor//Simplify;
chrtemp=chr//Simplify;
calctemp[i_,j_,k_]:= -  Sum[chrtemp[[k,i,α]]  schoutentemp[[α,j]],{α,dim}] - Sum[chrtemp[[k,j,β]]  schoutentemp[[i,β]],{β,dim}];
Table[calctemp[j,k,i]-calctemp[i,k,j] ,{i,dim},{j,dim},{k,dim}] 
];

DCottonTensor03:=Module[{cot03temp,chrtemp,calctemp,i,j,k,m,α,β,γ,a1,a2,a3,a4},
cot03temp=CottonTensor03//Simplify;
chrtemp=chr//Simplify;
calctemp[i_,j_,k_,m_]:=-  Sum[chrtemp[[m,i,α]]  cot03temp[[α,j,k]],{α,dim}] - Sum[chrtemp[[m,j,β]]  cot03temp[[i,β,k]],{β,dim}]-  Sum[chrtemp[[m,k,γ]]  cot03temp[[i,j,γ]],{γ,dim}];
Table[calctemp[a1,a2,a3,a4 ],{a1,dim},{a2,dim},{a3,dim},{a4,dim}] 
];

BachTensor02:=Module[{inv,weyltemp,ric02temp,dercotton03temp,i,j,k,l,α,β},
inv=Inverse[g]//Simplify;
weyltemp=W//Simplify;
ric02temp=ric02//Simplify;
dercotton03temp=DCottonTensor03//Simplify;
Table[(1/(dim-2)) (Sum[inv[[k,α]]  dercotton03temp[[k,i,j,α]],{α,dim},{k,dim}]+Sum[  ric02temp[[k,l]]  Sum[inv[[k,α]] inv[[l,β]]weyltemp[[i,α,j,β]],{α,1,dim},{β,1,dim}],{k,dim},{l,dim}]),{i,1,dim},{j,1,dim}]  
];

Section 4. Semi-direct extensions with Riemannian Lie group (ℝ^3)

(* Setting the Lie brackets and the left-invariant metric *)
dim=4;
ref={e1,e2,e3,e4};
cor=(({{0, 0, 0,  η1 e1-γ1 e2-γ2 e3}, {\[Placeholder], 0, 0, γ1 e1+η2 e2-γ3 e3}, {\[Placeholder], \[Placeholder], 0,  γ2 e1+γ3 e2+η3 e3}, {\[Placeholder], \[Placeholder], \[Placeholder], 0}}));
g=(({{1, 0, 0, 0}, {\[Placeholder], 1, 0, 0}, {\[Placeholder], \[Placeholder], 1, 0}, {\[Placeholder], \[Placeholder], \[Placeholder], -1}}));

memory;

(* We compute the expression of the Bach tensor *)
BachTMatrix=BachTensor02//Simplify;
BachTMatrix//MatrixForm

(* Since the Bach tensor is a symmetric (0,2)-tensor field we consider the following components *)
BachT={
BachTMatrix[[1,1]],BachTMatrix[[1,2]],BachTMatrix[[1,3]],BachTMatrix[[1,4]],
BachTMatrix[[2,2]],BachTMatrix[[2,3]],BachTMatrix[[2,4]],
BachTMatrix[[3,3]],BachTMatrix[[3,4]],
BachTMatrix[[4,4]]};
BachT//MatrixForm

Proof of Theorem 4.3.

Checking the components of the Bach tensor

B11=(1/6)(((η1)^4)-((η2)^4)-((η3)^4)-(η2+η3) ((η1)^3)+(3 η1+η3) ((η2)^3)+(3 η1+η2) ((η3)^3)-2 (((η2)^2)+((η3)^2)) ((η1)^2)
	-(7 η1-3 η2-3 η3) η1 η2 η3+9 (((γ1)^2)+((γ2)^2)) ((η1)^2)-3 (5 ((γ1)^2)+((γ3)^2)) ((η2)^2)-3 (5 ((γ2)^2)+((γ3)^2)) ((η3)^2)
	+6(((γ1)^2)+2 ((γ2)^2)) η1 η2+6 (2 ((γ1)^2)+((γ2)^2)) η1 η3-6 (2 ((γ1)^2)+2 ((γ2)^2)-((γ3)^2)) η2 η3 
	-18 γ1 γ2 γ3 (η2-η3));
B12=(1/6)(2 γ1 (((η1)^3)-((η2)^3))-2 γ1 (9 η2+5 η3) ((η1)^2)+2 γ1 (9 η1+5 η3) ((η2)^2)-γ1 (η1-η2) ((η3)^2)
	+6 γ2 γ3 (η1+η2-2 η3) (η1+η2+η3)+3 γ1 (4 ((γ1)^2)+4 ((γ2)^2)+((γ3)^2)) η1
	-3 γ1 (4 ((γ1)^2)+((γ2)^2)+4 ((γ3)^2)) η2-9 γ1 (((γ2)^2)-((γ3)^2)) η3);
B13=(1/6)(2 γ2 (((η1)^3)-((η3)^3))-2 γ2 (5 η2+9 η3) ((η1)^2)-γ2 (η1-η3) ((η2)^2)+2 γ2 (9 η1+5 η2) ((η3)^2)
	-6 γ1 γ3 (η1-2 η2+η3) (η1+η2+η3)+3 γ2 (4 ((γ1)^2)+4 ((γ2)^2)+((γ3)^2)) η1
	-9 γ2 (((γ1)^2)-((γ3)^2)) η2-3 γ2 (((γ1)^2)+4 ((γ2)^2)+4 ((γ3)^2)) η3);
B22=((-1)/6)(((η1)^4)-((η2)^4)+((η3)^4)-(3 η2+η3) ((η1)^3)+(η1+η3) ((η2)^3)-(η1+3 η2) ((η3)^3)+2 (((η1)^2)+((η3)^2)) ((η2)^2)
	-(3 η1-7 η2+3 η3) η1 η2 η3+3 (5 ((γ1)^2)+((γ2)^2)) ((η1)^2)-9 (((γ1)^2)+((γ3)^2)) ((η2)^2)+3 (((γ2)^2)+5 ((γ3)^2)) ((η3)^2)
	-6 (((γ1)^2)+2 ((γ3)^2)) η1 η2+6 (2 ((γ1)^2)-((γ2)^2)+2 ((γ3)^2)) η1 η3-6 (2 ((γ1)^2)+((γ3)^2)) η2 η3
	-18 γ1 γ2 γ3 (η1-η3));
B23=(1/6)(2 γ3 (((η2)^3)-((η3)^3))-γ3 (η2-η3) ((η1)^2)-2 γ3 (5 η1+9 η3) ((η2)^2)+2 γ3 (5 η1+9 η2) ((η3)^2)
	-6 γ1 γ2 (2 η1-η2-η3) (η1+η2+η3)-9 γ3 (((γ1)^2)-((γ2)^2)) η1
	+3 γ3 (4 ((γ1)^2)+((γ2)^2)+4 ((γ3)^2)) η2-3 γ3 (((γ1)^2)+4 ((γ2)^2)+4 ((γ3)^2)) η3);
B33=((-1)/6)(((η1)^4)+((η2)^4)-((η3)^4)-(η2+3 η3) ((η1)^3)-(η1+3 η3) ((η2)^3)+(η1+η2) ((η3)^3)+2 (((η1)^2)+((η2)^2)) ((η3)^2)
	-(3 η1+3 η2-7 η3) η1 η2 η3+3 (((γ1)^2)+5 ((γ2)^2)) ((η1)^2)+3 (((γ1)^2)+5 ((γ3)^2)) ((η2)^2)-9 (((γ2)^2)+((γ3)^2)) ((η3)^2)
	-6 (((γ1)^2)-2 ((γ2)^2)-2 ((γ3)^2)) η1 η2-6 (((γ2)^2)+2 ((γ3)^2)) η1 η3-6 (2 ((γ2)^2)+((γ3)^2)) η2 η3
	+18 γ1 γ2 γ3 (η1-η2));
B44=B11+B22+B33;

BachT=={B11,B12,B13,0,B22,B23,0,B33,0,B44}//Simplify

Case 2: (Subscript[η, 1])=(Subscript[η, 2])=1≠(Subscript[η, 3]), (Subscript[η, 3])≠0

(* We set the variables *)
varsI={γ1,γ2,γ3,η1,η2,η3};

(* Ideal *)
idealI=Join[BachT,{η1-1,η2-1}];

(* Polynomials g1, g2 and g3 *)
g1=η3(η3-1)(η3-4) ;
g2=((γ2)^2) (η3-1);
g3=((γ3)^2) (η3-1);

(* We compute a Gröbner basis for the ideal *)
gbI=GroebnerBasis[idealI,varsI,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI]];
gbI//Factor//MatrixForm

(* The polynomials g1, g2 and g3 are in the basis *)
g1==gbI[[1]]//Simplify
g2==gbI[[6]]//Simplify
g3==gbI[[4]]//Simplify

(* We double check that each of these polynomials belongs to the Gröbner basis *)
PolynomialReduce[g1,gbI,varsI,MonomialOrder->Lexicographic]
PolynomialReduce[g2,gbI,varsI,MonomialOrder->Lexicographic]
PolynomialReduce[g3,gbI,varsI,MonomialOrder->Lexicographic]

Case 3: (Subscript[η, 1])=1≠(Subscript[η, 2])≠(Subscript[η, 3]), (Subscript[η, 3])≠1, (Subscript[η, 2])(Subscript[η, 3])≠0

Initial basis

(* We set the variables *)
varsI={γ3,γ1,γ2,η1,η2,η3};

(* Ideal *)
idealI=Join[BachT,{η1-1}];

(* Polynomial g *)
polyp1=8 ((η3)^6)+12 ((η3)^5)-45 ((η3)^4)-100 ((η3)^3)-45 ((η3)^2)+12 η3+8;
polyp2=128 ((η3)^6)-48 ((η3)^5)-2985 ((η3)^4)-6340 ((η3)^3)-2985 ((η3)^2)-48 η3+128;
polyq1=(11723776 ((η3)^12)+57769728 ((η3)^11)+106929840 ((η3)^10)+192248648 ((η3)^9)+816486921 ((η3)^8) 
	+2162941668 ((η3)^7)+2962757238 ((η3)^6)+2162941668 ((η3)^5)+816486921 ((η3)^4)
	+192248648 ((η3)^3)+106929840 ((η3)^2)+57769728 η3
	+11723776)(2651014144 ((η3)^12)
	+24455976192 ((η3)^11)+98890309104 ((η3)^10)+43771506920 ((η3)^9) 
	-371896453695 ((η3)^8)+396884502852 ((η3)^7)+2252411448966 ((η3)^6) 
	+396884502852 ((η3)^5)-371896453695 ((η3)^4)+43771506920 ((η3)^3)+98890309104 ((η3)^2) 
	+24455976192 η3+2651014144) (((η3)^2)+η3+4);
polyq2=(η3-4)(7 ((η3)^2)-14 η3-23)(7 ((η3)^2)+7 η3-20);

polyg=γ1 γ2 polyp1 polyp2 η3 ( η3-1) polyq1 polyq2;

(* We compute a Gröbner basis for the ideal (time in a home computer: about 25 seconds) *)
gbI=GroebnerBasis[idealI,varsI,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI]];
gbI//Factor//MatrixForm

(* The polynomial g is in the basis *)
polyg==gbI[[10]]//Simplify

(* We double check that this polynomial belongs to the Gröbner basis *)
PolynomialReduce[polyg,gbI,varsI,MonomialOrder->Lexicographic]

Reduction of the case (Subscript[q, 2])=0 (the previous basis is needed)

(* We set the variables *)
varsIq2={η1,η2,η3,γ1,γ3,γ2};

(* Ideal *)
idealIq2=Join[gbI,{polyq2}];

(* Polynomial gq2 *)
polygq2=((γ2)^2) (14 ((γ2)^2)+33) (1421 ((γ2)^4)+28623 ((γ2)^2)+45600);

(* We compute a Gröbner basis for the ideal *)
gbIq2=GroebnerBasis[idealIq2,varsIq2,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbIq2]];
gbIq2//Factor//MatrixForm

(* The polynomial gq2 is in the basis *)
polygq2==gbIq2[[1]]//Simplify

(* We double check that this polynomial belongs to the Gröbner basis *)
PolynomialReduce[polygq2,gbIq2,varsIq2,MonomialOrder->Lexicographic]

Case 3.1: (Subscript[γ, 1])=0

First basis

(* We set the variables *)
varsI1={γ1,γ2,η1,η3,η2,γ3};

(* Ideal *)
idealI1=Join[BachT,{η1-1,γ1}];

(* Polynomial g1 *)
g1=γ2 γ3 (η2-1) (464 ((γ3)^4)+2175 ((γ3)^2)+1824) (2116 ((γ3)^4)+4884 ((γ3)^2)+1089) (80 ((γ3)^4)+95 ((γ3)^2)+32);

(* We compute a Gröbner basis for the ideal *)
gbI1=GroebnerBasis[idealI1,varsI1,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI1]];
gbI1//Factor//MatrixForm

(* The polynomial g1 is in the basis *)
g1==gbI1[[11]]//Simplify

(* We double check that this polynomial belongs to the Gröbner basis *)
PolynomialReduce[g1,gbI1,varsI1,MonomialOrder->Lexicographic]

Second basis

(* We set the variables *)
varsI1p={η1,η2,η3,γ1,γ2,γ3};

(* Ideal *)
idealI1p=Join[BachT,{γ1,γ2}];

(* Polynomial g1p *)
g1p=((γ3)^2) (η2-η3) (((η3)^4)+25 ((γ3)^4)+2 ((η3)^2) ((γ3)^2)) (9 ((η3)^4)+49 ((γ3)^4)+138 ((η3)^2) ((γ3)^2)) (81 ((η3)^4)+((γ3)^4)+90 ((η3)^2) ((γ3)^2));

(* We compute a Gröbner basis for the ideal *)
gbI1p=GroebnerBasis[idealI1p,varsI1p,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI1p]];
gbI1p//Factor//MatrixForm

(* The polynomial g1p is in the basis *)
g1p==gbI1p[[3]]//Simplify

(* We double check that this polynomial belongs to the Gröbner basis *)
PolynomialReduce[g1p,gbI1p,varsI1p,MonomialOrder->Lexicographic]

Case 3.2: (Subscript[γ, 1])≠0, (Subscript[p, 1])=0

(* We need de Initial basis (time in a home computer: about 25 seconds) *)
varsI={γ3,γ1,γ2,η1,η2,η3};
idealI=Join[BachT,{η1-1}];
gbI=GroebnerBasis[idealI,varsI,MonomialOrder->Lexicographic];

(* We set the variables *)
varsI2={γ1t,γ1,γ2,γ3,η1,η2,η3};

(* Ideal *)
polyp1=8 ((η3)^6)+12 ((η3)^5)-45 ((η3)^4)-100 ((η3)^3)-45 ((η3)^2)+12 η3+8;
idealI2=Join[gbI,{γ1 γ1t-1,polyp1}];

(* Polynomials g21 and g22 *)
g21=15 ((γ2)^2)-16 ((η3)^5)-24 ((η3)^4)+106 ((η3)^3)+196 ((η3)^2)+21 η3-28;
g22=15 ((γ3)^2)-88 ((η3)^5)-60 ((η3)^4)+535 ((η3)^3)+670 ((η3)^2)-15 η3-112;

(* We compute a Gröbner basis for the ideal *)
gbI2=GroebnerBasis[idealI2,varsI2,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI2]];
gbI2//Factor//MatrixForm

(* The polynomials g21 and g22 are in the basis *)
g21==gbI2[[5]]//Simplify
g22==gbI2[[4]]//Simplify

(* We double check that each of these polynomials belongs to the Gröbner basis *)
PolynomialReduce[g21,gbI2,varsI2,MonomialOrder->Lexicographic]
PolynomialReduce[g22,gbI2,varsI2,MonomialOrder->Lexicographic]

{(Subscript[p, 1])=0, (Subscript[g, 21])=0, (Subscript[g, 22])=0} has no real solution

(* We set the variables *)
varsI2p={η3,γ2,γ3};

(* Ideal *)
polyp1=8 ((η3)^6)+12 ((η3)^5)-45 ((η3)^4)-100 ((η3)^3)-45 ((η3)^2)+12 η3+8;
g21=15 ((γ2)^2)-16 ((η3)^5)-24 ((η3)^4)+106 ((η3)^3)+196 ((η3)^2)+21 η3-28;
g22=15 ((γ3)^2)-88 ((η3)^5)-60 ((η3)^4)+535 ((η3)^3)+670 ((η3)^2)-15 η3-112;
idealI2p={polyp1,g21,g22};

(* Polynomial g2p *)
g2p=216 ((γ2)^2) ((γ3)^4)+70 ((γ2)^4)+352 ((γ3)^4)+709 ((γ2)^2) ((γ3)^2)+234 ((γ2)^2)+1158 ((γ3)^2)+288;

(* We compute a Gröbner basis for the ideal *)
gbI2p=GroebnerBasis[idealI2p,varsI2p,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI2p]];
gbI2p//Factor//MatrixForm

(* The polynomial g2p is in the basis *)
g2p==gbI2p[[2]]//Simplify

(* We double check that this polynomial belongs to the Gröbner basis *)
PolynomialReduce[g2p,gbI2p,varsI2p,MonomialOrder->Lexicographic]

Case 3.3: (Subscript[γ, 1])≠0, (Subscript[p, 2])=0

(* We need de Initial basis (time in a home computer: about 25 seconds) *)
varsI={γ3,γ1,γ2,η1,η2,η3};
idealI=Join[BachT,{η1-1}];
gbI=GroebnerBasis[idealI,varsI,MonomialOrder->Lexicographic];

(* We set the variables *)
varsI3={γ1t,γ1,γ2,γ3,η1,η2,η3};

(* Ideal *)
polyp2=128 ((η3)^6)-48 ((η3)^5)-2985 ((η3)^4)-6340 ((η3)^3)-2985 ((η3)^2)-48 η3+128;
idealI3=Join[gbI,{γ1 γ1t-1,polyp2}];

(* Polynomials g31 and g32 *)
g31=180 ((γ2)^2)-1664 ((η3)^5)+624 ((η3)^4)+40469 ((η3)^3)+77204 ((η3)^2)+12969 η3-4592;
g32=3420 ((γ3)^2)-7040 ((η3)^5)+7248 ((η3)^4)+156623 ((η3)^3)+254948 ((η3)^2)+40023 η3-12512;

(* We compute a Gröbner basis for the ideal *)
gbI3=GroebnerBasis[idealI3,varsI3,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI3]];
gbI3//Factor//MatrixForm

(* The polynomials g31 and g32 are in the basis *)
g31==gbI3[[5]]//Simplify
g32==gbI3[[4]]//Simplify

(* We double check that each of these polynomials belongs to the Gröbner basis *)
PolynomialReduce[g31,gbI3,varsI3,MonomialOrder->Lexicographic]
PolynomialReduce[g32,gbI3,varsI3,MonomialOrder->Lexicographic]

{(Subscript[p, 2])=0, (Subscript[g, 31])=0, (Subscript[g, 32])=0} has no real solution

(* We set the variables *)
varsI3p={η3,γ2,γ3};

(* Ideal *)
polyp2=128 ((η3)^6)-48 ((η3)^5)-2985 ((η3)^4)-6340 ((η3)^3)-2985 ((η3)^2)-48 η3+128;
g31=180 ((γ2)^2)-1664 ((η3)^5)+624 ((η3)^4)+40469 ((η3)^3)+77204 ((η3)^2)+12969 η3-4592;
g32=3420 ((γ3)^2)-7040 ((η3)^5)+7248 ((η3)^4)+156623 ((η3)^3)+254948 ((η3)^2)+40023 η3-12512;
idealI3p={polyp2,g31,g32};

(* Polynomial g3p *)
g3p=384 ((γ2)^2) ((γ3)^4)+62 ((γ2)^4)+1616 ((γ3)^4)+485 ((γ2)^2) ((γ3)^2)+162 ((γ2)^2)+3582 ((γ3)^2)+1728;

(* We compute a Gröbner basis for the ideal *)
gbI3p=GroebnerBasis[idealI3p,varsI3p,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI3p]];
gbI3p//Factor//MatrixForm

(* The polynomial g3p is in the basis *)
g3p==gbI3p[[2]]//Simplify

(* We double check that this polynomial belongs to the Gröbner basis *)
PolynomialReduce[g3p,gbI3p,varsI3p,MonomialOrder->Lexicographic]

Section 5. Semi-direct extensions with Lorentzian Lie group (ℝ^3)

5.1. The self-adjoint part of the derivation (Subscript[D, sad]) is diagonalizable

(* Setting the Lie brackets and the left-invariant metric *)
dim=4;
ref={e1,e2,e3,e4};
cor=(({{0, 0, 0, η1 e1-γ1 e2+γ2 e3}, {\[Placeholder], 0, 0, γ1 e1+η2 e2+γ3 e3}, {\[Placeholder], \[Placeholder], 0, γ2 e1+γ3 e2+η3 e3}, {\[Placeholder], \[Placeholder], \[Placeholder], 0}}));
g=(({{1, 0, 0, 0}, {\[Placeholder], 1, 0, 0}, {\[Placeholder], \[Placeholder], -1, 0}, {\[Placeholder], \[Placeholder], \[Placeholder], 1}}));

memory;

(* We compute the expression of the Bach tensor *)
BachTMatrix=BachTensor02//Simplify;
BachTMatrix//MatrixForm

(* Since the Bach tensor is a symmetric (0,2)-tensor field we consider the following components *)
BachT={
BachTMatrix[[1,1]],BachTMatrix[[1,2]],BachTMatrix[[1,3]],BachTMatrix[[1,4]],
BachTMatrix[[2,2]],BachTMatrix[[2,3]],BachTMatrix[[2,4]],
BachTMatrix[[3,3]],BachTMatrix[[3,4]],
BachTMatrix[[4,4]]};
BachT//MatrixForm

Checking the components of the Bach tensor

B11=(1/6)(((η1)^4)-((η2)^4)-((η3)^4)-(η2+η3) ((η1)^3)+(3 η1+η3) ((η2)^3)+(3 η1+η2) ((η3)^3)-2 (((η2)^2)+((η3)^2)) ((η1)^2)
	-(7 η1-3 η2-3 η3) η1 η2 η3+9 (((γ1)^2)-((γ2)^2)) ((η1)^2)-3 (5 ((γ1)^2)-((γ3)^2)) ((η2)^2)+3 (5 ((γ2)^2)+((γ3)^2)) ((η3)^2)
	+6 (((γ1)^2)-2 ((γ2)^2)) η1 η2+6 (2 ((γ1)^2)-((γ2)^2)) η1 η3-6 (2 ((γ1)^2)-2 ((γ2)^2)+((γ3)^2)) η2 η3
	+18 γ1 γ2 γ3 (η2-η3));
B12=(1/6)(2 γ1 (((η1)^3)-((η2)^3))-2 γ1 (9 η2+5 η3) ((η1)^2)+2 γ1 (9 η1+5 η3) ((η2)^2)-γ1 (η1-η2) ((η3)^2)
	-6 γ2 γ3 (η1+η2-2 η3) (η1+η2+η3)+3 γ1 (4 ((γ1)^2)-4 ((γ2)^2)-((γ3)^2)) η1
	-3 γ1 (4 ((γ1)^2)-((γ2)^2)-4 ((γ3)^2)) η2+9 γ1 (((γ2)^2)-((γ3)^2)) η3);
B13=(1/6)(2 γ2 (((η1)^3)-((η3)^3))-2 γ2 (5 η2+9 η3) ((η1)^2)-γ2 (η1-η3) ((η2)^2)+2 γ2 (9 η1+5 η2) ((η3)^2)
	-6 γ1 γ3 (η1-2 η2+η3) (η1+η2+η3)+3 γ2 (4 ((γ1)^2)-4 ((γ2)^2)-((γ3)^2)) η1
	-9 γ2 (((γ1)^2)+((γ3)^2)) η2-3 γ2 (((γ1)^2)-4 ((γ2)^2)-4 ((γ3)^2)) η3);
B22=((-1)/6)(((η1)^4)-((η2)^4)+((η3)^4)-(3 η2+η3) ((η1)^3)+(η1+η3) ((η2)^3)-(η1+3 η2) ((η3)^3)+2 (((η1)^2)+((η3)^2)) ((η2)^2)
	-(3 η1-7 η2+3 η3) η1 η2 η3+3 (5 ((γ1)^2)-((γ2)^2)) ((η1)^2)-9 (((γ1)^2)-((γ3)^2)) ((η2)^2)-3 (((γ2)^2)+5 ((γ3)^2)) ((η3)^2)
	-6 (((γ1)^2)-2 ((γ3)^2)) η1 η2+6 (2 ((γ1)^2)+((γ2)^2)-2 ((γ3)^2)) η1 η3-6 (2 ((γ1)^2)-((γ3)^2)) η2 η3
	+18 γ1 γ2 γ3 (η1-η3));
B23=(1/6)(2 γ3 (((η2)^3)-((η3)^3))-γ3 (η2-η3) ((η1)^2)-2 γ3 (5 η1+9 η3) ((η2)^2)+2 γ3 (5 η1+9 η2) ((η3)^2)
	-6 γ1 γ2 (2 η1-η2-η3) (η1+η2+η3)-9 γ3 (((γ1)^2)+((γ2)^2)) η1
	+3 γ3 (4 ((γ1)^2)-((γ2)^2)-4 ((γ3)^2)) η2-3 γ3 (((γ1)^2)-4 ((γ2)^2)-4 ((γ3)^2)) η3);
B33=(1/6)(((η1)^4)+((η2)^4)-((η3)^4)-(η2+3 η3) ((η1)^3)-(η1+3 η3) ((η2)^3)+(η1+η2) ((η3)^3)+2 (((η1)^2)+((η2)^2)) ((η3)^2)
	-(3 η1+3 η2-7 η3) η1 η2 η3+3 (((γ1)^2)-5 ((γ2)^2)) ((η1)^2)+3 (((γ1)^2)-5 ((γ3)^2)) ((η2)^2)+9 (((γ2)^2)+((γ3)^2)) ((η3)^2)
	-6 (((γ1)^2)+2 ((γ2)^2)+2 ((γ3)^2)) η1 η2+6 (((γ2)^2)+2 ((γ3)^2)) η1 η3+6 (2 ((γ2)^2)+((γ3)^2)) η2 η3
	-18 γ1 γ2 γ3 (η1-η2));
B44=-B11-B22+B33;

BachT=={B11,B12,B13,0,B22,B23,0,B33,0,B44}//Simplify

Proof of Theorem 5.3.

Case 1: (Subscript[η, 1])=0

First basis

(* We set the variables *)
varsI1={γ1,γ2,γ3,η1,η3,η2};

(* Ideal *)
idealI1=Join[BachT,{η1}];

(* Polynomial g1 *)
g1=((η2)^8) (((η2-η3))^2);

(* We compute a Gröbner basis for the ideal *)
gbI1=GroebnerBasis[idealI1,varsI1,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI1]];
gbI1//Factor//MatrixForm

(* The polynomial g1 is in the basis *)
g1==gbI1[[1]]//Simplify

(* We double check that this polynomial belongs to the Gröbner basis *)
PolynomialReduce[g1,gbI1,varsI1,MonomialOrder->Lexicographic]

Second basis (if (Subscript[η, 2])=0)

(* We set the variables *)
varsI1p={γ1,γ2,γ3,η1,η3,η2};

(* Ideal *)
idealI1p=Join[BachT,{η1,η2}];

(* Polynomial g1p *)
g1p=((η3)^6);

(* We compute a Gröbner basis for the ideal *)
gbI1p=GroebnerBasis[idealI1p,varsI1p,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI1p]];
gbI1p//Factor//MatrixForm

(* The polynomial g1p is in the basis *)
g1p==gbI1p[[2]]//Simplify

(* We double check that this polynomial belongs to the Gröbner basis *)
PolynomialReduce[g1p,gbI1p,varsI1p,MonomialOrder->Lexicographic]

Third basis (if (Subscript[η, 3])=(Subscript[η, 2])≠0)

(* We set the variables *)
varsI1pp={γ1,γ2,γ3,η1,η3,η2};

(* Ideal *)
idealI1pp=Join[BachT,{η1,η2-η3}];

(* Polynomials g1pp and g2pp *)
g1pp=((γ2)^2) ((η2)^4);
g2pp=(((γ1)^2)-((γ2)^2))((η2)^2);

(* We compute a Gröbner basis for the ideal *)
gbI1pp=GroebnerBasis[idealI1pp,varsI1pp,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI1pp]];
gbI1pp//Factor//MatrixForm

(* The polynomials g1pp and g2pp are in the basis *)
g1pp==gbI1pp[[4]]//Simplify
g2pp==gbI1pp[[16]]//Simplify

(* We double check that each of these polynomials belongs to the Gröbner basis *)
PolynomialReduce[g1pp,gbI1pp,varsI1pp,MonomialOrder->Lexicographic]
PolynomialReduce[g2pp,gbI1pp,varsI1pp,MonomialOrder->Lexicographic]

Case 2: (Subscript[η, 1])≠0, (Subscript[η, 3])=0

(* We set the variables *)
varsI2={γ1,γ2,γ3,η2,η3,η1};

(* Ideal *)
idealI2=Join[BachT,{η3}];

(* Polynomial g2 *)
g2=((η1)^8) (((η1-η2))^2);

(* We compute a Gröbner basis for the ideal *)
gbI2=GroebnerBasis[idealI2,varsI2,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI2]];
gbI2//Factor//MatrixForm

(* The polynomial g2 is in the basis *)
g2==gbI2[[2]]//Simplify

(* We double check that this polynomial belongs to the Gröbner basis *)
PolynomialReduce[g2,gbI2,varsI2,MonomialOrder->Lexicographic]

Proof of Theorem 5.4.

Case 1: (Subscript[η, 1])=(Subscript[η, 2])=1≠(Subscript[η, 3]), (Subscript[η, 3])≠0

(* We set the variables *)
varsI1={γ1,γ2,γ3,η1,η2,η3};

(* Ideal *)
idealI1=Join[BachT,{η1-1,η2-1}];

(* Polynomials g11, g12 and g13 *)
g11=η3(η3-1)(η3-4);
g12=((γ2)^2) (η3-1);
g13=((γ3)^2) (η3-1);

(* We compute a Gröbner basis for the ideal *)
gbI1=GroebnerBasis[idealI1,varsI1,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI1]];
gbI1//Factor//MatrixForm

(* The polynomials g11, g12 and g13 are in the basis *)
g11==gbI1[[1]]//Simplify
g12==gbI1[[6]]//Simplify
g13==gbI1[[4]]//Simplify

(* We double check that each of these polynomials belongs to the Gröbner basis *)
PolynomialReduce[g11,gbI1,varsI1,MonomialOrder->Lexicographic]
PolynomialReduce[g12,gbI1,varsI1,MonomialOrder->Lexicographic]
PolynomialReduce[g13,gbI1,varsI1,MonomialOrder->Lexicographic]

Case 2: (Subscript[η, 1])=(Subscript[η, 3])=1≠(Subscript[η, 2]), (Subscript[η, 2])≠0

(* We set the variables *)
varsI2={γ2,γ1,γ3,η1,η2,η3};

(* Ideal *)
idealI2=Join[BachT,{η1-1,η3-1}];

(* Polynomials g21, g22 and g23 *)
g21= η2(η2-1)(η2-4);
g22=((γ1)^2) (η2-1);
g23=((γ3)^2) (η2-1);

(* We compute a Gröbner basis for the ideal *)
gbI2=GroebnerBasis[idealI2,varsI2,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI2]];
gbI2//Factor//MatrixForm

(* The polynomials g21, g22 and g23 are in the basis *)
g21==gbI2[[2]]//Simplify
g22==gbI2[[6]]//Simplify
g23==gbI2[[4]]//Simplify

(* We double check that each of these polynomials belongs to the Gröbner basis *)
PolynomialReduce[g21,gbI2,varsI2,MonomialOrder->Lexicographic]
PolynomialReduce[g22,gbI2,varsI2,MonomialOrder->Lexicographic]
PolynomialReduce[g23,gbI2,varsI2,MonomialOrder->Lexicographic]

Proof of Theorem 5.5.

Initial basis

(* We set the variables *)
varsI={γ3,γ1,γ2,η1,η2,η3};

(* Ideal *)
idealI=Join[BachT,{η1-1}];

(* Polynomial g *)
polyp1=8 ((η3)^6)+12 ((η3)^5)-45 ((η3)^4)-100 ((η3)^3)-45 ((η3)^2)+12 η3+8;
polyp2=128 ((η3)^6)-48 ((η3)^5)-2985 ((η3)^4)-6340 ((η3)^3)-2985 ((η3)^2)-48 η3+128;
polyq1=(11723776 ((η3)^12)+57769728 ((η3)^11)+106929840 ((η3)^10)+192248648 ((η3)^9)+816486921 ((η3)^8) 
	+2162941668 ((η3)^7)+2962757238 ((η3)^6)+2162941668 ((η3)^5)+816486921 ((η3)^4)
	+192248648 ((η3)^3)+106929840 ((η3)^2)+57769728 η3
	+11723776)(2651014144 ((η3)^12)
	+24455976192 ((η3)^11)+98890309104 ((η3)^10)+43771506920 ((η3)^9) 
	-371896453695 ((η3)^8)+396884502852 ((η3)^7)+2252411448966 ((η3)^6) 
	+396884502852 ((η3)^5)-371896453695 ((η3)^4)+43771506920 ((η3)^3)+98890309104 ((η3)^2) 
	+24455976192 η3+2651014144) (((η3)^2)+η3+4);
polyq2=(η3-4)(7 ((η3)^2)-14 η3-23)(7 ((η3)^2)+7 η3-20);

polyg=γ1 γ2 polyp1 polyp2 η3 (η3-1) polyq1 polyq2;

(* We compute a Gröbner basis for the ideal (time in a home computer: about 25 seconds) *)
gbI=GroebnerBasis[idealI,varsI,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI]];
gbI//Factor//MatrixForm

(* The polynomial g is in the basis *)
polyg==gbI[[10]]//Simplify

(* We double check that this polynomial belongs to the Gröbner basis *)
PolynomialReduce[polyg,gbI,varsI,MonomialOrder->Lexicographic]

Reduction of the case (Subscript[q, 2])=0 (the previous basis is needed)

(* We set the variables *)
varsIq2={η1,η2,η3,γ3,γ2,γ1};

(* Ideal *)
idealIq2=Join[gbI,{polyq2}];

(* Polynomial gq2 *)
polygq2=((γ2)^2) (14 ((γ1)^2)+33) (1421 ((γ1)^4)+28623 ((γ1)^2)+45600);

(* We compute a Gröbner basis for the ideal *)
gbIq2=GroebnerBasis[idealIq2,varsIq2,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbIq2]];
gbIq2//Factor//MatrixForm

(* The polynomial gq2 is in the basis *)
polygq2==gbIq2[[2]]//Simplify

(* We double check that this polynomial belongs to the Gröbner basis *)
PolynomialReduce[polygq2,gbIq2,varsIq2,MonomialOrder->Lexicographic]

Case 1: (Subscript[γ, 1])=0

Initial basis

(* We set the variables *)
varsI1={γ1,γ2,η1,η3,η2,γ3};

(* Ideal *)
idealI1=Join[BachT,{η1-1,γ1}];

(* Polynomial g1 *)
g1=γ2 γ3 (η2-1) (464 ((γ3)^4)-2175 ((γ3)^2)+1824) (2116 ((γ3)^4)-4884 ((γ3)^2)+1089) (80 ((γ3)^4)-95 ((γ3)^2)+32);

(* We compute a Gröbner basis for the ideal *)
gbI1=GroebnerBasis[idealI1,varsI1,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI1]];
gbI1//Factor//MatrixForm

(* The polynomial g1 is in the basis *)
g1==gbI1[[11]]//Simplify

(* We double check that this polynomial belongs to the Gröbner basis *)
PolynomialReduce[g1,gbI1,varsI1,MonomialOrder->Lexicographic]

Case 1.1: (Subscript[γ, 2])=0 - Initial basis

(* We set the variables *)
varsI11={η23t,γ1,γ2,η1,η3,η2,γ3};

(* Ideal *)
idealI1=Join[BachT,{η1-1,γ1}];

idealI11=Join[idealI1,{(η2-η3)η23t-1,γ2}];

(* Polynomial g11 *)
g11=((γ3)^2) (8 ((γ3)^2)-3) (8 ((γ3)^2)-9) (8 ((γ3)^2)+1);

(* We compute a Gröbner basis for the ideal *)
gbI11=GroebnerBasis[idealI11,varsI11,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI11]];
gbI11//Factor//MatrixForm

(* The polynomial g11 is in the basis *)
g11==gbI11[[1]]//Simplify

(* We double check that this polynomial belongs to the Gröbner basis *)
PolynomialReduce[g11,gbI11,varsI11,MonomialOrder->Lexicographic]

Case 1.1: (Subscript[γ, 2])=0 - (ii) 8(Subscript[γ, 3]^2)-3=0

(* We set the variables *)
varsI11p={η23t,γ1,γ2,η1,η3,η2,γ3};

(* Ideal *)
idealI1=Join[BachT,{η1-1,γ1}];
idealI11=Join[idealI1,{(η2-η3)η23t-1,γ2}];

idealI11p=Join[idealI11,{8 ((γ3)^2)-3}];

(* Polynomials g111p and g112p *)
g111p=8 ((η2)^2)-16 η2-7;
g112p=η2+η3-2;

(* We compute a Gröbner basis for the ideal *)
gbI11p=GroebnerBasis[idealI11p,varsI11p,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI11p]];
gbI11p//Factor//MatrixForm

(* The polynomials g111p and g112p are in the basis *)
g111p==gbI11p[[2]]//Simplify
g112p==gbI11p[[3]]//Simplify

(* We double check that each of these polynomials belongs to the Gröbner basis *)
PolynomialReduce[g111p,gbI11p,varsI11p,MonomialOrder->Lexicographic]
PolynomialReduce[g112p,gbI11p,varsI11p,MonomialOrder->Lexicographic]

Case 1.1: (Subscript[γ, 2])=0 - (iii) 8(Subscript[γ, 3]^2)-9=0

(* We set the variables *)
varsI11pp={η23t,γ1,γ2,η1,η3,η2,γ3};

(* Ideal *)
idealI1=Join[BachT,{η1-1,γ1}];
idealI11=Join[idealI1,{(η2-η3)η23t-1,γ2}];

idealI11pp=Join[idealI11,{8 ((γ3)^2)-9}];

(* Polynomials g111pp and g112pp *)
g111pp=8 ((η2)^2)+8 η2-1;
g112pp=η2+η3+1;

(* We compute a Gröbner basis for the ideal *)
gbI11pp=GroebnerBasis[idealI11pp,varsI11pp,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI11pp]];
gbI11pp//Factor//MatrixForm

(* The polynomials g111pp and g112pp are in the basis *)
g111pp==gbI11pp[[2]]//Simplify
g112pp==gbI11pp[[3]]//Simplify

(* We double check that each of these polynomials belongs to the Gröbner basis *)
PolynomialReduce[g111pp,gbI11pp,varsI11pp,MonomialOrder->Lexicographic]
PolynomialReduce[g112pp,gbI11pp,varsI11pp,MonomialOrder->Lexicographic]

Case 1.2: (Subscript[γ, 2])≠0, 464 ((γ3)^4)-2175 ((γ3)^2)+1824=0

(* We set the variables *)
varsI12={η23t,γ3,γ1,γ2,η1,η3,η2};

(* Ideal *)
idealI1=Join[BachT,{η1-1,γ1}];

idealI12=Join[idealI1,{(η2-η3)η23t-1,464 ((γ3)^4)-2175 ((γ3)^2)+1824}];

(* Polynomials g121, g122, g123 and g124 *)
g121=20 ((η2)^2)+47 η2+20;
g122=η2+η3+1;
g123=145 ((γ2)^2)+1316 η2+685;
g124=116 ((γ3)^2)-235 η2-548;

(* We compute a Gröbner basis for the ideal *)
gbI12=GroebnerBasis[idealI12,varsI12,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI12]];
gbI12//Factor//MatrixForm

(* The polynomials g121, g122, g123 and g124 are in the basis *)
g121==gbI12[[1]]//Simplify
g122==gbI12[[2]]//Simplify
g123==gbI12[[4]]//Simplify
g124==gbI12[[6]]//Simplify

(* We double check that each of these polynomials belongs to the Gröbner basis *)
PolynomialReduce[g121,gbI12,varsI12,MonomialOrder->Lexicographic]
PolynomialReduce[g122,gbI12,varsI12,MonomialOrder->Lexicographic]
PolynomialReduce[g123,gbI12,varsI12,MonomialOrder->Lexicographic]
PolynomialReduce[g124,gbI12,varsI12,MonomialOrder->Lexicographic]

Case 1.3: (Subscript[γ, 2])≠0, 2116 ((γ3)^4)-4884 ((γ3)^2)+1089=0

(* We set the variables *)
varsI13={η23t,γ3,γ1,γ2,η1,η3,η2};

(* Ideal *)
idealI1=Join[BachT,{η1-1,γ1}];

idealI13=Join[idealI1,{(η2-η3)η23t-1,2116 ((γ3)^4)-4884 ((γ3)^2)+1089}];

(* Polynomials g131, g132, g133 and g134 *)
g131=23((η2)^2)+74 η2+23;
g132=-η2+2η3-1;
g133=46 ((γ2)^2)+33 η2;
g134=46 ((γ3)^2)+33 η2;

(* We compute a Gröbner basis for the ideal *)
gbI13=GroebnerBasis[idealI13,varsI13,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI13]];
gbI13//Factor//MatrixForm

(* The polynomials g131, g132, g133 and g134 are in the basis *)
g131==gbI13[[1]]//Simplify
g132==gbI13[[2]]//Simplify
g133==gbI13[[4]]//Simplify
g134==gbI13[[6]]//Simplify

(* We double check that each of these polynomials belongs to the Gröbner basis *)
PolynomialReduce[g131,gbI13,varsI13,MonomialOrder->Lexicographic]
PolynomialReduce[g132,gbI13,varsI13,MonomialOrder->Lexicographic]
PolynomialReduce[g133,gbI13,varsI13,MonomialOrder->Lexicographic]
PolynomialReduce[g134,gbI13,varsI13,MonomialOrder->Lexicographic]

Case 2: (Subscript[γ, 1])≠0, (Subscript[γ, 2])=0

(* We set the variables *)
varsI2={γ3,η3,η1,η2,γ2,γ1};

(* Ideal *)
idealI2=Join[BachT,{η1-1,γ2}];

(* Polynomial g2 *)
g2=((γ1)^2) (η2-1) (((γ1)^4)+90 ((γ1)^2)+81) (5 ((γ1)^4)+5 ((γ1)^2)+2) (25 ((γ1)^4)+2 ((γ1)^2)+1) (49 ((γ1)^4)+138 ((γ1)^2)+9) (725 ((γ1)^4)+8613 ((γ1)^2)+2850) (2116 ((γ1)^4)+4884 ((γ1)^2)+1089);

(* We compute a Gröbner basis for the ideal *)
gbI2=GroebnerBasis[idealI2,varsI2,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI2]];
gbI2//Factor//MatrixForm

(* The polynomial g2 is in the basis *)
g2==gbI2[[2]]//Simplify

(* We double check that this polynomial belongs to the Gröbner basis *)
PolynomialReduce[g2,gbI2,varsI2,MonomialOrder->Lexicographic]

Case 3: (Subscript[γ, 1])≠0, (Subscript[p, 1])=0

(* We need de Initial basis (time in a home computer: about 25 seconds) *)
varsI={γ3,γ1,γ2,η1,η2,η3};
idealI=Join[BachT,{η1-1}];
gbI=GroebnerBasis[idealI,varsI,MonomialOrder->Lexicographic];

(* We set the variables *)
varsI3={γ1t,γ1,γ2,γ3,η1,η2,η3};

(* Ideal *)
polyp1=8 ((η3)^6)+12 ((η3)^5)-45 ((η3)^4)-100 ((η3)^3)-45 ((η3)^2)+12 η3+8;
idealI3=Join[gbI,{γ1 γ1t-1,polyp1}];

(* Polynomials g31, g32, g33 and g34 *)
g31=15 ((γ2)^2)+16 ((η3)^5)+24 ((η3)^4)-106 ((η3)^3)-196 ((η3)^2)-21 η3+28;
g32=15 ((γ3)^2)+88 ((η3)^5)+60 ((η3)^4)-535 ((η3)^3)-670 ((η3)^2)+15 η3+112;
g33=10 η2-56 ((η3)^5)-36 ((η3)^4)+339 ((η3)^3)+414 ((η3)^2)-η3-70;
g34=40 γ1+γ2 γ3 (408 ((η3)^5)+452 ((η3)^4)-2535 ((η3)^3)-4072 ((η3)^2)-359 η3+856);

(* We compute a Gröbner basis for the ideal *)
gbI3=GroebnerBasis[idealI3,varsI3,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI3]];
gbI3//Factor//MatrixForm

(* The polynomials g31, g32, g33 and g34 are in the basis *)
g31==gbI3[[5]]//Simplify
g32==gbI3[[4]]//Simplify
g33==gbI3[[2]]//Simplify
g34==gbI3[[6]]//Simplify

(* We double check that each of these polynomials belongs to the Gröbner basis *)
PolynomialReduce[g31,gbI3,varsI3,MonomialOrder->Lexicographic]
PolynomialReduce[g32,gbI3,varsI3,MonomialOrder->Lexicographic]
PolynomialReduce[g33,gbI3,varsI3,MonomialOrder->Lexicographic]
PolynomialReduce[g34,gbI3,varsI3,MonomialOrder->Lexicographic]

Case 4: (Subscript[γ, 1])≠0, (Subscript[p, 2])=0

(* We need de Initial basis (time in a home computer: about 25 seconds) *)
varsI={γ3,γ1,γ2,η1,η2,η3};
idealI=Join[BachT,{η1-1}];
gbI=GroebnerBasis[idealI,varsI,MonomialOrder->Lexicographic];

(* We set the variables *)
varsI4={γ1t,γ1,γ2,γ3,η1,η2,η3};

(* Ideal *)
polyp2=128 ((η3)^6)-48 ((η3)^5)-2985 ((η3)^4)-6340 ((η3)^3)-2985 ((η3)^2)-48 η3+128;
idealI4=Join[gbI,{γ1 γ1t-1,polyp2}];

(* Polynomials g41, g42, g43 and g44 *)
g41=180 ((γ2)^2)+1664 ((η3)^5)-624 ((η3)^4)-40469 ((η3)^3)-77204 ((η3)^2)-12969 η3+4592;
g42=3420 ((γ3)^2)+7040 ((η3)^5)-7248 ((η3)^4)-156623 ((η3)^3)-254948 ((η3)^2)-40023 η3+12512;
g43=760 η2-5504 ((η3)^5)+5136 ((η3)^4)+124131 ((η3)^3)+207636 ((η3)^2)+33731 η3-12040;
g44=82080 γ1-γ2 γ3 (477824 ((η3)^5)-328688 ((η3)^4)-11131997 ((η3)^3)-
				19885316 ((η3)^2)-3537085 η3+1274912);

(* We compute a Gröbner basis for the ideal *)
gbI4=GroebnerBasis[idealI4,varsI4,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI4]];
gbI4//Factor//MatrixForm

(* The polynomials g41, g42, g43 and g44 are in the basis *)
g41==gbI4[[5]]//Simplify
g42==gbI4[[4]]//Simplify
g43==gbI4[[2]]//Simplify
g44==gbI4[[6]]//Simplify

(* We double check that each of these polynomials belongs to the Gröbner basis *)
PolynomialReduce[g41,gbI4,varsI4,MonomialOrder->Lexicographic]
PolynomialReduce[g42,gbI4,varsI4,MonomialOrder->Lexicographic]
PolynomialReduce[g43,gbI4,varsI4,MonomialOrder->Lexicographic]
PolynomialReduce[g44,gbI4,varsI4,MonomialOrder->Lexicographic]

5.2. The self-adjoint part of the derivation (Subscript[D, sad]) has complex eigenvalues

(* Setting the Lie brackets and the left-invariant metric *)
dim=4;
ref={e1,e2,e3,e4};
cor=(({{0, 0, 0, η e1-γ1 e2+γ2 e3}, {\[Placeholder], 0, 0, γ1 e1+δ e2+(γ3-ν) e3}, {\[Placeholder], \[Placeholder], 0, γ2 e1+(γ3+ν) e2+δ e3}, {\[Placeholder], \[Placeholder], \[Placeholder], 0}}));
g=(({{1, 0, 0, 0}, {\[Placeholder], 1, 0, 0}, {\[Placeholder], \[Placeholder], -1, 0}, {\[Placeholder], \[Placeholder], \[Placeholder], 1}}));

memory;

(* We compute the expression of the Bach tensor *)
BachTMatrix=BachTensor02//Simplify;
BachTMatrix//MatrixForm

(* Since the Bach tensor is a symmetric (0,2)-tensor field we consider the following components *)
BachT={
BachTMatrix[[1,1]],BachTMatrix[[1,2]],BachTMatrix[[1,3]],BachTMatrix[[1,4]],
BachTMatrix[[2,2]],BachTMatrix[[2,3]],BachTMatrix[[2,4]],
BachTMatrix[[3,3]],BachTMatrix[[3,4]],
BachTMatrix[[4,4]]};
BachT//MatrixForm

Proof of Theorem 5.9.

Checking the components of the Bach tensor

B11=(1/6)((η-δ) (η+3 δ) ((η^2)-4 η δ+9 ((γ1)^2)-9 ((γ2)^2))-3 (ν^2) ((η^2)-4 (δ^2)+4 η δ)
	+12 γ1 γ2 ν (η+5 δ)-ν(18 γ3 (((γ1)^2)+((γ2)^2))-3 ν (((γ1)^2)-((γ2)^2)-4 ((γ3)^2))+4 (ν^3)));
B12=(1/6)((η-δ)(γ1 (2 (η^2)-9 (δ^2)-26 η δ)-6 γ2 γ3 (η+2 δ)+3 γ1 (4 ((γ1)^2)-4 ((γ2)^2)-((γ3)^2)))
	+γ2 ν (8 (η^2)-3 (δ^2)-38 η δ)+γ1 ν((18 γ3-17 ν) η+(36 γ3+17 ν) δ)
	+γ2 ν (21 ((γ1)^2)-3 ((γ2)^2)-21 ((γ3)^2)-11 (ν^2)));
B13=(1/6)((η-δ)(γ2 (2 (η^2)-9 (δ^2)-26 η δ)-6 γ1 γ3 (η+2 δ)+3 γ2 (4 ((γ1)^2)-4 ((γ2)^2)-((γ3)^2)))
	-γ1 ν (8 (η^2)-3 (δ^2)-38 η δ)-γ2 ν((36 γ3-17 ν) δ+(18 γ3+17 ν) η)
	-γ1 ν (3 ((γ1)^2)-21 ((γ2)^2)-21 ((γ3)^2)-11 (ν^2)));
B22=((-1)/6)((η-δ)((η^3)-3 (η^2) δ-4 η (δ^2)+3(5 ((γ1)^2)-((γ2)^2)) η+3 (7 ((γ1)^2)+((γ2)^2)) δ+18 γ1 γ2 γ3)
	+ν((2 γ3-5 ν) (η^2)+4 (6 γ3+ν) (δ^2)+4 (10 γ3+ν) η δ)+12 γ1 γ2 ν (2 η+δ)
	+ν(24 ((γ3)^3)-24 ((γ1)^2) γ3+6 ((γ2)^2) γ3-3 ν (((γ1)^2)-((γ2)^2)-4 ((γ3)^2))+40 (ν^2) γ3+4 (ν^3)));
B23=((-1)/6)(3 (η-δ)(4 γ1 γ2 (η+2 δ)+3 γ3(((γ1)^2)+((γ2)^2)))
	-2 (η+2 δ) ν ((η^2)-4 η δ+3 ((γ1)^2)-3 ((γ2)^2)-12 ((γ3)^2)-4 (ν^2)));
B44=(1/6)((((η-δ))^2) ((η^2)-4 η δ+9 ((γ1)^2)-9 ((γ2)^2))-(ν^2) (7 (η^2)+4 (δ^2)-20 η δ)+36 γ1 γ2 ν (η-δ)
	-3 (ν^2) (3 ((γ1)^2)-3 ((γ2)^2)-12 ((γ3)^2)-4 (ν^2)));
B33=B11+B22+B44;

BachT=={B11,B12,B13,0,B22,B23,0,B33,0,B44}//Simplify

Initial basis

(* We set the variables *)
varsI={γ1,γ2,γ3,δ,η,ν};

(* Ideal *)
idealI=Join[BachT,{ν-1}];

(* Polynomial g *)
polyg=γ2 γ3 ((η^2)+4) (9 (η^2)+4) (9 (η^4)-12 (η^2)+64) (9 (η^4)+2676 (η^2)+1600)(225 (η^4)+
	2220 (η^2)+2116) (75 (η^6)+600 (η^4)+1200 (η^2)+256) (6075 (η^6)+
	81000 (η^4)+74160 (η^2)+4096) (150921225 (η^12)+
	352923480 (η^10)+5525229024 (η^8)-1113219072 (η^6) +
	7682072832 (η^4)+37029888 (η^2)+750321664) (41280080625 (η^12)+
	135680265000 (η^10)+3185958636000 (η^8)+
	719476262400 (η^6)+802494470400 (η^4)-671994193920 (η^2)+169664905216);

(* We compute a Gröbner basis for the ideal (time in a home computer: about 10 seconds) *)
gbI=GroebnerBasis[idealI,varsI,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI]];
gbI//Factor//MatrixForm

(* The polynomial g is in the basis *)
polyg==gbI[[9]]//Simplify

(* We double check that this polynomial belongs to the Gröbner basis *)
PolynomialReduce[polyg,gbI,varsI,MonomialOrder->Lexicographic]

Case 1: ν=1, (Subscript[γ, 2])=0

(* We set the variables *)
varsI1={δ,η,γ1,γ2,γ3,ν};

(* Ideal *)
idealI=Join[BachT,{ν-1}];
idealI1=Join[idealI,{γ2}];

(* Polynomials g11 and g12 *)
g11=((γ3)^2) (((γ3)^2)+3) (3 ((γ3)^2)+1) (5 ((γ3)^2)+1);
g12=γ1 (30 (η^2)+7);

(* We compute a Gröbner basis for the ideal *)
gbI1=GroebnerBasis[idealI1,varsI1,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI1]];
gbI1//Factor//MatrixForm

(* The polynomials g11 and g12 are in the basis *)
g11==gbI1[[2]]//Simplify
g12==gbI1[[7]]//Simplify

(* We double check that each of these polynomials belongs to the Gröbner basis *)
PolynomialReduce[g11,gbI1,varsI1,MonomialOrder->Lexicographic]
PolynomialReduce[g12,gbI1,varsI1,MonomialOrder->Lexicographic]

Case 2: ν=1, (Subscript[γ, 2])≠0, (Subscript[γ, 3])=0

First basis

(* We set the variables *)
varsI2={γ1,γ2t,γ2,γ3,δ,η,ν};

(* Ideal *)
idealI=Join[BachT,{ν-1}];
idealI2=Join[idealI,{γ2 γ2t-1,γ3}];

(* Polynomial g2 *)
g2=(15 (η^2)-4) ((η^2)+4) (9 (η^2)+4) (30 (η^2)+7) (87 (η^2)+28) (225 (η^4)+1416 (η^2)+16) (225 (η^4)+
2760 (η^2)+784) (441 (η^4)-696 (η^2)+400);

(* We compute a Gröbner basis for the ideal *)
gbI2=GroebnerBasis[idealI2,varsI2,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI2]];
gbI2//Factor//MatrixForm

(* The polynomial g2 is in the basis *)
g2==gbI2[[2]]//Simplify

(* We double check that this polynomial belongs to the Gröbner basis *)
PolynomialReduce[g2,gbI2,varsI2,MonomialOrder->Lexicographic]

Second basis

(* We set the variables *)
varsI2p={γ1,γ2t,γ2,γ3,δ,η,ν};

(* Ideal *)
idealI=Join[BachT,{ν-1}];
idealI2=Join[idealI,{γ2 γ2t-1,γ3}];
idealI2p=Join[idealI2,{15 (η^2)-4}];

(* Polynomials g21p, g22p and g23p *)
g21p=η+2 δ;
g22p=375 ((γ2)^4)+500 ((γ2)^2)-4;
g23p=4 γ1+25 γ2(3 ((γ2)^2)+4) η;

(* We compute a Gröbner basis for the ideal *)
gbI2p=GroebnerBasis[idealI2p,varsI2p,MonomialOrder->Lexicographic];
Print["Number of polynomials in the basis: ",Length[gbI2p]];
gbI2p//Factor//MatrixForm

(* The polynomials g21p, g22p and g23p are in the basis *)
g21p==gbI2p[[3]]//Simplify
g22p==gbI2p[[5]]//Simplify
g23p==gbI2p[[7]]//Simplify

(* We double check that each of these polynomials belongs to the Gröbner basis *)
PolynomialReduce[g21p,gbI2p,varsI2p,MonomialOrder->Lexicographic]
PolynomialReduce[g22p,gbI2p,varsI2p,MonomialOrder->Lexicographic]
PolynomialReduce[g23p,gbI2p,varsI2p,MonomialOrder->Lexicographic]

5.3. The self-adjoint part of the derivation (Subscript[D, sad]) has a double root

(* Setting the Lie brackets and the left-invariant metric *)
dim=4;
ref={u1,u2,u3,u4};
cor=(({{0, 0, 0, (η1+γ1) u1+ϵ u2-γ3 u3}, {\[Placeholder], 0, 0, (η1-γ1) u2-γ2 u3}, {\[Placeholder], \[Placeholder], 0,  γ2 u1+γ3 u2+η2 u3}, {\[Placeholder], \[Placeholder], \[Placeholder], 0}}));
g=(({{0, 1, 0, 0}, {\[Placeholder], 0, 0, 0}, {\[Placeholder], \[Placeholder], 1, 0}, {\[Placeholder], \[Placeholder], \[Placeholder], 1}}));

memory;

(* We compute the expression of the Bach tensor *)
BachTMatrix=Simplify[BachTensor02,{(ϵ^2)==1}];
BachTMatrix//MatrixForm

(* Since the Bach tensor is a symmetric (0,2)-tensor field we consider the following components *)
BachT={
BachTMatrix[[1,1]],BachTMatrix[[1,2]],BachTMatrix[[1,3]],BachTMatrix[[1,4]],
BachTMatrix[[2,2]],BachTMatrix[[2,3]],BachTMatrix[[2,4]],
BachTMatrix[[3,3]],BachTMatrix[[3,4]],
BachTMatrix[[4,4]]};
BachT//MatrixForm

Proof of Theorem 5.13.

Checking the components of the Bach tensor

B11=(1/3)(3 ((γ3)^2) (4 η1+2 η2+3 γ1) (η1-η2)-ϵ(η2 (4 η1-η2)+12 (((γ1)^2)-γ2 γ3)) (2 η1+η2)
	-ϵ γ1 (12 ((η1)^2)+((η2)^2)+20 η1 η2+12 ((γ1)^2)-30 γ2 γ3)-3 ϵ γ3 (4 η1+2 η2+5 γ1) γ2);
B12=((-1)/6)(η2 (η1-η2) (4 η1-η2) (η1+η2)-3 ϵ (2 η1+4 η2+3 γ1) ((γ2)^2)-18 γ3 (((η1)^2)-((η2)^2)) γ2);
B13=((-1)/6)(γ3 (η1-η2)(9 ((η1)^2)-2 ((η2)^2)+26 η1 η2+6 γ1 (2 η1+η2)+3 ((γ1)^2))-15 ϵ γ3 ((γ2)^2)
	+ϵ((3 ((η1)^2)-8 ((η2)^2)+38 η1 η2)-12 (2 ϵ ((γ3)^2)-3 γ1) η1+6 (4 ϵ ((γ3)^2)+3 γ1) η2+21 ((γ1)^2)) γ2);
B22=((γ2)^2) (η1-η2) (4 η1+2 η2-3 γ1);
B23=((-1)/6)(γ2((η1-η2) (9 ((η1)^2)-2 ((η2)^2)+26 η1 η2-6 γ1 (2 η1+η2)+3 (((γ1)^2)-8 γ2 γ3))-9 ϵ ((γ2)^2)));
B44=((-1)/6)(η2 (((η1-η2))^2) (4 η1-η2)-18 (η1-η2)(ϵ γ2+γ3 (η1-η2)) γ2);
B33=-2B12-B44;

Simplify[BachT=={B11,B12,B13,0,B22,B23,0,B33,0,B44},{(ϵ^2)==1}]

5.4. The self-adjoint part of the derivation (Subscript[D, sad]) has a triple root

(* Setting the Lie brackets and the left-invariant metric *)
dim=4;
ref={u1,u2,u3,u4};
cor=(({{0, 0, 0, (η+γ1) u1-γ3 u3}, {\[Placeholder], 0, 0, (η-γ1) u2-(γ2-1) u3}, {\[Placeholder], \[Placeholder], 0, (γ2+1) u1+γ3 u2+η u3}, {\[Placeholder], \[Placeholder], \[Placeholder], 0}}));
g=(({{0, 1, 0, 0}, {\[Placeholder], 0, 0, 0}, {\[Placeholder], \[Placeholder], 1, 0}, {\[Placeholder], \[Placeholder], \[Placeholder], 1}}));

memory;

(* We compute the expression of the Bach tensor *)
BachTMatrix=BachTensor02//Simplify;
BachTMatrix//MatrixForm

(* Since the Bach tensor is a symmetric (0,2)-tensor field we consider the following components *)
BachT={
BachTMatrix[[1,1]],BachTMatrix[[1,2]],BachTMatrix[[1,3]],BachTMatrix[[1,4]],
BachTMatrix[[2,2]],BachTMatrix[[2,3]],BachTMatrix[[2,4]],
BachTMatrix[[3,3]],BachTMatrix[[3,4]],
BachTMatrix[[4,4]]};
BachT//MatrixForm

Proof of Theorem 5.17.

Checking the components of the Bach tensor

B11=-(3 ((γ3)^3));
B12=(1/2)(γ3(11 (η^2)-6 γ1 η+((γ1)^2)-γ3(8 γ2+3)));
B13=9 ((γ3)^2) η;
B22=11 γ2 (η^2)-18 γ1 γ2 η+7 ((γ1)^2) γ2-γ3(5 ((γ2)^2)+3);
B23=((-1)/2)(6 (η^3)-11 γ1 (η^2)+6 (((γ1)^2)-5 γ2 γ3) η-γ1(((γ1)^2)-14 γ2 γ3));
B33=((-1)/2)(γ3(22 (η^2)-12 γ1 η+2 ((γ1)^2)-γ3 (16 γ2-3)));
B44=(1/2)(9 ((γ3)^2));

BachT=={B11,B12,B13,0,B22,B23,0,B33,0,B44}//Simplify

Section 6. Semi-direct extensions with degenerate Lie group (ℝ^3)

(* Setting the Lie brackets and the left-invariant metric *)
dim=4;
ref={u1,u2,u3,u4};
cor=(({{0, 0, 0,  γ1 u1+γ2 u2+γ3 u3}, {\[Placeholder], 0, 0, γ4 u1+γ5 u2+γ6 u3}, {\[Placeholder], \[Placeholder], 0, γ7 u1+γ8 u2+γ9 u3}, {\[Placeholder], \[Placeholder], \[Placeholder], 0}}));
g=(({{1, 0, 0, 0}, {\[Placeholder], 1, 0, 0}, {\[Placeholder], \[Placeholder], 0, 1}, {\[Placeholder], \[Placeholder], \[Placeholder], 0}}));

memory;

(* We compute the expression of the Bach tensor *)
BachTMatrix=BachTensor02//Simplify;
BachTMatrix//MatrixForm

(* Since the Bach tensor is a symmetric (0,2)-tensor field we consider the following components *)
BachT={
BachTMatrix[[1,1]],BachTMatrix[[1,2]],BachTMatrix[[1,3]],BachTMatrix[[1,4]],
BachTMatrix[[2,2]],BachTMatrix[[2,3]],BachTMatrix[[2,4]],
BachTMatrix[[3,3]],BachTMatrix[[3,4]],
BachTMatrix[[4,4]]};
BachT//MatrixForm

Checking the components of the Bach tensor

B11=((-1)/6)((((γ7)^2)+((γ8)^2))(5 ((γ7)^2)+((γ8)^2)));
B12=((-1)/3)(2 γ7 γ8 (((γ7)^2)+((γ8)^2)));
B14=(1/3)(2((2 γ1+γ5) γ7+γ2 γ8) (((γ7)^2)+((γ8)^2)));
B22=((-1)/6)((((γ7)^2)+((γ8)^2)) (((γ7)^2)+5 ((γ8)^2)));
B24=(1/3)(2(γ4 γ7+(γ1+2 γ5) γ8) (((γ7)^2)+((γ8)^2)));
B34=(1/2)((((((γ7)^2)+((γ8)^2)))^2));
B44=((-1)/12)(4 γ3 ((γ7)^3)+4 γ6 ((γ8)^3)+4 γ6 ((γ7)^2) γ8+4 γ3 γ7 ((γ8)^2)+18 (γ1 γ2+γ4 γ5) γ7 γ8
	+(16 ((γ1)^2)-((γ2)^2)+8 ((γ4)^2)+7 ((γ5)^2)+12 γ1 γ5+4 γ2 γ4-4 γ1 γ9-4 γ5 γ9) ((γ7)^2)
	+(7 ((γ1)^2)+8 ((γ2)^2)-((γ4)^2)+16 ((γ5)^2)+12 γ1 γ5+4 γ2 γ4-4 γ1 γ9-4 γ5 γ9) ((γ8)^2));

BachT=={B11,B12,0,B14,B22,0,B24,0,B34,B44}//Simplify