needsPackage("LieAlgebraRepresentations");
debug LieAlgebraRepresentations;


matApprox = (A,B,e) -> (
    max apply(flatten entries(A-B), x->abs(x) ) < e
);




end


load "Proposition_A.4_statement.m2";


-- Example 1
g = simpleLieAlgebra("A",3);
lambda = {2,0,0};
Vlambda = irreducibleLieAlgebraModule(lambda,g);
N = dim Vlambda;
Llambda = GTrepresentationMatrices(Vlambda);
lambdastar = starInvolution(lambda,g);
Vlambdastar = irreducibleLieAlgebraModule(lambdastar,g);
Llambdastar = GTrepresentationMatrices(Vlambdastar);
Blambda = gtPatterns("A",dynkinToPartition("A",lambda));
Blambdastar = gtPatterns("A",dynkinToPartition("A",lambdastar));
sigma = apply(N, i -> first select(N, k -> Blambdastar_k==dualGTPattern(Blambda_i)));
Blambdastar= apply(sigma, j -> Blambdastar_j);
Llambdastar = apply(Llambdastar, M -> M_sigma^sigma);
MC = apply(Blambda, p -> molevCoefficient(gtp(p)));
MCstar = apply(#Blambdastar, i -> molevCoefficient(gtp(Blambdastar_i)));
signs = apply(#Blambdastar, i -> (-1)^(1+level((gtp(Blambdastar_i))#"weight","A",#lambdastar)));
Q1 = diagonalMatrix apply(#MC, i-> sqrt(1/(MC_i)));
Q2 = diagonalMatrix apply(#MCstar, i-> signs_i*sqrt(1/(MCstar_i)));
Q2plus = diagonalMatrix apply(#MCstar, i-> sqrt(1/(MCstar_i)));

-- Check the main statement
all(dim(g), j -> matApprox((inverse(Q2)*(Llambdastar_j)*Q2),-transpose(inverse(Q1)*(Llambda_j)*Q1),10^-12))


-- Check the statement about Bstarplus
LAB = lieAlgebraBasis(g);
LOI = LAB#"LoweringOperatorIndices";
ROI = LAB#"RaisingOperatorIndices";
m=g#"LieAlgebraRank";
all(m, j -> matApprox((inverse(Q2plus)*(Llambdastar_(LOI_j))*Q2plus),transpose(inverse(Q1)*(Llambda_(LOI_j))*Q1),10^-12))
all(m, j -> matApprox((inverse(Q2plus)*(Llambdastar_(ROI_j))*Q2plus),transpose(inverse(Q1)*(Llambda_(ROI_j))*Q1),10^-12))





-- Example 2
g = simpleLieAlgebra("A",3);
lambda = {3,2,1};
Vlambda = irreducibleLieAlgebraModule(lambda,g);
N = dim Vlambda;
Llambda = GTrepresentationMatrices(Vlambda);
lambdastar = starInvolution(lambda,g);
Vlambdastar = irreducibleLieAlgebraModule(lambdastar,g);
Llambdastar = GTrepresentationMatrices(Vlambdastar);
Blambda = gtPatterns("A",dynkinToPartition("A",lambda));
Blambdastar = gtPatterns("A",dynkinToPartition("A",lambdastar));
sigma = apply(N, i -> first select(N, k -> Blambdastar_k==dualGTPattern(Blambda_i)));
Blambdastar= apply(sigma, j -> Blambdastar_j);
Llambdastar = apply(Llambdastar, M -> M_sigma^sigma);
MC = apply(Blambda, p -> molevCoefficient(gtp(p)));
MCstar = apply(#Blambdastar, i -> molevCoefficient(gtp(Blambdastar_i)));
signs = apply(#Blambdastar, i -> (-1)^(1+level((gtp(Blambdastar_i))#"weight","A",#lambdastar)));
Q1 = diagonalMatrix apply(#MC, i-> sqrt(1/(MC_i)));
Q2 = diagonalMatrix apply(#MCstar, i-> signs_i*sqrt(1/(MCstar_i)));
Q2plus = diagonalMatrix apply(#MCstar, i-> sqrt(1/(MCstar_i)));

-- Check the main statement
all(dim(g), j -> matApprox((inverse(Q2)*(Llambdastar_j)*Q2),-transpose(inverse(Q1)*(Llambda_j)*Q1),10^-12))


-- Check the statement about Bstarplus
LAB = lieAlgebraBasis(g);
LOI = LAB#"LoweringOperatorIndices";
ROI = LAB#"RaisingOperatorIndices";
m=g#"LieAlgebraRank";
all(m, j -> matApprox((inverse(Q2plus)*(Llambdastar_(LOI_j))*Q2plus),transpose(inverse(Q1)*(Llambda_(LOI_j))*Q1),10^-12))
all(m, j -> matApprox((inverse(Q2plus)*(Llambdastar_(ROI_j))*Q2plus),transpose(inverse(Q1)*(Llambda_(ROI_j))*Q1),10^-12))




-- Example 3
g = simpleLieAlgebra("A",5);
lambda = {1,0,0,0,1};
Vlambda = irreducibleLieAlgebraModule(lambda,g);
N = dim Vlambda;
Llambda = GTrepresentationMatrices(Vlambda);
lambdastar = starInvolution(lambda,g);
Vlambdastar = irreducibleLieAlgebraModule(lambdastar,g);
Llambdastar = GTrepresentationMatrices(Vlambdastar);
Blambda = gtPatterns("A",dynkinToPartition("A",lambda));
Blambdastar = gtPatterns("A",dynkinToPartition("A",lambdastar));
sigma = apply(N, i -> first select(N, k -> Blambdastar_k==dualGTPattern(Blambda_i)));
Blambdastar= apply(sigma, j -> Blambdastar_j);
Llambdastar = apply(Llambdastar, M -> M_sigma^sigma);
MC = apply(Blambda, p -> molevCoefficient(gtp(p)));
MCstar = apply(#Blambdastar, i -> molevCoefficient(gtp(Blambdastar_i)));
signs = apply(#Blambdastar, i -> (-1)^(1+level((gtp(Blambdastar_i))#"weight","A",#lambdastar)));
Q1 = diagonalMatrix apply(#MC, i-> sqrt(1/(MC_i)));
Q2 = diagonalMatrix apply(#MCstar, i-> signs_i*sqrt(1/(MCstar_i)));
Q2plus = diagonalMatrix apply(#MCstar, i-> sqrt(1/(MCstar_i)));

-- Check the main statement
all(dim(g), j -> matApprox((inverse(Q2)*(Llambdastar_j)*Q2),-transpose(inverse(Q1)*(Llambda_j)*Q1),10^-12))


-- Check the statement about Bstarplus
LAB = lieAlgebraBasis(g);
LOI = LAB#"LoweringOperatorIndices";
ROI = LAB#"RaisingOperatorIndices";
m=g#"LieAlgebraRank";
all(m, j -> matApprox((inverse(Q2plus)*(Llambdastar_(LOI_j))*Q2plus),transpose(inverse(Q1)*(Llambda_(LOI_j))*Q1),10^-12))
all(m, j -> matApprox((inverse(Q2plus)*(Llambdastar_(ROI_j))*Q2plus),transpose(inverse(Q1)*(Llambda_(ROI_j))*Q1),10^-12))






-- Example 4
g = simpleLieAlgebra("A",6);
lambda = {2,0,0,0,1,0};
Vlambda = irreducibleLieAlgebraModule(lambda,g);
N = dim Vlambda;
Llambda = GTrepresentationMatrices(Vlambda);
lambdastar = starInvolution(lambda,g);
Vlambdastar = irreducibleLieAlgebraModule(lambdastar,g);
Llambdastar = GTrepresentationMatrices(Vlambdastar);
Blambda = gtPatterns("A",dynkinToPartition("A",lambda));
Blambdastar = gtPatterns("A",dynkinToPartition("A",lambdastar));
sigma = apply(N, i -> first select(N, k -> Blambdastar_k==dualGTPattern(Blambda_i)));
Blambdastar= apply(sigma, j -> Blambdastar_j);
Llambdastar = apply(Llambdastar, M -> M_sigma^sigma);
MC = apply(Blambda, p -> molevCoefficient(gtp(p)));
MCstar = apply(#Blambdastar, i -> molevCoefficient(gtp(Blambdastar_i)));
signs = apply(#Blambdastar, i -> (-1)^(1+level((gtp(Blambdastar_i))#"weight","A",#lambdastar)));
Q1 = diagonalMatrix apply(#MC, i-> sqrt(1/(MC_i)));
Q2 = diagonalMatrix apply(#MCstar, i-> signs_i*sqrt(1/(MCstar_i)));
Q2plus = diagonalMatrix apply(#MCstar, i-> sqrt(1/(MCstar_i)));

-- Check the main statement
all(dim(g), j -> matApprox((inverse(Q2)*(Llambdastar_j)*Q2),-transpose(inverse(Q1)*(Llambda_j)*Q1),10^-12))
-- true

-- Check the statement about Bstarplus
LAB = lieAlgebraBasis(g);
LOI = LAB#"LoweringOperatorIndices";
ROI = LAB#"RaisingOperatorIndices";
m=g#"LieAlgebraRank";
all(m, j -> matApprox((inverse(Q2plus)*(Llambdastar_(LOI_j))*Q2plus),transpose(inverse(Q1)*(Llambda_(LOI_j))*Q1),10^-12))
all(m, j -> matApprox((inverse(Q2plus)*(Llambdastar_(ROI_j))*Q2plus),transpose(inverse(Q1)*(Llambda_(ROI_j))*Q1),10^-12))
-- true