// first group

R<u>:= PolynomialRing(Integers());
v:= u^2-2;
K<t>:= NumberField(v);
print "\n";

z:= ZeroMatrix(K,7,7);

a:= z;     // alpha
a[1,6]:=t; a[3,1]:=-t; a[5,4]:=1; a[7,2]:=-1;

b:= z;     // beta
b[2,3]:=1; b[6,5]:=-1;

c:= z;     // alpha + beta
c[1,5]:=t; c[2,1]:=-t; c[6,4]:=-1; c[7,3]:=1;
print "c = ", c, "\n";

d:= z;     // 2 alpha + beta
d[1,4]:=t; d[7,1]:=-t; d[2,6]:=-1; d[3,5]:=1;
print "d = ", d, "\n";

e:= z;     // 3 alpha + beta
e[3,4]:=1; e[7,6]:=-1;

f:= z;     // 3 alpha + 2 beta
f[2,4]:=1; f[7,5]:=-1;


// second group

print "ab-ba = -c", a*b-b*a eq -c;     // prints ab-ba = -c true
print "ac-ca = -2d", a*c-c*a eq (-2)*d;
print "ad-da = -3e", a*d-d*a eq (-3)*e;
print "ae-ea = 0", a*e-e*a eq z;
print "af-fa = 0", a*f-f*a eq z, "\n";

print "bc-cb = 0", b*c-c*b eq z;
print "bd-db = 0", b*d-d*b eq z;
print "be-eb = f", b*e-e*b eq f;
print "bf-fb = 0", b*f-f*b eq z, "\n";

print "cd-dc = -3f", c*d-d*c eq (-3)*f;
print "ce-ec = 0", c*e-e*c eq z;
print "cf-fc = 0", c*f-f*c eq z, "\n";

print "de-ed = 0", d*e-e*d eq z;
print "df-fd = 0", d*f-f*d eq z, "\n";

print "ef-fe = 0", e*f-f*e eq z, "\n";

print "fa is not 0", f*a ne z;     // prints fa is not 0 true
print "fb = 0", f*b eq z, "\n";


// third group

h1:=z;
h1[2,2]:=1; h1[5,5]:=-1;
h1[3,3]:=-1; h1[6,6]:=1;

print "h1 a - a h1 = -a", h1*a - a*h1 eq -a;     // alpha = -1
print "h1 b - b h1 = 2b", h1*b - b*h1 eq 2*b;     // beta = 2
print " ";

h2:=z;
h2[3,3]:=1; h2[6,6]:=-1;
h2[4,4]:=-1; h2[7,7]:=1;

print "h2 a - a h2 = a", h2*a - a*h2 eq a;     // alpha = 1
print "h2 b - b h2 = -b", h2*b - b*h2 eq -b;    // beta = -1
print " ";

h3:=z;
h3[2,2]:=-1; h3[5,5]:=1;
h3[4,4]:=1; h3[7,7]:=-1;

print "h3 a - a h3 = 0", h3*a - a*h3 eq z;     // alpha = 0
print "h3 b - b h3 = -b", h3*b - b*h3 eq -b;    // beta = -1
print " ";


//  fourth group

V:= VectorSpace(K,7);

VNN:= V*c + V*d + V*e + V*f;     // V [N,N]
print "dim V [N,N] = ", Dimension(VNN);     // dim = 5
print "basis of V [N,N]", Basis(VNN), "\n";      // v1, v3, v4, v5, v6

VA:= V*a + VNN;      // V N_a
print "dim V N_a = ", Dimension(VA);    // dim = 6

W:= VA;
WA:= W*a  + W*c + W*d + W*e + W*f;     // W N_a
print "dim W N_a = ", Dimension(WA);     // dim = 4

print "dim (W N_a + W h1) = ", Dimension(WA + W*h1);     // dim = 6
print "dim (W N_a + W h2) = ", Dimension(WA + W*h2);     // dim = 5
print "dim (W N_a + W h3) = ", Dimension(WA + W*h3), "\n";     // dim = 5

VB:= V*b + VNN;     // V N_b
print "dim V N_b = ", Dimension(VB);     // dim = 5
print "basis of V N_b", Basis(VB), "\n";     // v1, v3, v4, v5, v6

W1:= VB + V*h1;     // V N_b + V h1
print "basis of W1", Basis(W1);     // v1, v2, v3, v4, v5, v6
W1B:= W1*b + W1*c + W1*d + W1*e + W1*f;     // W1 N_b
print "dim W1 N_b = ", Dimension(W1B), "\n";     // dim = 5

W2:= VB + V*h2;     // V N_b + V h2
print "basis of W2", Basis(W2);     // v1, v3, v4, v5, v6, v7
W2B:= W2*b + W2*c + W2*d + W2*e + W2*f;     // W2 N_b
print "dim W2 N_b = ", Dimension(W2B), "\n";     // dim = 5

VBB:= VB*b + VB*c + VB*d + VB*e + VB*f;     // V N_b N_b
print "dim V N_b N_b = ", Dimension(VBB);     // dim = 2
print "basis of V N_b N_b", Basis(VBB);     // v4, v5