#the K field with the parameters we will need
K.<a,b,c,A,B,C,x1,x2,x3,x4,x5,x6,x7,x8,x9,x10,x11,x12,x13,x14,x15,x16,x17,x18,x19,x20,x21,x22,x23,x24,x25,x26,x27,x28,x29,x30,x31,x32,x33,x34,x35,x36,x37,x38,x39,x40,x41,x42,x43,x44,x45,x46,x47,x48,x49,x50,x51,x52,x53,x54,x55,x56,x57,x58,x59,x60,x61,x62,x63,x64,y1,y2,y3,y4,y5,y6,y7,y8,y9,y10,y11,y12,y13,y14,y15,y16,y17,y18,y19,y20,y21,y22,y23,y24,y25,y26,y27,y28,y29,y30,y31,y32,y33,y34,y35,y36,y37,y38,y39,y40,y41,y42,y43,y44,y45,y46,y47,y48,y49,y50,y51,y52,y53,y54,y55,y56,y57,y58,y59,y60,y61,y62,y63,y64,z1,z2,z3,z4,z5,z6,z7,z8,z9,z10,z11,z12,z13,z14,z15,z16,z17,z18,z19,z20,z21,z22,z23,z24,z25,z26,z27,z28,z29,z30,z31,z32,z33,z34,z35,z36,z37,z38,z39,z40,z41,z42,z43,z44,z45,z46,z47,z48,z49,z50,z51,z52,z53,z54,z55,z56,z57,z58,z59,z60,z61,z62,z63,z64,u1,u2,u3,u4,u5,u6,u7,u8,u9,u10,u11,u12,u13,u14,u15,u16,u17,u18,u19,u20,u21,u22,u23,u24,u25,u26,u27,u28,u29,u30,u31,u32,u33,u34,u35,u36,u37,u38,u39,u40,u41,u42,u43,u44,u45,u46,u47,u48,u49,u50,u51,u52,u53,u54,u55,u56,u57,u58,u59,u60,u61,u62,u63,u64,s1,s2,s3,s4,s5,s6,s7,s8,s9,s10,s11,s12,s13,s14,r1,r2,r3,r4,r5,r6,r7,r8,r9,r10,r11,r12,r13,r14>=QQ[]
K=K.fraction_field()

#definition of two parent classes of octonion algebras
octpar1=octonion_parent(K,a,b,c)
octpar2=octonion_parent(K,A,B,C)

#definition of the parent class of bioctonions
OxOpar=OxOparent(octpar1,octpar2)

#we create 4 symbolic elements
x=OxOpar([x1,x2,x3,x4,x5,x6,x7,x8,x9,x10,x11,x12,x13,x14,x15,x16,x17,x18,x19,x20,x21,x22,x23,x24,x25,x26,x27,x28,x29,x30,x31,x32,x33,x34,x35,x36,x37,x38,x39,x40,x41,x42,x43,x44,x45,x46,x47,x48,x49,x50,x51,x52,x53,x54,x55,x56,x57,x58,x59,x60,x61,x62,x63,x64])
y=OxOpar([y1,y2,y3,y4,y5,y6,y7,y8,y9,y10,y11,y12,y13,y14,y15,y16,y17,y18,y19,y20,y21,y22,y23,y24,y25,y26,y27,y28,y29,y30,y31,y32,y33,y34,y35,y36,y37,y38,y39,y40,y41,y42,y43,y44,y45,y46,y47,y48,y49,y50,y51,y52,y53,y54,y55,y56,y57,y58,y59,y60,y61,y62,y63,y64])
z=OxOpar([z1,z2,z3,z4,z5,z6,z7,z8,z9,z10,z11,z12,z13,z14,z15,z16,z17,z18,z19,z20,z21,z22,z23,z24,z25,z26,z27,z28,z29,z30,z31,z32,z33,z34,z35,z36,z37,z38,z39,z40,z41,z42,z43,z44,z45,z46,z47,z48,z49,z50,z51,z52,z53,z54,z55,z56,z57,z58,z59,z60,z61,z62,z63,z64])
u=OxOpar([u1,u2,u3,u4,u5,u6,u7,u8,u9,u10,u11,u12,u13,u14,u15,u16,u17,u18,u19,u20,u21,u22,u23,u24,u25,u26,u27,u28,u29,u30,u31,u32,u33,u34,u35,u36,u37,u38,u39,u40,u41,u42,u43,u44,u45,u46,u47,u48,u49,u50,u51,u52,u53,u54,u55,u56,u57,u58,u59,u60,u61,u62,u63,u64])
#we create a symbolic skew elements
s=OxOpar.embed_in_OxO(octpar1([0,s1,s2,s3,s4,s5,s6,s7]),octpar2.one())+OxOpar.embed_in_OxO(octpar1.one(),octpar2([0,s8,s9,s10,s11,s12,s13,s14]))
r=OxOpar.embed_in_OxO(octpar1([0,r1,r2,r3,r4,r5,r6,r7]),octpar2.one())+OxOpar.embed_in_OxO(octpar1.one(),octpar2([0,r8,r9,r10,r11,r12,r13,r14]))

print "Let O1 be the octonion algebra (a,b,c) and O2 be the octonion algebra (A,B,C). O1xO2 is the tensor product of these algebras. x,y,z,u are symbolic elements of O1xO2. s,r are symbolic skew elements of O1xO2."

#test if the standard involution is an involution
print("si(x*y)==si(y)*si(x)? %s" %((x*y).si()==y.si()*x.si()))

#tests on skew element by Lemma 3.13
print ("verification of Lemma 3.13.(i): %s" %(s*(r*s)==(s*r)*s))
print ("verification of Lemma 3.13.(ii): %s" %((s*(r*s))*x==(s*(r*(s*x)))))
print ("verification of Lemma 3.13.(iii): %s" %(s*(s.inverse()*x)==x))

#to test if OxO is structurable, we introduce the V-operator and the operator T_x=V_{x,1}
def V(a,b,c):
    return a.V(b,c)
def T(x,z):
    return x*z+z*x-z*(x.si())

#testing if identity (3.2) holds for x,y,z,u as above asks to much memory. We can only do this for random values of x,y,z,u
x_rand=OxOpar.random()
y_rand=OxOpar.random()
z_rand=OxOpar.random()
u_rand=OxOpar.random()

print("\nWe verify if the multiplication and involution we defined on OxO, gives indeed rise to a structurable algebra. We can only verify this for x,y,z and u random elemts in OxO with coefficients in Q. We verify if formula (3.2) holds.")
print T(z_rand,V(x_rand,y_rand,u_rand))-V(x_rand,y_rand,T(z_rand,u_rand))==V(T(z_rand,x_rand),y_rand,u_rand)-V(x_rand,T(z_rand.si(),y_rand),u_rand)
#print T(z,V(x,y,u))-V(x,y,T(z,u))==V(T(z,x),y,u)-V(x,T(z.si(),y),u)