version 1.1, 2009/10/09 07:15:44 |
version 1.9, 2009/11/12 01:39:54 |
Line 28 localf weyl_minipolym, weyl_minipoly, weyl_nf, weyl_nf |
|
Line 28 localf weyl_minipolym, weyl_minipoly, weyl_nf, weyl_nf |
|
localf weyl_nf_quo, weyl_nf_mod, b_subst, v_factorial, w_tdeg$ |
localf weyl_nf_quo, weyl_nf_mod, b_subst, v_factorial, w_tdeg$ |
localf replace_vars_f, replace_vars_v, replace_var$ |
localf replace_vars_f, replace_vars_v, replace_var$ |
localf action_on_gfs, action_on_gfs_1$ |
localf action_on_gfs, action_on_gfs_1$ |
|
localf nd_gb_candidate$ |
|
localf in_gb_oaku$ |
|
|
/* stratification */ |
/* stratification */ |
|
|
|
|
|
|
def in_ww(F) |
def in_ww(F) |
{ |
{ |
|
F = ptozp(F); |
V = vars(F); |
V = vars(F); |
N = length(V); |
N = length(V); |
D = newvect(N); |
D = newvect(N); |
Line 203 def in_ww_main(F,VW1,VW2) |
|
Line 206 def in_ww_main(F,VW1,VW2) |
|
VDVH = append(VDV,[TMP_H]); |
VDVH = append(VDV,[TMP_H]); |
FH = map(dp_dtop,map(dp_homo,map(dp_ptod,F,VDV)),VDVH); |
FH = map(dp_dtop,map(dp_homo,map(dp_ptod,F,VDV)),VDVH); |
|
|
/* compute a groebner basis of FH w.r.t. MWH */ |
/* |
|
* FH is a GB w.r.t. any term order s.t. LT(FH)=[t,dx1,...,dxn] |
|
* Compute a groebner basis of FH w.r.t. MWH by modular change of |
|
* ordering. |
|
* Since F is Z-coef, LC(FH)=[1,...,1] and we can use any prime p |
|
* for trace algorithm. |
|
*/ |
/* dp_gr_flags(["Top",1,"NoRA",1]); */ |
/* dp_gr_flags(["Top",1,"NoRA",1]); */ |
#if 0 |
for ( I = 0, HC=[]; I <= N; I++ ) HC = cons(1,HC); |
GH = dp_weyl_gr_main(FH,VDVH,0,0,11); |
GH = nd_gb_candidate(FH,VDVH,11,0,HC,1); |
#else |
|
#if 0 |
|
GH = dp_weyl_gr_main(FH,VDVH,0,-1,11); |
|
#else |
|
GH = nd_weyl_gr_trace(FH,VDVH,0,-1,11); |
|
#endif |
|
#endif |
|
/* dp_gr_flags(["Top",0,"NoRA",0]); */ |
/* dp_gr_flags(["Top",0,"NoRA",0]); */ |
|
|
/* dehomigenize GH */ |
/* dehomigenize GH */ |
Line 232 def in_ww_main(F,VW1,VW2) |
|
Line 234 def in_ww_main(F,VW1,VW2) |
|
|
|
AllData = [G,GIN,VDV,W,T,WtV]; |
AllData = [G,GIN,VDV,W,T,WtV]; |
if ( VW2 ) { |
if ( VW2 ) { |
|
/* take LC(GIN) w.r.t. DRL */ |
|
dp_set_weight(WtV); dp_ord(0); |
|
HC = map(dp_hc,map(dp_ptod,GIN,VDV)); |
V2 = VW2[0]; DV2 = VW2[1]; WtV2 = VW2[2]; |
V2 = VW2[0]; DV2 = VW2[1]; WtV2 = VW2[2]; |
VDV2 = append(V2,DV2); |
VDV2 = append(V2,DV2); |
dp_set_weight(WtV2); |
dp_set_weight(WtV2); |
GIN2 = nd_weyl_gr_trace(GIN,VDV2,0,-1,0); |
GIN2 = nd_gb_candidate(GIN,VDV2,0,0,HC,1); |
InData = [GIN2,VDV2,V2,DV2,WtV2]; |
InData = [GIN2,VDV2,V2,DV2,WtV2]; |
} else { |
} else { |
if ( 0 ) { |
if ( 0 ) { |
|
|
{ |
{ |
if ( member(s,vars(F)) ) |
if ( member(s,vars(F)) ) |
error("ann : the variable 's' is reserved."); |
error("ann : the variable 's' is reserved."); |
|
F = ptozp(F); |
V = vars(F); |
V = vars(F); |
N = length(V); |
N = length(V); |
D = newvect(N); |
D = newvect(N); |
|
if ( type(Wt=getopt(weight)) == -1 ) |
|
for ( I = N-1, Wt = []; I >= 0; I-- ) Wt = append([V[I],1],Wt); |
|
|
for ( I = 0; I < N; I++ ) |
Wt1 = vector(N); |
D[I] = [deg(F,V[I]),V[I]]; |
for ( I = 0, F1 =F; I < N; I++ ) { |
|
VI = Wt[2*I]; WI = Wt[2*I+1]; |
|
for ( J = 0; J < N; J++ ) |
|
if ( VI == V[J] ) break; |
|
F1 = subst(F1,VI,VI^WI); |
|
} |
|
for ( I = 0; I < N; I++ ) D[I] = [deg(F1,V[I]),V[I]]; |
qsort(D,compare_first); |
qsort(D,compare_first); |
for ( V = [], I = N-1; I >= 0; I-- ) |
for ( V = [], I = 0; I < N; I++ ) V = cons(D[I][1],V); |
V = cons(D[I][1],V); |
V = reverse(V); |
|
for ( I = 0; I < N; I++ ) { |
|
VI = Wt[2*I]; WI = Wt[2*I+1]; |
|
for ( J = 0; J < N; J++ ) if ( VI == V[J] ) break; |
|
Wt1[J] = WI; |
|
} |
|
Wt = vtol(Wt1); |
|
|
for ( I = N-1, DV = []; I >= 0; I-- ) |
for ( I = N-1, DV = []; I >= 0; I-- ) |
DV = cons(strtov("d"+rtostr(V[I])),DV); |
DV = cons(strtov("d"+rtostr(V[I])),DV); |
|
|
B = cons(DV[I]+TMP_Y1*diff(F,V[I])*TMP_DT,B); |
B = cons(DV[I]+TMP_Y1*diff(F,V[I])*TMP_DT,B); |
} |
} |
|
|
/* homogenized (heuristics) */ |
Tdeg = w_tdeg(F,V,Wt); |
dp_nelim(2); |
/* y1*y2-1, t-y1*f, dx1+y1*df/dx1*dt ... */ |
G0 = nd_weyl_gr(B,append(W,DW),0,[[0,2],[0,length(W)*2-2]]); |
/* weight for [y1,y2,t, x1,...,xn, dy1,dy2, dt,dx1,...,dxn, h] */ |
|
/* 0 1 2 3 N3-1 N3 N3+1 N3+2 2*N3 */ |
|
/* 1 1 D+1 w1 wn 1 1 1 D D 1 */ |
|
N3 = N+3; |
|
WtV = newvect(2*N3+1); |
|
WtV[0] = WtV[1] = 1; WtV[2] = Tdeg+1; |
|
for ( I = 3; I < N3; I++ ) WtV[I] = Wt[I-3]; |
|
for ( ; I <= N3+2; I++ ) WtV[I] = 1; |
|
for ( ; I < 2*N3; I++ ) WtV[I] = Tdeg; |
|
WtV[2*N3] = 1; |
|
|
|
/* B is already a GB => modular change of ordering can be applied */ |
|
/* any prime is available => HC=[1] */ |
|
dp_set_weight(WtV); |
|
G0 = nd_gb_candidate(B,append(W,DW),[[0,2],[0,length(W)*2-2]],0,[1],1); |
|
dp_set_weight(0); |
G1 = []; |
G1 = []; |
for ( T = G0; T != []; T = cdr(T) ) { |
for ( T = G0; T != []; T = cdr(T) ) { |
E = car(T); VL = vars(E); |
E = car(T); VL = vars(E); |
|
|
return G3; |
return G3; |
} |
} |
|
|
|
def in_gb_oaku(F) |
|
{ |
|
if ( member(s,vars(F)) ) |
|
error("ann : the variable 's' is reserved."); |
|
F = ptozp(F); |
|
V = vars(F); |
|
N = length(V); |
|
D = newvect(N); |
|
if ( type(Wt=getopt(weight)) == -1 ) |
|
for ( I = N-1, Wt = []; I >= 0; I-- ) Wt = append([V[I],1],Wt); |
|
|
|
Wt1 = vector(N); |
|
for ( I = 0, F1 =F; I < N; I++ ) { |
|
VI = Wt[2*I]; WI = Wt[2*I+1]; |
|
for ( J = 0; J < N; J++ ) |
|
if ( VI == V[J] ) break; |
|
F1 = subst(F1,VI,VI^WI); |
|
} |
|
for ( I = 0; I < N; I++ ) D[I] = [deg(F1,V[I]),V[I]]; |
|
qsort(D,compare_first); |
|
for ( V = [], I = 0; I < N; I++ ) V = cons(D[I][1],V); |
|
V = reverse(V); |
|
for ( I = 0; I < N; I++ ) { |
|
VI = Wt[2*I]; WI = Wt[2*I+1]; |
|
for ( J = 0; J < N; J++ ) if ( VI == V[J] ) break; |
|
Wt1[J] = WI; |
|
} |
|
Wt = vtol(Wt1); |
|
|
|
for ( I = N-1, DV = []; I >= 0; I-- ) |
|
DV = cons(strtov("d"+rtostr(V[I])),DV); |
|
|
|
W = append([TMP_Y1,TMP_Y2,TMP_T],V); |
|
DW = append([TMP_DY1,TMP_DY2,TMP_DT],DV); |
|
|
|
B = [TMP_T-TMP_Y1*F]; |
|
for ( I = 0; I < N; I++ ) { |
|
B = cons(DV[I]+TMP_Y1*diff(F,V[I])*TMP_DT,B); |
|
} |
|
|
|
Tdeg = w_tdeg(F,V,Wt); |
|
/* y1*y2-1, t-y1*f, dx1+y1*df/dx1*dt ... */ |
|
/* weight for [y1,y2,t, x1,...,xn, dy1,dy2, dt,dx1,...,dxn, h] */ |
|
/* 0 1 2 3 N3-1 N3 N3+1 N3+2 2*N3 */ |
|
/* 1 1 D+1 1 1 1 1 1 D D 1 */ |
|
N3 = N+3; |
|
WtV = newvect(2*N3+1); |
|
WtV[0] = WtV[1] = 1; WtV[2] = Tdeg+1; |
|
for ( I = 3; I <= N3+2; I++ ) WtV[I] = 1; |
|
for ( ; I < 2*N3; I++ ) WtV[I] = Tdeg; |
|
WtV[2*N3] = 1; |
|
|
|
/* B is already a GB => modular change of ordering can be applied */ |
|
/* any prime is available => HC=[1] */ |
|
dp_set_weight(WtV); |
|
G0 = nd_gb_candidate(B,append(W,DW),[[0,2],[0,length(W)*2-2]],0,[1],1); |
|
dp_set_weight(0); |
|
G1 = map(subst,G0,TMP_Y1,1); |
|
return [G1,append(V,DV)]; |
|
} |
|
|
/* F = [F0,F1,...] */ |
/* F = [F0,F1,...] */ |
|
|
def ann_n(F) |
def ann_n(F) |
|
|
VA = append(U,append(W,V)); |
VA = append(U,append(W,V)); |
DVA = append(DU,append(DW,DV)); |
DVA = append(DU,append(DW,DV)); |
VDV = append(VA,DVA); |
VDV = append(VA,DVA); |
|
#if 0 |
G0 = nd_weyl_gr(B,VDV,0,[[0,2*L],[0,length(VDV)-2*L]]); |
G0 = nd_weyl_gr(B,VDV,0,[[0,2*L],[0,length(VDV)-2*L]]); |
|
#else |
|
G0 = nd_gb_candidate(B,VDV,[[0,2*L],[0,length(VDV)-2*L]],0,[1],1); |
|
#endif |
G1 = []; |
G1 = []; |
for ( T = G0; T != []; T = cdr(T) ) { |
for ( T = G0; T != []; T = cdr(T) ) { |
E = car(T); VL = vars(E); |
E = car(T); VL = vars(E); |
Line 1304 def bf_strat_stage2(L) |
|
Line 1404 def bf_strat_stage2(L) |
|
M = elim_mat(VDV,DVRV); |
M = elim_mat(VDV,DVRV); |
for ( K = 0; K < N; K++ ) |
for ( K = 0; K < N; K++ ) |
M[1][K] = W[K]; |
M[1][K] = W[K]; |
dp_ord(0); H1 = map(dp_ht,map(dp_ptod,G1,VDV)); |
dp_ord(0); D1 = map(dp_ptod,G1,VDV); |
|
H1 = map(dp_ht,D1); HC1 = map(dp_hc,D1); |
dp_ord(M); H2 = map(dp_ht,map(dp_ptod,G1,VDV)); |
dp_ord(M); H2 = map(dp_ht,map(dp_ptod,G1,VDV)); |
if ( H1 == H2 ) |
if ( H1 == H2 ) |
G2 = G1; |
G2 = G1; |
else |
else |
G2 = nd_weyl_gr_trace(G1,VDV,0,-1,M); |
G2 = nd_gb_candidate(G1,VDV,M,0,HC1,1); |
G1 = elimination(G2,DVRV); |
G1 = elimination(G2,DVRV); |
} |
} |
T1 = time(); |
T1 = time(); |
Line 1332 def bf_strat_stage3(L) |
|
Line 1433 def bf_strat_stage3(L) |
|
QQ = L[0]; V0 = L[2]; B = L[3]; BF = L[4]; W0 = L[5]; |
QQ = L[0]; V0 = L[2]; B = L[3]; BF = L[4]; W0 = L[5]; |
NF = length(BF); |
NF = length(BF); |
Data = vector(NF); |
Data = vector(NF); |
|
W1 = W0? cons(1,append(W0,[1])) : 0; |
for ( I = J = 0; I < NF; I++ ) { |
for ( I = J = 0; I < NF; I++ ) { |
DI = tower_in_p(QQ,B,BF[I],V0,W0); |
DI = tower_in_p(QQ,B,BF[I],V0,W0); |
NDI = length(DI); |
NDI = length(DI); |
|
dp_set_weight(W1); |
for ( K = 0; K < J; K++ ) { |
for ( K = 0; K < J; K++ ) { |
if ( length(DK=Data[K]) == NDI ) { |
if ( length(DK=Data[K]) == NDI ) { |
for ( L = 0; L < NDI; L++ ) { |
for ( L = 0; L < NDI; L++ ) { |
Line 1345 def bf_strat_stage3(L) |
|
Line 1448 def bf_strat_stage3(L) |
|
if ( L == NDI ) break; |
if ( L == NDI ) break; |
} |
} |
} |
} |
|
dp_set_weight(0); |
if ( K < J ) { |
if ( K < J ) { |
for ( L = 0, T = []; L < NDI; L++ ) |
for ( L = 0, T = []; L < NDI; L++ ) { |
|
#if 0 |
|
NewId = DK[L][1]; |
|
#else |
|
NewId = ideal_intersection(DK[L][1],DI[L][1],V0,0); |
|
#endif |
T = cons([[DK[L][0][0]*DI[L][0][0],DK[L][0][1]], |
T = cons([[DK[L][0][0]*DI[L][0][0],DK[L][0][1]], |
DK[L][1],DK[L][2]],T); |
NewId,DK[L][2]],T); |
|
} |
Data[K] = reverse(T); |
Data[K] = reverse(T); |
} else |
} else |
Data[J++] = DI; |
Data[J++] = DI; |
Line 1357 def bf_strat_stage3(L) |
|
Line 1467 def bf_strat_stage3(L) |
|
for ( I = 0; I < J; I++ ) |
for ( I = 0; I < J; I++ ) |
Data1[I] = Data[I]; |
Data1[I] = Data[I]; |
T1 = time(); |
T1 = time(); |
Str = stratify_bf(Data1,V0); |
Str = stratify_bf(Data1,V0,W0); |
T2 = time(); |
T2 = time(); |
print(["tower",(T1[0]+T1[1])-(T0[0]+T0[1])]); |
print(["tower",(T1[0]+T1[1])-(T0[0]+T0[1])]); |
print(["strat",(T2[0]+T2[1])-(T1[0]+T1[1])]); |
print(["strat",(T2[0]+T2[1])-(T1[0]+T1[1])]); |
Line 1371 def bf_strat_stage3(L) |
|
Line 1481 def bf_strat_stage3(L) |
|
|
|
def bf_local(F,P) |
def bf_local(F,P) |
{ |
{ |
|
/* F -> F/Fcont */ |
|
F1 = ptozp(F); Fcont = sdiv(F,F1); F = F1; |
if ( type(Heu=getopt(heuristic)) == -1 ) Heu = 0; |
if ( type(Heu=getopt(heuristic)) == -1 ) Heu = 0; |
if ( type(Vord=getopt(vord)) == -1 || type(Vord) != 4 ) Vord = 0; |
if ( type(Vord=getopt(vord)) == -1 || type(Vord) != 4 ) Vord = 0; |
if ( type(Wt=getopt(weight)) == -1 ) Wt = 0; |
if ( type(Wt=getopt(weight)) == -1 ) Wt = 0; |
|
if ( type(Op=getopt(op)) == -1 ) Op = 0; |
L = in_ww(F|weight=Wt,heuristic=Heu,vord=Vord); |
L = in_ww(F|weight=Wt,heuristic=Heu,vord=Vord); |
InData = L[0]; AllData = L[1]; VData = L[2]; |
InData = L[0]; AllData = L[1]; VData = L[2]; |
G = InData[0]; VDV = InData[1]; |
G = InData[0]; VDV = InData[1]; |
Line 1407 def bf_local(F,P) |
|
Line 1520 def bf_local(F,P) |
|
break; |
break; |
if ( List == [] ) error("bf_local : inconsitent intersection"); |
if ( List == [] ) error("bf_local : inconsitent intersection"); |
Ax = car(List); |
Ax = car(List); |
for ( BPT = 1, List = BP; List != []; List = cdr(List) ) |
LB = []; |
|
for ( BPT = 1, List = BP; List != []; List = cdr(List) ) { |
BPT *= car(List)[0]^car(List)[1]; |
BPT *= car(List)[0]^car(List)[1]; |
|
LB = cons([subst(car(List)[0],s,-s-1),car(List)[1]],LB); |
|
} |
|
LB = reverse(LB); |
|
if ( !Op ) return LB; |
|
|
BPT = weyl_subst(BPT,T*DT,VDV); |
BPT = weyl_subst(BPT,T*DT,VDV); |
|
|
/* computation using G0,GIN0,VDV0 */ |
/* computation using G0,GIN0,VDV0 */ |
Line 1428 def bf_local(F,P) |
|
Line 1547 def bf_local(F,P) |
|
CR = conv_tdt(R,F,V0,DV0,T,DT); |
CR = conv_tdt(R,F,V0,DV0,T,DT); |
|
|
dp_set_weight(0); |
dp_set_weight(0); |
return [BP,Ax,CR]; |
Cont = cont(CR); CR /= Cont; |
|
Cont *= dn(Fcont); Den *= nm(Fcont); |
|
Gcd = igcd(Den,Cont); |
|
return [LB,(Den/Gcd)*Ax,(Cont/Gcd)*CR]; |
} |
} |
|
|
/* t^(l+k)*dt^l (k>l) -> (s-k)(s-k-1)...(s-(k+l-1))t^k */ |
/* t^(l+k)*dt^l (k>l) -> (s-k)(s-k-1)...(s-(k+l-1))t^k */ |
Line 1453 def conv_tdt(P,F,V0,DV0,T,DT) |
|
Line 1575 def conv_tdt(P,F,V0,DV0,T,DT) |
|
return subst(R,T,s); |
return subst(R,T,s); |
} |
} |
|
|
def merge_tower(Str,Tower,V) |
/* W1=[W,1], W2=[1,W,1] */ |
|
|
|
def merge_tower(Str,Tower,V,W1,W2) |
{ |
{ |
Prev = car(Tower); T = cdr(Tower); |
Prev = car(Tower); T = cdr(Tower); |
Str1 = []; |
Str1 = []; |
Line 1468 def merge_tower(Str,Tower,V) |
|
Line 1592 def merge_tower(Str,Tower,V) |
|
U = []; |
U = []; |
for ( T = Str; T != []; T = cdr(T) ) { |
for ( T = Str; T != []; T = cdr(T) ) { |
for ( S = Str1; S != []; S = cdr(S) ) { |
for ( S = Str1; S != []; S = cdr(S) ) { |
Int = int_cons(T[0],S[0],V); |
Int = int_cons(T[0],S[0],V,W1,W2); |
if ( Int[0] != [1] ) |
if ( Int[0] != [1] ) |
U = cons(append(Int,[append(T[0][2],S[0][2])]),U); |
U = cons(append(Int,[append(T[0][2],S[0][2])]),U); |
} |
} |
Line 1477 def merge_tower(Str,Tower,V) |
|
Line 1601 def merge_tower(Str,Tower,V) |
|
return U; |
return U; |
} |
} |
|
|
def stratify_bf(Data,V) |
def stratify_bf(Data,V,W) |
{ |
{ |
N = length(Data); |
N = length(Data); |
R = []; |
R = []; |
|
if ( W ) { |
|
W1 = append(W,[1]); |
|
W2 = cons(1,W1); |
|
} else |
|
W1 = W2 = 0; |
for ( I = 0; I < N; I++ ) |
for ( I = 0; I < N; I++ ) |
R = merge_tower(R,Data[I],V); |
R = merge_tower(R,Data[I],V,W1,W2); |
return R; |
return R; |
} |
} |
|
|
Line 1588 def zero_inclusion(A,B,V) |
|
Line 1717 def zero_inclusion(A,B,V) |
|
NV = ttttt; |
NV = ttttt; |
for ( T = A; T != []; T = cdr(T) ) { |
for ( T = A; T != []; T = cdr(T) ) { |
F = car(T); |
F = car(T); |
G = nd_gr_trace(cons(NV*F-1,B),cons(NV,V),1,-1,0); |
G = nd_gr_trace(cons(NV*F-1,B),cons(NV,V),1,1,0); |
if ( type(car(G)) != 1 ) return 0; |
if ( type(car(G)) != 1 ) return 0; |
} |
} |
return 1; |
return 1; |
Line 1622 def elim_mat(V,W) |
|
Line 1751 def elim_mat(V,W) |
|
=(P cap R)-(Q cup S) |
=(P cap R)-(Q cup S) |
*/ |
*/ |
|
|
def int_cons(A,B,V) |
def int_cons(A,B,V,W1,W2) |
{ |
{ |
AZ = A[0]; AN = A[1]; |
AZ = A[0]; AN = A[1]; |
BZ = B[0]; BN = B[1]; |
BZ = B[0]; BN = B[1]; |
|
if ( W1 ) dp_set_weight(W1); |
CZ = nd_gr_trace(append(AZ,BZ),V,1,1,0); |
CZ = nd_gr_trace(append(AZ,BZ),V,1,1,0); |
|
if ( W2 ) dp_set_weight(W2); |
CN = ideal_intersection(AN,BN,V,0); |
CN = ideal_intersection(AN,BN,V,0); |
if ( zero_inclusion(CN,CZ,V) ) |
ZI = zero_inclusion(CN,CZ,V); |
|
dp_set_weight(0); |
|
if ( ZI ) |
return [[1],[]]; |
return [[1],[]]; |
else |
else |
return [CZ,CN]; |
return [CZ,CN]; |
Line 1641 def ideal_intersection(A,B,V,Ord) |
|
Line 1774 def ideal_intersection(A,B,V,Ord) |
|
cons(T,V),1,1,[[0,1],[Ord,length(V)]]); |
cons(T,V),1,1,[[0,1],[Ord,length(V)]]); |
return elimination(G,V); |
return elimination(G,V); |
} |
} |
|
|
|
def nd_gb_candidate(G,V,Ord,Homo,HC,Weyl) |
|
{ |
|
Ind = 0; |
|
N = length(HC); |
|
while ( 1 ) { |
|
P = lprime(Ind++); |
|
for ( I = 0; I < N; I++ ) |
|
if ( !(HC[I]%P) ) break; |
|
if ( I < N ) continue; |
|
if ( Weyl ) |
|
G = nd_weyl_gr_trace(G,V,Homo,-P,Ord); |
|
else |
|
G = nd_gr_trace(G,V,Homo,-P,Ord); |
|
if ( G ) return G; |
|
} |
|
} |
|
|
endmodule $ |
endmodule $ |
end$ |
end$ |
|
|