##### GENERAL ROOT SYSTEM PROCEDURES #####
my_rank := proc (R) local i, j, res;
res := 0;
for i from 2 to length(R) do
member(substring(R,i..i),[`0`,`1`,`2`,`3`,`4`,`5`,`6`,`7`,`8`,`9`],'j');
res := res * 10 + j - 1;
od:
res;
end proc:

my_cartan_matrix := proc (R) option remember; local Car, L, r, i, j;
L := substring (R, 1..1);
r := my_rank (R);
Car := matrix (r, r);
for i to r do for j to r do Car [i, j] := 0; od; od;
for i to r do Car [i, i] := 2; od;
for i from 2 to r do Car [i - 1, i] := -1; od;
for i to r - 1 do Car [i + 1, i] := -1; od;
if L = 'A' then
elif L = 'B' then Car [2, 1] := -2;
elif L = 'C' then Car [1, 2] := -2;
elif L = 'D' then Car [1, 2] := 0; Car [1, 3] := -1;
Car [2, 1] := 0; Car [3, 1] := -1;
elif L = 'E' then Car [1, 2] := 0; Car [1, 3] := -1;
Car [2, 1] := 0; Car [2, 3] := 0; Car [2, 4] := -1;
Car [3, 1] := -1; Car [3, 2] := 0; Car [4, 2] := -1;
elif R = F4 then Car [3, 2] := -2;
elif R = G2 then Car [2, 1] := -3;
else ERROR (cat(R,` is an unknown root system`)); fi;
Car;
end proc:

my_cartan_inverse := proc (R) option remember;
linalg:-inverse(my_cartan_matrix(R)); end proc:

pos_roots_number := proc (R) local L, r;
L:=substring(R,1..1); r:=my_rank(R);
if L='A' then return r*(r+1)/2; fi;
if (L='B') or (L='C') then return r^2; fi;
if (L='D') then return r*(r-1); fi;
if R=G2 then return 6; fi;
if R=F4 then return 24; fi;
if R=E6 then return 36; fi;
if R=E7 then return 63; fi;
if R=E8 then return 120; fi;
ERROR (cat(R,` is an unknown root system`));
end proc:

my_root_length := proc (R, k) local L;
L:=substring(R,1..1);
if (R=G2) and (k=2) then return 3; fi;
if (L='B') and (k>1) or (L='C') and (k=1) or (R=F4) and (k>2) then return 2; fi;
return 1;
end proc:

my_opp_inv := proc (R) local L, r, i;
L:=substring(R,1..1); r:=my_rank(R);
if L='A' then return [seq(r+1-i,i=1..r)]; fi;
if (L='D') and (r mod 2=1) then return [2,1,$3..r]; fi;
if R=E6 then return [6,2,5,4,3,1]; fi;
[$1..r];
end proc:

##### NUMERICAL REPRESENTATION AND REDUCTION MODULO W_P #####
my_root_coords := proc (R, p) local C, i, k;
C:=my_cartan_inverse(R);
[seq(add(coeff(p,omega[i])*C[i,k],i=1..my_rank(R)),k=1..my_rank(R))];
end proc:

reflect_weights := proc (R) option remember; local Car,i,j;
Car:=my_cartan_matrix(R);
[seq(omega[i]-add(Car[i,j]*omega[j],j=1..my_rank(R)),i=1..my_rank(R))];
end proc:

numer_repr := proc (Car, r, w) local i, t;
t:=add(omega[i],i=1..r);
for i to nops(w) do t:=subs(omega[w[-i]]=Car[w[-i]],t) od; t;
end proc:

weyl_length := proc (Car, r, w) local res, i, j, t;
t:=add(omega[i],i=1..r); res := 0;
for i to nops(w) do if coeff(t,omega[w[i]])>0 then res:=res+1 else
res:=res-1; fi;
t := subs(omega[w[i]]=Car[w[i]],t); od; res;
end proc:

is_minimal_repr := proc (J, Car, r, w) local j, t;
t:=numer_repr(Car,r,ListTools:-Reverse(w));
for j in J do if coeff(t,omega[j])<0 then return false; fi; od;
return true; end proc:

reduce_left := proc (J, Car, r, q) local i, flag, t, res;
t:=q; res:=[]; flag:=true;
while flag do
flag:=false;
for i in J do if coeff(t,omega[i])<0 then res:=[op(res),i];
t:=subs(omega[i]=Car[i],t); flag:=true; break; fi; od; od; [res,t];
end proc:

