Annotation of OpenXM_contrib2/asir2000/lib/weight, Revision 1.2
1.2 ! kimura 1: load("solve")$
! 2:
! 3: def nonposdegchk(Res){
! 4:
! 5: for(I=0;I<length(Res);I++)
! 6: if(Res[I][1]<=0)
! 7: return 0$
! 8:
! 9: return 1$
! 10: }
! 11:
! 12: def resvars(Res,Vars){
! 13:
! 14: ResVars=newvect(length(Vars),Vars)$
! 15:
! 16: for(I=0;I<length(Res);I++){
! 17:
! 18: for(J=0;J<size(ResVars)[0];J++)
! 19: if(Res[I][0]==ResVars[J])
! 20: break$
! 21:
! 22: ResVars[J]=Res[I][1]$
! 23: }
! 24:
! 25: return(ResVars)$
! 26: }
! 27:
! 28: def makeret(Res,Vars,B){
! 29:
! 30: VarsNum=length(Vars)$
! 31:
! 32: ResMat=newvect(VarsNum)$
! 33: for(I=0;I<VarsNum;I++)
! 34: ResMat[I]=newvect(2)$
! 35:
! 36: for(I=0;I<VarsNum;I++){
! 37: ResMat[I][0]=Vars[I]$
! 38: ResMat[I][1]=Vars[I]$
! 39: }
! 40:
! 41: for(F=0,L=1,D=0,M=1,I=0;I<length(Res);I++){
! 42:
! 43: for(J=0;J<size(ResMat)[0];J++)
! 44: if(Res[I][0]==ResMat[J][0])
! 45: break$
! 46:
! 47: if(J<VarsNum){
! 48: K=Res[I][1]$
! 49: ResMat[J][1]=K$
! 50:
! 51: if(F==0 && type(K)==1){
! 52: if(B==2){
! 53: L=ilcm(L,dn(K))$
! 54: D=igcd(D,nm(K))$
! 55: }
! 56: else{
! 57: if(K<M)
! 58: M=K$
! 59: }
! 60: }
! 61: else
! 62: F=1$
! 63:
! 64: }
! 65: }
! 66:
! 67: if(F==0)
! 68: if(B==2)
! 69: for(I=0;I<VarsNum;I++)
! 70: ResMat[I][1]=ResMat[I][1]*L/D$
! 71: else
! 72: for(I=0;I<VarsNum;I++)
! 73: ResMat[I][1]=ResMat[I][1]/M*1.0$
! 74:
! 75: for(I=0;I<VarsNum;I++)
! 76: for(J=0;J<length(Vars);J++)
! 77: ResMat[I][1]=subst(ResMat[I][1],Vars[J],
! 78: strtov(rtostr(Vars[J])+"_deg"))$
! 79:
! 80: ResMat=map(vtol,ResMat)$
! 81: return(vtol(ResMat))$
! 82:
! 83: }
! 84:
! 85: def afo(A,B){
! 86:
! 87: for(I=0;I<size(A)[0];I++){
! 88: if(A[I]<B[I])
! 89: return 1$
! 90:
! 91: if(A[I]>B[I])
! 92: return -1$
! 93: }
! 94:
! 95: return 0$
! 96: }
! 97:
! 98: def weight(PolyList,Vars){
! 99:
! 100: dp_ord(2)$
! 101:
! 102: PolyListNum=length(PolyList)$
! 103:
! 104: ExpMat=[]$
! 105: for(I=0;I<PolyListNum;I++)
! 106: for(Poly=dp_ptod(PolyList[I],Vars);Poly!=0;Poly=dp_rest(Poly))
! 107: ExpMat=cons(dp_etov(dp_ht(Poly)),ExpMat)$
! 108:
! 109: ExpMat=reverse(ExpMat)$
! 110: ExpMat=newvect(length(ExpMat),ExpMat)$
! 111:
! 112: ExpMatRowNum=size(ExpMat)[0]$
! 113: ExpMatColNum=size(ExpMat[0])[0]$
! 114: ExtMatColNum=ExpMatColNum+PolyListNum$
! 115:
! 116: OneMat=newvect(PolyListNum+1,[0])$
! 117: for(I=0,SUM=0;I<PolyListNum;I++){
! 118: SUM+=nmono(PolyList[I])$
! 119: OneMat[I+1]=SUM$
! 120: }
! 121:
! 122: RevOneMat=newvect(ExpMatRowNum)$
! 123: for(I=0;I<PolyListNum;I++)
! 124: for(J=OneMat[I];J<OneMat[I+1];J++)
! 125: RevOneMat[J]=I$
! 126:
! 127: NormMat=newmat(ExpMatColNum,ExtMatColNum)$
! 128:
! 129: for(I=0;I<ExpMatColNum;I++)
! 130: for(J=0;J<ExpMatColNum;J++)
! 131: for(K=0;K<ExpMatRowNum;K++)
! 132: NormMat[I][J]+=ExpMat[K][I]*ExpMat[K][J]$
! 133:
! 134: for(I=0;I<ExpMatColNum;I++)
! 135: for(J=0;J<PolyListNum-1;J++)
! 136: for(K=OneMat[J];K<OneMat[J+1];K++)
! 137: NormMat[I][J+ExpMatColNum]-=ExpMat[K][I]$
! 138:
! 139: for(I=0;I<ExpMatColNum;I++)
! 140: for(J=OneMat[PolyListNum-1];J<OneMat[PolyListNum];J++)
! 141: NormMat[I][ExtMatColNum-1]+=ExpMat[J][I]$
! 142:
! 143: NormMat2=newmat(PolyListNum-1,ExpMatColNum+1)$
! 144:
! 145: for(I=0;I<PolyListNum-1;I++)
! 146: for(J=0;J<ExpMatColNum;J++)
! 147: for(K=OneMat[I];K<OneMat[I+1];K++)
! 148: NormMat2[I][J]-=ExpMat[K][J]$
! 149:
! 150: for(I=0;I<PolyListNum-1;I++)
! 151: NormMat2[I][ExpMatColNum]=OneMat[I+1]-OneMat[I]$
! 152:
! 153: ExtVars=Vars$
! 154: for(I=0;I<PolyListNum-1;I++)
! 155: ExtVars=append(ExtVars,[uc()])$
! 156:
! 157: SolveList=[]$
! 158: for(I=0;I<ExpMatColNum;I++){
! 159: TMP=0$
! 160: for(J=0;J<ExtMatColNum-1;J++)
! 161: TMP+=NormMat[I][J]*ExtVars[J]$
! 162:
! 163: TMP-=NormMat[I][ExtMatColNum-1]$
! 164: SolveList=cons(TMP,SolveList)$
! 165: }
! 166:
! 167: for(I=0;I<PolyListNum-1;I++){
! 168: TMP=0$
! 169: for(J=0;J<ExpMatColNum;J++)
! 170: TMP+=NormMat2[I][J]*ExtVars[J]$
! 171:
! 172: TMP+=NormMat2[I][ExpMatColNum]*ExtVars[I+ExpMatColNum]$
! 173:
! 174: SolveList=cons(TMP,SolveList)$
! 175: }
! 176:
! 177: ReaVars=vars(SolveList)$
! 178: Res=solve(SolveList,reverse(ExtVars))$
! 179:
! 180: RET=[]$
! 181: if(nonposdegchk(Res)){
! 182:
! 183: ResVars=resvars(Res,ExtVars)$
! 184: ResVars=append(vtol(ResVars),[1])$
! 185:
! 186: for(I=0;I<ExpMatRowNum;I++){
! 187: TMP=0$
! 188: for(J=0;J<ExpMatColNum;J++)
! 189: TMP+=ExpMat[I][J]*ResVars[J]$
! 190:
! 191: TMP-=ResVars[RevOneMat[I]+ExpMatColNum]$
! 192:
! 193: if(TMP!=0)
! 194: break$
! 195: }
! 196:
! 197: if(I==ExpMatRowNum){
! 198: print("complitely homogenized")$
! 199: return([makeret(Res,Vars,2)])$
! 200: }
! 201: else
! 202: RET=cons(makeret(Res,Vars,1),RET)$
! 203: }
! 204:
! 205: ExpMat=qsort(ExpMat,afo)$
! 206: ExpMat2=[]$
! 207: for(I=0;I<size(ExpMat)[0];I++)
! 208: if(car(ExpMat2)!=ExpMat[I])
! 209: ExpMat2=cons(ExpMat[I],ExpMat2)$
! 210:
! 211: ExpMat=newvect(length(ExpMat2),ExpMat2)$
! 212: ExpMatRowNum=size(ExpMat)[0]$
! 213: ExpMatColNum=size(ExpMat[0])[0]$
! 214:
! 215: NormMat=newmat(ExpMatColNum,ExpMatColNum+1)$
! 216:
! 217: for(I=0;I<ExpMatColNum;I++)
! 218: for(J=0;J<ExpMatColNum;J++)
! 219: for(K=0;K<ExpMatRowNum;K++)
! 220: NormMat[I][J]+=ExpMat[K][I]*ExpMat[K][J]$
! 221:
! 222: for(I=0;I<ExpMatColNum;I++)
! 223: for(J=0;J<ExpMatRowNum;J++)
! 224: NormMat[I][ExpMatColNum]+=ExpMat[J][I]$
! 225:
! 226: SolveList=[]$
! 227: for(I=0;I<ExpMatColNum;I++){
! 228: TMP=0$
! 229: for(J=0;J<ExpMatColNum;J++)
! 230: TMP+=NormMat[I][J]*Vars[J]$
! 231:
! 232: TMP-=NormMat[I][ExpMatColNum]$
! 233: SolveList=cons(TMP,SolveList)$
! 234: }
! 235:
! 236: Res=solve(SolveList,Vars)$
! 237:
! 238: if(nonposdegchk(Res))
! 239: RET=cons(makeret(Res,Vars,1),RET)$
! 240:
! 241: return reverse(RET)$
! 242: }
! 243:
! 244: end$
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>