version 1.8, 2003/11/20 10:41:12 |
version 1.17, 2004/01/06 17:23:23 |
|
|
load("solve")$ |
load("solve")$ |
load("gr")$ |
load("gr")$ |
|
|
def junban(A,B){ |
def nonzerovec(A){ |
|
|
|
for(I=0;I<size(A)[0];I++) |
|
if(A[I]!=0) |
|
return 1$ |
|
|
|
return 0$ |
|
} |
|
|
|
def junban(A,B){ |
return (A<B ? 1:(A>B ? -1:0))$ |
return (A<B ? 1:(A>B ? -1:0))$ |
} |
} |
|
|
def worder(A,B){ |
def worder(A,B){ |
|
|
return (A[0]<B[0] ? 1:(A[0]>B[0] ? -1:0))$ |
return (A[0]<B[0] ? 1:(A[0]>B[0] ? -1:0))$ |
} |
} |
|
|
def wsort(A,B,C){ |
def bsort(A){ |
|
|
D=newvect(length(B))$ |
K=size(A)[0]-1$ |
for(I=0;I<length(B);I++) |
while(K>=0){ |
D[I]=[A[I],B[I],C[I]]$ |
J=-1$ |
|
for(I=1;I<=K;I++) |
|
if(A[I-1][0]<A[I][0]){ |
|
J=I-1$ |
|
X=A[J]$ |
|
A[J]=A[I]$ |
|
A[I]=X$ |
|
} |
|
K=J$ |
|
} |
|
return A$ |
|
} |
|
|
D=qsort(D,worder)$ |
def perm(I,P,TMP){ |
E=[]$ |
|
for(I=0;I<length(B);I++) |
if(I>0){ |
E=cons(D[I][1],E)$ |
TMP=perm(I-1,P,TMP)$ |
E=reverse(E)$ |
for(J=I-1;J>=0;J--){ |
F=[]$ |
T=P[I]$ |
for(I=0;I<length(B);I++) |
P[I]=P[J]$ |
F=cons(D[I][2],F)$ |
P[J]=T$ |
F=reverse(F)$ |
TMP=perm(I-1,P,TMP)$ |
|
T=P[I]$ |
|
P[I]=P[J]$ |
|
P[J]=T$ |
|
} |
|
|
|
return TMP$ |
|
} |
|
else{ |
|
for(TMP0=[],K=0;K<size(P)[0];K++) |
|
TMP0=cons(P[K],TMP0)$ |
|
|
|
TMP=cons(TMP0,TMP)$ |
|
return TMP$ |
|
} |
|
} |
|
|
|
def marge(A,B){ |
|
|
|
RET=[]$ |
|
for(I=0;I<length(A);I++) |
|
for(J=0;J<length(B);J++) |
|
RET=cons(append(A[I],B[J]),RET)$ |
|
|
|
return RET$ |
|
} |
|
|
|
def wsort(A,B,C,FLAG){ |
|
|
|
if(FLAG==0){ |
|
D=newvect(length(B))$ |
|
for(I=0;I<length(B);I++) |
|
D[I]=[A[I],B[I],C[I]]$ |
|
|
|
D=bsort(D)$ |
|
E=[]$ |
|
for(I=0;I<length(B);I++) |
|
E=cons(D[I][1],E)$ |
|
E=reverse(E)$ |
|
F=[]$ |
|
for(I=0;I<length(B);I++) |
|
F=cons(D[I][2],F)$ |
|
F=reverse(F)$ |
|
|
return [E,F]$ |
return [[E,F]]$ |
|
} |
|
else{ |
|
D=newvect(length(B))$ |
|
for(I=0;I<length(B);I++) |
|
D[I]=[A[I],B[I],C[I]]$ |
|
|
|
D=qsort(D,worder)$ |
|
D0=[]$ |
|
|
|
for(I=0,J=0,TMP=[],X=0;I<size(D)[0];I++){ |
|
if(X==D[I][0]) |
|
TMP=cons(cdr(D[I]),TMP)$ |
|
else{ |
|
D0=cons(TMP,D0)$ |
|
TMP=[]$ |
|
TMP=cons(cdr(D[I]),TMP)$ |
|
X=car(D[I])$ |
|
} |
|
} |
|
D0=cdr(reverse(cons(TMP,D0)))$ |
|
D0=map(ltov,D0)$ |
|
for(I=0,TMP=[[]];I<length(D0);I++){ |
|
TMP0=perm(length(D0[I])-1,D0[I],[])$ |
|
TMP=marge(TMP,TMP0)$ |
|
} |
|
|
|
RET=[]$ |
|
for(I=0;I<length(TMP);I++){ |
|
TMP0=[]$ |
|
TMP1=[]$ |
|
for(J=0;J<length(TMP[I]);J++){ |
|
TMP0=cons(TMP[I][J][0],TMP0)$ |
|
TMP1=cons(TMP[I][J][1],TMP1)$ |
|
} |
|
TMP0=reverse(TMP0)$ |
|
TMP1=reverse(TMP1)$ |
|
|
|
RET=cons([TMP0,TMP1],RET)$ |
|
} |
|
|
|
return RET$ |
|
} |
} |
} |
|
|
def derase(A){ |
def derase(A){ |
|
Line 163 def getgcd(A,B){ |
|
for(I=0;I<VarsNumA;I++){ |
for(I=0;I<VarsNumA;I++){ |
|
|
for(J=0;J<VarsNumB;J++) |
for(J=0;J<VarsNumB;J++) |
if(C[J]==A[I][0]) |
if(B[J]==A[I][0]) |
break$ |
break$ |
|
|
if(J<VarsNumB) |
if(J<VarsNumB) |
|
Line 198 def getgcd(A,B){ |
|
return RET$ |
return RET$ |
} |
} |
|
|
def resvars(Res,Vars){ |
def makeret(Res,Vars,FLAG){ |
|
|
ResVars=newvect(length(Vars),Vars)$ |
|
for(I=0;I<length(Res);I++){ |
|
|
|
for(J=0;J<size(ResVars)[0];J++) |
|
if(Res[I][0]==ResVars[J]) |
|
break$ |
|
|
|
if(J<size(ResVars)[0]) |
|
ResVars[J]=Res[I][1]$ |
|
} |
|
return(ResVars)$ |
|
} |
|
|
|
def makeret(Res,Vars){ |
|
|
|
ResNum=length(Res)$ |
ResNum=length(Res)$ |
VarsNum=length(Vars)$ |
VarsNum=length(Vars)$ |
|
|
ResVec=newvect(ResNum)$ |
ResVec=newvect(ResNum)$ |
for(M=0,I=0;I<ResNum;I++){ |
|
if(member(Res[I][0],Vars)){ |
|
ResVec[I]=Res[I][1]$ |
|
|
|
if(type(ResVec[I])==1){ |
for(M=0,I=0;I<ResNum;I++){ |
if(M==0) |
if(member(Res[I][0],Vars)){ |
M=ResVec[I]$ |
ResVec[I]=Res[I][1]$ |
else |
|
if(ResVec[I]<M) |
if(FLAG && type(ResVec[I])==1){ |
M=ResVec[I]$ |
if(M==0) |
} |
M=ResVec[I]$ |
} |
else |
} |
if(ResVec[I]<M) |
|
M=ResVec[I]$ |
|
} |
|
} |
|
} |
|
|
if(M!=0) |
if(M!=0) |
ResVec=ResVec/M; |
ResVec=ResVec/M; |
|
|
RET=newvect(VarsNum,Vars)$ |
RET=newvect(VarsNum,Vars)$ |
|
|
for(I=0;I<ResNum;I++){ |
for(I=0;I<ResNum;I++){ |
for(J=0;J<VarsNum;J++) |
for(J=0;J<VarsNum;J++) |
Line 145 def makeret(Res,Vars){ |
|
Line 233 def makeret(Res,Vars){ |
|
RET[J]=ResVec[I]$ |
RET[J]=ResVec[I]$ |
} |
} |
|
|
|
|
|
for(J=0;J<length(Vars);J++) |
|
RET=map(subst,RET,Vars[J], |
|
strtov(rtostr(Vars[J])+"_deg"))$ |
|
|
for(I=0;I<VarsNum;I++) |
for(I=0;I<VarsNum;I++) |
if(type(RET[I])!=1) |
if(type(RET[I])!=1) |
return [1,RET]$ |
return [1,RET]$ |
Line 257 def checktd(PolyList,Vars,ResVars){ |
|
Line 350 def checktd(PolyList,Vars,ResVars){ |
|
return 1$ |
return 1$ |
} |
} |
|
|
def qcheck(PolyList,Vars){ |
def qcheck(PolyList,Vars,FLAG){ |
|
|
|
RET=[]$ |
Res=qcheckmain(PolyList,Vars)$ |
Res=qcheckmain(PolyList,Vars)$ |
VarsNum=length(Vars)$ |
VarsNum=length(Vars)$ |
|
|
Line 287 def qcheck(PolyList,Vars){ |
|
Line 381 def qcheck(PolyList,Vars){ |
|
|
|
if(nonposdegchk(Res)){ |
if(nonposdegchk(Res)){ |
|
|
ResVars=resvars(Res,Vars)$ |
ResVars=makeret(Res,Vars,0)$ |
|
|
if(checktd(PolyList,Vars,ResVars)==1){ |
if(checktd(PolyList,Vars,ResVars[1])==1){ |
|
if(ResVars[0]==0){ |
for(J=0;J<length(Vars);J++) |
RET=append(RET,wsort(ResVars[1],Vars, |
ResVars=map(subst,ResVars,Vars[J], |
ResVars[1],FLAG))$ |
strtov(rtostr(Vars[J])+"_deg"))$ |
return RET$ |
|
} |
return [wsort(ResVars,Vars,ResVars)]$ |
else{ |
|
RET=append(RET,[[Vars,vtol(ResVars[1])]])$ |
|
return RET$ |
|
} |
} |
} |
else |
else |
return []$ |
return []$ |
Line 305 def qcheck(PolyList,Vars){ |
|
Line 402 def qcheck(PolyList,Vars){ |
|
|
|
} |
} |
|
|
def leastsq(NormMat,ExpMat,Vars){ |
def leastsq(NormMat,ExpMat,Vars,FLAG){ |
|
|
RET=[]$ |
RET=[]$ |
|
|
Line 352 def leastsq(NormMat,ExpMat,Vars){ |
|
Line 449 def leastsq(NormMat,ExpMat,Vars){ |
|
Res=getgcd(Res,Rea)$ |
Res=getgcd(Res,Rea)$ |
|
|
if(nonposdegchk(Res)){ |
if(nonposdegchk(Res)){ |
TMP1=makeret(Res,Vars)$ |
|
|
TMP1=makeret(Res,Vars,1)$ |
|
|
if(TMP1[0]==0){ |
if(TMP1[0]==0){ |
TMP=roundret(TMP1[1]*1.0)$ |
TMP=roundret(TMP1[1])$ |
if(TMP!=[]) |
if(TMP!=[]) |
RET=cons(wsort(TMP1[1],Vars,TMP),RET)$ |
RET=append(RET,wsort(TMP1[1],Vars,TMP,FLAG))$ |
|
|
RET=cons(wsort(TMP1[1],Vars, |
RET=append(RET,wsort(TMP1[1],Vars, |
map(drint,TMP1[1]*1.0)),RET)$ |
map(drint,TMP1[1]*1.0),FLAG))$ |
|
|
return RET$ |
return RET$ |
} |
} |
else{ |
else{ |
RET=cons(wsort(TMP1[1],Vars,TMP1[1]*1.0),RET)$ |
RET=append(RET,[[Vars,vtol(TMP1[1]*1.0)]])$ |
return RET$ |
return RET$ |
} |
} |
} |
} |
Line 373 def leastsq(NormMat,ExpMat,Vars){ |
|
Line 472 def leastsq(NormMat,ExpMat,Vars){ |
|
|
|
} |
} |
|
|
def weightr(ExpMat,Vars,PolyListNum,OneMat){ |
def weightr(ExpMat,Vars,PolyListNum,OneMat,FLAG){ |
|
|
RET=[]$ |
RET=[]$ |
|
|
Line 458 def weightr(ExpMat,Vars,PolyListNum,OneMat){ |
|
Line 557 def weightr(ExpMat,Vars,PolyListNum,OneMat){ |
|
Res=solve(SolveList,SolVars)$ |
Res=solve(SolveList,SolVars)$ |
Res=cons([ExtVars[F],1],Res)$ |
Res=cons([ExtVars[F],1],Res)$ |
|
|
Rea=cons(ExtVars[F],Rea)$ |
TMP=[]$ |
Res=getgcd(Res,Rea)$ |
for(I=0;I<length(Rea);I++) |
|
if(member(Rea[I],Vars)) |
|
TMP=cons(Rea[I],TMP)$ |
|
|
|
if(member(ExtVars[F],Vars)) |
|
TMP=cons(ExtVars[F],TMP)$ |
|
|
|
Res=getgcd(Res,TMP)$ |
|
|
if(nonposdegchk(Res)){ |
if(nonposdegchk(Res)){ |
|
|
TMP1=makeret(Res,Vars)$ |
TMP1=makeret(Res,Vars,1)$ |
|
|
if(TMP1[0]==0){ |
if(TMP1[0]==0){ |
TMP=roundret(TMP1[1]*1.0)$ |
TMP=roundret(TMP1[1])$ |
if(TMP!=[]) |
if(TMP!=[]) |
RET=cons(wsort(TMP1[1],Vars,TMP),RET)$ |
RET=append(RET,wsort(TMP1[1],Vars, |
|
TMP,FLAG))$ |
|
|
RET=cons(wsort(TMP1[1],Vars, |
RET=append(RET,wsort(TMP1[1],Vars, |
map(drint,TMP1[1]*1.0)),RET)$ |
map(drint,TMP1[1]*1.0),FLAG))$ |
} |
} |
else{ |
else{ |
RET=cons(wsort(TMP1[1],Vars,TMP1[1]*1.0),RET)$ |
RET=append(RET,[[Vars,vtol(TMP1[1]*1.0)]])$ |
} |
} |
} |
} |
|
|
Line 482 def weightr(ExpMat,Vars,PolyListNum,OneMat){ |
|
Line 590 def weightr(ExpMat,Vars,PolyListNum,OneMat){ |
|
return [NormMat,RET]$ |
return [NormMat,RET]$ |
} |
} |
|
|
def weight(PolyList,Vars,FLAG){ |
def weight(PolyList,Vars,FLAG1,FLAG2){ |
|
|
Vars0=vars(PolyList)$ |
Vars0=vars(PolyList)$ |
Vars1=[]$ |
Vars1=[]$ |
Line 494 def weight(PolyList,Vars,FLAG){ |
|
Line 602 def weight(PolyList,Vars,FLAG){ |
|
|
|
RET=[]$ |
RET=[]$ |
|
|
/* first */ |
TMP=qcheck(PolyList,Vars,FLAG2)$ |
|
|
TMP=qcheck(PolyList,Vars)$ |
|
|
|
if(TMP!=[]){ |
if(TMP!=[]){ |
RET=append(RET,TMP)$ |
RET=append(RET,TMP)$ |
return cons(1,RET)$ |
return cons(1,RET)$ |
Line 507 def weight(PolyList,Vars,FLAG){ |
|
Line 613 def weight(PolyList,Vars,FLAG){ |
|
|
|
PolyListNum=length(PolyList)$ |
PolyListNum=length(PolyList)$ |
|
|
OneMat=newvect(PolyListNum+1,[0])$ |
if(FLAG1){ |
ExpMat=[]$ |
|
for(I=0;I<PolyListNum;I++){ |
OneMat=newvect(PolyListNum+1,[0])$ |
for(Poly=dp_ptod(PolyList[I],Vars); |
ExpMat=[]$ |
Poly!=0;Poly=dp_rest(Poly)){ |
for(I=0;I<PolyListNum;I++){ |
ExpMat=cons(dp_etov(dp_ht(Poly)),ExpMat)$ |
for(Poly=dp_ptod(PolyList[I],Vars); |
|
Poly!=0;Poly=dp_rest(Poly)){ |
|
ExpMat=cons(dp_etov(dp_ht(Poly)),ExpMat)$ |
|
} |
|
OneMat[I+1]=length(ExpMat)$ |
} |
} |
OneMat[I+1]=length(ExpMat)$ |
|
} |
|
|
|
ExpMat=reverse(ExpMat)$ |
ExpMat=reverse(ExpMat)$ |
ExpMat=newvect(length(ExpMat),ExpMat)$ |
ExpMat=newvect(length(ExpMat),ExpMat)$ |
|
|
/* second */ |
TMP=weightr(ExpMat,Vars,PolyListNum,OneMat,FLAG2)$ |
|
|
if(FLAG){ |
|
TMP=weightr(ExpMat,Vars,PolyListNum,OneMat)$ |
|
RET=append(RET,TMP[1])$ |
RET=append(RET,TMP[1])$ |
|
RET=append(RET,leastsq(TMP[0],ExpMat,Vars,FLAG2))$ |
} |
} |
|
else{ |
|
ExpMat=[]$ |
|
for(I=0;I<PolyListNum;I++){ |
|
for(Poly=dp_ptod(PolyList[I],Vars); |
|
Poly!=0;Poly=dp_rest(Poly)){ |
|
if(nonzerovec(TMP=dp_etov(dp_ht(Poly)))) |
|
ExpMat=cons(TMP,ExpMat)$ |
|
} |
|
} |
|
|
/* third */ |
ExpMat=reverse(ExpMat)$ |
|
ExpMat=newvect(length(ExpMat),ExpMat)$ |
|
|
if(FLAG) |
RET=append(RET,leastsq(0,ExpMat,Vars,FLAG2))$ |
RET=append(RET,leastsq(TMP[0],ExpMat,Vars))$ |
} |
else |
|
RET=append(RET,leastsq(0,ExpMat,Vars))$ |
|
|
|
/* forth */ |
|
|
|
ExpMat=qsort(ExpMat,junban)$ |
ExpMat=qsort(ExpMat,junban)$ |
|
|
ExpMat2=[]$ |
ExpMat2=[]$ |
Line 545 def weight(PolyList,Vars,FLAG){ |
|
Line 657 def weight(PolyList,Vars,FLAG){ |
|
|
|
if(size(ExpMat)[0]!=length(ExpMat2)){ |
if(size(ExpMat)[0]!=length(ExpMat2)){ |
ExpMat=newvect(length(ExpMat2),ExpMat2)$ |
ExpMat=newvect(length(ExpMat2),ExpMat2)$ |
RET=append(RET,leastsq(0,ExpMat,Vars))$ |
RET=append(RET,leastsq(0,ExpMat,Vars,FLAG2))$ |
} |
} |
|
|
RET=derase(RET)$ |
RET=derase(RET)$ |