my_right_coset := proc (J, Car, r, w) local t;
t:=reduce_left(J,Car,r,numer_repr (Car, r, w));
[t[1],reduce_left([$1..r], Car, r, t[2])[1]];
end proc:

right_coset := proc (J, R, w);
my_right_coset(J,reflect_weights(R),my_rank(R),w);
end proc:

my_reduce := proc (Car, r, w); my_right_coset([],Car,r,w)[2]; end proc:

my_left_coset := proc (J, Car, r, w) local t;
t:=my_right_coset (J, Car, r, ListTools:-Reverse(w));
[my_reduce(Car,r,ListTools:-Reverse(t[2])),ListTools:-Reverse(t[1])];
end proc:

left_coset := proc (J, R, w)
my_left_coset(J,reflect_weights(R),my_rank(R),w);
end proc:

my_double_coset := proc (JL, JR, Car, r, w) local t;
t:=my_left_coset (JR, Car, r, w);[op(my_right_coset(JL, Car, r, t[1])),t[2]];
end proc:

double_coset := proc (JL, JR, R, w)
my_double_coset(JL,JR,reflect_weights(R),my_rank(R),w);
end proc:

my_longest_elt := proc (R) option remember; local r,i;
r := my_rank(R);
reduce_left([$1..r],reflect_weights(R),r,-add(omega[i],i=1..r))[1];
end proc:

my_longest_elt2 := proc (J, R) option remember;
my_right_coset(J,reflect_weights(R),my_rank(R),my_longest_elt(R))[1];
end proc:

my_pos_roots := proc (R) option remember;
local Car, r, w, j, i, t, res1, res2, res3, res4;
Car := reflect_weights (R); r := my_rank (R);
w := my_longest_elt (R);
res1 := []; res2 := []; res3 := []; res4 := [];
for i to nops(w) do
t := omega[w[i]]-Car[w[i]];
for j from i-1 to 1 by -1 do
t := subs (omega[w[j]]=Car[w[j]], t);
od:
res1 := [op(res1), t];
res2 := [op(res2),
numer_repr(Car, r, [op(w[1..i]),op(ListTools:-Reverse(w[1..i-1]))])];
res3 := [op(res3), my_root_coords (R, t)];
res4 := [op(res4), my_root_length (R, w[i])];
od:
[res1, res2, res3, res4];
end proc:

##### CHOW GENERATORS #####
list_of_generators := proc (J, R, k) option remember; local Car, X, L,
N, S, T, r, i, j, m, num, flag; r:=my_rank(R);
if k=0 then return [[[]],[add(omega[i],i=1..r)-add(omega[i],i in J)],[],0]; fi;
Car:=reflect_weights(R);  X:=list_of_generators(J,R,k-1); L:=[]; N :=[]; T:=[];
if k>1 then S:={op(list_of_generators(J,R,k-2)[2])}; else S:={}; fi;
for j to nops(X[1]) do
flag:=true;
for i to r do
num := subs(omega[i]=Car[i],X[2][j]);
if (num<>X[2][j]) and not (num in S) then
m:=ListTools:-Search(num,N);
if m=0 then
N:=[op(N),num]; L:=[op(L),my_reduce(Car,r,[i,op(X[1][j])])];
if flag then T:=[op(T),[i,nops(N)]]; flag:=false; fi;
else
if flag then T:=[op(T),[i,m]]; flag:=false; fi;
fi;
fi;
od;
od;
[L,N,T,X[4]+nops(X[1])];
end proc:

chow_generators := proc (J,R,k);
list_of_generators(J,R,k)[1];
end proc:

chow_dual := proc (J, R, w) local w0;
w0 := my_longest_elt(R);
my_reduce(reflect_weights(R),my_rank(R),[op(w0),op(w),op(my_longest_elt2(J,R))]);
end proc:

chow_dim := proc (J, R) local t;
pos_roots_number(R)-add(pos_roots_number(t[2]),t in levi_part2(J,R));
end proc:

##### MULTIPLICATION ALGORITHM #####
my_indets := proc (p) local res, j, s; res:=[];
for s in indets(p) do if op(0,s)=Z then
res:=[op(res),[seq(op(j,s),j=1..nops(s))]]; fi; od; res;
end proc:

my_weight_eval := proc (R) option remember; local C, i, k;
C:=my_cartan_inverse(R);
[seq(add(C[i,k],k=1..my_rank(R)),i=1..my_rank(R))];
end proc:

