version 1.23, 2003/04/20 11:59:57 |
version 1.26, 2003/10/20 00:58:47 |
|
|
* DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE, |
* DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE, |
* PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE. |
* PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE. |
* |
* |
* $OpenXM: OpenXM_contrib2/asir2000/lib/bfct,v 1.22 2003/04/20 08:54:28 noro Exp $ |
* $OpenXM: OpenXM_contrib2/asir2000/lib/bfct,v 1.25 2003/04/28 03:02:52 noro Exp $ |
*/ |
*/ |
/* requires 'primdec' */ |
/* requires 'primdec' */ |
|
|
#define TMP_S ssssssss |
#define TMP_S ssssssss |
#define TMP_T tttttttt |
#define TMP_DS dssssssss |
|
#define TMP_T dtttttttt |
|
#define TMP_DT tttttttt |
#define TMP_Y1 yyyyyyyy1 |
#define TMP_Y1 yyyyyyyy1 |
|
#define TMP_DY1 dyyyyyyyy1 |
#define TMP_Y2 yyyyyyyy2 |
#define TMP_Y2 yyyyyyyy2 |
|
#define TMP_DY2 dyyyyyyyy2 |
|
|
extern LIBRARY_GR_LOADED$ |
if (!module_definedp("gr")) load("gr") $$ |
extern LIBRARY_PRIMDEC_LOADED$ |
if (!module_definedp("primdec")) load("primdec") $$ |
|
module bfct $ |
|
/* Empty for now. It will be used in a future. */ |
|
endmodule $ |
|
|
if(!LIBRARY_GR_LOADED) load("gr"); else ; LIBRARY_GR_LOADED = 1$ |
|
if(!LIBRARY_PRIMDEC_LOADED) load("primdec"); else ; LIBRARY_PRIMDEC_LOADED = 1$ |
|
|
|
/* toplevel */ |
/* toplevel */ |
|
|
def bfunction(F) |
def bfunction(F) |
{ |
{ |
/* XXX */ |
|
F = replace_vars_f(F); |
|
|
|
V = vars(F); |
V = vars(F); |
N = length(V); |
N = length(V); |
D = newvect(N); |
D = newvect(N); |
|
|
|
|
def ann(F) |
def ann(F) |
{ |
{ |
/* XXX */ |
if ( member(s,vars(F)) ) |
F = replace_vars_f(F); |
error("ann : the variable 's' is reserved."); |
|
|
V = vars(F); |
V = vars(F); |
N = length(V); |
N = length(V); |
D = newvect(N); |
D = newvect(N); |
|
|
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); |
|
|
W = append([y1,y2,t],V); |
W = append([TMP_Y1,TMP_Y2,TMP_T],V); |
DW = append([dy1,dy2,dt],DV); |
DW = append([TMP_DY1,TMP_DY2,TMP_DT],DV); |
|
|
B = [1-y1*y2,t-y1*F]; |
B = [1-TMP_Y1*TMP_Y2,TMP_T-TMP_Y1*F]; |
for ( I = 0; I < N; I++ ) { |
for ( I = 0; I < N; I++ ) { |
B = cons(DV[I]+y1*diff(F,V[I])*dt,B); |
B = cons(DV[I]+TMP_Y1*diff(F,V[I])*TMP_DT,B); |
} |
} |
|
|
/* homogenized (heuristics) */ |
/* homogenized (heuristics) */ |
|
|
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); |
if ( !member(y1,VL) && !member(y2,VL) ) |
if ( !member(TMP_Y1,VL) && !member(TMP_Y2,VL) ) |
G1 = cons(E,G1); |
G1 = cons(E,G1); |
} |
} |
G2 = map(psi,G1,t,dt); |
G2 = map(psi,G1,TMP_T,TMP_DT); |
G3 = map(subst,G2,t,-1-s); |
G3 = map(subst,G2,TMP_T,-1-s); |
return G3; |
return G3; |
} |
} |
|
|
|
|
|
|
def ann0(F) |
def ann0(F) |
{ |
{ |
/* XXX */ |
F = subst(F,s,TMP_S); |
F = replace_vars_f(F); |
Ann = ann(F); |
|
Bf = bfunction(F); |
|
|
V = vars(F); |
|
N = length(V); |
|
D = newvect(N); |
|
|
|
for ( I = 0; I < N; I++ ) |
|
D[I] = [deg(F,V[I]),V[I]]; |
|
qsort(D,compare_first); |
|
for ( V = [], I = 0; I < N; I++ ) |
|
V = cons(D[I][1],V); |
|
|
|
for ( I = N-1, DV = []; I >= 0; I-- ) |
|
DV = cons(strtov("d"+rtostr(V[I])),DV); |
|
|
|
/* XXX : heuristics */ |
|
W = append([y1,y2,t],reverse(V)); |
|
DW = append([dy1,dy2,dt],reverse(DV)); |
|
WDW = append(W,DW); |
|
|
|
B = [1-y1*y2,t-y1*F]; |
|
for ( I = 0; I < N; I++ ) { |
|
B = cons(DV[I]+y1*diff(F,V[I])*dt,B); |
|
} |
|
|
|
/* homogenized (heuristics) */ |
|
dp_nelim(2); |
|
G0 = dp_weyl_gr_main(B,WDW,1,0,6); |
|
G1 = []; |
|
for ( T = G0; T != []; T = cdr(T) ) { |
|
E = car(T); VL = vars(E); |
|
if ( !member(y1,VL) && !member(y2,VL) ) |
|
G1 = cons(E,G1); |
|
} |
|
G2 = map(psi,G1,t,dt); |
|
G3 = map(subst,G2,t,-1-s); |
|
|
|
/* G3 = J_f(s) */ |
|
|
|
V1 = cons(s,V); DV1 = cons(ds,DV); V1DV1 = append(V1,DV1); |
|
G4 = dp_weyl_gr_main(cons(F,G3),V1DV1,0,1,0); |
|
Bf = weyl_minipoly(G4,V1DV1,0,s); |
|
|
|
FList = cdr(fctr(Bf)); |
FList = cdr(fctr(Bf)); |
for ( T = FList, Min = 0; T != []; T = cdr(T) ) { |
for ( T = FList, Min = 0; T != []; T = cdr(T) ) { |
LF = car(car(T)); |
LF = car(car(T)); |
|
|
if ( dn(Root) == 1 && Root < Min ) |
if ( dn(Root) == 1 && Root < Min ) |
Min = Root; |
Min = Root; |
} |
} |
return [Min,map(subst,G3,s,Min)]; |
return [Min,map(subst,Ann,s,Min,TMP_S,s,TMP_DS,ds)]; |
} |
} |
|
|
def indicial1(F,V) |
|
{ |
|
W = append([y1,t],V); |
|
N = length(V); |
|
B = [t-y1*F]; |
|
for ( I = N-1, DV = []; I >= 0; I-- ) |
|
DV = cons(strtov("d"+rtostr(V[I])),DV); |
|
DW = append([dy1,dt],DV); |
|
for ( I = 0; I < N; I++ ) { |
|
B = cons(DV[I]+y1*diff(F,V[I])*dt,B); |
|
} |
|
dp_nelim(1); |
|
|
|
/* homogenized (heuristics) */ |
|
G0 = dp_weyl_gr_main(B,append(W,DW),1,0,6); |
|
G1 = map(subst,G0,y1,1); |
|
G2 = map(psi,G1,t,dt); |
|
G3 = map(subst,G2,t,-s-1); |
|
return G3; |
|
} |
|
|
|
def psi(F,T,DT) |
def psi(F,T,DT) |
{ |
{ |
D = dp_ptod(F,[T,DT]); |
D = dp_ptod(F,[T,DT]); |
|
|
return Minipoly; |
return Minipoly; |
} |
} |
|
|
|
/* called from bfct() only */ |
|
|
|
def indicial1(F,V) |
|
{ |
|
W = append([y1,t],V); |
|
N = length(V); |
|
B = [t-y1*F]; |
|
for ( I = N-1, DV = []; I >= 0; I-- ) |
|
DV = cons(strtov("d"+rtostr(V[I])),DV); |
|
DW = append([dy1,dt],DV); |
|
for ( I = 0; I < N; I++ ) { |
|
B = cons(DV[I]+y1*diff(F,V[I])*dt,B); |
|
} |
|
dp_nelim(1); |
|
|
|
/* homogenized (heuristics) */ |
|
G0 = dp_weyl_gr_main(B,append(W,DW),1,0,6); |
|
G1 = map(subst,G0,y1,1); |
|
G2 = map(psi,G1,t,dt); |
|
G3 = map(subst,G2,t,-s-1); |
|
return G3; |
|
} |
|
|
/* b-function computation via generic_bfct() (experimental) */ |
/* b-function computation via generic_bfct() (experimental) */ |
|
|
def bfct_via_gbfct(F) |
def bfct_via_gbfct(F) |
{ |
{ |
/* XXX */ |
|
F = replace_vars_f(F); |
|
|
|
V = vars(F); |
V = vars(F); |
N = length(V); |
N = length(V); |
D = newvect(N); |
D = newvect(N); |
Line 445 def bfct_via_gbfct(F) |
|
Line 404 def bfct_via_gbfct(F) |
|
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 = [t-F]; |
B = [TMP_T-F]; |
for ( I = 0; I < N; I++ ) { |
for ( I = 0; I < N; I++ ) { |
B = cons(DV[I]+diff(F,V[I])*dt,B); |
B = cons(DV[I]+diff(F,V[I])*TMP_DT,B); |
} |
} |
V1 = cons(t,V); DV1 = cons(dt,DV); |
V1 = cons(TMP_T,V); DV1 = cons(TMP_DT,DV); |
W = newvect(N+1); |
W = newvect(N+1); |
W[0] = 1; |
W[0] = 1; |
R = generic_bfct(B,V1,DV1,W); |
R = generic_bfct(B,V1,DV1,W); |
Line 461 def bfct_via_gbfct(F) |
|
Line 420 def bfct_via_gbfct(F) |
|
|
|
def bfct_via_gbfct_weight(F,V) |
def bfct_via_gbfct_weight(F,V) |
{ |
{ |
/* XXX */ |
|
F = replace_vars_f(F); |
|
V = replace_vars_v(V); |
|
|
|
N = length(V); |
N = length(V); |
D = newvect(N); |
D = newvect(N); |
Wt = getopt(weight); |
Wt = getopt(weight); |
Line 486 def bfct_via_gbfct_weight(F,V) |
|
Line 441 def bfct_via_gbfct_weight(F,V) |
|
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 = [t-F]; |
B = [TMP_T-F]; |
for ( I = 0; I < N; I++ ) { |
for ( I = 0; I < N; I++ ) { |
B = cons(DV[I]+diff(F,V[I])*dt,B); |
B = cons(DV[I]+diff(F,V[I])*TMP_DT,B); |
} |
} |
V1 = cons(t,V); DV1 = cons(dt,DV); |
V1 = cons(TMP_T,V); DV1 = cons(TMP_DT,DV); |
W = newvect(N+1); |
W = newvect(N+1); |
W[0] = 1; |
W[0] = 1; |
R = generic_bfct_1(B,V1,DV1,W); |
R = generic_bfct_1(B,V1,DV1,W); |
Line 502 def bfct_via_gbfct_weight(F,V) |
|
Line 457 def bfct_via_gbfct_weight(F,V) |
|
|
|
def bfct_via_gbfct_weight_1(F,V) |
def bfct_via_gbfct_weight_1(F,V) |
{ |
{ |
/* XXX */ |
|
F = replace_vars_f(F); |
|
V = replace_vars_v(V); |
|
|
|
N = length(V); |
N = length(V); |
D = newvect(N); |
D = newvect(N); |
Wt = getopt(weight); |
Wt = getopt(weight); |
Line 527 def bfct_via_gbfct_weight_1(F,V) |
|
Line 478 def bfct_via_gbfct_weight_1(F,V) |
|
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 = [t-F]; |
B = [TMP_T-F]; |
for ( I = 0; I < N; I++ ) { |
for ( I = 0; I < N; I++ ) { |
B = cons(DV[I]+diff(F,V[I])*dt,B); |
B = cons(DV[I]+diff(F,V[I])*TMP_DT,B); |
} |
} |
V1 = append(V,[t]); DV1 = append(DV,[dt]); |
V1 = append(V,[TMP_T]); DV1 = append(DV,[TMP_DT]); |
W = newvect(N+1); |
W = newvect(N+1); |
W[N] = 1; |
W[N] = 1; |
R = generic_bfct_1(B,V1,DV1,W); |
R = generic_bfct_1(B,V1,DV1,W); |
Line 541 def bfct_via_gbfct_weight_1(F,V) |
|
Line 492 def bfct_via_gbfct_weight_1(F,V) |
|
|
|
def bfct_via_gbfct_weight_2(F,V) |
def bfct_via_gbfct_weight_2(F,V) |
{ |
{ |
/* XXX */ |
|
F = replace_vars_f(F); |
|
V = replace_vars_v(V); |
|
|
|
N = length(V); |
N = length(V); |
D = newvect(N); |
D = newvect(N); |
Wt = getopt(weight); |
Wt = getopt(weight); |
Line 581 def bfct_via_gbfct_weight_2(F,V) |
|
Line 528 def bfct_via_gbfct_weight_2(F,V) |
|
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 = [t-F]; |
B = [TMP_T-F]; |
for ( I = 0; I < N; I++ ) { |
for ( I = 0; I < N; I++ ) { |
B = cons(DV[I]+diff(F,V[I])*dt,B); |
B = cons(DV[I]+diff(F,V[I])*TMP_DT,B); |
} |
} |
V1 = cons(t,V); DV1 = cons(dt,DV); |
V1 = cons(TMP_T,V); DV1 = cons(TMP_DT,DV); |
V2 = append(V,[t]); DV2 = append(DV,[dt]); |
V2 = append(V,[TMP_T]); DV2 = append(DV,[TMP_DT]); |
W = newvect(N+1,[1]); |
W = newvect(N+1,[1]); |
dp_weyl_set_weight(W); |
dp_weyl_set_weight(W); |
|
|
Line 637 def bfct_via_gbfct_weight_2(F,V) |
|
Line 584 def bfct_via_gbfct_weight_2(F,V) |
|
return subst(R,s,-s-1); |
return subst(R,s,-s-1); |
} |
} |
|
|
|
/* minimal polynomial of s; modular computation */ |
|
|
def weyl_minipolym(G,V,O,M,V0) |
def weyl_minipolym(G,V,O,M,V0) |
{ |
{ |
N = length(V); |
N = length(V); |
Line 678 def weyl_minipolym(G,V,O,M,V0) |
|
Line 627 def weyl_minipolym(G,V,O,M,V0) |
|
} |
} |
} |
} |
|
|
|
/* minimal polynomial of s over Q */ |
|
|
def weyl_minipoly(G0,V0,O0,P) |
def weyl_minipoly(G0,V0,O0,P) |
{ |
{ |
HM = hmlist(G0,V0,O0); |
HM = hmlist(G0,V0,O0); |
|
|
S = append(S,[F]); |
S = append(S,[F]); |
return S; |
return S; |
} |
} |
|
|
def member(A,L) { |
|
for ( ; L != []; L = cdr(L) ) |
|
if ( A == car(L) ) |
|
return 1; |
|
return 0; |
|
} |
|
|
|
def intersection(A,B) |
def intersection(A,B) |
{ |
{ |