===================================================================
RCS file: /home/cvs/OpenXM_contrib2/asir2000/lib/gr,v
retrieving revision 1.16
retrieving revision 1.22
diff -u -p -r1.16 -r1.22
--- OpenXM_contrib2/asir2000/lib/gr 2002/09/03 08:12:25 1.16
+++ OpenXM_contrib2/asir2000/lib/gr 2006/07/24 06:36:01 1.22
@@ -45,8 +45,13 @@
* DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE,
* PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE.
*
- * $OpenXM: OpenXM_contrib2/asir2000/lib/gr,v 1.15 2002/06/12 08:19:04 noro Exp $
+ * $OpenXM: OpenXM_contrib2/asir2000/lib/gr,v 1.21 2005/08/02 07:21:48 noro Exp $
*/
+
+module gr $
+ /* Empty for now. It will be used in a future. */
+endmodule $
+
extern INIT_COUNT,ITOR_FAIL$
extern REMOTE_MATRIX,REMOTE_NF,REMOTE_VARS$
@@ -126,6 +131,8 @@ def tolex_tl(G0,V,O,W,H)
def tolex(G0,V,O,W)
{
+ Procs = getopt(procs);
+
TM = TE = TNF = 0;
N = length(V); HM = hmlist(G0,V,O); ZD = zero_dim(HM,V,O);
if ( ZD )
@@ -157,7 +164,10 @@ def tolex(G0,V,O,W)
T0 = time()[0]; NF = gennf(G0,TL,V,O,W[N-1],ZD)[0];
TNF += time()[0] - T0;
T0 = time()[0];
- R = tolex_main(V,O,NF,GM,M,MB);
+ if ( type(Procs) != -1 )
+ R = tolex_d_main(V,O,NF,GM,M,MB,Procs);
+ else
+ R = tolex_main(V,O,NF,GM,M,MB);
TE += time()[0] - T0;
if ( R ) {
if ( dp_gr_print() )
@@ -357,6 +367,85 @@ def tolex_main(V,O,NF,GM,M,MB)
return SL;
}
+def tolex_d_main(V,O,NF,GM,M,MB,Procs)
+{
+ map(ox_reset,Procs);
+ /* register data in servers */
+ map(ox_cmo_rpc,Procs,"register_data_for_find_base",NF,V,O,MB,M);
+ /* discard return value in stack */
+ map(ox_pop_cmo,Procs);
+ Free = Procs;
+ Busy = [];
+ T = GM;
+ SL = [];
+ while ( T != [] || Busy != [] ){
+ if ( Free == [] || T == [] ) {
+ /* someone is working; wait for data */
+ Ready = ox_select(Busy);
+ Busy = setminus(Busy,Ready);
+ Free = append(Ready,Free);
+ for ( ; Ready != []; Ready = cdr(Ready) )
+ SL = cons(ox_get(car(Ready)),SL);
+ } else {
+ P = car(Free);
+ Free = cdr(Free);
+ Busy = cons(P,Busy);
+ Template = car(T);
+ T = cdr(T);
+ ox_cmo_rpc(P,"find_base",Template);
+ ox_push_cmd(P,262); /* 262 = OX_popCMO */
+ }
+ }
+ return SL;
+}
+
+struct find_base_data { NF,V,O,MB,M,PosDim,DV }$
+extern Find_base$
+
+def register_data_for_find_base(NF,V,O,MB,M)
+{
+ Find_base = newstruct(find_base_data);
+ Find_base->NF = NF;
+ Find_base->V = V;
+ Find_base->O = O;
+ Find_base->M = M;
+ Find_base->MB = MB;
+
+ if ( MB ) {
+ Find_base->PosDim = 0;
+ DIM = length(MB);
+ Find_base->DV = newvect(DIM);
+ } else
+ Find_base->PosDim = 1;
+}
+
+def find_base(S) {
+ NF = Find_base->NF;
+ V = Find_base->V;
+ O = Find_base->O;
+ MB = Find_base->MB;
+ M = Find_base->M;
+ PosDim = Find_base->PosDim;
+ DV = Find_base->DV;
+
+ S = p_terms(S,V,2);
+ if ( PosDim ) {
+ MB = gather_nf_terms(S,NF,V,O);
+ DV = newvect(length(MB));
+ }
+ dp_ord(O); RHS = termstomat(NF,map(dp_ptod,cdr(S),V),MB,M);
+ dp_ord(O); NHT = nf_tab_gsl(dp_ptod(car(S),V),NF);
+ dptov(NHT[0],DV,MB);
+ dp_ord(O); B = hen_ttob_gsl([DV,NHT[1]],RHS,cdr(S),M);
+ if ( !B )
+ return 0;
+ Len = length(S);
+ for ( U = B[1]*car(S), I = 1; I < Len; I++ )
+ U += B[0][I-1]*S[I];
+ R = ptozp(U);
+ return R;
+}
+
/*
* NF = [Pairs,DN]
* Pairs = [[NF1,T1],[NF2,T2],...]
@@ -421,6 +510,15 @@ def minipoly(G0,V,O,P,V0)
def gennf(G,TL,V,O,V0,FLAG)
{
+ F = dp_gr_flags();
+ for ( T = F; T != []; T = cdr(T) ) {
+ Key = car(T); T = cdr(T);
+ if ( Key == "Demand" ) {
+ Dir = car(T); break;
+ }
+ }
+ if ( Dir )
+ return gennf_demand(G,TL,V,O,V0,FLAG,Dir);
N = length(V); Len = length(G); dp_ord(O); PS = newvect(Len);
for ( I = 0, T = G, HL = []; T != []; T = cdr(T), I++ ) {
PS[I] = dp_ptod(car(T),V); HL = cons(dp_ht(PS[I]),HL);
@@ -472,6 +570,78 @@ def gennf(G,TL,V,O,V0,FLAG)
return [[map(adj_dn,H,LCM),LCM],PS,GI];
}
+def gennf_demand(G,TL,V,O,V0,FLAG,Dir)
+{
+ N = length(V); Len = length(G); dp_ord(O); PS = newvect(Len);
+ NTL = length(TL);
+ for ( I = 0, T = G, HL = []; T != []; T = cdr(T), I++ ) {
+ PS[I] = dp_ptod(car(T),V); HL = cons(dp_ht(PS[I]),HL);
+ }
+ for ( I = 0, DTL = []; TL != []; TL = cdr(TL) )
+ DTL = cons(dp_ptod(car(TL),V),DTL);
+ for ( I = Len - 1, GI = []; I >= 0; I-- )
+ GI = cons(I,GI);
+
+ USE_TAB = (FLAG != 0);
+ if ( USE_TAB ) {
+ T0 = time()[0];
+ MB = dp_mbase(HL); DIM = length(MB);
+ U = dp_ptod(V0,V);
+ UTAB = newvect(DIM);
+ for ( I = 0; I < DIM; I++ ) {
+ UTAB[I] = [MB[I],remove_cont(dp_true_nf(GI,U*MB[I],PS,1))];
+ if ( dp_gr_print() )
+ print(".",2);
+ }
+ if ( dp_gr_print() )
+ print("");
+ TTAB = time()[0]-T0;
+ }
+
+ T0 = time()[0];
+ for ( LCM = 1, Index = 0, H = []; DTL != []; Index++ ) {
+ if ( dp_gr_print() )
+ print(".",2);
+ T = car(DTL); DTL = cdr(DTL);
+ if ( L = search_redble(T,H) ) {
+ L = nf_load(Dir,L[0]);
+ DD = dp_subd(T,L[1]);
+ if ( USE_TAB && (DD == U) ) {
+ NF = nf_tab(L[0],UTAB);
+ NF = [NF[0],dp_hc(L[1])*NF[1]*T];
+ } else
+ NF = nf(GI,L[0]*dp_subd(T,L[1]),dp_hc(L[1])*T,PS);
+ } else
+ NF = nf(GI,T,T,PS);
+ NF = remove_cont(NF);
+ nf_save(NF,Dir,Index);
+ H = cons([Index,NF[1]],H);
+ LCM = ilcm(LCM,dp_hc(NF[1]));
+ }
+ TNF = time()[0]-T0;
+ if ( dp_gr_print() )
+ print("gennf(TAB="+rtostr(TTAB)+" NF="+rtostr(TNF)+")");
+
+ for ( I = 0; I < NTL; I++ ) {
+ NF = nf_load(Dir,I);
+ NF = adj_dn(NF,LCM);
+ nf_save(NF,Dir,I);
+ }
+ for ( H = [], I = NTL-1; I >= 0; I-- )
+ H = cons(nf_load(Dir,I),H);
+ return [[H,LCM],PS,GI];
+}
+
+def nf_load(Dir,I)
+{
+ return bload(Dir+"/nf"+rtostr(I));
+}
+
+def nf_save(NF,Dir,I)
+{
+ bsave(NF,Dir+"/nf"+rtostr(I));
+}
+
def adj_dn(P,D)
{
return [(idiv(D,dp_hc(P[1])))*P[0],dp_ht(P[1])];
@@ -993,9 +1163,15 @@ def gb_comp(A,B)
LB = length(B);
if ( LA != LB )
return 0;
- A1 = qsort(newvect(LA,A));
- B1 = qsort(newvect(LB,B));
+ A = newvect(LA,A);
+ B = newvect(LB,B);
for ( I = 0; I < LA; I++ )
+ A[I] *= headsgn(A[I]);
+ for ( I = 0; I < LB; I++ )
+ B[I] *= headsgn(B[I]);
+ A1 = qsort(A);
+ B1 = qsort(B);
+ for ( I = 0; I < LA; I++ )
if ( A1[I] != B1[I] && A1[I] != -B1[I] )
break;
return I == LA ? 1 : 0;
@@ -1492,6 +1668,105 @@ def check_trace(NF,NFIndex,HL)
error("check_trace");
}
+/*
+ * Trace = [Input,[[j1,[[c,i,m,d],...]],[j2,[[...],...]],...]]
+ * if c != 0
+ * g = 0
+ * g = (c*g + m*gi)/d
+ * ...
+ * finally fj = g
+ */
+
+def show_trace(Trace,V)
+{
+ Input = Trace[0];
+ for ( I = 0, T = Input; T != []; T = cdr(T), I++ ) {
+ print("F"+rtostr(I)+"=",0);
+ print(dp_dtop(car(T),V));
+ }
+ Trace = cdr(Trace);
+ for ( T = Trace; T != []; T = cdr(T) ) {
+ HL = car(T);
+ J = car(HL); HL = HL[1];
+ L = length(HL);
+ print("F"+rtostr(J)+"=",0);
+ for ( I = 0; I < L; I++ ) print("(",0);
+ for ( First = 1, S = HL; S != []; S = cdr(S) ) {
+ H = car(S);
+
+ Coeff = H[0];
+ Index = H[1];
+ Monomial = H[2];
+ Denominator = H[3];
+ if ( First ) {
+ if ( Monomial != 1 ) {
+ print("(",0);
+ print(type(Monomial)==9?dp_dtop(Monomial,V):Monomial,0);
+ print(")*",0);
+ }
+ print("F"+rtostr(Index)+")",0);
+ } else {
+ if ( Coeff != 1 ) {
+ print("*(",0); print(Coeff,0); print(")",0);
+ }
+ print("+",0);
+ if ( Monomial != 1 ) {
+ print("(",0);
+ print(type(Monomial)==9?dp_dtop(Monomial,V):Monomial,0);
+ print(")*",0);
+ }
+ print("F"+rtostr(Index)+")",0);
+ if ( Denominator != 1 ) {
+ print("/",0); print(Denominator,0);
+ }
+ }
+ if ( First ) First = 0;
+ }
+ print("");
+ }
+}
+
+def generating_relation(Trace,V)
+{
+ Trace = cdr(Trace);
+ Tab = [];
+ for ( T = Trace; T != []; T = cdr(T) ) {
+ HL = car(T);
+ J = car(HL); HL = HL[1];
+ L = length(HL);
+ LHS = strtov("f"+rtostr(J));
+ Dn = 1;
+ for ( First = 1, S = HL; S != []; S = cdr(S) ) {
+ H = car(S);
+
+ Coeff = H[0];
+ Index = H[1];
+ Monomial = type(H[2])==9?dp_dtop(H[2],V):H[2];
+ Denominator = H[3];
+ F = strtov("f"+rtostr(Index));
+ for ( Z = Tab; Z != []; Z = cdr(Z) )
+ if ( Z[0][0] == F ) break;
+ if ( Z != [] ) Value = Z[0][1];
+ else Value = [F,1];
+ if ( First ) {
+ RHS = Monomial*Value[0];
+ Dn *= Value[1];
+ } else {
+ RHS = RHS*Coeff*Value[1]+Dn*Value[0]*Monomial;
+ Dn = Value[1]*Dn*Denominator;
+ }
+ VVVV = tttttttt;
+ P = ptozp(Dn*VVVV+RHS);
+ RHS = coef(P,0,VVVV);
+ Dn = coef(P,1,VVVV);
+ if ( First ) First = 0;
+ }
+ Tab = cons([LHS,[RHS,Dn]],Tab);
+ }
+ return Tab;
+}
+
+end$
/*
* realloc NFArray so that it can hold * an element as NFArray[Ind].
*/