generate_common_table := proc (J, R) option remember; local r, C, Car, L, N, n,
B, B1, B2, B3, B4, B5, X, S, T, k, i, j, s, t, m, ch, x, y, l, inv;
r:=my_rank(R); Car:=reflect_weights(R); C:=my_weight_eval(R);
inv := my_opp_inv(R); L:=[]; n:=0;
while true do
X:=list_of_generators(J, R, n);
L:=[op(L), X];
if nops(X[1])=0 then break; fi;
n:=n+1;
od;
n:=n-1;
N:=X[4];
B:=array(1..N);
for k from 1 to n+1 do
for i to nops(L[k][1]) do
B[L[k][4]+i]:=L[k][1][i];
od;
od;
B1:=array(1..N);
B1[N]:=[seq(Omega[j]=omega[j],j=1..r)];
for k from n to 1 by -1 do
for i to nops(L[k][1]) do
m:=inv[L[k+1][3][i][1]];
B1[L[k][4]+i]:=subs(omega[m]=Car[m],B1[L[k+1][4]+L[k+1][3][i][2]]);
od;
od;
B2:=array(1..N);
for i to N do B2[i]:=[subs(seq(omega[j]=C[j],j=1..r),B1[i])]; od;
S:=unipotent_rad(J,R);
ch:=product(S[j],j=1..nops(S));
ch:=subs(seq(omega[j]=Omega[j],j=1..r),ch);
B3:=array(1..N);
for i to N do B3[i]:=subs(op(B2[i]),ch); od;
B4:=array(1..N);
B4[N]:=S;
T:=my_pos_roots(R)[1];
for k from n to 1 by -1 do
for i to nops(L[k][1]) do
B4[L[k][4]+i]:={};
m:=inv[L[k+1][3][i][1]];
S:=B4[L[k+1][4]+L[k+1][3][i][2]];
for s in T do
y:=subs(omega[m]=Car[m],s);
x:=my_root_coords(R,y);
for l to r do if x[l]<>0 then break; fi; od;
if (x[l]>0) and (y in S) or (x[l]<0) and not ((-y) in S) then
B4[L[k][4]+i]:=B4[L[k][4]+i] union {s}; fi;
od;
od;
od;
B5:=array(1..N);
for i to N do
B5[i]:=subs(seq(omega[j]=C[j],j=1..r),product(B4[i][t],t=1..nops(B4[i])));
od;
[B,B1,B2,B3,B4,B5,N,n];
end proc:

degree_map := proc (J, R, p) local L, i;
L := generate_common_table (J, R);
add(subs(op(L[3][i]),p)/L[4][i],i=1..L[7]);
end proc:

find_expansion := proc (J, R, p) local A, x, q, L, X, i, j, k, d, res, res2, flag;
d:=my_degree(my_rank(R),p);
L := generate_common_table (J, R);
X := list_of_generators (J, R, d);
A:=array(1..X[4]+nops(X[1]));
for i to X[4]+nops(X[1]) do A[i]:=subs(op(L[3][i]),p); od;
res:=0; k:=1;
while (k<=X[4]) do
flag:=true;
for i from k to X[4] do if A[i]<>0 then flag:=false; break; fi; od;
if flag then break; fi;
k:=i+1;
x:=get_generator(J, R, i);
q:=A[i]/subs(op(L[3][i]),x);
res:=res+q*Z[op(L[1][i])];
for j from i+1 to X[4]+nops(X[1]) do A[j]:=A[j]-subs(op(L[3][j]),x)*q; od;
od;
res2:=0;
for i from X[4]+1 to X[4]+nops(X[1]) do
res2:=res2+A[i]/L[6][i]*Z[op(L[1][i])]
od;
[res+res2,res2];
end proc:

find_expansion2 := proc (J, R, p) local A, x, q, L, X, i, j, k, d, m, res, flag;
d:=my_degree(my_rank(R),p);
L := generate_common_table (J, R);
X := list_of_generators (J, R, d);
A:=array(1..X[4]+nops(X[1]));
for i to X[4]+nops(X[1]) do A[i]:=expand(subs(op(L[2][i]),p)); od;
res:=0; k:=1;
while (k<=X[4]) do
flag:=true;
for i from k to X[4] do if A[i]<>0 then flag:=false; break; fi; od;
if flag then break; fi;
k:=i+1;
x:=get_generator2(J, R, i);
q:=simplify(A[i]/subs(op(L[2][i]),x));
res:=res+q*Z[op(L[1][i])];
for j from i+1 to X[4]+nops(X[1]) do A[j]:=A[j]-expand(subs(op(L[2][j]),x)*q); od;
od;
for i from X[4]+1 to X[4]+nops(X[1]) do
res:=res+simplify(A[i]/product(L[5][i][m],m=1..nops(L[5][i])))*Z[op(L[1][i])]
od;
res;
end proc:

