/* set of differential operators for Stiefel manifold */ import("yang.rr")$ def dsf(P,R) { S=[]; V=[]; DV=[]; for (I=1; I<=P; I++) { for (J=1; J<=R; J++) { V = cons(util_v(x,[I,J]),V); DV = cons(util_v(dx,[I,J]),DV); } } /* type 1 */ for (I=1; I<=R ; I++) { for (J=I; J<=R ; J++) { if (I==J) L = -1; else L = 0; for (K=1; K<=P; K++) { L += util_v(dx,[K,I])*util_v(dx,[K,J]); } S = cons(L,S); } } /* type 2 */ for (I=1; I<=R ; I++) { for (J=I+1; J<=R ; J++) { L=0; for (K=1; K<=P; K++) { L += -util_v(x,[K,J])*util_v(dx,[K,I]) +util_v(x,[K,I])*util_v(dx,[K,J]); } S = cons(L,S); } } for (I=1; I<=P ; I++) { for (J=I+1; J<=P ; J++) { L=0; for (K=1; K<=R; K++) { L += -util_v(x,[J,K])*util_v(dx,[I,K]) +util_v(x,[I,K])*util_v(dx,[J,K]); } S = cons(L,S); } } return [reverse(S),V,DV]; } def dblock(N) { M=newmat(2+N*2,2*N); for (I=N; I<2*N; I++) M[0][I]=1; for (I=2*N-1; I>=N; I--) M[1+(2*N-1)-I][I]= -1; for (I=0; I=0; I--) M[N+1+1+(N-1)-I][I] = -1; return M; } /* cf. misc-2009/10/fish/sei13.rr */ def bynd() { V=reverse([x_1_1,x_1_2,x_2_1,x_2_2,x_3_1,x_3_2]); Dv=reverse([dx_1_1,dx_1_2,dx_2_1,dx_2_2,dx_3_1,dx_3_2]); Ls=dsf(3,2)[0]; dp_gr_print(1); G=nd_weyl_gr(Ls,append(V,Dv),0,dblock(length(V))); bsave(G,"tmp-dsf32.ab"); return 0; } def mypf() { /* P=3,R=2 */ P=dsf(3,2); V=P[1]; Dv=P[2]; Ls = bload("tmp-dsf32-reduced.ab"); /* Ls = bload("tmp-dsf32.ab"); */ yang.define_ring(V); L=yang.util_pd_to_euler(Ls,V); L=map(nm,L); L=map(dp_ptod,L,Dv); yang.verbose(); /* G=yang.buchberger(L); */ G=L; print(yang.stdmon(G)); S0=yang.constant(1); S1=yang.operator(x_1_1); S2=yang.operator(x_1_2); S3=S1*S1; Base=[S0,S1,S2,S3]; Pf = yang.pfaffian(Base,G); bsave(Pf,"tmp-pf-sf32.ab"); return Pf; } /* 余分な元の消去 */ def mysupp1(P) { Dv=reverse([dx_1_1,dx_1_2,dx_2_1,dx_2_2,dx_3_1,dx_3_2]); P = dp_ptod(P,Dv); Q=0; while (P!=0) { Q += dp_ht(P); P = dp_rest(P); } return(Q); } def mypolysort(F) { return poly_sort(F | v=reverse([d11,d12,d21,d22,d31,d32])); } def mysupp() { /* G=bload("tmp-dsf32.ab"); */ G=bload("tmp-dsf32-reduced.ab"); Dvc=reverse([d11,d12,d21,d22,d31,d32]); Rule=[[x_1_1,3],[x_1_2,2], [x_2_1,1],[x_2_2,1/5], [x_3_1,1/6],[x_3_2,-1/3]]; G0=base_replace(G,Rule); G0=map(mysupp1,G0); G0=map(dp_dtop,G0,Dvc); G1=map(mypolysort,G0); for (I=0; I