get_generator := proc (J, R, i) option remember; local L, p, j, k;
L := generate_common_table (J, R);
p := subs(seq(omega[j]=Omega[j],j=1..my_rank(R)),
c_func_inv(J, R, Z[op(L[1][i])]));
p-subs(seq(Z[op(L[1][k])]=get_generator(J,R,k),
k=1..(list_of_generators (J, R, nops(L[1][i])))[4]),
(find_expansion(J,R,p))[1]-Z[op(L[1][i])]
);
end proc:

get_generator2 := proc (J, R, i) option remember; local L, p, j, k;
L := generate_common_table (J, R);
p := subs(seq(omega[j]=Omega[j],j=1..my_rank(R)),
c_func_inv(J, R, Z[op(L[1][i])]));
expand(p-subs(seq(Z[op(L[1][k])]=get_generator2(J,R,k),
k=1..(list_of_generators (J, R, nops(L[1][i])))[4]),
find_expansion2(J,R,p)-Z[op(L[1][i])]
));
end proc:

steenrod1 := proc (J, R, u, p, k) option remember; local X, X2, L, A, x,
r, t, i, j, m, q, res;
r:=my_rank(R);
L:=generate_common_table(J,R);
X:=list_of_generators(J,R,nops(u));
X2:=list_of_generators(J,R,nops(u)+k*(p-1));
m:=ListTools:-Search(u,X[1]);
if m=0 then ERROR(`A wrong generator`); fi;
m:=X[4]+m;
A:=array(m..X2[4]+nops(X2[1]));
x:=get_generator2(J,R,m);
for i from m to X2[4]+nops(X2[1]) do
A[i]:=expand(coeff(subs(seq(omega[s]=omega[s]+t*omega[s]^p,s=1..r),
expand(subs(op(L[2][i]),x))),
t^k)) mod p;
od;
for i from m to X2[4] do if A[i]<>0 then
x:=get_generator2(J,R,i);
if not (Divide(A[i],expand(subs(op(L[2][i]),x)),'q') mod p) then ERROR(`Internal error in Steenrod`); fi;
for j from i+1 to X2[4]+nops(X2[1]) do A[j]:=(A[j]-expand(subs(op(L[2][j]),x)*q)) mod p; od;
fi; od;
res:=0;
for i from X2[4]+1 to X2[4]+nops(X2[1]) do
if not Divide(A[i],product(L[5][i][s],s=1..nops(L[5][i])),'q') mod p
then ERROR(`Internal error in Steenrod`); fi;
res:=res+q*Z[op(L[1][i])];
od;
res;
end proc:

my_degree := proc (r, p) local i; degree(p,{seq(omega[i],i=1..r),seq(Omega[i],i=1..r)});
end proc:

c_func := proc (J, R, p) local i;
(find_expansion(J,R,subs(seq(omega[i]=Omega[i],i=1..my_rank(R)),p)))[2];
end proc:

steenrod := proc (J, R, p, k, q) local L,l; L:=my_indets(q);
if nops(L)=0 then if k=0 then return q else return 0; fi; fi;
subs(seq(Z[op(l)]=steenrod1(J,R,l,p,k),l in L),q)
end proc:

##### PIERI FORMULA #####

piericoef := proc (R, k, w, u) option remember;
local L, i;
L:=my_pos_roots(R);
i:=ListTools:-Search(numer_repr(reflect_weights(R),my_rank(R),
[op(ListTools:-Reverse(u)),op(w)]),L[2]);
if i=0 then return 0; fi;
L[3][i][k]*my_root_length(R,k)/L[4][i];
end proc:

pieri1 := proc (J, R, k, u) option remember; local w;
add(piericoef (R, k, w, u) * Z[op(w)], w in chow_generators(J, R, nops(u)+1));
end proc:

pieri2 := proc (J, R, k, p) local L,l; L:=my_indets(p);
if nops(L)=0 then return Z[k]*p; fi;
subs(seq(Z[op(l)]=pieri1(J, R, k, l), l in L), p)
end proc:

pieri := proc (J, R, k, p, n) local i,t;
if nargs=4 then return pieri2 (J, R, k, p); fi;
t:=p; for i to n do t:=pieri2(J, R, k, t); od; t;
end proc:

##### LEVI PART COMPUTATION #####
root_sys_name := proc (type, rk); if rk=1 then return A1; fi;
if substring(type,1..1)='F' then return cat(A,rk); fi;
return cat(type,rk);
end proc:

levi_part := proc (J, Ord, R) local L, i, k, r;
if Ord=[] then return []; fi;
k:=0; for i to nops(Ord) do if not (Ord[i] in {op(J)}) then k:=i;
break; fi; od;
if k=0 then return [[Ord, R]]; fi;
r:=my_rank(R); L:=substring(R,1..1); 
if (L='A') or (L='B') or (L='C') or (R=G2) or ((L='D') and (k>4)) or
((R=F4) and ((k=2) or (k=3))) or ((L='E') and (k>6)) then 
return [op(levi_part(J,Ord[1..k-1],root_sys_name(L,k-1))),
op(levi_part(J,Ord[k+1..r],root_sys_name(A,r-k)))]; fi;
if L='D' then if k=4 then return
[op(levi_part(J,Ord[5..r],root_sys_name(A,r-4))),
op(levi_part(J,[Ord[1],Ord[3],Ord[2]],A3))]; fi;
if k=3 then return [op(levi_part(J,Ord[4..r],root_sys_name(A,r-3))),
op(levi_part(J,[Ord[1]],A1)),op(levi_part(J,[Ord[2]],A1))]; fi;
return levi_part(J,[Ord[3-k],op(Ord[3..r])],root_sys_name(A,r-1)); fi;
if R=F4 then if k=1 then return levi_part(J,Ord[2..4],B3); fi;
if k=4 then return levi_part(J,[Ord[3],Ord[2],Ord[1]],C3); fi;
fi;
if L<>'E' then ERROR (cat(R,` is an unknown root system`)); fi;
if k=1 then return levi_part(J,Ord[2..r],root_sys_name(D,r-1)); fi;
if k=2 then return
levi_part(J,[Ord[1],op(Ord[3..r])],root_sys_name(A,r-1)); fi;
if k=3 then return [op(levi_part(J,[Ord[1]],A1)),
op(levi_part(J,[Ord[2],op(Ord[4..r])],root_sys_name(A,r-2)))]; fi;
if k=4 then return
[op(levi_part(J,[Ord[1],Ord[3]],A2)),op(levi_part(J,[Ord[2]],A1)),
op(levi_part(J,Ord[5..r],root_sys_name(A,r-4)))]; fi;
if k=5 then return [op(levi_part(J,[Ord[1],Ord[3],Ord[4],Ord[2]],A4)),
op(levi_part(J,Ord[6..r],root_sys_name(A,r-5)))]; fi;
if k=6 then return [op(levi_part(J,[Ord[2],Ord[5],Ord[4],Ord[3],Ord[1]],D5)),
op(levi_part(J,Ord[7..r],root_sys_name(A,r-6)))]; fi;
end proc:

levi_part2 := proc (J, R) option remember;
levi_part(J,[$1..my_rank(R)],R); end proc:

##### FUNDAMENTAL INVARIANTS #####
deg_fundam_invariant := proc (R, k) option remember; local L;
L:=substring(R,1..1);
if L='A' then return k+1; fi;
if (L='B') or (L='C') then return 2*k; fi;
if (L='D') then if k<my_rank(R) then return 2*k; fi; return k; fi;
if R=G2 then return [2,6][k]; fi;
if R=F4 then return [2,6,8,12][k]; fi;
if R=E6 then return [2,5,6,8,9,12][k]; fi;
if R=E7 then return [2,6,8,10,12,14,18][k]; fi;
if R=E8 then return [2,8,12,14,18,20,24,30][k]; fi;
end proc:

fundam_invariant := proc (R, Ord, R2, k) option remember;
local Car, res, r, L, x, deg, hyp, cond, i, j, m;
r:=my_rank(R); L:=substring(R,1..1); deg:=deg_fundam_invariant(R,k);
cond:={};
if L='A' then res:=add(x[j]^deg,j=1..r+1);
hyp:=[seq(x[i+1]-x[i],i=1..r)]; cond:={add(x[i],i=1..r+1)=0}; fi;
if L='B' then res:=add(x[j]^deg,j=1..r);
hyp:=[x[1],seq(x[i+1]-x[i],i=1..r)]; fi;
if L='C' then res:=add(x[j]^deg,j=1..r);
hyp:=[2*x[1],seq(x[i+1]-x[i],i=1..r)]; fi;
if L='D' then if k<r then res:=add(x[j]^deg,j=1..r); else
res:=mul(x[j],j=1..r); fi;
hyp:=[x[1]+x[2],seq(x[i+1]-x[i],i=1..r)]; fi;
if R=G2 then if k=1 then res:=x[1]^2+x[2]^2; else
res:=add((x[1]*cos(Pi/3*j)+x[2]*sin(Pi/3*j))^6,j=1..6); fi;
hyp:=[x[1],-3*x[1]/2+sqrt(3)/2*x[2]]; fi;
if R=F4 then res:=0; for j to 4 do for i to j-1 do
res:=res+(x[i]+x[j])^deg+(x[i]-x[j])^deg; od; od;
hyp:=[(x[1]-x[2]-x[3]-x[4])/2,x[4],x[3]-x[4],x[2]-x[3]]; fi;
if R=E6 then res:=add((-x[m]+x[8])^deg+(-x[m]-x[8])^deg,m=1..6);
for j to 6 do for i to j-1 do res:=res+(x[i]+x[j])^deg; od; od;
hyp:=[x[2]-x[1],(x[1]+x[2]+x[3]+x[8]-x[4]-x[5]-x[6]-x[7])/2,x[3]-x[2],
x[4]-x[3],x[5]-x[4],x[6]-x[5]];
cond:={x[1]+x[2]+x[3]+x[4]+x[5]+x[6]=0,x[7]+x[8]=0}; fi;
if R=E7 then res:=0; for j to 8 do for i to j-1 do
res:=res+(x[i]+x[j])^deg; od; od;
hyp:=[x[2]-x[1],(x[1]+x[2]+x[3]+x[8]-x[4]-x[5]-x[6]-x[7])/2,x[3]-x[2],
x[4]-x[3],x[5]-x[4],x[6]-x[5],x[7]-x[6]];
cond:={x[1]+x[2]+x[3]+x[4]+x[5]+x[6]+x[7]+x[8]=0}; fi;
if R=E8 then res:=0; for j to 9 do for i to j-1 do
res:=res+(x[i]-x[j])^deg; od; od;
for m to 9 do for j to m-1 do for i to j-1 do
res:=res+(x[i]+x[j]+x[m])^deg; od; od; od;
hyp:=[x[2]-x[1],(2*(x[1]+x[2]+x[3])-x[4]-x[5]-x[6]-x[7]-x[8]-x[9])/3,
x[3]-x[2],x[4]-x[3],x[5]-x[4],x[6]-x[5],x[7]-x[6],x[8]-x[7]];
cond:={x[1]+x[2]+x[3]+x[4]+x[5]+x[6]+x[7]+x[8]+x[9]=0};
fi;
Car:=my_cartan_matrix(R2);
subs(op(solve({seq(hyp[i]=add(Car[Ord[i],t]*omega[t],t=1..my_rank(R2)),i=1..r)}
union cond,{seq(x[i],i=1..r+2)})),res);
end proc:

##### CHOW RING MULTIPLICATION #####
generators := proc (J, R, n) option remember;
local G, G1, L, L1, S, LP, B, M, x, i, j, k, l, rk, d, N;
if n=0 then return [[1],[[1]],[[1]]]; fi;
G1:=chow_generators(J,R,n-1);
G:=chow_generators(J,R,n);
LP:=levi_part2(J,R);
L:=generators(J,R,n-1); B:=[];
M:=matrix(nops(G),nops(G)); rk:=0;
for i to nops(G) do for j to nops(G) do M[i,j]:=0; od; od;
S:={$1..my_rank(R)} minus {op(J)};
for i to nops(S) do for j to nops(G1) do
x:=pieri(J,R,S[i],add(L[2][j,k]*Z[op(G1[k])],k=1..nops(G1)));
for k to nops(G) do M[rk+1,k]:=coeff(x,Z[op(G[k])]) od;
if linalg:-rank(M)>rk then rk:=rk+1; B:=[op(B),L[1][j]*omega[S[i]]];
fi;
if rk=nops(G) then break; fi;
od; if rk=nops(G) then break; fi; od;
if rk<nops(G) then
for i to nops(LP) do
for j to my_rank(LP[i][2]) do
d:=deg_fundam_invariant(LP[i][2],j);
if d<=n then
L1:=generators(J,R,n-d);
for l to nops(L1[1]) do if my_degree(my_rank(R),L1[1][l])=0 then
x:=chow_expand(J,R,P[i,j]*L1[1][l]);
for k to nops(G) do M[rk+1,k]:=coeff(x,Z[op(G[k])]) od;
if linalg:-rank(M)>rk then rk:=rk+1; B:=[op(B),P[i,j]*L1[1][l]];
if rk=nops(G) then return [B,M,linalg:-inverse(M)]; fi; fi; fi; od;
fi; od; od; fi;
[B,M,linalg:-inverse(M)];
end proc:

add_to_mult:=proc (J,R,p) local v, L, SUBS, j, k;
SUBS:=[];
for v in my_indets(p) do
L:=generators(J,R,nops(v));
k:=ListTools:-Search(my_reduce(reflect_weights(R),my_rank(R),v),
chow_generators(J,R,nops(v)));
if k=0 then ERROR(`A wrong generator`); fi;
SUBS:=[op(SUBS),Z[op(v)]=add(L[3][k,j]*L[1][j],j=1..nops(L[1]))];
od;
subs(op(SUBS),p);
end proc:

c_func_inv:=proc (J,R,p) local V,LP,S,i,j;
LP:=levi_part2(J,R); S:={};
for i to nops(LP) do
for j to my_rank(LP[i][2]) do
S:=S union {P[i,j]=fundam_invariant(LP[i][2],LP[i][1],R,j)};
od; od;
subs(op(S),add_to_mult(J,R,p));
end proc:

chow_expand:=proc (J,R,p);
c_func(J,R,c_func_inv(J,R,p));
end proc:

##### CHERN CLASSES #####
unipotent_rad := proc (J, R) option remember; local j, S, i, res;
res:={}; S:=my_pos_roots(R);
for i to nops (S[1]) do
for j in {$ 1..my_rank(R)} minus {op(J)} do
if S[3][i][j]>0 then res:=res union {S[1][i]}; break; fi; od; od; res;
end proc:

relative_root := proc (J, Pat, R) option remember; local j, S, i, flag, res;
res:={}; S:=my_pos_roots(R);
for i to nops (S[1]) do flag := true;
for j in {$ 1..my_rank(R)} minus {op(J)} do
if S[3][i][j]<>Pat[j] then flag := false; break; fi; od;
if flag then res:=res union {S[1][i]}; fi; od; res;
end proc:

sym_funct := proc (T, R, k, p) local x, y, j;
x:=array(0..k);
x[0]:=1; for j to k do x[j]:=0; od;
for y in T do for j from k by -1 to 1 do
x[j]:=x[j]+expand(x[j-1]*y); if nargs>3 then x[j]:=x[j] mod p; fi; od;
od;
x[k];
end proc:

chern_class := proc (J, R, k, p) option remember; local res;
res:=c_func(J,R,sym_funct(unipotent_rad(J,R),args[2..nargs]));
if nargs>3 then res := res mod p; fi;
res;
end proc:

weyl_orbit := proc (J, R, lambda) local Car, T, S, t, i, flag;
Car := reflect_weights(R); T := {lambda}; flag := false;
do
S := T;
for t in T do for i in J do S := S union {subs(omega[i]=Car[i], t)}; od; od;
if S = T then break; fi;
T := S;
od;
T;
end proc:

steinberg_weight := proc (R, w) option rememeber; local Car, t, x, i, j, k, flag, r, res;
r := my_rank (R); Car := reflect_weights (R); res := 0;
for i to r do t := omega[i]-Car[i]; for j to nops(w) do t := subs (omega[w[j]] = Car[w[j]], t); od;
x := my_root_coords (R, t); flag := true; for k to r do if x[k] < 0 then flag := false; break; fi; od;
if flag then t := omega[i]; for j to nops(w) do t := subs (omega[w[j]] = Car[w[j]], t); od;
res := res + t;
fi;
od;
res;
end proc:

##### DOUBLE PARABOLIC FILTRATION #####

prodbase := proc (J,K,R,w) local flag, i, j, x, y, res, K1, Car;
Car := reflect_weights (R);
K1:={$1..my_rank(R)} minus {op(K)}; res:=[];
for i in J do
y:=omega[i]-Car[i]; for j to nops(w) do y:=subs(omega[w[j]]=Car[w[j]],y); od; x:=my_root_coords(R,y);
flag:=true;
for j in K1 do if x[j]>0 then flag:=false; break; fi; od;
if flag then res:=[op(res), i]; fi;
od;
res;
end proc:

my_dcoset_reps := proc (J,K,R) option remember; local i, w, S;
S := {};
for i from 0 to chow_dim(K,R) do
for w in chow_generators(K,R,i) do
S := S union {right_coset(J,R,w)[2]};
od; od;
S;
end proc:

prodbases := proc (J,K,R) option remember; local w, res;
res := [];
for w in my_dcoset_reps(J,K,R) do res := [op(res), [prodbase(J,K,R,w),w]]; od;
res;
end proc:

prodimagecor := proc (J,K,R,w,u,v) local M;
M := prodbase(J,K,R,w);
my_pushforward(M,K,R,w,chow_expand(M,R,Z[op(u)]*Z[op(v)]));
end proc:

prodimage := proc(J,K,R,w,u) local v,i,t,res;
res := [];
for i from 0 to chow_dim(J,R) do
for v in chow_generators(J,R,i) do
t := prodimagecor (J,K,R,w,u,v);
if t<>0 then res := [op(res), [Z[op(chow_dual(J,R,v))],t]]; fi;
od;
od;
res;
end proc:

my_pushforward1 := proc(J,K,R,w,v) local x, y;
x := chow_dual(J,R,v);
y := left_coset(K,R,[op(x),op(w)])[1];
if nops(y)=nops(x)+nops(w) then Z[op(chow_dual(K,R,y))] else 0; fi;
end proc:

my_pushforward := proc (J, K, R, w, p) local L,l; L:=my_indets(p);
if nops(L)=0 then return my_pushforward(J,K,R,w,Z[])*p; fi;
subs(seq(Z[op(l)]=my_pushforward1(J, K, R, w, l), l in L), p)
end proc:

pushforward := proc(J,K,R,p);
my_pushforward(J,K,R,[],p);
end proc:

##### LIE ALGEBRA #####

number_of_root := proc (R, p) local L, t;
L:=my_pos_roots(R)[1];
t:=ListTools:-Search(p,L);
if t>0 then return t; fi;
-ListTools:-Search(-p,L);
end proc:

number_to_root := proc (R, i) local L;
L:=my_pos_roots(R)[1];
if i>0 then L[i] else -L[-i]; fi;
end proc:

sum_of_roots := proc (R, i, j) option remember;
number_of_root(R,number_to_root(R,i)+number_to_root(R,j));
end proc:

structure_sign := proc (R, i, j) option remember; local Car, r, q1, q2, k, l, t;
r:=my_rank(R); Car:=my_cartan_matrix(R);
q1:=my_root_coords(R,number_to_root(R,i));
q2:=my_root_coords(R,number_to_root(R,j));
t:=0;
for k to r do for l to k-1 do t:=t+q1[k]*q2[l]*Car[k,l]; od;
t:=t+q1[k]*q2[k]; od;
if t mod 2=0 then t:=1 else t:=-1; fi;
signum(i)*signum(j)*signum(sum_of_roots(R,i,j))*t;
end proc:

lie_bracket := proc (R, p1, p2) local N, L, r, res, i, j, k, p, q, t;
N:=pos_roots_number(R); r:=my_rank(R);
res:=0;
for i from -N to N do
if i<>0 then
t:=coeff(p1,e[i]);
if t<>0 then
p:=number_to_root(R,i); q:=my_root_coords(R,p);
for k to r do
res:=res-t*coeff(p,omega[k])*coeff(p2,h[k])*e[i];
od;
for k to r do
res:=res+t*coeff(p2,e[-i])*q[k]*h[k];
od;
for j from -N to N do
if (j<>0) and (j<>-i) then
res:=res+t*coeff(p2,e[j])*structure_sign(R,i,j)*e[sum_of_roots(R,i,j)];
fi;
od;
fi;
fi;
od;
for j from -N to N do
if j<>0 then
t:=coeff(p2,e[j]);
if t<>0 then
p:=number_to_root(R,j);
for k to r do
res:=res+t*coeff(p,omega[k])*coeff(p1,h[k])*e[j];
od;
fi;
fi;
od;
res;
end proc:

lprint(`Chow ring package v. 3.0 loaded`):