[BACK]Return to os_muldif.rr CVS log [TXT][DIR] Up to [local] / OpenXM / src / asir-contrib / packages / src

Diff for /OpenXM/src/asir-contrib/packages/src/os_muldif.rr between version 1.8 and 1.58

version 1.8, 2017/05/10 02:37:37 version 1.58, 2020/02/25 02:47:35
Line 1 
Line 1 
 /* $OpenXM: OpenXM/src/asir-contrib/packages/src/os_muldif.rr,v 1.7 2017/05/04 11:12:16 takayama Exp $ */  /* $OpenXM: OpenXM/src/asir-contrib/packages/src/os_muldif.rr,v 1.57 2020/02/21 05:36:17 takayama Exp $ */
 /* The latest version will be at ftp://akagi.ms.u-tokyo.ac.jp/pub/math/muldif  /* The latest version will be at https://www.ms.u-tokyo.ac.jp/~oshima/index-j.html
  scp os_muldif.[dp]* ${USER}@lemon.math.kobe-u.ac.jp:/home/web/OpenXM/Current/doc/other-docs   scp os_muldif.[dp]* ${USER}@lemon.math.kobe-u.ac.jp:/home/web/OpenXM/Current/doc/other-docs
 */  */
 #define USEMODULE 1  #define USEMODULE 1
 /* #undef USEMODULE */  /* #undef USEMODULE */
   
 /*             os_muldif.rr (Library for Risa/Asir)  /*             os_muldif.rr (Library for Risa/Asir)
  *          Toshio Oshima (Nov. 2007 - Apr. 2017)   *          Toshio Oshima (Nov. 2007 - Feb. 2020)
  *   *
  *   For polynomials and differential operators with coefficients   *   For polynomials and differential operators with coefficients
  *   in rational funtions (See os_muldif.pdf)   *   in rational funtions (See os_muldif.pdf)
Line 22  static Muldif.rr$
Line 22  static Muldif.rr$
 static TeXEq$  static TeXEq$
 static TeXLim$  static TeXLim$
 static DIROUT$  static DIROUT$
   static DIROUTD$
 static DVIOUTL$  static DVIOUTL$
 static DVIOUTA$  static DVIOUTA$
 static DVIOUTB$  static DVIOUTB$
Line 42  static Canvas$
Line 43  static Canvas$
 static ID_PLOT$  static ID_PLOT$
 static Rand$  static Rand$
 static LQS$  static LQS$
   static SVORG$
 localf spType2$  localf spType2$
 localf erno$  localf erno$
 localf chkfun$  localf chkfun$
Line 51  localf makenewv$
Line 53  localf makenewv$
 localf vweyl$  localf vweyl$
 localf mycat$  localf mycat$
 localf mycat0$  localf mycat0$
   localf fcat$
 localf findin$  localf findin$
 localf countin$  localf countin$
 localf mycoef$  localf mycoef$
 localf mydiff$  localf mydiff$
 localf myediff$  localf myediff$
   localf mypdiff$
   localf pTaylor$
   localf pwTaylor$
 localf m2l$  localf m2l$
 localf m2ll$  localf m2ll$
 localf mydeg$  localf mydeg$
Line 114  localf lchange$
Line 120  localf lchange$
 localf llsize$  localf llsize$
 localf llbase$  localf llbase$
 localf lsort$  localf lsort$
   localf rsort$
   localf lpair$
 localf lmax$  localf lmax$
 localf lmin$  localf lmin$
 localf lgcd$  localf lgcd$
Line 142  localf getel$
Line 150  localf getel$
 localf ptol$  localf ptol$
 localf rmul$  localf rmul$
 localf mtransbys$  localf mtransbys$
   localf trcolor$
 localf drawopt$  localf drawopt$
 localf execdraw$  localf execdraw$
 localf execproc$  localf execproc$
Line 166  localf myasin$
Line 175  localf myasin$
 localf myacos$  localf myacos$
 localf myatan$  localf myatan$
 localf mylog$  localf mylog$
   localf nlog$
 localf mypow$  localf mypow$
   localf scale$
 localf arg$  localf arg$
 localf sqrt$  localf sqrt$
 localf gamma$  localf gamma$
Line 178  localf eta$
Line 189  localf eta$
 localf jell$  localf jell$
 localf frac$  localf frac$
 localf erfc$  localf erfc$
   localf orthpoly$
   localf schurpoly$
 localf fouriers$  localf fouriers$
 localf todf$  localf todf$
 localf f2df$  localf f2df$
Line 235  localf expower$
Line 248  localf expower$
 localf seriesHG$  localf seriesHG$
 localf seriesMc$  localf seriesMc$
 localf seriesTaylor$  localf seriesTaylor$
   localf mulpolyMod$
   localf solveEq$
   localf baseODE$
   localf taylorODE$
 localf evalred$  localf evalred$
 localf toeul$  localf toeul$
 localf fromeul$  localf fromeul$
Line 248  localf expat$
Line 265  localf expat$
 localf polbyroot$  localf polbyroot$
 localf polbyvalue$  localf polbyvalue$
 localf pcoef$  localf pcoef$
   localf pmaj$
 localf prehombf$  localf prehombf$
 localf prehombfold$  localf prehombfold$
 localf sub3e$  localf sub3e$
Line 269  localf okuboetos$
Line 287  localf okuboetos$
 localf heun$  localf heun$
 localf fspt$  localf fspt$
 localf abs$  localf abs$
   localf sgn$
 localf calc$  localf calc$
 localf isint$  localf isint$
 localf israt$  localf israt$
Line 284  localf iscoef$
Line 303  localf iscoef$
 localf iscombox$  localf iscombox$
 localf sproot$  localf sproot$
 localf spgen$  localf spgen$
   localf spbasic$
 localf chkspt$  localf chkspt$
 localf cterm$  localf cterm$
 localf terms$  localf terms$
Line 293  localf cutgrs$
Line 313  localf cutgrs$
 localf mcgrs$  localf mcgrs$
 localf mc2grs$  localf mc2grs$
 localf mcmgrs$  localf mcmgrs$
   localf spslm$
 localf anal2sp$  localf anal2sp$
 localf delopt$  localf delopt$
 localf str_char$  localf str_char$
Line 323  localf divmattex$
Line 344  localf divmattex$
 localf dviout0$  localf dviout0$
 localf myhelp$  localf myhelp$
 localf isMs$  localf isMs$
   localf getline$
 localf showbyshell$  localf showbyshell$
 localf readcsv$  localf readcsv$
 localf tocsv$  localf tocsv$
Line 338  localf texsp$
Line 360  localf texsp$
 localf getbygrs$  localf getbygrs$
 localf mcop$  localf mcop$
 localf shiftop$  localf shiftop$
   localf shiftPfaff;
 localf conf1sp$  localf conf1sp$
   localf confexp$
   localf confspt$
   localf mcvm$
   localf s2csp$
   localf partspt$
 localf pgen$  localf pgen$
 localf diagm$  localf diagm$
 localf mgen$  localf mgen$
Line 355  localf fimag$
Line 383  localf fimag$
 localf trig2exp$  localf trig2exp$
 localf intpoly$  localf intpoly$
 localf integrate$  localf integrate$
   localf rungeKutta$
 localf simplog$  localf simplog$
 localf fshorter$  localf fshorter$
 localf isshortneg$  localf isshortneg$
Line 371  localf primroot$
Line 400  localf primroot$
 localf varargs$  localf varargs$
 localf ptype$  localf ptype$
 localf pfargs$  localf pfargs$
   localf regress$
 localf average$  localf average$
   localf tobig$
 localf sint$  localf sint$
 localf frac2n$  localf frac2n$
   localf openGlib$
 localf xyproc$  localf xyproc$
 localf xypos$  localf xypos$
 localf xyput$  localf xyput$
Line 394  localf periodicf$
Line 426  localf periodicf$
 localf cmpf$  localf cmpf$
 localf areabezier$  localf areabezier$
 localf saveproc$  localf saveproc$
   localf xyplot$
 localf xygraph$  localf xygraph$
 localf xy2graph$  localf xy2graph$
   localf addIL$
   localf xy2curve$
   localf xygrid$
 localf xyarrow$  localf xyarrow$
 localf xyarrows$  localf xyarrows$
 localf xyang$  localf xyang$
 localf xyoval$  localf xyoval$
   localf xypoch$
 localf ptcommon$  localf ptcommon$
 localf ptcopy$  localf ptcopy$
 localf ptaffine$  localf ptaffine$
Line 423  extern Muldif.rr$
Line 460  extern Muldif.rr$
 extern TeXEq$  extern TeXEq$
 extern TeXLim$  extern TeXLim$
 extern DIROUT$  extern DIROUT$
   extern DIROUTD$
 extern DVIOUTL$  extern DVIOUTL$
 extern DVIOUTA$  extern DVIOUTA$
 extern DVIOUTB$  extern DVIOUTB$
Line 444  extern Canvas$
Line 482  extern Canvas$
 extern ID_PLOT$  extern ID_PLOT$
 extern Rand$  extern Rand$
 extern LQS$  extern LQS$
   extern SV=SVORG$
 #endif  #endif
 static S_Fc,S_Dc,S_Ic,S_Ec,S_EC,S_Lc$  static S_Fc,S_Dc,S_Ic,S_Ec,S_EC,S_Lc$
 static S_FDot;  static S_FDot$
 extern AMSTeX$  extern AMSTeX$
 Muldif.rr="00170510"$  extern Glib_math_coordinate$
   extern Glib_canvas_x$
   extern Glib_canvas_y$
   Muldif.rr="00200223"$
 AMSTeX=1$  AMSTeX=1$
 TeXEq=5$  TeXEq=5$
 TeXLim=80$  TeXLim=80$
Line 467  LCOPT=["red","green","blue","yellow","cyan","magenta",
Line 509  LCOPT=["red","green","blue","yellow","cyan","magenta",
 COLOPT=[0xff,0xff00,0xff0000,0xffff,0xffff00,0xff00ff,0,0xffffff,0xc0c0c0]$  COLOPT=[0xff,0xff00,0xff0000,0xffff,0xffff00,0xff00ff,0,0xffffff,0xc0c0c0]$
 LPOPT=["above","below","left","right"]$  LPOPT=["above","below","left","right"]$
 LFOPT=["very thin","thin","dotted","dashed"]$  LFOPT=["very thin","thin","dotted","dashed"]$
   SVORG=["x","y","z","w","u","v","p","q","r","s"]$
 Canvas=[400,400]$  Canvas=[400,400]$
 LQS=[[1,0]]$  LQS=[[1,0]]$
   
Line 547  def makenewv(L)
Line 590  def makenewv(L)
         if((V=getopt(var))<2) V="z_";          if((V=getopt(var))<2) V="z_";
         else if(isvar(V)) V=rtostr(V);          else if(isvar(V)) V=rtostr(V);
         if(type(N=getopt(num))!=1) N=0;          if(type(N=getopt(num))!=1) N=0;
         Var=vars(L);          Var=varargs(L|all=2);
         for(Va=Var;Va!=[];Va=cdr(Va))  
                 if(vtype(car(Va))==2) Var=append(vars(args(car(Va))),Var);  
         for(XX=[],I=J=0;;I++){          for(XX=[],I=J=0;;I++){
                 X=strtov(V+rtostr(I));                  X=strtov(V+rtostr(I));
                 if(findin(X,Var)<0){                  if(findin(X,Var)<0){
Line 608  def mycat(L)
Line 649  def mycat(L)
                 Do = 1;                  Do = 1;
         }          }
         if(CR) print("");          if(CR) print("");
           else print("",2);
 }  }
   
   def fcat(S,X)
   {
           if(type(S)!=7){
                   if(type(DIROUTD)!=7){
                           DIROUTD=str_subst(DIROUT,["%HOME%","%ASIRROOT%","\\"],
                                   [getenv("HOME"),get_rootdir(),"/"])+"/";
                           if(isMs()) DIROUTD=str_subst(DIROUTD,"/","\\"|sjis=1);
                   }
                   if(S==-1) return;
                   T="fcat";
                   if(S>=2&&S<=9) T+=rtostr(S);
                   T=DIROUTD+T+".txt";
                   if(S==-1) return T;
                   if(S!=0&&access(T)) remove_file(T);
                   S=T;
           }
           R=output(S);
           print(X);
           output();
           if(getopt(exe)==1) shell("\""+S+"\"");
           return R;
   }
   
 def mycat0(L,T)  def mycat0(L,T)
 {  {
         Opt = getopt(delim);          Opt = getopt(delim);
         Del = (type(Opt) >= 0)?Opt:"";          Del = (type(Opt) >= 0)?Opt:"";
           if(type(L)!=4) L=[L];
         while(L != []){          while(L != []){
                 if(Do==1)                  if(Do==1)
                         print(Del,0);                          print(Del,0);
Line 622  def mycat0(L,T)
Line 688  def mycat0(L,T)
                 Do = 1;                  Do = 1;
         }          }
         if(T) print("");          if(T) print("");
           else print("",2);
 }  }
   
 def findin(M,L)  def findin(M,L)
Line 639  def findin(M,L)
Line 706  def findin(M,L)
   
 def countin(S,M,L)  def countin(S,M,L)
 {  {
         if(((Step=getopt(step))==1)||Step==-1){          Step=getopt(step);
           if(type(Step)==1){
                   N=(Step>0)?Step:-Step;
                 if(type(L)==5) L=vtol(L);                  if(type(L)==5) L=vtol(L);
                 L=qsort(L);                  L=qsort(L);
                 while(car(L)<S&&L!=[]) L=cdr(L);                  while(car(L)<S&&L!=[]) L=cdr(L);
                 S+=M;                  S+=M;
                 for(R=[],C=0;L!=[];){                  for(R=[],C=I=0;L!=[];){
                         if(car(L)<S||(Step==1&&car(L)==S)){                          if(car(L)<S||(Step>0&&car(L)==S)){
                                 C++;                                  C++;
                                 L=cdr(L);                                  L=cdr(L);
                         }else{                          }else{
                                 R=cons(C,R);C=0;S+=M;                                  R=cons(C,R);C=0;S+=M;
                                   if(N>1&&++I>=N) break;
                         }                          }
                 }                  }
                 if(C>0) R=cons(C,R);                  if(C>0) R=cons(C,R);
                   if(N>1&&(N-=length(R))>0) while(N-->0) R=cons(0,R);
                 return reverse(R);                  return reverse(R);
         }          }
         if(type(L)==4){          if(type(L)==4){
Line 705  def mydiff(P,X)
Line 776  def mydiff(P,X)
                 for(;X!=[];X=cdr(X)) P=mydiff(P,car(X));                  for(;X!=[];X=cdr(X)) P=mydiff(P,car(X));
                 return P;                  return P;
         }          }
         if(deg(dn(P),X) == 0)          if(ptype(dn(P),X)<2)
                 return red(diff(nm(P),X)/dn(P));                  return red(diff(nm(P),X)/dn(P));
         return red(diff(P,X));          return red(diff(P,X));
 }  }
Line 727  def myediff(P,X)
Line 798  def myediff(P,X)
         return red(X*diff(P,X));          return red(X*diff(P,X));
 }  }
   
   def mypdiff(P,L)
   {
           if(type(P)>3) return map(os_md.mypdiff,P,L);
           for(Q=0;L!=[];L=cdr(L)){
                   Q+=mydiff(P,car(L))*L[1];
                   L=cdr(L);
           }
           return red(Q);
   }
   
   def pTaylor(S,X,N)
   {
           if(!isvar(T=getopt(time))) T=t;
           if(type(S)<4) S=[S];
           if(type(X)<4) X=[X];
           if(findin(T,varargs(S|all=2))>=0){
                   S=cons(z_z,S);X=cons(z_z,X);FT=1;
           }else FT=0;
           LS=length(S);
           FR=(getopt(raw)==1)?1:0;
           if(!FR) R=newvect(LS);
           else R=R1=[];
           for(L=[],I=0,TS=S,TX=X;I<LS;I++,TS=cdr(TS),TX=cdr(TX)){
                   if(!FR) R[I]=car(TX)+car(TS)*T;
                   else{
                           R=cons(car(TX),R);R1=cons(car(TS),R1);
                   }
                   L=cons(car(TS),cons(car(TX),L));
           }
           L=reverse(L);
           if(FR) R=[reverse(R1),reverse(R)];
           for(K=M=1;N>1;N--){
                   S=mypdiff(S,L);
                   K*=++M;
                   for(TS=S,I=0,R1=[];TS!=[];TS=cdr(TS),I++){
                           if(!FR) R[I]+=car(TS)*t^M/K;
                           else R1=cons(car(TS)/K,R1);
                   }
                   if(FR) R=cons(reverse(R1),R);
           }
           if(FT){
                   if(!FR){
                           S=newvect(LS-1);
                           for(I=1;I<LS;I++) S[I-1]=R[I];
                   }else{
                           for(S=[];R!=[];R=cdr(R)){
                                   S=cons(cdr(car(R)),S);
                           }
                           R=S;
                   }
                   R=subst(S,z_z,0);
           }
           return (FR&&!FT)?reverse(R):R;
   }
   
 def m2l(M)  def m2l(M)
 {  {
         if(type(M) < 4)          if(type(M) < 4)
Line 749  def m2l(M)
Line 875  def m2l(M)
   
 def mydeg(P,X)  def mydeg(P,X)
 {  {
         if(type(P) < 3)          if(type(P) < 3 && type(X)==2)
                 return deg(P,X);                  return deg(P,X);
         II = -1;          II=(type(X)==4)?-100000:-1;
         Opt = getopt(opt);          Opt = getopt(opt);
         if(type(P) >= 4){          if(type(P) >= 4){
                 S=(type(P) == 6)?size(P)[0]:0;                  S=(type(P) == 6)?size(P)[0]:0;
                 P = m2l(P);                  P = m2l(P);
                 for(I = 0, Deg = -3; P != []; P = cdr(P), I++){                  for(I = 0, Deg = -100000; P != []; P = cdr(P), I++){
                         if( (DT = mydeg(car(P),X)) == -2)                          if( (DT = mydeg(car(P),X)) == -2&&type(X)!=4)
                                 return -2;                                  return -2;
                         if(DT > Deg){                          if(DT > Deg){
                                 Deg = DT;                                  Deg = DT;
Line 767  def mydeg(P,X)
Line 893  def mydeg(P,X)
                 return (Opt==1)?([Deg,(S==0)?II:[idiv(II,S),irem(II,S)]]):Deg;                  return (Opt==1)?([Deg,(S==0)?II:[idiv(II,S),irem(II,S)]]):Deg;
         }          }
         P = red(P);          P = red(P);
         if(deg(dn(P),X) == 0)          if(type(X)==2){
                 return deg(nm(P),X);                  if(deg(dn(P),X) == 0)
                           return deg(nm(P),X);
           }else{
                   P=nm(red(P));
                   for(D=-100000,I=deg(P,X[1]);I>=0;I--){
                           if(TP=mycoef(P,I,X[1])){
                                   TD=mydeg(TP,X[0])-I;
                                   if(D<TD) D=TD;
                           }
                   }
                   return D;
           }
         return -2;          return -2;
 }  }
   
Line 845  def mulsubst(F,L)
Line 982  def mulsubst(F,L)
         if(N == 0)          if(N == 0)
                 return F;                  return F;
         if(type(L[0])!=4)       L=[L];          if(type(L[0])!=4)       L=[L];
           if(getopt(lpair)==1||(type(L[0])==4&&length(L[0])>2)) L=lpair(L[0],L[1]);
         if(getopt(inv)==1){          if(getopt(inv)==1){
                 for(R=[];L!=[];L=cdr(L)) R=cons([car(L)[1],car(L)[0]],R);                  for(R=[];L!=[];L=cdr(L)) R=cons([car(L)[1],car(L)[0]],R);
                 L=reverse(R);                  L=reverse(R);
Line 1161  def vprod(V1,V2)
Line 1299  def vprod(V1,V2)
 def dnorm(V)  def dnorm(V)
 {  {
         if(type(V)<2) return dabs(V);          if(type(V)<2) return dabs(V);
           if((M=getopt(max))==1||M==2){
                   if(type(V)==5) V=vtol(V);
                   for(S=0;V!=[];V=cdr(V)){
                           if(M==2) S+=dabs(car(V));
                           else{
                                   if((T=dabs(car(V)))>S) S=T;
                           }
                   }
                   return S;
           }
         R=0;          R=0;
         if(type(V)!=4)          if(type(V)!=4)
                 for (I = length(V)-1; I >= 0; I--) R+= V[I]^2;                  for (I = length(V)-1; I >= 0; I--) R+= real(V[I])^2+imag(V[I])^2;
         else{          else{
                 if(type(V[0])>3){                  if(type(V[0])>3){
                         V=ltov(V[0])-ltov(V[1]);                          V=ltov(V[0])-ltov(V[1]);
                         return dnorm(V);                          return dnorm(V);
                 }                  }
                 for(;V!=[]; V=cdr(V))   R+=car(V)^2;                  for(;V!=[]; V=cdr(V))   R+=real(car(V))^2+imag(car(V))^2;
         }          }
         return dsqrt(R);          return dsqrt(R);
 }  }
Line 1219  def mulseries(V1,V2)
Line 1367  def mulseries(V1,V2)
         return VV;          return VV;
 }  }
   
   def scale(L)
   {
           T=F=0;LS=1;
           Pr=getopt(prec);
           Inv=getopt(inv);
           Log10=dlog(10);
           if(type(L)==7){
                   V=findin(L,["CI","DI","CIF","CIF'","DIF","DIF'","SI","TI1","TI2","STI"]);
                   if(V>=0){
                           L=["C","D","CF","CF'","DF","DF'","S","T1","T2","ST"];
                           Inv=1;L=L[V];
                   }
                   V=findin(L,["C","A","K","CF","CF'","S","T1","T2","ST","LL0","LL1","LL2","LL3","LL00",
                           "LL01","LL02","LL03"])+1;
                   if(V==0) V=findin(L,["D","B","K","DF","DF'"])+1;
                   if(V>0) L=V;
           }
           if(type(OL=L)!=4){
                   if(L==2){
                           L=(Pr==0)?
                             [[[1,2,1/20],[2,5,1/10],[5,10,1/5], [10,20,1/2],[20,50,1],[50,100,2]],
                      [[1,2,1/10],[2,5,1/2], [10,20,1],[20,50,5]],
                      [[1,2,1/2],[2,10,1], [10,20,5],[20,100,10]]]:
                              [[[1,2,1/50],[2,5,1/20],[5,10,1/10], [10,20,1/5],[20,50,1/2],[50,100,1]],
                                   [[1,5,1/10],[5,10,1/2], [10,20,1],[50,100,5]],
                                   [[1,5,1/2],[5,10,1], [10,50,5],[50,100,10]]];
                           LS=2;M2=[[1,10,1],[10,100,10]];
                   }else if(L==3){
                           L=(Pr==0)?
                             [[[1,2,1/20],[2,5,1/10],[5,10,1/5], [10,20,1/2],[20,50,1],[50,100,2],
                     [100,200,5],[200,500,10],[500,1000,20]],
                           [[1,2,1/10],[2,5,1/2], [10,20,1],[20,50,5], [100,200,10],[200,500,50]],
                           [[1,2,1/2],[2,10,1], [10,20,5],[20,100,10], [100,200,50],[200,1000,100]]]:
                             [[[1,2,1/50],[2,5,1/20],[5,10,1/10],[10,20,1/5],[20,50,1/2],[50,100,1],
                                 [100,200,2],[200,500,5],[500,1000,10]],
                                   [[1,5,1/10],[5,10,1/2], [10,50,1],[50,100,5], [100,500,10],[500,1000,50]],
                                   [[1,5,1/2],[5,10,1],[10,50,5],[50,100,10], [100,500,50],[500,1000,100]]];
                           LS=3;M2=[[1,5,1],[10,50,10],[100,500,100],[500,1000,500]];
                   }else if(L>9&&L<18){
                           if(L<18){       /* LL0 - LL3, LL00 - LL03 */
                                   if(L==10){
                                           L=[ [[1.001,1.002,0.00001],[1.002,1.005,0.00002],[1.005,1.0105,0.00005]],
                                                   [[1.001,1.002,0.00005],[1.002,1.005,0.0001], [1.005,1.0105,0.0001]],
                                                   [[1.001,1.002,0.0001],[1.002,1.005,0.0005], [1.005,1.0105,0.0005]]];
                                           M2=[1.001,1.0015,1.002,1.003,1.004,1.005,1.006,1.007,1.008,1.009,1.01];
                                   }
                                   if(L==11){
                                           L=[ [[1.01,1.02,0.0001],[1.02,1.05,0.0002],[1.05,1.105,0.0005]],
                                                   [[1.01,1.02,0.0005],[1.02,1.05,0.001], [1.05,1.105,0.001]],
                                                   [[1.01,1.02,0.001],[1.02,1.05,0.005], [1.05,1.105,0.005]]];
                                           M2=[1.01,1.015,1.02,1.03,1.04,1.05,1.06,1.07,1.08,1.09,1.10];
                                   }else if(L==12){
                                           L=[ [[1.105,1.2,0.001],[1.2,1.4,0.002],[1.4,1.8,0.005],[1.8,2.5,0.01],
                                                 [2.5,2.72,0.02]],
                                                   [[1.105,1.2,0.005],[1.2,1.4,0.01],[1.4,1.8,0.01],[1.8,2.5,0.05],
                                                 [2.5,2.72,0.1]],
                                                   [[1.105,1.2,0.01],[1.2,1.4,0.05],[1.4,1.8,0.05],[1.8,2.5,0.1],
                                                 [2.5,2.72,0.1]]];
                                           M2=[1.11,1.15,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2.0,2.2,2.5];
                                   }else if(L==13){
                                           L=[ [[2.72,4,0.02],[4,6,0.05],[6,10,0.1],[10,15,0.2],[15,30,0.5],[30,50,1],
                                                    [50,100,2],[100,200,5],[200,400,10],[400,500,20],[500,1000,50],
                                                    [1000,2000,100],[2000,5000,200],[5000,10000,500],[10000,22000,1000]],
                                                   [[2.7,4,0.1],[4,6,0.1],[6,10,0.5],[10,15,1],[15,30,1],[30,50,5],
                                                    [50,100,10],[100,200,10],[200,400,50],[400,500,100],[500,1000,100],
                                                    [1000,2000,500],[2000,5000,1000],[5000,10000,1000],[10000,22000,5000]],
                                                   [[3,4,0.5],[4,6,0.5],[6,10,1],[10,15,5],[15,30,5],[30,50,10],
                                                    [50,100,50],[100,200,50],[200,400,100],[400,500,100],[500,1000,500],
                                                    [1000,2000,1000],[2000,5000,3000],[5000,10000,5000],[10000,22000,10000]]];
                                           M2=[3,4,5,6,7,8,9,10,15,20,30,40,50,100,200,500,1000,2000,5000,10000,20000];
                                   }else if(L==14){
                                           L=[ [[0.998,0.999,0.00001],[0.995,0.998,0.00002],[0.99,0.995,0.00005]],
                                                   [[0.998,0.999,0.00005],[0.995,0.998,0.0001],[0.99,0.995,0.0001]],
                                                   [[0.998,0.999,0.0001],[0.995,0.998,0.0005],[0.99,0.995,0.0005]]];
                                           M2=[0.999,0.9985,0.998,0.997,0.996,0.995,0.994,0.993,0.992,0.991,0.99];
                                   }else if(L==15){
                                           L=[ [[0.98,0.9901,0.0001],[0.95,0.98,0.0002],[0.905,0.95,0.0005]],
                                                   [[0.98,0.99,0.0005],[0.95,0.98,0.001], [0.905,0.95,0.001]],
                                                   [[0.98,0.99,0.001],[0.95,0.98,0.005], [0.91,0.95,0.005]]];
                                           M2=[0.99,0.985,0.98,0.97,0.96,0.95,0.94,0.93,0.92,0.91];
                                   }else if(L==16){
                                           L=[ [[0.8,0.906,0.001],[0.6,0.8,0.002],[0.37,0.6,0.005]],
                                                   [[0.8,0.906,0.005],[0.6,0.8,0.01],[0.37,0.6,0.01]],
                                                   [[0.8,0.9,0.01],[0.6,0.8,0.05],[0.4,0.6,0.05]]];
                                           M2=[0.9,0.85,0.8,0.75,0.7,0.65,0.6,0.55,0.5,0.45,0.4];
                                   }else{
                                           L=[ [[0.05,0.37,0.002],[0.02,0.05,0.001],[0.01,0.02,0.0005],
                                                    [0.005,0.01,0.0002],[0.001,0.005,0.0001],
                                                    [0.0005,0.001,0.00002],[0.0001,0.0005,0.00001],[0.00005,0.0001,0.000002]],
                                                   [[0.05,0.37,0.01],[0.02,0.05,0.002],[0.01,0.02,0.001],
                                                    [0.005,0.01,0.001],[0.001,0.005,0.0002],
                                                    [0.0005,0.001,0.0001],[0.0001,0.0005,0.00002],[0.00005,0.0001,0.00001]],
                                                   [[0.05,0.37,0.05],[0.02,0.05,0.01],[0.01,0.02,0.005],
                                                    [0.005,0.01,0.005],[0.002,0.005,0.001],
                                                    [0.0005,0.001,0.0005],[0.0001,0.0005,0.0001],[0.00005,0.0001,0.00005]]];
                                           M2=[0.3,0.2,0.1,0.05,0.03,0.02,0.01,0.005,0.002,0.001,0.0005,0.0002,0.0001];
                                   }
                           }
                   }else{
                           if(L==6){       /* S */
                                   L=[ [[6-3/12,15,1/12],[15,30,1/6],[30,50,1/3],[50,70,1/2],[70,80,1],[80,90,5]],
                                           [[6-1/6,15,1/6],[15,30,1/2],[30,70,1],[70,80,5],[80,90,10]],
                                           [[6,15,1/2],[15,30,1],[30,70,5],[70,90,10]] ];
                                   M2=[6,7,8,9,10,15,20,30,40,50,60,70,90];
                           }else if(L==7){ /* T1 */
                                   F=log(tan(x*3.1416/180))/Log10+1;
                                   L=[ [[6-1/3,15,1/12],[15,45,1/6]],
                                           [[6-1/3,15,1/6],[15,45,1/2]],
                                           [[6,45,1]] ];
                                   M2=[6,7,8,9,10,15,20,30,40,45];
                           }else if(L==8){ /* T2 */
                                   L=[ [[45,75,1/6],[75,84+1/6,1/12]],
                                           [[45,75,1],[75,84+1/6,1/6]],
                                           [[45,84,1]] ];
                                   M2=[45,50,60,70,75,80,81,82,83,84];
                           }else if(L==9){ /* ST */
                                   L=[ [[35/60,1,1/120],[1,2,1/60],[2,5+9/12,1/30]],
                                           [[35/60,1,1/60],[1,2,1/6],[2,5+9/12,1/6]],
                                           [[40/60,1,1/6],[1,2,1/2],[2,5+9/12,1]] ];
                                   M2=[1,2,3,4,5];
                           }else{
                                   M2=(L==4||L==5)?[[1,2,1/2],[2,9,1]]:[[1,2,1/2],[2,10,1]];
                                   L=(Pr==0)?
                                   [ [[1,2,1/50],[2,5,1/20],[5,10,1/10]],
                                     [[1,5,1/10],[5,10,1/2]],
                                     [[1,5,1/2],[5,10,1]] ]:
                                   [[[1,2,1/100],[2,5,1/50],[5,10,1/20]],
                                 [[1,2,1/20],[2,10,1/10]],
                                 [[1,2,1/10],[2,10,1/2]] ];
                           }
                   }
           }else if(type(L[0])!=4){
                   L=[L];
                   if(length(L)!=3||L[0]+L[2]>L[1]) T=L;
           }
           if(T==0){
                   if(type(L[0][0])!=4) L=[L];
                   for(R=[];L!=[];L=cdr(L)){
                           for(RR=[],LT=car(L);LT!=[];LT=cdr(LT))
                                   for(I=car(LT)[0];I<=car(LT)[1];I+=car(LT)[2]) RR=cons(I,RR);
                   RR=lsort(RR,[],1);
                           R=cons(RR,R);
                   }
                   R=reverse(R);
                   for(T=[];R!=[];R=cdr(R)){
                           if(length(R)>1) T=cons(lsort(R[0],R[1],"setminus"),T);
                           else T=cons(R[0],T);
                   }
           }
           V0=dlog(10);
           S0=S1=1;D0=D1=0;
           SC=getopt(scale);
           if(type(SC)==4){
                   S0=SC[0];S1=SC[1];
           }else if(type(SC)==1){
                   S0=SC;S1=0;
           }else return T;
           if(type(D=getopt(shift))==4){
                   D0=D[0];D1=D[1];
           }else if(type(D)<2&&type(D)>=0){
                   D0=0;D1=D;
           };
           if(Inv==1){
                   D0+=S0;S0=-S0;
           }
           if(type(TF=getopt(f))>1) F=TF;
           if(F) F=f2df(F);
           if(type(I=getopt(ol))==1&&OL>3) OL=I;
           for(M=M0=[],I=length(T);T!=[];T=cdr(T),I--){
                   for(S=car(T);S!=[];S=cdr(S)){
                           VS=car(S);
                           if(F) V=myfdeval(F,car(S));
                           else if(OL==4) V=frac(dlog(VS)/Log10+0.5);
                           else if(OL==5) V=frac(dlog(VS*3.1416)/Log10);
                           else if(OL>5&&OL<10){
                                   VS=VS*3.1416/180;
                                   if(OL==6) V=dlog(dsin(VS))/Log10+1;
                                   else if(OL==9) V=dlog(VS)/Log10+2;
                                   else V=dlog(dtan(VS))/Log10+8-OL;
                           }
                           else if(OL>9&&OL<14) V=dlog(dlog(VS))/Log10+13-OL;
                           else if(OL>13&&OL<18) V=dlog(-dlog(VS))/Log10+17-OL;
                           else V=dlog(VS)/Log10/LS;
                           V*=S0;
                           if(S1!=0){
                                   M=cons([V+D0,D1],M);
                                   M=cons([V+D0,((length(SC)>2)?SC[I]:(I*S1))+D1],M);
                                   M=cons(0,M);
                           }else M0=cons(V+D0,M0);
                   }
                   if(S1==0) M=cons(reverse(M0),M);
           }
           if(S1!=0) M=cdr(M);
           if(S1==0||getopt(TeX)!=1) return M;
           M=reverse(M);
           if(type(U=getopt(line))==4){
                   if(Inv==1) U=[U[0]+S0,U[1]+S0];
                   M=cons([U[0]+D0,D1],cons([U[1]+D0,D1],cons(0,M)));
           }
           if((VT=getopt(vert))==1){
                   for(N=[];M!=[];M=cdr(M)){
                           if(type(TM=car(M))==4) N=cons([TM[1],TM[0]],N);
                           else N=cons(TM,N);
                   }
                   M=reverse(N);
           }
           if(type(Col=getopt(col))<1) S=xylines(M);
           else S=xylines(M|opt=Col);
           if(type(Mes=getopt(mes))==4){
                   if(length(Mes)==1&&type(M2)==4) Mes=cons(car(Mes),M2);
                   S3=car(Mes);
                   if(type(S3)==4){
                           Col=S3[1];
                           S3=car(S3);
                   }else Col=0;
                   V=car(scale(cdr(Mes)));
                   if(!F) Mes=scale(cdr(Mes)|scale=[S0/LS,0],shift=[D0,D1],ol=OL);
                   else Mes=scale(cdr(Mes)|f=F,scale=[S0,0],shift=[D0,D1]);
                   for(M=car(Mes);M!=[];M=cdr(M),V=cdr(V)){
                           TV=deval(car(V));
                           if(Col!=0) TV=[Col,TV];
                           S+=(VT==1)?xyput([S3+D1,car(M),TV]):xyput([car(M),S3+D1,TV]);
                   }
           }
           if(type(Mes=getopt(mes2))==4){
                   if(type(car(Mes))!=4) Mes=[Mes];
                   for(;Mes!=[];Mes=cdr(Mes)){
                           TM=car(Mes);
                           if(!F) V=scale([car(TM)]|scale=[S0/LS,0],shift=[D0,D1],ol=OL);
                           else V=scale([car(TM)]|f=F,scale=[S0,0],shift=[D0,D1]);
                           V=car(car(V));
                           TM=cdr(TM);
                           if(type(Col=car(TM))==4){
                                   C0=Col[0];C1=Col[1];
                                   if(length(Col)==3){
                                           S+=(VT==1)?xyline([D1+C0,V],[D1+C1,V]|opt=Col[2])
                                                   :xyline([V,D1+C0],[V,D1+C1]|opt=Col[2]);
                                   }else S+=(VT==1)?xyline([D1+C0,V],[D1+C1,V]):xyline([V,D1+C0],[V,D1+C1]);
                           }
                           if(type(TM[1]<2)){
                                   TM=cdr(TM);
                                   S3=car(TM);
                           }
                           S+=(VT==1)?xyput([S3+D1,V,TM[1]]):xyput([V,S3+D1,TM[1]]);
                   }
           }
           return S;
   }
   
 def pluspower(P,V,N,M)  def pluspower(P,V,N,M)
 {  {
         RR = 1;          RR = 1;
Line 1393  def mtoupper(MM, F)
Line 1790  def mtoupper(MM, F)
         if(type(St = getopt(step))!=1) St=0;          if(type(St = getopt(step))!=1) St=0;
         Opt = getopt(opt);          Opt = getopt(opt);
         if(type(Opt)!=1) Opt=0;          if(type(Opt)!=1) Opt=0;
           if(type(Main=getopt(main))!=1) Main=0;
         TeX=getopt(dviout);          TeX=getopt(dviout);
         if(type(Tab=getopt(tab))!=1 && Tab!=0) Tab=2;          if(type(Tab=getopt(tab))!=1 && Tab!=0) Tab=2;
         Line="\\text{line}";          Line="\\text{line}";
Line 1423  def mtoupper(MM, F)
Line 1821  def mtoupper(MM, F)
                         Top+=(TeX)?"\\ ":" ";                          Top+=(TeX)?"\\ ":" ";
         }          }
         PC=IF=1;          PC=IF=1;
           if(Opt>3){
                   for(P=[1],K=0;K<Size[1]-F;K++){
                           for(J=0;J<Size[0];J++)
                                   if(type(dn(M[J][K]))==2) P=cons(dn(M[J][K]),P);
                   }
                   PC=llcm(P|poly=1);
           }
         for(K = JJ = 0; K < Size[1] - F; K++){          for(K = JJ = 0; K < Size[1] - F; K++){
                 for(J = JJ; J < Size[0]; J++){                  for(J = JJ; J < Size[0]; J++){
                         if(M[J][K] != 0){               /* search simpler element */                          if(M[J][K] != 0){               /* search simpler element */
Line 1503  def mtoupper(MM, F)
Line 1908  def mtoupper(MM, F)
                                                                 KRC=-KRC;Sgn=1;                                                                  KRC=-KRC;Sgn=1;
                                                         }else                                                          }else
                                                                 Sgn=0;                                                                  Sgn=0;
                                                         if(St){                                                          if(St&&!Main){
                                                                 if(TeX){                                                                  if(TeX){
                                                                         if(KRC==1)                                                                          if(KRC==1)
                                                                                 Lout=cons([Top+"\\xrightarrow{", Line,KJ0+1,TeXs[Sgn],                                                                                  Lout=cons([Top+"\\xrightarrow{", Line,KJ0+1,TeXs[Sgn],
Line 1528  def mtoupper(MM, F)
Line 1933  def mtoupper(MM, F)
                                 }                                  }
                         /* a parameter Var */                          /* a parameter Var */
                                 Var=0;                                  Var=0;
   /* mycat(["start",J,K]); */
                                 if(St && Opt>4 && length(Var=vars(nm(M[J][K])))==1){                                  if(St && Opt>4 && length(Var=vars(nm(M[J][K])))==1){
                                         J0=J;Jv=mydeg(nm(M[J0][K]),car(Var));                                          J0=J;Jv=mydeg(nm(M[J0][K]),car(Var));
                                         for(I=JJ;I<Size[0]; I++){                                          for(I=JJ;I<Size[0]; I++){
Line 1537  def mtoupper(MM, F)
Line 1943  def mtoupper(MM, F)
                                                 }                                                  }
                                                 if(length(T)>1) continue;                                                  if(length(T)>1) continue;
                                                 if(mydeg(MIK,T[0])<Jv){                                                  if(mydeg(MIK,T[0])<Jv){
                                                         J0=I;Jv=mydeg(MIK);Var=T;       /* search minimal degree */                                                          J0=I;Jv=mydeg(MIK,T[0]);Var=T;  /* search minimal degree */
                                                 }                                                  }
                                         }                                          }
                                         if(length(Var)==1){                                          if(length(Var)==1){
                                                 Var=car(Var);                                                  Var=car(Var);
                                                 Q=nm(M[J0][K]);                                                  Q=nm(M[J0][K]);
   /* mycat(["min",Q,M[J0][K],"J0=",J0,"J=",J,"JJ=",JJ,K,M]); */
   J=J0;
                                                 for(I=JJ; I<Size[0]; I++){                                                  for(I=JJ; I<Size[0]; I++){
                                                         if(I==J0 || mydeg(nm(M[I][K]),Var)<0) continue;                                                          if(I==J0 || mydeg(nm(M[I][K]),Var)<0) continue;
                                                         T=rpdiv(nm(M[I][K]),Q,Var);                                                          T=rpdiv(nm(M[I][K]),Q,Var);
Line 1553  def mtoupper(MM, F)
Line 1961  def mtoupper(MM, F)
                                 if(type(Var)==2){ /* 1 variable */                                  if(type(Var)==2){ /* 1 variable */
                                         if(I==Size[0]){                                          if(I==Size[0]){
                                                 for(QF=0,Q0=1,QR=getroot(Q,Var|mult=1);QR!=[];QR=cdr(QR)){                                                  for(QF=0,Q0=1,QR=getroot(Q,Var|mult=1);QR!=[];QR=cdr(QR)){
   /* mycat(["root",Q,QR,PC]); */
                                                         if(deg(T=QR[0][1],Var)>0){                                                          if(deg(T=QR[0][1],Var)>0){
                                                                 QF=1;Q0*=T; continue;                                                                  QF=1;Q0*=T; continue;
                                                         }                                                          }
Line 1563  def mtoupper(MM, F)
Line 1972  def mtoupper(MM, F)
                                                                 if(TeX){                                                                  if(TeX){
                                                                         Lout=cons(["\\hspace{",Tab*(St-3)-1,"mm}\\text{If }",                                                                          Lout=cons(["\\hspace{",Tab*(St-3)-1,"mm}\\text{If }",
                                                                                 Var,"=",T,","] ,Lout);                                                                                  Var,"=",T,","] ,Lout);
                                                                         Lout=append(mtoupper(M0,F|step=St+1,opt=Opt,dviout=-2,tab=Tab),Lout);                                                                          Lout=append(mtoupper(M0,F|step=St+1,opt=Opt,dviout=-2,tab=Tab,main=Main),Lout);
                                                                 }else{                                                                  }else{
                                                                         mycat([str_times(" ",St-1)+"If",Var,"=",T,","]);                                                                          mycat([str_times(" ",St-1)+"If",Var,"=",T,","]);
                                                                         mtoupper(M0,F|step=St+1,opt=Opt);                                                                          mtoupper(M0,F|step=St+1,opt=Opt,main=Main);
                                                                 }                                                                  }
                                                         }                                                          }
                                                 }                                                  }
Line 1583  def mtoupper(MM, F)
Line 1992  def mtoupper(MM, F)
                                                 KRC=-red((T[2]*dn(M[J0][K]))/(T[1]*dn(M[I][K])));                                                  KRC=-red((T[2]*dn(M[J0][K]))/(T[1]*dn(M[I][K])));
                                                 for(II=K;II<Size[1];II++)                                                  for(II=K;II<Size[1];II++)
                                                         M[I][II]=radd(M[I][II],rmul(M[J0][II],KRC));                                                          M[I][II]=radd(M[I][II],rmul(M[J0][II],KRC));
                                                 if(TeX)                                                  if(!Main){
                                                         Lout=cons([Top+"\\xrightarrow{", Line,I+1,"\\ +=\\ ",Line,                                                          if(TeX)
                                                                 J0+1,"\\times\\left(",KRC,"\\right)}",dupmat(M)],Lout);                                                                  Lout=cons([Top+"\\xrightarrow{", Line,I+1,"\\ +=\\ ",Line,
                                                 else                                                                          J0+1,"\\times\\left(",KRC,"\\right)}",dupmat(M)],Lout);
                                                         mycat([Top+"line",I+1,"+=",Line,J0+1," * (",KRC,")\n",M,"\n"]);                                                          else
                                                                   mycat([Top+"line",I+1,"+=",Line,J0+1," * (",KRC,")\n",M,"\n"]);
                                                   }
                                                 J=JJ-1;                                                  J=JJ-1;
                                                 continue;                                                  continue;
                                         }                                          }
Line 1622  def mtoupper(MM, F)
Line 2033  def mtoupper(MM, F)
                                                                         if(TeX){                                                                          if(TeX){
                                                                                 Lout=cons(["\\hspace{",Tab*(St-3)-1,"mm}\\text{If }",                                                                                  Lout=cons(["\\hspace{",Tab*(St-3)-1,"mm}\\text{If }",
                                                                                         X,"=",T,","] ,Lout);                                                                                          X,"=",T,","] ,Lout);
                                                                                 Lout=append(mtoupper(M0,F|step=St+1,opt=Opt,dviout=-2,tab=Tab),                                                                                  Lout=append(mtoupper(M0,F|step=St+1,opt=Opt,dviout=-2,tab=Tab,main=Main),
                                                                                         Lout);                                                                                          Lout);
                                                                         }else{                                                                          }else{
                                                                                 mycat([str_times(" ",St-1)+"If",X,"=",T,","]);                                                                                  mycat([str_times(" ",St-1)+"If",X,"=",T,","]);
                                                                                 mtoupper(M0,F|step=St+1,opt=Opt);                                                                                  mtoupper(M0,F|step=St+1,opt=Opt,main=Main);
                                                                         }                                                                          }
                                                                         break;                                                                          break;
                                                                 }                                                                  }
Line 1653  def mtoupper(MM, F)
Line 2064  def mtoupper(MM, F)
                                                                                         if(TeX){                                                                                          if(TeX){
                                                                                                 Lout=cons(["\\hspace{",Tab*(St-3)-1,"mm}\\text{If }",                                                                                                  Lout=cons(["\\hspace{",Tab*(St-3)-1,"mm}\\text{If }",
                                                                                                         X0,"=",T0,","] ,Lout);                                                                                                          X0,"=",T0,","] ,Lout);
                                                                                                 Lout=append(mtoupper(M0,F|step=St+1,opt=Opt,dviout=-2,tab=Tab),                                                                                                  Lout=append(mtoupper(M0,F|step=St+1,opt=Opt,dviout=-2,tab=Tab,main=Main),
                                                                                                         Lout);                                                                                                          Lout);
                                                                                         }else{                                                                                          }else{
                                                                                                 mycat([str_times(" ",St-1)+"If",X0,"=",T0,","]);                                                                                                  mycat([str_times(" ",St-1)+"If",X0,"=",T0,","]);
                                                                                                 mtoupper(M0,F|step=St+1,opt=Opt);                                                                                                  mtoupper(M0,F|step=St+1,opt=Opt,main=Main);
                                                                                         }                                                                                          }
                                                                                 }                                                                                  }
   
Line 1701  def mtoupper(MM, F)
Line 2112  def mtoupper(MM, F)
                                                 for(I = K+1; I < Size[1]; I++)                                                  for(I = K+1; I < Size[1]; I++)
                                                         M[J][I] = radd(M[J][I],rmul(M[JJ][I],Mul));                                                          M[J][I] = radd(M[J][I],rmul(M[JJ][I],Mul));
                                                 M[J][K] = 0;                                                  M[J][K] = 0;
                                                 if(St){                                                  if(St&&!Main){
                                                         if(Mul<0){                                                          if(Mul<0){
                                                                 Mul=-Mul;Sgn=0;                                                                  Mul=-Mul;Sgn=0;
                                                         }else   Sgn=1;                                                          }else   Sgn=1;
Line 1930  def vgen(V,W,S)
Line 2341  def vgen(V,W,S)
 def mmc(M,X)  def mmc(M,X)
 {  {
         Mt=getopt(mult);          Mt=getopt(mult);
         if(type(M)==7) M=os_md.s2sp(M);          if(type(M)==7) M=s2sp(M);
         if(type(M)!=4||type(M[0])!=6) return 0;          if(type(M)!=4) return 0;
           if(type(M[0])<=3){
                   for(RR=[];M!=[];M=cdr(M)) RR=cons(mat([car(M)]),RR);
                   M=reverse(RR);
           }
         if(type(M[0])!=6){                      /* spectre type -> GRS */          if(type(M[0])!=6){                      /* spectre type -> GRS */
                 G=s2sp(M|std=1);                  G=s2sp(M|std=1);
                 L=length(G);                  L=length(G);
Line 3040  def llbase(VV,L)
Line 3455  def llbase(VV,L)
         return V;          return V;
 }  }
   
   def rsort(L,T,K)
   {
           for(R=[];L!=[];L=cdr(L))
                   R=cons((type(car(L))==4)?rsort(car(L),T-1,K):car(L),R);
           if(T>0||iand(T,iand(K,2)/2)) return reverse(R);
           R=qsort(R);
           return (iand(K,1))? reverse(R):R;
   }
   
   
 def lsort(L1,L2,T)  def lsort(L1,L2,T)
 {  {
           C1=getopt(c1);C2=getopt(c2);
         if(type(T)==4){          if(type(T)==4){
                 K=T;                  K=T;
                 T=K[0];                  if(length(T)>0){
                 K=cdr(K);                          T=K[0];
                           K=cdr(K);
                   }else T=0;
         }else K=0;          }else K=0;
         if(type(T)==7)          if(type(TT=T)==7)
                 T = findin(T,["cup","setminus","cap","reduce"]);                  T = findin(T,["cup","setminus","cap","reduce","sum","subst"]);
         if(K){          if(type(L2)==7&&T<0)
                 C1=getopt(c1);C2=getopt(c2);                  T=findin(TT,["put","get","sub"]);
                 KN=K[0];          if(K){           /* [[..],..] */
                 if(L2==[]){ /* sort or deduce duplication */                  if(K!=[]) KN=K[0];
                         if((T!=0&&T!=1)||length(K)!=1) return L1;                  if(L2==[]||L2=="sort"){ /* sort or deduce duplication */
                           if((T!=0&&T!=3)||length(K)!=1) return L1;
                         if(KN<0){                          if(KN<0){
                                 KN=-KN-1;                                  KN=-KN-1;
                                 F=-1;                                  F=-1;
                         }else F=1;                          }else F=1;
                         L1=msort(L1,[F,0,KN]);                          L1=msort(L1,[F,0,KN]);
                         if(T==1){                          if(T==3){
                                 R=[car(L1)];L1=cdr(L1);                                  R=[car(L1)];L1=cdr(L1);
                                 for(;L1!=[];L1=cdr(L1)){                                  for(;L1!=[];L1=cdr(L1)){
                                         if(car(L1)[KN]!=car(R)[KN]) R=cons(car(L1),R);                                          if(car(L1)[KN]!=car(R)[KN]) R=cons(car(L1),R);
Line 3067  def lsort(L1,L2,T)
Line 3496  def lsort(L1,L2,T)
                                 L1=reverse(R);                                  L1=reverse(R);
                         }                          }
                         return L1;                          return L1;
                 }else if(L2==0&&type(C1)==4){                  }else if((L2==0||L2=="col")&&type(C1)==4){
                         if(T==0||T==1){ /* extract or delete columns */                          if(T==0||T==1){ /* extract or delete columns */
                                 for(R=[];L1!=[];L1=cdr(L1)){                                  for(R=[];L1!=[];L1=cdr(L1)){
                                         if(T==1&&L=[0]){        /* delete top column */                                          if(T==1&&C1==[0]){      /* delete top column */
                                                 R=cons(cdr(car(L1)),R);                                                  R=cons(cdr(car(L1)),R);
                                                 continue;                                                  continue;
                                         }                                          }
                                         LT=car(L1);                                          LT=car(L1);RT=[];
                                         if(T==0){                                          if(T==0){
                                                 for(CT=C1;CT!=[];CT=cdr(CT)){                                                  for(CT=C1;CT!=[];CT=cdr(CT)) RT=cons(LT[car(CT)],RT);
                                                         for(RT=[],CT=C1;CT!=[];CT=cdr(CT))  
                                                                 RT=cons(LT[car(CT)],R);  
                                                 }  
                                         }else{                                          }else{
                                                 for(I=0,RT=[];LT!=[];I++,LT=cdr(LT))                                                  for(I=0;LT!=[];I++,LT=cdr(LT))
                                                         if(findin(I,C1)<0) RT=cons(car(LT),RT);                                                          if(findin(I,C1)<0) RT=cons(car(LT),RT);
                                                 RT=reverse(RT);                                                  RT=reverse(RT);
                                         }                                          }
Line 3089  def lsort(L1,L2,T)
Line 3515  def lsort(L1,L2,T)
                                 }                                  }
                                 return reverse(R);                                  return reverse(R);
                         }                          }
                 }else{                  }else if(type(L2)==1||type(L2)==7){
                           if(L2==1||L2=="num"){
                                   if(T==4) T=3;
                                   I=(length(K)<2)?(-1):K[1];
                                   if(T==0||T==1||T==2||T==3){
                                           S=F=CT=0;R=[];
                                           if(K==[] || type((S=K[0]))==1 || S==0){
                                                   if(T==0||T==1||T==2){
                                                           for(J;L1!=[];L1=cdr(L1),J++){
                                                                   if(T==0) R=cons(cons(J+S,car(L1)),R);
                                                                   else if(T==1){
                                                                           for( ;C1!=[]; C1=cdr(C1))
                                                                                   R=cons(L1[car(C1)],R);
                                                                   }else{
                                                                           if(findin(J,C1)<0) R=cons(car(L1),R);
                                                                   }
                                                           }
                                                           return reverse(R);
                                                   }else if(T==3) return length(L1);
                                           }else{
                                                   if(type(S)==2&&vtype(S)>2) F=1;
                                                   else if(type(S)==4) F=2;
                                                   else if(S=="+") F=3;
                                                   else return L1;
                                           }
                                           for(R=[];L1!=[];L1=cdr(L1)){
                                                   L1T=car(L1);
                                                   if(F==1) V=call(S,(I<0)?L1T:L1T[I]);
                                                   else if(F==2) V=calc((I<0)?L1T:L1T[I],S);
                                                   else if(F==3){
                                                           for(C=C1,V=0;C!=[];C=cdr(C))
                                                                   if(type(X=L1T[car(C)])==1) V+=X;
                                                   }
                                                   if(T==0) R=cons(cons(V,L1T),R);
                                                   else if(T==1){
                                                           if(V) R=cons(L1T,R);
                                                   }else if(T==2){
                                                           if(!V) R=cons(L1T,R);
                                                   }else if(T==3){
                                                           if(F==3) CT+=V;
                                                           else if(V) CT++;
                                                   }
                                           }
                                           return (T==3)?CT:reverse(R);
                                   }else if(TT=="col"){
                                           J=(length(K)>0)?car(K):0;
                                           I=length(car(L1))+J;
                                           for(V=[];I>J;)
                                                   V=cons(--I,V);
                                           return cons(V,L1);
                                   }
                           }else if(L2=="transpose") return mtranspose(L1);
                           else if(L2=="subst"||L2=="adjust"){
                                   Null=(!K)?"":car(K);
                                   if(L2=="adjust") C1=[];
                                   R=lv2m(L1|null="");
                                   for(;C1!=[];C1=cdr(C1)) R[car(C1)[0]][car(C1)[1]]=car(C1)[2];
                                   return m2ll(R);
                           }
                           return L1;
                   }else{           /* [[..],..], [[..],..] */
                           if(type(L2[0])<4){
                                   for(R=[];L2!=[];L2=cdr(L2)) R=cons([car(L2)],R);
                                   L2=reverse(R);
                           }
                           if(TT=="sum")  T=3;
                           if(TT=="over") T=4;
                           if(findin(T,[0,1,2,3,4,5])<0) return L1;
                           if(T==4||T==5){
                                   if(type(C1)<2) C1=[C1];
                                   if(type(C2)<2) C2=[C2];
                           }
                         if(type(car(L2))!=4){                          if(type(car(L2))!=4){
                                 for(R=[];L2!=[];L2=cdr(L2)) R=cons([car(L2)],R);                                  for(R=[];L2!=[];L2=cdr(L2)) R=cons([car(L2)],R);
                                 R=reverse(R);                                  R=reverse(R);
                                 if(length(K)==1) K=[K[0],0];                                  if(length(K)==1) K=[K[0],0];
                         }  
                         for(R=[],I=0;L1!=[];I++,L1=cdr(L1))  
                                 R=cons(cons(I,car(L1)),R);  
                         KN++;  
                         KR=K[1];  
                         K0=K[0];K1=K[1];  
                         L1=msort(R,[],[1,0,K0]);  
                         if(type(C2)==4){  
                                 L2=msort(L2,0,0|c1=cons(K1,C2)); /* extract columns */  
                                 C2=0;                                  C2=0;
                                 K1=0;  
                         }                          }
                         L2=msort(L2,[],[1,0,K1]);                          L1=lsort(L1,"num",["put",0]);           /* insert number */
                         S=size(L2);                          K0=(length(K)>0)?K[0]+1:1;
                         for(R0=[];S>0;S--) R0=cons("",R0);                          K1=(length(K)>1)?K1=K[1]:0;
                         R0=[R0];                          L1=lsort(L1,"sort",[0,K0]);
                         if(T==0||T==1||T==3){ /* cup or setminus or cap */                          if(T<4&&type(C2)==4&&length(L2[0])>1){
                                 for(R=[];L1!=[];L1=cdr(L1)){                                  L2=lsort(L2,"col",["put"]|c1=cons(K1,C2)); /* add key and extract columns */
                                         while(L2!=[]&&car(L1)[K0]>car(L2)[K1]) L2=cdr(L2);                                  C2=0;K1=0;
                                         if(L2==[]||car(L1)[K0]<car(L2)[K1]){                          }
                                                 if(T1==0||T1==1)                          L2=lsort(L2,"sort",[0,K1]);
                                                                 R=cons((T1==0)?append(car(L1),R0):car(L1),R);                          for(R0=[],S=S1=length(L1[0]);S>0;S--) R0=cons("",R0);
                                         }else if(T==1||T1==3) R=cons([car(L1),car(L2)],R);                          for(R1=[],S=length(L2[0]);S>0;S--) R1=cons("",R1);
                           if(!K1&&T!=3) R1=cdr(R1);
                           for(R=[];L1!=[];L1=cdr(L1)){
                                   while(L2!=[]&&car(L1)[K0]>car(L2)[K1]){
                                           if(T==3) R=cons(append(R0,car(L2)),R);
                                           L2=cdr(L2);
                                 }                                  }
                                   if(L2==[]||car(L1)[K0]<car(L2)[K1]){
                                           if(T!=2) R=cons((T==1||T>3||R1==[])?car(L1):append(car(L1),R1),R);
                                   }else if(T==0||T==2||T==3){
                                           if(R0==[]) R=append(car(L1),R);
                                           else R=cons(append(car(L1),(!K1&&T!=3)?cdr(car(L2)):car(L2)),R);
                                           L2=cdr(L2);
                                   }else if(T==4||T==5){
                                           V1=ltov(car(L1));V2=ltov(car(L2));
                                           for(D1=C1,D2=C2;D1!=[];D1=cdr(D1),D2=cdr(D2))
                                                   if((I=V2[car(D2)])!=""||T==4) V1[car(D1)+1]=I;
                                           R=cons(vtol(V1),R);
                                   }
                         }                          }
                         R=msort(R,[],[1,0,K0]);                          if(T==3){
                         R=msort(R,0,[1]|c1=(T1!=0&&T!=3)?[0]:[0,length(L1)+K[1]]);                                  while(L2!=[]){
                         if(type(C1)!=4&&type(C2)!=4) return R;                                          R=cons(append(R0,car(L2)),R);
                         C=[];S0=size(L1);                                          L2=cdr(L2);
                         if(type(C1)==4)  
                                 for(;C1!=[];C1=cdr(C1)) C=cons(car(C1),C);  
                         else for(I=0;I<S0; I++) C=cons(I,C);  
 /*  
                         if(type(C2)==4){  
                                 for(F=0,I=0;C1!=[];I++,C1=cdr(C1)){  
                                         if(I==I) continue;  
                                         if(car(C1)==K1) continue;  
                                         F=(car(C1)>K1)?S0:(S0-1);  
                                         C=cons(car(C1)+F,C);  
                                 }                                  }
                         }else                          }
 */                          R=lsort(R,"sort",["put",0]);    /* original order */
                         for(I=0;I<S-1;I++) C=cons(I+S0,C);                          D=(((T==0||T==2)&&!K1)||T==3)?[0]:[0,S1+K1];
                           R=lsort(R,0,[1]|c1=D); /* delete */
                           if(type(C1)!=4||T==1||T==4||T==5) return R;
                           C=[];S0=size(L1[0]);
                           for(;C1!=[];C1=cdr(C1)) C=cons(car(C1),C);
                           for(I=0;I<S0-S1;I++) C=cons(I+S1,C);
                         C=reverse(C);                          C=reverse(C);
                         return msort(L,0,[0]|c1=C);                          return lsort(R,"col",[1]|c1=C);
                 }                  }
         }          }
         if(L2 == []){          if(L2 == []){           /* [...] */
                 if(T == 2) return L2;                  if(T==8||TT=="count") return [length(L1),length(lsort(L1,[],1))];
                 if(T == 3)      return [L1,L2];                  if(T==7||TT=="cut"){
                           K=length(L1);
                           if(C1<0) C1=K+C1;
                           for(R=[],I=0;I<C1&&L1!=[];I++,L1=cdr(L1))
                                   R=cons(car(L1),R);
                           for(S=[];L1!=[];L1=cdr(L1))
                                   S=cons(car(L1),S);
                           return [reverse(R),reverse(S)];
                   }
                   if(T==2) return L2;
                   if(T==3) return [L1,L2];
                 L1 = ltov(L1); qsort(L1);                  L1 = ltov(L1); qsort(L1);
                 if(T != 1)                  if(T != 1)
                         return vtol(L1);                          return vtol(L1);
Line 3155  def lsort(L1,L2,T)
Line 3664  def lsort(L1,L2,T)
                 }                  }
                 return L3;                  return L3;
         }          }
           if(T==8||TT=="count"){
                   K=length(lsort(L1,L2,3)[0]);
                   R=[length(L2),length(L1)];
                   L1 = lsort(L1,[],1);
                   L2 = lsort(L2,[],1);
                   R=append([length(L2),length(L1)],R);
                   R=cons(length(lsort(L1,L2,2)),R);
                   return reverse(cons(K,R));
           }
           if((T==9||TT=="cons")&&type(car(L1))==4){
                   if(type(L2)!=4) L2=[L2];
                   for(R=[];L1!=[];L1=cdr(L1)){
                           R=cons(cons(car(L2),car(L1)),R);
                           if(length(L2)>1) L2=cdr(L2);
                   }
                   return reverse(R);
           }
           if(T==10||TT=="cmp"){
                   if(length(L1)!=length(L2)){
                           mycat("Different length!");
                           return 1;
                   }
                   R=[];
                   if(type(car(L1))==4){
                           for(U=[],I=0;L1!=[];I++,L1=cdr(L1),L2=cdr(L2)){
                                   if(length(S=car(L1))!=length(T=car(L2))){
                                           mycat(["Different size : line ",I]);
                                           return 0;
                                   }
                                   for(J=0;S!=[];S=cdr(S),T=cdr(T),J++)
                                           if(car(S)!=car(T)) U=cons([[I,J],car(S),car(T)],U);
                           }
                           if(U!=[]) R=cons(reverse(U),R);
                   }else{
                           for(I=0;L1!=[];L1=cdr(L1),L2=cdr(L2),I++)
                                   if(car(L1)!=car(L2)) R=cons([I,car(L1),car(L2)],R);
                   }
                   return reverse(R);
           }
           if(T==11||TT=="append"){
                   if(type(car(L1))!=4) return append(L1,L2);
                   for(R=[];L1!=[];L1=cdr(L1),L2=cdr(L2))
                           R=cons(append(car(L1),car(L2)),R);
                   return reverse(R);
           }
         if(T == 1 || T == 2){          if(T == 1 || T == 2){
                 L1 = lsort(L1,[],1);                  L1 = lsort(L1,[],1);
                 L2 = lsort(L2,[],1);                  L2 = lsort(L2,[],1);
Line 3260  def msort(L,S)
Line 3814  def msort(L,S)
         return qsort(L,os_md.mqsub);          return qsort(L,os_md.mqsub);
 }  }
   
   def lpair(A,B)
   {
           if(B==0){
                   for(S=T=[];A!=[];A=cdr(A)){
                           S=cons(car(A)[0],S);T=cons(car(A)[1],T);
                   }
                   return [reverse(S),reverse(T)];
           }else{
                   for(R=[];A!=[];A=cdr(A),B=cdr(B))
                           R=cons([car(A),car(B)],R);
                   return reverse(R);
           }
   }
   
 def lmax(L)  def lmax(L)
 {  {
         if(type(L)==4){          if(type(L)==4){
Line 3299  def lgcd(L)
Line 3867  def lgcd(L)
         return [];          return [];
 }  }
   
   #if0
 def llcm(L)  def llcm(L)
 {  {
           if(type(L)==5||type(L)==6) L=m2l(L);
           if(type(L)<4) L=[L];
         if(type(L)==4){          if(type(L)==4){
                 F=getopt(poly);                  F=getopt(poly);
                 V=car(L);                  V=car(L);
Line 3314  def llcm(L)
Line 3885  def llcm(L)
                 if(F!=1&&V<0) V=-V;                  if(F!=1&&V<0) V=-V;
                 return V;                  return V;
         }          }
         else if(type(L)==5||type(L)==6)  
                 return llcm(m2l(L)|option_list=getopt());  
         return [];          return [];
 }  }
   #else
   def llcm(R)
   {
           if(type(R)==4||type(R)==5) R=m2l(R);
           if(type(R)<4) R=[R];
           if(type(R)!=4) return 0;
           V=getopt(poly);
           if(type(V)<1){
                   for(L=R;L!=[];L=cdr(L)){
                           if(type(car(L))>1){
                                   V=1; break;
                           }
                   }
           }
           if(getopt(dn)!=1){
                   for(L=[];R!=[];R=cdr(R)) if(R!=0) L=cons(1/car(R),L);
                   R=L;
           }
           P=1;
           if(type(V)<1){
                   for(;R!=[];R=cdr(R)){
                           if(!(TL=car(R))) continue;
                           else P=ilcm(P,dn(TL));
                   }
                   return P;
           }
           for(;R!=[];R=cdr(R)){
                   if(!car(R)) continue;
                   D=dn(red(car(R)));
                   N=red(P/D);
                   if(type(V)<2){
                           if(type(N)!=3) continue;
                           P*=dn(N);
                           continue;
                   }
                   if(ptype(N,V)>2){
                           L=fctr(dn(N));
                           for(;L!=[];L=cdr(L)){
                                   if(ptype(car(L)[0],V)<2) continue;
                                   P*=car(L)[0]^car(L)[1];
                           }
                   }
           }
           return P;
   }
   #fi
   
 def ldev(L,S)  def ldev(L,S)
 {  {
Line 3393  def lnsol(VV,L)
Line 4008  def lnsol(VV,L)
   
 def ladd(X,Y,M)  def ladd(X,Y,M)
 {  {
         if(type(X)==4) X=ltov(X);          if(Y==0){
                   Y=X[1];X=X[0];
           }
         if(type(Y)==4) Y=ltov(Y);          if(type(Y)==4) Y=ltov(Y);
           if(type(X)==4) X=ltov(X);
         return vtol(X+M*Y);          return vtol(X+M*Y);
 }  }
   
 def mrot(X)  def mrot(X)
 {  {
           if(type(X)==4){
                   if(getopt(deg)==1)
                           X=[deval(@pi*X[0]/180),deval(@pi*X[1]/180),deval(@pi*X[2]/180)];
                   if(getopt(conj)==1)
                           return mrot([-X[2],-X[1],0])*mrot([X[0],X[1],X[2]]);
                   if(X[1]==0){
                           X=[X[0]+X[2],0,0];
                           if(X[0]==0) return diagm(3,[1]);
                   }
                   if(X[0]!=0){
                           M=mat([dcos(X[0]),-dsin(X[0]),0],[dsin(X[0]),dcos(X[0]),0],[0,0,1]);
                           if(X[1]==0) return M;
                   }
                   N=mat([dcos(X[1]),0,-dsin(X[1])],[0,1,0],[dsin(X[1]),0,dcos(X[1])]);
                   if(X[0]!=0) N=M*N;
                   if(X[2]==0) return N;
                   return N*mrot([X[2],0,0]);
           }
         if(getopt(deg)==1) X=@pi*X/180;          if(getopt(deg)==1) X=@pi*X/180;
         X=deval(X);          X=deval(X);
         return mat([dcos(X),dsin(X)],[-dsin(X),dcos(X)]);          return mat([dcos(X),-dsin(X)],[dsin(X),dcos(X)]);
 }  }
   
 def m2v(M)  def m2v(M)
Line 4084  def mtransbys(FN,F,LL)
Line 4720  def mtransbys(FN,F,LL)
         return call(FN, cons(F,LL)|option_list=Opt);          return call(FN, cons(F,LL)|option_list=Opt);
 }  }
   
   def trcolor(S)
   {
           if(type(S)!=7) return S;
           return ((I=findin(S,LCOPT))>=0)?COLOPT[I]:0;
   }
   
 def drawopt(S,T)  def drawopt(S,T)
 {  {
         if(type(S)!=7) return -1;          if(type(S)!=7) return -1;
Line 4115  def drawopt(S,T)
Line 4757  def drawopt(S,T)
         return -1;          return -1;
 }  }
   
   def openGlib(W)
   {
           extern Glib_canvas_x;
           extern Glib_canvas_y;
           extern Glib_math_coordinate;
   
           if(W==0){
                   glib_clear();
                   return;
           }
           if(type(W)==4&&length(W)==2){
                   Glib_canvas_x=W[0];
                   Glib_canvax_y=W[1];
           }
           Glib_math_coordinate=1;
           if(getopt(null)!=1) return glib_open();
   }
   
 def execdraw(L,P)  def execdraw(L,P)
 {  {
         if((Proc=getopt(proc))!=1) Proc=0;          if((Proc=getopt(proc))!=1) Proc=0;
Line 4367  def execdraw(L,P)
Line 5027  def execdraw(L,P)
                                                 LOut=cons(T[2],Out);                                                  LOut=cons(T[2],Out);
                                         }                                          }
                                 }                                  }
                           }else if(T[0]==6){      /* plot */
                                   F++;
                                   if((T1=findin(T[1],LCOPT))>-1) T1=COLOPT(T1);
                                   else if(type(T1)!=1 && T1!=0) T1=0xffffff;
                                   for(T2=ptaffine(M,T[2]|option_list=Org);T2!=[];T2=cdr(T2))
                                           draw_obj(Id,Ind,[rint(car(T2)[0]),rint(car(T2)[1])],T1);
                         }else if(Proc==1&&type(T[0])==2){                          }else if(Proc==1&&type(T[0])==2){
                                 if(length(T)<3) call(T[0],T[1]);                                  if(length(T)<3) call(T[0],T[1]);
                                 else call(T[0],T[1]|option_list=T[2]);                                  else call(T[0],T[1]|option_list=T[2]);
Line 4446  def execdraw(L,P)
Line 5112  def execdraw(L,P)
                                         if(P[0]==2)     dviout(T[2]|option_list=T[1]);                                          if(P[0]==2)     dviout(T[2]|option_list=T[1]);
                                         else LOut=cons(T[2],Out);                                          else LOut=cons(T[2],Out);
                                 }                                  }
                           }else if(T[0]==6){      /* plot */
                                   F++;
                                   if(type(T[1])==7) T1=[T[1],"."];
                                   else T1=".";
                                   for(T2=ptaffine(M,T[2]|option_list=Org);T2!=[];T2=cdr(T2))
                                           str_tb(xypos([car(T2)[0],car(T2)[1],T1]),Out);
                         }else if(T[0]==-2)                          }else if(T[0]==-2)
                                 str_tb(["%",T[1],"\n"],Out);                                  str_tb(["%",T[1],"\n"],Out);
                         else if(Proc==1&&type(T[0])==2){                           else if(Proc==1&&type(T[0])==2){
                                 if(length(T)<3) call(T[0],T[1]);                                  if(length(T)<3) call(T[0],T[1]);
                                 else call(T[0],T[1]|option_list=T[2]);                                  else call(T[0],T[1]|option_list=T[2]);
                         }                          }
Line 4496  def myswap(P,L)
Line 5168  def myswap(P,L)
 def mysubst(P,L)  def mysubst(P,L)
 {  {
         if(P==0) return 0;          if(P==0) return 0;
           if(getopt(lpair)==1||(type(L[0])==4&&length(L[0])>2)) L=lpair(L[0],L[1]);
         Inv=getopt(inv);          Inv=getopt(inv);
         if(type(L[0]) == 4){          if(type(L[0]) == 4){
                 while((L0 = car(L))!=[]){                  while((L0 = car(L))!=[]){
Line 4614  def mmulbys(FN,P,F,L)
Line 5287  def mmulbys(FN,P,F,L)
   
 def appldo(P,F,L)  def appldo(P,F,L)
 {  {
           if(getopt(Pfaff)==1){
                   L = vweyl(L);
                   X = L[0]; DX = L[1];
                   for(I=mydeg(P,DX);I>0;I--){
                           if(!(TP=mycoef(P,D,DX))) continue;
                           P=red(P+TP*(muldo(D^(I-1),F,L)-D^I));
                   }
                   return P;
           }
         if(type(F) <= 3){          if(type(F) <= 3){
                 if(type(L) == 4 && type(L[0]) == 4)                  if(type(L) == 4 && type(L[0]) == 4)
                         return applpdo(P,F,L);                          return applpdo(P,F,L);
Line 4691  def muldo(P,Q,L)
Line 5373  def muldo(P,Q,L)
 def jacobian(F,X)  def jacobian(F,X)
 {  {
         F=ltov(F);X=ltov(X);          F=ltov(F);X=ltov(X);
         N=length(F);          N=length(F);L=length(X);
         M=newmat(N,N);          M=newmat(N,L);
         for(I=0;I<N;I++)          for(I=0;I<N;I++)
                 for(J=0;J<N;J++) M[I][J]=red(diff(F[I],X[J]));                  for(J=0;J<L;J++) M[I][J]=red(diff(F[I],X[J]));
         if(getopt(mat)==1) return M;          if(N!=L||getopt(mat)==1) return M;
         return mydet(M);          return mydet(M);
 }  }
   
Line 4778  def mce(P,L,V,R)
Line 5460  def mce(P,L,V,R)
 {  {
         L = vweyl(L);          L = vweyl(L);
         X = L[0]; DX = L[1];          X = L[0]; DX = L[1];
         P = sftexp(laplace1(P,L),L,V,R);          P = sftexp(laplace1(P,L),L,V,R|option_list=getopt());
         return laplace(P,L);          return laplace(P,L);
 }  }
   
 def mc(P,L,R)  def mc(P,L,R)
 {  {
         return mce(P,L,0,R);          return mce(P,L,0,R|option_list=getopt());
 }  }
   
 def rede(P,L)  def rede(P,L)
Line 4961  def mulpdo(P,Q,L);
Line 5643  def mulpdo(P,Q,L);
   
 def transpdosub(P,LL,K)  def transpdosub(P,LL,K)
 {  {
           if(type(P)>3) return
   #ifdef  USEMODULE
                   mtransbys(os_md.transpdosub,P,[LL,K]);
   #else
                   mtransbys(transpdosub,P,[LL,K]);
   #endif
         Len = length(K)-1;          Len = length(K)-1;
         if(Len < 0 || P == 0)          if(Len < 0 || P == 0)
                 return P;                  return P;
Line 4986  def transpdosub(P,LL,K)
Line 5674  def transpdosub(P,LL,K)
   
 def transpdo(P,LL,K)  def transpdo(P,LL,K)
 {  {
         if(type(K[0]) < 4)  
                 K = [K];  
         Len = length(K)-1;          Len = length(K)-1;
         K1=K2=[];          K1=K2=[];
         if(type(LL)!=4) LL=[LL];          if(type(LL)!=4) LL=[LL];
         if(type(LL[0])!=4) LL=[LL];          if(type(LL[0])!=4) LL=[LL];
           if(type(car(K)) < 4 && length(LL)!=length(K)) K = [K];
         if(getopt(ex)==1){          if(getopt(ex)==1){
                 for(LT=LL, KT=K; KT!=[]; LT=cdr(LT), KT=cdr(KT)){                  for(LT=LL, KT=K; KT!=[]; LT=cdr(LT), KT=cdr(KT)){
                         L = vweyl(LL[J]);                          L = vweyl(LL[J]);
Line 5000  def transpdo(P,LL,K)
Line 5687  def transpdo(P,LL,K)
                 }                  }
                 K2=append(K1,K2);                  K2=append(K1,K2);
         }else{          }else{
                   if(length(LL)==length(K) && type(car(K))!=4){
                           for(DV=V=TL=[],J=length(LL)-1;J>=0;J--){
                                   TL=cons(vweyl(LL[J]),TL);
                                   V=cons(car(TL)[0],V);
                                   DV=cons(car(TL)[1],DV);
                           }
                           LL=TL;
                           if(type(RK=solveEq(K,V|inv=1))!=4) return TK;
                           if(!isint(Inv=getopt(inv))) Inv=0;
                           if(iand(Inv,1)){J=K;K=RK;RK=J;}
                           M=jacobian(RK,V|mat=1);
                           M=mulsubst(M,[V,K]|lpair=1);
                           RK=vtol(M*ltov(DV));
                           if(Inv>1) return RK;
                           K=lpair(K,RK);
                   }
                 for(J = length(K)-1; J >= 0; J--){                  for(J = length(K)-1; J >= 0; J--){
                         L = vweyl(LL[J]);                          L = vweyl(LL[J]);
                         if(L[0] != K[J][0])                          if(L[0]!= K[J][0]) K1=cons([L[0],K[J][0]],K1);
                                 K1 = cons([L[0],K[J][0]],K1);  
                         K2 = cons(K[J][1],K2);                          K2 = cons(K[J][1],K2);
                 }                  }
                 P = mulsubst(P, K1);                  P = mulsubst(P, K1);
Line 5056  def texbegin(T,S)
Line 5758  def texbegin(T,S)
 {  {
         if(type(Opt=getopt(opt))==7) Opt="["+Opt+"]\n";          if(type(Opt=getopt(opt))==7) Opt="["+Opt+"]\n";
         else Opt="\n";          else Opt="\n";
         return "\\begin{"+T+"}"+Opt+S+"%\n\\end{"+T+"}\n";          U=(str_chr(S,str_len(S)-1,"\n")<0)?"%\n":"";
           return "\\begin{"+T+"}"+Opt+S+U+"\\end{"+T+"}\n";
 }  }
   
 def mygcd(P,Q,L)  def mygcd(P,Q,L)
Line 5548  def pol2sft(F,A)
Line 6251  def pol2sft(F,A)
   
 def binom(P,N)  def binom(P,N)
 {  {
         if(type(N)!=1 || N<0) return 1;          if(type(N)!=1 || N<=0) return 1;
         for(S=1;N>0;N--,P-=1)   S*=P/N;          for(S=1;N>0;N--,P-=1)   S*=P/N;
         return red(S);          return red(S);
 }  }
Line 5564  def expower(P,R,N)
Line 6267  def expower(P,R,N)
   
 def seriesHG(A,B,X,N)  def seriesHG(A,B,X,N)
 {  {
           if(N==0) return 1;
         if(type(N)!=1 || N<0) return 0;          if(type(N)!=1 || N<0) return 0;
         if(type(X)<4){          if(type(X)<4){
                 for(K=0,S=S0=1;K<N;K++){                  for(K=0,S=S0=1;K<N;K++){
Line 5598  def evalred(F)
Line 6302  def evalred(F)
                  Opt=[];                   Opt=[];
         }else if(length(Opt)==2 && type(Opt[0])!=4)     Opt=[Opt];          }else if(length(Opt)==2 && type(Opt[0])!=4)     Opt=[Opt];
         for(;;){          for(;;){
                 G=mysubst(F,[[sin(0),0],[tan(0),0],[asin(0),0],[atan(0),0],[sinh(0),0],[tanh(0),0],                  G=mysubst(F,[[tan(0),0],[asin(0),0],[atan(0),0],[sinh(0),0],[tanh(0),0],
                         [log(1),0],[cos(0),1],[cosh(0),1],[exp(0),1]]);                          [log(1),0],[cosh(0),1],[exp(0),1]]);
                 for(Rep=Opt; Rep!=[]; Rep=cdr(Rep))                  for(Rep=Opt; Rep!=[]; Rep=cdr(Rep))
                         G=subst(G,car(Rep)[0],car(Rep)[1]);                          G=subst(G,car(Rep)[0],car(Rep)[1]);
                 Var=vars(G);                  Var=vars(G);
                 for(V=Var; V!=[]; V=cdr(V)){                  for(V=Var; V!=[]; V=cdr(V)){
                         if(functor(car(V))!=pow || (P=args(car(V))[0])!=1) continue;                          if(!(VV=args(CV=car(V)))) continue;
                         G=subst(G,car(V),1);                          if((functor(CV)==sin||functor(CV)==cos)){
                                   P=2*red(VV[0]/@pi);
                                   if(functor(CV)==sin) P=1-P;
                                   if(isint(P)){
                                           if(iand(P,1)) G=subst(G,CV,0);
                                           else if(!iand(P,3)) G=subst(G,CV,1);
                                           else G=subst(G,CV,-1);
                                           continue;
                                   }
                                   if(isint(P*=3/2)){
                                           if(iand(P,3)==1) G=subst(G,CV,1/2);
                                           else G=subst(G,CV,-1/2);
                                   }
                           }
                           for(;VV!=[];VV=cdr(VV))
                                   if(car(VV)!=(TV=evalred(car(VV)))) G=subst(G,car(VV),TV);
                           if(functor(CV)!=pow || (args(CV)[0])!=1) continue;
                           G=subst(G,CV,1);
                 }                  }
                 if(G==F) return F;                  if(G==F) return F;
                 F=G;                  F=G;
Line 5680  def seriesTaylor(F,N,V)
Line 6401  def seriesTaylor(F,N,V)
         return F;          return F;
 }  }
   
   def mulpolyMod(P,Q,X,N)
   {
           Red=(type(P)>2||type(Q)>2)?1:0;
           for(I=R=0;I<=N;I++){
                   P0=mycoef(P,I,X);
                   for(J=0;J<=N-I;J++){
                           R+=P0*mycoef(Q,J,X)*X^(I+J);
                           if(Red) R=red(R);
                   }
           }
           return R;
   }
   
   def solveEq(L,V)
   {
           Inv=0;K=length(V);
           H=(getopt(h)==1)?1:0;
           if(getopt(inv)==1){
                   if(K!=length(L)) return -5;
                   Inv=1;
                   VN=makenewv(vars(L)|num=K);
                   for(TL=[],I=K-1;I>=0;I--) TL=cons(VN[I]-L[I],TL);
                   S=solveEq(TL,V|h=H);
                   if(type(S)!=4) return S;
                   return mysubst(S,[VN,V]|lpair=1);
           }
           for(TL=[];L!=[];L=cdr(L)) TL=cons(nm(red(car(L))),TL);
           S=gr(TL,reverse(V),2);
           if(length(S)!=K) return -1;
           for(R=[],I=F=0;I<K;I++){
                   TS=S[I];
                   VI=lsort(vars(TS),V,2);
                   if(length(VI)!=1) return -2;
                   if((VI=car(VI))!=V[I]) return -3;
                   if(mydeg(TS,VI)!=1){
                           F=1;R=cons([VI,TS],R);
                   }else R=cons(-red(mycoef(TS,0,VI)/mycoef(TS,1,VI)),R);
           }
           R=reverse(R);
           if(!F||H==1) return R;
           return -4;
   }
   
   /* Opt: f, var, ord, to, in, TeX */
   def baseODE(L)
   {
           SV=SVORG;
           if(type(TeX=getopt(TeX))!=1) TeX=0;
           if(type(F=getopt(f))!=1) F=0;
           if(isint(In=getopt(in))!=1) In=0;
           if(type(Ord=getopt(ord))!=1&&Ord!=0) Ord=2;
           if(Ord>3){
                   Ord-=4; Hgr=1;
           }else Hgr=0;
           if(type(car(L))==4&&type(L[1])==7){
                   Tt=L[1];L=car(L);
           }
           M=N=length(L);  SV=SVORG;
           if(type(Var=getopt(var))==4&&(In>0||length(Var)==N)){
                   SV=Var;
                   M=length(SV);
                   if(type(car(SV))==2){
                           for(R=[];SV!=[];SV=cdr(SV)) R=cons(rtostr(car(SV)),R);
                           SV=reverse(R);
                   }
           }else{
                   if(N>10){
                           R=[];
                           for(K=M-1;K>9;K++) R=cons(SV[floor(K/10)-1]+SV[K%10],R);
                           SV=append(SV,R);
                   }
                   for(Var=[],I=M-1;I>=0;I--) Var=cons(makev([SV[I]]),Var);
           }
           if(type(To=getopt(to))<2||type(To)>4) To=0;
           else if(!isvar(To)){
                   if(type(To)!=4){
                           To=red(To);
                           for(K=0;K<length(Var);K++){
                                   I=mydeg(nm(To),Var[K]);J=mydeg(dn(To),Var[K]);
                                   if(I+J>0&&I<2&&J<2) break;
                           }
                           if(K==length(Var)) return -9;
                           J=To;
                           for(To=[],I=length(Var)-1;I>=0;I--)
                                   if(I!=K) To=cons(Var[I],To);
                           To=cons(J,To);
                   }
                   if(type(To)==4){
                           if(type(car(To))==4){
                                   R=1;To=car(To);
                           }else R=0;
                           if(type(IL=solveEq(To,Var|inv=1))!=4) return IL;
                           if(R==1){
                                   R=To;To=IL;IL=R;
                           }
                           L=mulsubst(L,[Var,IL]|lpair=1);
                           if(!In){   /* X_i'=\sum_j(\p_{x_j}X_i)*x_j' */
                                   for(TL=[],I=M-1;I>=0;I--){
                                           P=To[I];Q=mydiff(P,t);
                                           for(J=0;J<M;J++) Q=red(Q+mydiff(P,Var[J])*L[J]);
                                           TL=cons(Q,TL);
                                   }
                                   L=TL;
                           }else{  /* x_i'=\sum_j(\p_{X_j}x_i)*X_j' */
                                   for(I=M-1;I>=0;I--){
                                           P=IL[I];Q=mydiff(P,t);
                                           for(J=0;J<M;J++){
                                                   V=makev([SV[J],1]);
                                                   Q=red(Q+mydiff(P,V)*V);
                                           }
                                           L=mysubst(L,[makev([SV[I],1]),TL[I]]);
                                   }
                                   for(TL=L,L=[],I=M-1;I>=0;I--) L=cons(num(TL[I]),L);
                           }
                   }
           }
           if(F==-3) return [Var,L];
           for(I=0;I<M;I++) L=subst(L,Var[I],makev([SV[I],0]));
           if(TeX){
                   for(TL=L,I=0;I<M;I++)
                           TL=subst(TL,makev([SV[I],0]),Var[I]);
                   for(I=0;I<N;I++){
                           if(I) S0+=",\\\\\n";
                           if(In) S0+=" "+my_tex_form(TL[I])+"=0";
                           else S0+=" "+SV[I]+"'\\!\\!\\! &= "+my_tex_form(TL[I]);
                   }
                   S0+=".\n";
                   S0=texbegin("cases", S0);
                   S0=texbegin("align",S0);
                   if(type(Tt)==7) S0=Tt+"\n"+S0;
                   if(F<0){
                           if(TeX==2) dviout(S0);
                           return S0;
                   }
           }
           for(I=0,TL=[];L!=[];L=cdr(L),I++){
                   T=car(L);
                   if(!In) T=makev([SV[I],1])-T;
                   TL=cons(nm(red(T)),TL);
           }
           if(isvar(To)){
                   T=rtostr(To);
                   IT=findin(T,SV);
                   if(IT>=0 && IT<M){
                           R=[SV[IT]];
                           for(J=0;SV!=[];SV=cdr(SV),J++){
                                   if(J==IT) continue;
                                   R=cons(car(SV),R);
                           }
                           SV=reverse(R);
                   }else{
                           IT=0;
                           mycat(["Cannot find variable", T, "!\n"]);
                   }
           }
           for(S=1;S<M;S++){
                   L=append(TL,L);
                   TL=reverse(TL);
                   for(RL=[];TL!=[];TL=cdr(TL)){
                           if(In==0&&S==N-1&&length(TL)!=N-IT) continue;
                           T=car(TL);R=mydiff(V,t);
                           for(I=0;I<M;I++){
                                   for(J=0;J<=S;J++){
                                           V=makev([SV[I],J]|num=1);
                                           if((DR=mydiff(T,V))!=0) R+=DR*makev([SV[I],J+1]|num=1);
                                   }
                           }
                           RL=cons(R,RL);
                   }
                   TL=RL;
           }
           L=append(TL,L);
           for(I=0;I<M;I++) L=subst(L,makev([SV[I],0]),Var[I]);
           for(V=VV=[],I=0;I<M;I++){
                   for(J=0;J<M;J++) V=cons(J?makev([SV[I],J]):makev([SV[I]]),V);
                   if(!I||In) V=cons(makev([SV[0],M]),V);
                   if(F==-2){
                           VV=cons(V,VV);
                           V=[];
                   }
           }
           if(F>=0&&!chkfun("gr",0)){
                   mycat("load(\"gr\"); /* <- do! */\n");
                   F=-1;
           }
           if(F==-2) return [VV,L];
           if(F<0) return [V,L];
            LL=(Hgr==1)?hgr(L,V,Ord):gr(L,V,Ord);
           if(F==2) return [V,L,LL];
           if(Ord==2) P=LL[0];
           else{
                   P=LL[length(LL)-1];
                   for(RV=reverse(V), I=0;I<M+1;I++) RV=cdr(RV);
                   if(lsort(vars(P),RV,2)!=[]){
                           LL=tolex_tl(LL,V,Ord,V,2);P=LL[0];
                   }
           }
           V0=makev([car(SV),M]);
           CP=mycoef(P,mydeg(P,V0),V0);
           if(cmpsimple(-CP,CP)<0) P=-P;
           if(TeX){
                   for(V0=[makev([car(SV)])],I=1;I<=M;I++) V0=cons(makev([car(SV),I]),V0);
                   T="&\\!\\!\\!"+fctrtos(P|var=VV,dic=1,TeX=3);
                   S=((F==1)?(Tt+"\n"):S0)+texbegin("align*",texbegin("split",T));
                   if(TeX==2) dviout(S);
                   return S;
           }
           return (F==1)? P:[P,V,L,LL];
   }
   
   def taylorODE(D){
           Dif=(getopt(dif)==1)?1:0;
           if(D==0) return Dif?f:f_00;
           if(type(T=getopt(runge))!=1||ntype(T)!=0) T=0;
           if(type(F=getopt(f))!=7&&type(F)<2) F="f_";
           if(type(D)!=1||ntype(D)!=0||D<0||D>30) return 0;
           if(type(H=getopt(taylor))==4&&length(H)==2){
                   if(type(Lim=getopt(lim))==2) DD=D;
                   else if(type(Lim)==4){
                           DD=Lim[1];Lim=Lim[0];
                   }else Lim=0;
                   for(R=I=0;I<=D;I++){
                           if(I){
                                   if(Lim) H0=mulpolyMod(H0,H[0],Lim,DD);
                                   else H0*=H[0];
                           }else  H0=1;
                           if(type(F)!=7) G=I?mydiff(G,x):F;
                           for(J=0;J<=D-I;J++){
                                   if(J){
                                           if(Lim) H1=mulpolyMod(H1,H[1],Lim,DD);
                                           else H1*=H[1];
                                   }else H1=H0;
                                   if(type(F)==7) G=makev([F,I,J]);
                                   else if(J) G=mydiff(G,y);
                                   R+=G*H1/fac(I)/fac(J);
                           }
                   }
                   if(Lim) R=os_md.polcut(R,DD,Lim);
                   return R;
           }else{
                   if(type(H=getopt(series))>=0||getopt(list)==1){
                           if(type(F)!=7){
                                   for(PP=[F],I=1;I<D;I++)
                                           PP=cons(mydiff(car(PP),x)+mydiff(car(PP),y)*F,PP);
                                   if(type(H)<0) return PP;
                                   for(R=0,DD=D;DD>=1;DD--,PP=cdr(PP)) R+=car(PP)*H^DD/fac(DD);
                                   return red(R);
                           }
                           if(type(H)>=0) D--;
                           PP=taylorODE(D-1|list=1);
                           if(type(PP)!=4) PP=[PP];
                           P=car(PP);
                   }else P=taylorODE(D-1);
                   for(R=I=0;I<D;I++){
                           for(J=0;J<D-I;J++){
                                   Q=diff(P,makev([F,I,J]));
                                   if(Q!=0) R+=Q*(f_00*makev([F,I,J+1])+makev([F,I+1,J]));
                           }
                   }
                   if(getopt(list)==1){
                           R=cons(R,PP);
                           if(Dif!=1) return R;
                   }else if(type(H)>=0){
                           R=y+R*H^(D+1)/fac(D+1);
                           for(DD=D;DD>0;PP=cdr(PP),DD--) R+=car(PP)*H^(DD)/fac(DD);
                           if(T){
                                   if(T<0){
                                           Dif=0;TT=-T;
                                   }else TT=T;
                                   K=newvect(TT);K[0]=Dif?f:f_00;
                                   if(getopt(c1)==1) K[0]=taylorODE(D|taylor=[c_1*H,0]);
                                   for(I=1;I<TT;I++){
                                           for(S=J=0;J<I;J++) S+=makev(["a_",I+1,J+1])*K[J];
                                           K[I]=taylorODE(D|taylor=[makev(["c_",I+1])*H,S*H],lim=[H,D]);
                                   }
                                   for(S=I=0;I<TT;I++) S+=makev(["b_",I+1])*K[I];
                                   S=S*H+y;
                                   R=S-R;
                                   if(T<0){
                                           for(V=[H],I=0;I<=D;I++)
                                                   for(J=0;J<=D-I;J++) V=cons(makev([F,I,J]),V);
                                           return os_md.ptol(R,reverse(V)|opt=0);
                                   }
                           }else T=0;
                   }
           }
           if(Dif){
                   for(I=0;I<=D;I++){
                           for(J=0;J<=D;J++){
                                   if(I==0&&J==0){
                                           R=subst(R,f_00,f);
                                           continue;
                                   }
                                   V=makev([F,str_times("x",I),str_times("y",J)]);
                                   R=subst(R,makev([F,I,J]),V);
                           }
                   }
           }
           return R;
   }
   
 def toeul(F,L,V)  def toeul(F,L,V)
 {  {
         L = vweyl(L);          L = vweyl(L);
         X = L[0]; DX = L[1];          X = L[0]; DX = L[1];
         I = mydeg(F,DX);          I = mydeg(F,DX);
         if(V == "infty"){          if(getopt(raw)!=1){
                 for(II=I; II>=0; II--){                  for(II=I; II>=0; II--){
                         J = mydeg(P=mycoef(F,I,DX),X);                          J = mydeg(P=mycoef(F,II,DX),X);
                         if(II==I) S=II-J;                          if(II==I) S=II-J;
                         else if(P!=0 && II-J>S) S=II-J;                          else if(P!=0 && II-J>S) S=II-J;
                 }                  }
                 F *= X^S;                  F *= X^S;
                 R = 0;          }
                 for( ; I >= 0; I--)          if(V == "infty"){
                   for(R=0; I >= 0; I--)
                          R += red((mysubst(mycoef(F,I,DX),[X,1/X])*(x*DX)^I));                           R += red((mysubst(mycoef(F,I,DX),[X,1/X])*(x*DX)^I));
                 return(subst(pol2sft(R,DX),DX,-DX));                  return(subst(pol2sft(R,DX),DX,-DX));
         }          }
         F = subst(F,X,X+V);          for(R=0; I >= 0; I--)
         for(II=I; II>=0; II--){  
                 J = mymindeg(P=mycoef(F,II,DX),X);  
                 if(II==I) S=II-J;  
                 else if(P!=0 && II-J>S) S=II-J;  
         }  
         F *= X^S;  
         R = 0;  
         for( ; I >= 0; I--)  
                 R += (red(mycoef(F,I,DX)/X^I))*DX^I;                  R += (red(mycoef(F,I,DX)/X^I))*DX^I;
         return pol2sft(R,DX);          return pol2sft(R,DX);
 }  }
Line 5741  def fromeul(P,L,V)
Line 6756  def fromeul(P,L,V)
                 S = DX*(S*X + mydiff(S,DX));                  S = DX*(S*X + mydiff(S,DX));
                 R += mycoef(P,J,DX)*S;                  R += mycoef(P,J,DX)*S;
         }          }
         while(mycoef(R,0,X) == 0)          if(getopt(raw)!=1){
                 R = tdiv(R,X);                  while(mycoef(R,0,X) == 0)
                           R = tdiv(R,X);
           }
         if(V != "infty" && V != 0)          if(V != "infty" && V != 0)
                 R = mysubst(R,[X,X-V]);                  R = mysubst(R,[X,X-V]);
         return R;          return R;
Line 5751  def fromeul(P,L,V)
Line 6768  def fromeul(P,L,V)
 def sftexp(P,L,V,N)  def sftexp(P,L,V,N)
 {  {
         L = vweyl(L); DX = L[1];          L = vweyl(L); DX = L[1];
         P = mysubst(toeul(P,L,V),[DX,DX+N]);          P = mysubst(toeul(P,L,V|opt_list=getpt()),[DX,DX+N]);
         return fromeul(P,L,V);          return fromeul(P,L,V|opt_list=getopt());
 }  }
   
   
Line 6147  def expat(F,L,V)
Line 7164  def expat(F,L,V)
   
 def polbyroot(P,X)  def polbyroot(P,X)
 {  {
           if(isvar(V=getopt(var))&&length(P)>1&&isint(car(P))){
                   for(Q=[],I=car(P);I<=P[1];I++) Q=cons(makev([V,I]),Q);
                   P=Q;
           }
         R = 1;          R = 1;
         while(length(P)){          while(length(P)){
                 R *= X-car(P);                  R *= X-car(P);
Line 6335  def pcoef(P,L,Q)
Line 7356  def pcoef(P,L,Q)
         return Coef;          return Coef;
 }  }
   
   def pmaj(P)
   {
           if(type(P)==4){
                   Opt=delopt(getopt(),"var"|inv=1);
                   P=mtransbys(os_md.pmaj,P,[]|optilon_list=Opt);
                   if(Opt==[]) return P;
                   X=Opt[0][1];
                   D=mydeg(P,X);
                   for(S=0;D>=0;D--) S+=lmax(mycoef(P,D,X))*X^D;
                   return S;
           }
           V=vars(P);
           if(!(K=length(V))) return abs(P);
           for(R=0,D=deg(P,X=V[0]);D>=0;D--){
                   Q=coef(P,D,X);
                   if(Q!=0) R+=((type(Q)>1)?pmaj(Q):abs(Q))*X^D;
           }
           if(isvar(Y=getopt(var))) for(;V!=[];V=cdr(V)) R=subst(R,car(V),Y);
           return R;
   }
   
 def prehombf(P,Q)  def prehombf(P,Q)
 {  {
         if((Mem=getopt(mem))!=1 && Mem!=-1)          if((Mem=getopt(mem))!=1 && Mem!=-1)
Line 7167  def abs(X)
Line 8209  def abs(X)
         return X;          return X;
 }  }
   
   def sgn(X)
   {
           if(X==0) return 0;
           if(type(X)==1){
                   return (X>0)?1:-1;
           }
           if(type(X)==5) X=vtol(X);
           if(type(X)==4){
                   for(W=0,Y=X;Y!=[];Y=cdr(Y))
                           for(Z=cdr(Y);Z!=[];Z=cdr(Z))
                                   if(car(Y)>car(Z)) W++;
                   if(getopt(val)==1) return W;
                   return (iand(W,1))?-1:1;
           }
   }
   
 def calc(X,L)  def calc(X,L)
 {  {
         if(type(X)<4){          if(type(X)<4||type(X)==7){
                 if(type(L)==4){                  if(type(L)==4||type(L)==7){
                         V=L[1];                          V=L[1];
                         if((L0=L[0])=="+") X+=V;                          if(type(X)!=7){
                         else if(L0=="-")   X-=V;                                  if((L0=L[0])=="+") X+=V;
                         else if(L0=="*")   X*=V;                                  else if(L0=="-")   X-=V;
                         else if(L0=="/")   X/=V;                                  else if(L0=="*")   X*=V;
                         else if(L0=="^")   X^=V;                                  else if(L0=="/")   X/=V;
                         else if(L0==">")   X=(X>V);                                  else if(L0=="^")   X^=V;
                         else if(L0=="<")   X=(X<V);                          }
                         else if(L0=="=")   X=(X==V);                          if((L0=L[0])==">")      X=(X>V);
                           else if(L0=="<")        X=(X<V);
                           else if(L0=="=")        X=(X==V);
                         else if(L0==">=")   X=(X>=V);                          else if(L0==">=")   X=(X>=V);
                         else if(L0=="<=")   X=(X<=V);                          else if(L0=="<=")   X=(X<=V);
                         else if(L0=="!=")       X=(X!=V);                          else if(L0=="!=")       X=(X!=V);
                 }else if(type(L)==7){                  }else if(type(L)==7&&type(X)<4){
                         if(L=="neg") X=-X;                          if(L=="neg") X=-X;
                         else if(L=="abs") X=abs(X);                          else if(L=="abs") X=abs(X);
                         else if(L=="neg") X=-X;                          else if(L=="neg") X=-X;
Line 7198  def calc(X,L)
Line 8258  def calc(X,L)
         return X;          return X;
 }  }
   
   def tobig(X)
   {
           if((type(X)==1 && ntype(X)==3)||type(X)>3) return X;
           return eval(X*exp(0));
   }
   
 def isint(X)  def isint(X)
 {  {
         if(X==0||(type(X)==1 && ntype(X)==0 && dn(X)==1)) return 1;          if(X==0||(type(X)==1 && ntype(X)==0 && dn(X)==1)) return 1;
Line 7443  def spgen(MO)
Line 8509  def spgen(MO)
         if(F!=1&&F!=-1) F=0;          if(F!=1&&F!=-1) F=0;
         if(type(LP)==4){          if(type(LP)==4){
                 L0=LP[0]; L1=LP[1];                  L0=LP[0]; L1=LP[1];
           }else if(type(LP)==1){
                   L0=L1=LP;
         }else{          }else{
                 L0=0; L1=MO+1;                  L0=0; L1=MO+1;
         }          }
         if(MO<=0){          if(M0<=0){
                 MO=-MO;                  MO=-MO;
                 if(iand(MO,1)==1) return [];                  if(iand(MO,1)==1) return [];
                 if(MO>1){                  MO=MO/2;
                         if(isMs()==0) return [];                  B=spbasic(-2*MO,0|str=1);
                         Cmd="okubo "+rtostr(-MO);                  if(L1<3) L1=MO+4;
                         MO/=2;  
                         if(L1>0) Cmd=Cmd+"+"+rtostr(L0)+"-"+rtostr(L1);  
                         else L1=MO+4;  
                         Cmd=Cmd+" B";  
                         Id=getbyshell(Cmd);  
                         if(Id<0) return [];  
                         B=[];  
                         while((S=get_line(Id)) !=0){  
                                 P0=str_chr(S,1,":")+1;  
                                 if(P0>1){  
                                         P1=str_chr(S,P,"\n");  
                                         if(P1<0) P1=str_len(S);  
                                         B=cons(sub_str(S,P0,P1-1),B);  
                                 }  
                         }  
                 }else{  
                         MO/=2;  
                         if(L1<=1) L1=MO+4;  
 BB=[  
 ["11,11,11,11","111,111,111","1^4,1^4,22","1^6,222,33"],  
 ["11,11,11,11,11","1^4,1^4,211","211,22,22,22","1^6,2211,33",  
 "2211,222,222","22211,2^4,44","2^511,444,66","1^4,22,22,31",  
 "2^5,3331,55","1^5,1^5,32","1^8,332,44","111,111,21,21","1^5,221,221"],  
 ["11,11,11,11,11,11","1^4,1^4,1^4","1^4,22,22,22","111,111,111,21",  
 "1^6,21^4,33","21^4,222,222","221^4,2^4,44","2^41^4,444,66",  
 "1^5,1^5,311","1^8,3311,44","1^6,222,321","321,33,33,33",  
 "3321,333,333","33321,3^4,66","3^721,666,99","2^5,3322,55",  
 "1^6,1^6,42","222,33,33,42","1^a,442,55","1^6,33,33,51",  
 "222,222,33,51","1^9,333,54","2^7,554,77","1^5,2111,221",  
 "2^41,333,441","1^7,2221,43","211,211,22,22","2211,2211,222",  
 "22211,22211,44","1^4,211,22,31","2^411,3331,55","1^4,1^4,31,31",  
 "22,22,22,31,31","1^7,331,331","2221,2221,331","111,21,21,21,21"],  
 ["11,11,11,11,11,11,11","111,111,111,111","1^6,1^6,33",  
 "1^6,222,222","222,33,33,33","1^5,1^5,221",  
 "1^4,211,22,22","1^4,1^4,22,31","22,22,22,22,31",  
 "111,111,21,21,21","21^6,2^4,44","2221^6,444,66",  
 "1^6,222,3111","3111,33,33,33","33111,333,333",  
 "333111,3^4,66","3^5111,666,99","2^5,33211,55",  
 "1^8,3221,44","3222,333,333","33222,3^4,66",  
 "3^4222,666,99","1^6,1^6,411","222,33,33,411",  
 "1^a,4411,55","2^4,2^4,431","431,44,44,44",  
 "2^6,4431,66","4431,444,444","44431,4^4,88",  
 "4^531,888,cc","1^a,433,55","1^7,1^7,52",  
 "1^c,552,66","3^4,444,552","1^8,2^4,53",  
 "1^8,44,44,71","3^5,555,771","21^4,2211,222",  
 "221^4,22211,44","2221^4,3331,55","1^6,2211,321",  
 "2^411,3322,55","1^7,322,331","2211,33,33,42",  
 "3^42,4442,77","2211,222,33,51","3^51,5551,88",  
 "2^611,554,77","2221,2221,322","2^41,2^41,54",  
 "1^5,2111,2111","222111,333,441","1^7,22111,43",  
 "1^5,1^5,41,41","1^9,441,441","22111,2221,331",  
 "1^5,221,32,41","221,221,221,41","211,211,211,22",  
 "2211,2211,2211","1^4,211,211,31","211,22,22,31,31",  
 "1^4,22,31,31,31","1^5,32,32,32","221,221,32,32","21,21,21,21,21,21"],  
 ["11,11,11,11,11,11,11,11","1^4,1^4,22,22","1^8,2^4,44",  
 "1^6,2211,222","2211,33,33,33","111,111,111,21,21",  
 "1^5,1^5,2111","1^4,211,211,22","1^4,1^4,211,31",  
 "211,22,22,22,31","1^4,22,22,31,31","111,21,21,21,21,21",  
 "221^8,444,66","2^5,331^4,55","1^8,32111,44",  
 "32211,333,333","332211,3^4,66","3^42211,666,99",  
 "2^5,32221,55","1^7,1^7,511","1^c,5511,66",  
 "3^4,444,5511","541,55,55,55","5541,555,555",  
 "55541,5^4,aa","5^541,aaa,ff","1^8,1^8,62",  
 "1^a1^4,662,77","1^a,55,55,91","2^71,555,87",  
 "21^6,22211,44","221^6,3331,55","1^6,2211,3111",  
 "2^411,33211,55","1^7,3211,331","2211,33,33,411",  
 "3^42,44411,77","22211,2^4,431","2^511,4431,66",  
 "1^8,332,431","3^42,4433,77","1^8,22211,53",  
 "2221,2221,3211","221^5,333,441","1^7,21^5,43",  
 "1^b,443,65","21^5,2221,331","2^51,3332,65",  
 "21^4,21^4,222","221^4,221^4,44","1^6,21^4,321",  
 "2221^4,3322,55","21^4,33,33,42","21^4,222,33,51",  
 "2^51^4,554,77","2^4,3311,3311","3^411,4442,77",  
 "321,321,33,33","3321,3321,333","33321,33321,66",  
 "222,321,33,42","1^6,321,33,51","222,222,321,51",  
 "1^9,3321,54","1^7,322,322","3^422,5551,88",  
 "1^6,33,42,42","1^6,222,42,51","33,33,33,42,51",  
 "1^6,1^6,51,51","222,33,33,51,51","1^b,551,551",  
 "1^5,221,311,41","2^41,3321,441","22111,2221,322",  
 "2^51,443,551","222111,2^41,54","21^4,2211,2211",  
 "1^5,311,32,32","3331,3331,442","2211,2211,33,51",  
 "221,221,311,32","22111,22111,331","1^5,2111,32,41",  
 "2111,221,221,41","2111,221,32,32","211,211,211,211",  
 "211,211,22,31,31","1^4,211,31,31,31","22,22,31,31,31,31"],  
 ["11,11,11,11,11,11,11,11,11","1^5,1^5,1^5","2^5,2^5,55",  
 "111,111,111,111,21","2^41,333,333","1^4,1^4,211,22",  
 "211,22,22,22,22","1^8,22211,44","1^4,1^4,1^4,31",  
 "1^4,22,22,22,31","1^7,1^7,43","1^7,2221,331",  
 "2221,2221,2221","1^6,21^4,222","21^4,33,33,33",  
 "1^6,1^6,321","222,321,33,33","1^6,33,33,42",  
 "222,222,33,42","1^6,222,33,51","222,222,222,51",  
 "33,33,33,33,51","1^6,2211,2211","111,111,21,21,21,21",  
 "1^5,1^5,32,41","1^5,221,221,41","1^5,221,32,32",  
 "221,221,221,32","1^4,211,211,211","211,211,22,22,31",  
 "1^4,211,22,31,31","1^4,1^4,31,31,31","22,22,22,31,31,31",  
 "21,21,21,21,21,21,21","21^a,444,66","1^8,31^5,44",  
 "321^4,333,333","3321^4,3^4,66","3^421^4,666,99",  
 "2^5,322111,55","32^41,3^4,66","3332^41,666,99",  
 "1^8,1^8,611","2^4,44,44,611","1^d,6611,77",  
 "4^5,66611,aa","2^6,444,651","3^4,3^4,651",  
 "651,66,66,66","3^6,6651,99","6651,666,666",  
 "66651,6^4,cc","6^551,ccc,ii","2^8,655,88",  
 "1^9,1^9,72","1^g,772,88","1^c,444,75",  
 "2^6,3^4,75","1^c,66,66,b1","3^4,444,66,b1",  
 "3^7,777,ba","1^7,2221,4111","2^41,333,4311",  
 "1^9,2^41,63","21^8,3331,55","2^411,331^4,55",  
 "1^7,31^4,331","2^411,32221,55","22211,2^4,422",  
 "2^511,4422,66","1^8,332,422","2^5,3331,541",  
 "22211,44,44,62","2^411,2^5,64","2^711,664,88",  
 "1^a,3331,64","2221,2221,31^4","21^7,333,441",  
 "333,333,441,81","2^6111,555,87","21^6,221^4,44",  
 "221^6,3322,55","2^41^6,554,77","1^6,21^4,3111",  
 "3111,321,33,33","33111,3321,333","333111,33321,66",  
 "222,3111,33,42","1^6,3111,33,51","222,222,3111,51",  
 "1^9,33111,54","2221^4,33211,55","1^7,3211,322",  
 "3^4211,5551,88","2^4,3221,3311","333221,4442,77",  
 "3222,3321,333","33222,33321,66","1^9,3222,54",  
 "21^4,33,33,411","3^411,44411,77","222,321,33,411",  
 "1^6,33,411,42","1^6,222,411,51","33,33,33,411,51",  
 "221^4,2^4,431","2^41^4,4431,66","1^8,3311,431",  
 "3^411,4433,77","33321,444,552","1^8,221^4,53",  
 "3311,44,44,53","4^42,5553,99","2^4,3311,44,71",  
 "3^421,555,771","4^52,7771,bb","3^611,776,aa",  
 "2^41,33111,441","22111,2221,3211","2^41,3222,441",  
 "2^61,4441,76","3331,3331,4411","22211,22211,431",  
 "3331,3331,433","3^41,3^41,76","1^7,1^7,61,61",  
 "1^d,661,661","21^5,2221,322","221^5,2^41,54",  
 "2^51,33311,65","21^5,22111,331","3^41,4441,661",  
 "1^7,331,43,61","2221,2221,43,61","2221,331,331,61",  
 "21^4,21^4,2211","21^4,2211,33,51","22211,3311,3311",  
 "1^5,311,311,32","2211,321,33,42","2211,222,321,51",  
 "3322,3331,442","2211,222,42,42","2^411,442,442",  
 "1^6,2211,42,51","2211,33,33,51,51","221,221,311,311",  
 "1^5,2111,311,41","222111,3321,441","22111,22111,322",  
 "222111,222111,54","2111,221,311,32","2111,2111,221,41",  
 "1^5,221,41,41,41","2221,43,43,43","1^5,32,32,41,41",  
 "331,331,43,43","221,221,32,41,41","221,32,32,32,41",  
 "211,211,211,31,31","211,22,31,31,31,31","1^4,31,31,31,31,31"]];  
                         B=BB[MO];  
                 }  
                 if(St!=1){                  if(St!=1){
                         for(R=[]; B!=[]; B=cdr(B)){                          for(R=[]; B!=[]; B=cdr(B)){
                                 RT=F?s2sp(car(B)|std=F):s2sp(car(B));                                  RT= F?s2sp(car(B)|std=F): s2sp(car(B));
                                 if(length(RT)<L0 || length(RT)>L1) continue;                                  if(length(RT)<L0 || length(RT)>L1) continue;
                                 R=cons(RT,R);                                  R=cons(RT,R);
                         }                          }
Line 7696  BB=[
Line 8624  BB=[
         return LL;          return LL;
 }  }
   
   def spbasic(Idx,D)
   {
   /*
     D<=3|Idx|+6,  D<=|Idx|+2 (p>3),  p<=|Idx|/2+4
     Idx=2*D^2-(D^2-\sum m_{j,\nu}^2); \sum(D-m_{j,1})>=2*D;
     \sum (m_{j,1)-m_{j,\nu})*m_{j,\nu)
     0<=(2*D-\sum(D-m_{j,1})})*D=\sum_(m_{j,1}-m_{j,\mu})*m_{j,\nu} -|Idx|
     (-2,0)                                    13ŒÂ (9+3+?)
     (-4,0)                                    37ŒÂ (25+9+?)
     (-6,0) :  8.5sec  ?sec          0.05sec   69ŒÂ (46+17+?)
     (-8,0) : 97  sec  1sec          0.13sec  113ŒÂ (73+29+?)   <- (-2,0)
     (-10,0):          4sec          0.27sec  198ŒÂ (127+50+?)
   @(-12,0)          28sec   4.2sec 0.64sec  291ŒÂ (182+76+?)
     (-14,0)          27sec  10.2sec 1.31sec  415ŒÂ (249+115+?)
     (-16,0)                 34.0sec 2.47sec  647ŒÂ (395+172+?) <- (-4,0)
     (-18,0)                         4.42sec  883ŒÂ (521+243+?) <- (-2,0)
     (-20,0)                         8.17sec 1186ŒÂ (680+345+?)
   */
           Idx=-Idx;
           if((Str=getopt(str))!=1) Str=0;
           if(!isint(Idx)||!isint(Idx/2)||Idx<0||!isint(D)||D<0||D==1||D>3*Idx+6) return [];
           if(D==0){
                   for(R=[],D=3*Idx+6;D>=2;D--) R=append(spbasic(-Idx,D|str=Str),R);
                   return R;
           }
           if(!Idx){
                   R=0;
                   if(D==2) R="11,11,11,11";
                   if(D==3) R="111,111,111";
                   if(D==4) R="22,1111,1111";
                   if(D==6) R="33,222,111111";
                   if(!R) return [];
                   return [(Str==1)?R:s2sp(R)];
           }
           if(D>Idx+2){
                   L=3;
                   if(D==3*Idx+6){
                           R=[[D/2,D/2],[D/3,D/3,D/3],[D/6,D/6,D/6,D/6,D/6,D/6-1,1]];
                           return [(Str==1)?s2sp(R):R];
                   }
                   if(iand(D,1)&&(D-3)/2>Idx) return [];
           }else L=Idx/2+4;
           V=newvect(L);SV=newvect(L);
           for(S1=[],I=0;I<D;I++) S1=cons(1,S1);
           for(T=D-1;T>1;T--){
                   K=D%T;
                   if((T-K)*K<=Idx) break;
           }
           J=(T-K)*K;SJ=K^2+(D-K)*T;
           TV=K?[K]:[];
           for(I=(D-K)/T;I>0;I--) TV=cons(T,TV);
           for(I=0;I<L;I++){
                   SV[I]=2*D^2-(I+1)*(D^2-J)-Idx;
                   V[I]=TV;
           }
           if(SV[2]>0) return [];
           if(D>Idx+2 && V[0][0]+V[1][0]>=D && V[1][0]>1){
                   T=V[1][0]-1;K=D%T;TV=K?[K]:[];
                   for(I=(D-K)/T;I>0;I--) TV=cons(T,TV);
                   V[1]=V[2]=TV;
           }
           for(R=[];;){
                   if(D>Idx+2){
                           if(3*V[0][0]<D) break;
                           if(V[0][0]+V[1][0]>=D && (T=D-V[0][0]-1)>0){
                                   K=D%T;TV=K?[K]:[];
                                   for(I=(D-K)/T;I>0;I--) TV=cons(T,TV);
                                   V[1]=V[2]=TV;
                           }
                           S2=V[0][0]+V[1][0]+V[2][0]-D;
                           if(V[0][0]+2*V[1][0]<D ||(S2<0&&V[1][0]==1) ){
                                   V[0]=V[1]=V[2]=nextpart(V[0]);
                                   T=V[0][0];
                                   T=D-2*T;
                                   if(T==0){
                                           V[1]=[D/2-1,1];
                                           V[2]=S1;
                                   }else if(T>0){
                                           J=D%T;
                                           K=J?[J]:[];
                                           for(J=(D-J)/T;J>0;J--) K=cons(T,K);
                                           V[2]=K;
                                   }
                                   continue;
                           }
                           if(S2<0||V[2][0]<=S2){
                                   V[1]=V[2]=nextpart(V[1]);
                                   continue;
                           }else if(S2>0){
                                   T=V[2][0]-S2;J=D%T;
                                   K=J?[J]:[];
                                   for(J=(D-J)/T;J>0;J--) K=cons(T,K);
                                   V[2]=K;
                           }
                   }
                   for(S=-2*D,IL=0;IL<L;IL++){
                           S+=D-car(V[IL]);
                           if(S>=0) break;
                   }
                   if((I=IL)==L){  /* reducible i.e. IL=L && S<0 */
                           for(LL=L-1;LL>=0;LL--){
                                   if((K=car(V[LL]))+S>0){
                                           K+=S;
                                           for(TV=[],TD=D;TD>=K;TD-=K) TV=cons(K,TV);
                                           if(TD>0) V[LL]=append(TV,[TD]);
                                           else V[LL]=TV;
                                           break;
                                   }else{
                                           S+=K-1;
                                           V[LL]=S1;
                                   }
                           }
                           if(LL<0) break;
                           continue;
                   }
                   for(S0=K=0;K<=IL;K++){
                           ST=car(V[K]);J=V[K][length(V[K])-1];S0+=(ST-J)*J;
                           if(S0>Idx) break;
                   }
                   if(S0>Idx && car(V[K])!=1){
                           ST=car(V[K]);
                           S0-=(ST-J)*J;
                           for(ST--;ST>0;ST--){
                                   J=D%ST;
                                   if(S0+(ST-J)*J <= Idx) break;
                           }
                           V[K]=J?[J]:[];
                           for(J=D-J;J>0;J-=ST) V[K]=cons(ST,V[K]);
                           for(J=K+1;J<L;J++) V[J]=V[K];
                           continue;
                   }
   
                   for(K=SS=0;K<L&&SS<=Idx;K++){
                           ST=car(V[K]);
                           for(S0=0,TV=cdr(V[K]);TV!=[];TV=cdr(TV)) S0+=(ST-car(TV))*car(TV);
                           SS+=S0;
                   }
                   if(SS>Idx && K<=IL && K!=L){
                           SS0=Idx-SS+S0;
                           for(TV=car(V[K]);TV>1;TV--){
                                   U=D%TV;
                                   if((D-U)*U<=SS0) break;
                           }
                           if(TV==car(V[K])){
                                   K=K-1;
                                   V[K]=nextpart(V[K]); /* to be improves */
                           }else{
                                   V[K]=U?[U]:[];  /* to be improved */
                                   for(J=D-U;J>0;J-=TV) V[K]=cons(TV,V[K]);
                           }
                           for(J=K+1;J<L;J++) V[J]=V[K];
                           continue;
                   }
   
                   for(Ix=2*D^2+Idx,J=0;J<L;J++){
                           IxF=Ix;
                           for(Ix-=D^2,TV=V[J];TV!=[];TV=cdr(TV)) Ix+=car(TV)^2;
                           if(Ix<=0) break;
                   }
                   if(Ix==0&&(J>=I||IL==2)){
                           for(TR=[],K=J;K>=0;K--) TR=cons(V[K],TR);
                           R=cons((Str==1)?s2sp(TR):TR,R);
                   }
                   if(J>=0 && J<L && Ix<=0){
                           I=V[J][0];K=D%I;S0=(D-K)*I+K^2;
                           if(I>1&& IxF-D^2+S0<0){
                                   for(V[J]=[],K=D-I;K>0;K--) V[J]=cons(1,V[J]);
                                   V[J]=cons(I,V[J]);
                                   V[J]=nextpart(V[J]);
                                   for(I=J+1;I<L;I++) V[I]=V[J];
                                   continue;
                           }
                   }
                   if(J>=0 && J<L && Ix<=0 && car(V[J])>(U=V[J][length(V[J])-1])+1){
                           TV=reverse(V[J]);
                           for(S0=0,K=[];TV!=[];TV=cdr(TV),S0++){
                                   if((I=car(TV))<U+2||(length(TV)>1&&S0<2)){
                                           while(I-->0) K=cons(1,K);
                                   }else K=cons(car(TV),K);
                           }
                           V[I=J]=K;
                   }else{
                           if(J>=L) J=L-1;
                           for(I=J;I>=0&&length(V[I])==D;I--);
                           if(I<0) break;
                   }
                   V[I]=nextpart(V[I]);                    /* to be improved */
                   for(J=I+1;J<L;J++) V[J]=V[I];
           }
           return R;
   }
   
 def spType2(L)  def spType2(L)
 {  {
         C=0;R=[];          C=0;R=[];
Line 7732  def chkspt(M)
Line 8852  def chkspt(M)
         Opt= getopt(opt);          Opt= getopt(opt);
         Mat= getopt(mat);          Mat= getopt(mat);
         if(type(M)==7) M=s2sp(M);          if(type(M)==7) M=s2sp(M);
         if(type(Opt) >= 0){          if(type(Opt) >= 0&&Opt!="idx"){
                 if(type(Opt) == 7)                  if(type(Opt) == 7)
                         Opt = findin(Opt, ["sp","basic","construct","strip","short","long","sort","root"]);                          Opt = findin(Opt, ["sp","basic","construct","strip","short","long","sort","root"]);
                 if(Opt < 0){                  if(Opt < 0){
Line 7741  def chkspt(M)
Line 8861  def chkspt(M)
                 }                  }
                 return fspt(M,Opt);                  return fspt(M,Opt);
         }          }
         MR = fspt(M,1);  
         P  = length(M);          P  = length(M);
         OD = -1;          OD = -1;
         XM = newvect(P);          XM = newvect(P);
Line 7764  def chkspt(M)
Line 8883  def chkspt(M)
                 if(OD < 0)                  if(OD < 0)
                         OD = SM;                          OD = SM;
                 else if(OD != SM){                  else if(OD != SM){
                         print("irregal partitions");                          if(getopt(dumb)!=1) print("irregal partitions");
                         return 0;                          return -1;
                 }                  }
                 XM[I] = JM;                  XM[I] = JM;
         }          }
Line 7782  def chkspt(M)
Line 8901  def chkspt(M)
                 SM += MV;                  SM += MV;
         }          }
         SM -= (P-2)*OD;          SM -= (P-2)*OD;
           if(Opt=="idx") return SSM;
         if(SM > SMM && SM != 2*OD){          if(SM > SMM && SM != 2*OD){
                 print("not realizable");                  if(getopt(dumb)!=1) print("not realizable");
                 return -1;                  return 0;
         }          }
         if(JM==1 && Mat!=1)          if(JM==1 && Mat!=1)
                 Fu -= OD - SSM/2;                  Fu -= OD - SSM/2;
         return [P, OD, SSM, Fu, SM, XM, MR];          return [P, OD, SSM, Fu, SM, XM, fspt(M,1)];
 }  }
   
 def cterm(P)  def cterm(P)
Line 7903  def redgrs(M)
Line 9023  def redgrs(M)
                                 }                                  }
                                 L = cons([VM,EV], L);                                  L = cons([VM,EV], L);
 /*  /*
                                 if(R[2] >= 2){ */ /* digid */                                  if(R[2] >= 2){ */ /* rigid */
 /*          P = dx^(R[1]);  /*          P = dx^(R[1]);
                                 } */                                  } */
                         }                          }
Line 7930  def mcgrs(G, R)
Line 9050  def mcgrs(G, R)
 {  {
         NP = length(G);          NP = length(G);
         Mat = (getopt(mat)==1)?0:1;          Mat = (getopt(mat)==1)?0:1;
           if(Mat==0 && type(SM=getopt(slm))==4){
                   SM0=SM[0];SM1=anal2sp(SM[1],["*",-1]);
                   if(findin(0,SM0)>=0){
                           for(SM=[],I=length(G)-1;I>0;I--)
                                   if(findin(I,SM0)<0) SM=cons(I,SM);
                           SM=[SM,SM1];
                           G=mcgrs(G,R|mat=1,slm=SM);
                           return [G[0],anal2sp(G[1],["*",-1])];
                   }
           }else SM0=0;
         for(R = reverse(R) ; R != []; R = cdr(R)){          for(R = reverse(R) ; R != []; R = cdr(R)){
                 GN = [];                  GN = [];
                 L = length(G)-1;                  L = length(G)-1;
                 RT = car(R);                  RT = car(R);
                 if(type(RT) == 4){                  if(type(RT) == 4){
                         RT = reverse(RT); S = 0;                          if(length(RT)==L+1&&RT[0]!=0){
                         for(G = reverse(G); G != []; G = cdr(G), L--){                                  R=cons(cdr(RT),cdr(R));
                                 AD = car(RT); RT = cdr(RT);                                  R=cons(RT[0],R);
                                 if(L > 0)                                  R=cons(0,R);
                                   continue;
                           }               /* addition */
                           RT = reverse(RT); S = ADS = 0;
                           for(G = reverse(G); G != []; G = cdr(G), L--, RT=cdr(RT)){
                                   AD = car(RT);
                                   if(L > 0){
                                         S += AD;                                          S += AD;
                                 else                                          if(SM && findin(L,SM0)>=0) ADS+=AD;
                                   }else
                                         AD = -S;                                          AD = -S;
                                 for(GTN = [], GT = reverse(car(G)); GT != []; GT = cdr(GT))                                  for(GTN = [], GT = reverse(car(G)); GT != []; GT = cdr(GT))
                                         GTN = cons([car(GT)[0],car(GT)[1]+AD], GTN);                                          GTN = cons([car(GT)[0],car(GT)[1]+AD], GTN);
                                 GN = cons(GTN, GN);                                  GN = cons(GTN, GN);
                         }                          }
                         G = GN;                          G = GN;
                           if(SM0){
                                   for(ST=reverse(SM1),SM1=[]; ST!=[]; ST=cdr(ST))
                                           SM1 = cons([car(ST)[0],car(ST)[1]+ADS], SM1);
                           }
                         continue;                          continue;
                 }                  }
                 VP = newvec(L+1); GV = ltov(G);                  if(RT==0) continue;
                   VP = newvect(L+1); GV = ltov(G);        /* middle convolution */
                 for(I = S = OD = 0; I <= L; I++){                  for(I = S = OD = 0; I <= L; I++){
                         RTT = (I==0)?(Mat-RT):0;                          RTT = (I==0)?(Mat-RT):0;
                         VP[I] = -1;                          VP[I] = -1;
                         for(J = M = 0, GT = GV[I]; GT != []; GT = cdr(GT), J++){                          for(J = M = K = 0, GT = GV[I]; GT != []; GT = cdr(GT), J++){
                                 if(I == 0)                                  if(I == 0)
                                         OD += car(GT)[0];                                          OD += car(GT)[0];
                                 if(car(GT)[1] == RTT && car(GT)[0] > M){                                  if(car(GT)[1] == RTT && car(GT)[0] > M){
                                         S += car(GT)[0]-M;                                          S += car(GT)[0]-M;
                                           M=car(GT)[0];
                                         VP[I] = J;                                          VP[I] = J;
                                 }                                  }
                         }                          }
                         S -= (L-1)*OD;                  }
                         for(GN = [] ; L >= 0; L--){                  S -= (L-1)*OD;
                                 GT = GV[L];                  for(GN = []; L >= 0; L--){
                                 RTT = (L==0)?(-RT):RT;                          GT = GV[L];
                                 FTN = (VP[L] >= 0 || S == 0)?[]:[-S,(L==0)?(Mat-RT):0];                          RTT = (L==0)?(-RT):RT;
                                 for(J = 0; GT != []; GT = cdr(GT), J++){                          GTN = (VP[L]>=0 || S == 0)?[]:[[-S,(L==0)?(Mat-RT):0]];
                                         if(J != VP[L]){                          for(J = 0; GT != []; GT = cdr(GT), J++){
                                                 GTN = cons([car(GT)[0],car(GT)[1]+RTT], GTN);                                  if(J != VP[L]){
                                                 continue;                                          GTN = cons([car(GT)[0],car(GT)[1]+RTT], GTN);
                                         }                                          continue;
                                         K = car(GT)[0] - S;  
                                         if(K < 0){  
                                                 print("Not realizable");  
                                                 return;  
                                         }  
                                         GTN = cons([K,(L==0)?(Mat-RT):0], GTN);  
                                 }                                  }
                                 GN = cons(reverse(GTN), GN);                                  K = car(GT)[0] - S;
                                   if(K < 0){
                                           print("Not realizable");
                                           return;
                                   }
                                   if(K>0) GTN = cons([K,(L==0)?(Mat-RT):0], GTN);
                         }                          }
                           GN = cons(reverse(GTN), GN);
                 }                  }
                   if(SM0&&RT!=0){
                           for(M0=M1=-OD,L=length(G)-1;L>=0;L--){
                                   if(findin(L,SM0)>=0){
                                           M0+=OD;
                                           if(VP[L]>=0) M0-=GV[L][VP[L]][0];
                                   }else{
                                            M1+=OD;
                                           if(VP[L]>=0) M1-=GV[L][VP[L]][0];
                                   }
                           }
                           SM2=[];
                           if((Mx1=anal2sp(SM1,["max",1,-RT])[0])<0){
                                   if(M1>0) SM2=cons([M1,0],SM2);
                           }else M1+=car(SM1[Mx1]);
                           if((Mx0=anal2sp(SM1,["max",1,0])[0])<0){
                                   if(M0>0) SM2=cons([M0,RT],SM2);
                           }else M0+=car(SM1[Mx0]);
                           for(J=0;SM1!=[];J++,SM1=cdr(SM1)){
                                   if(J==Mx0){
                                           if(M0>0) SM2=cons([M0,-RT],SM2);
                                   }else if(J==Mx1){
                                           if(M1>0) SM2=cons([M1,0],SM2);
                                   }else SM2=cons([car(SM1)[0],car(SM1)[1]+RT],SM2);
                           }
                           SM1=reverse(SM2);
                   }
                 G = cutgrs(GN);                  G = cutgrs(GN);
         }          }
         return G;          return SM0?[G,SM1]:G;
 }  }
   
   def spslm(M,TT)
   {
           R=getbygrs(M,1|mat=1);
           if(type(R)!=4||type(R[0])!=4||type(S=R[0][1])!=4){
                   errno(0);return0;
           }
           if(S[1]!=[[1,0]]){
                   print("Not rigid!");return0;
           }
           if((F=S[0][0][1])!=0){
                   for(V=vars(F);V!=[];V=cdr(V)){
                           if(mydeg(F,car(V))==1){
                                   T=lsol(F,car(V));
                                   break;
                           }
                   }
                   if(V==[]){
                           print("Violate Fuchs condition!");
                           return0;
                   }
           }
           for(P=[];R!=[];R=cdr(R))
                   P=cons(car(R)[0],P);
           if(F!=0){
                   S=mysubst(S,[car(V),T]);P=mysubst(P,[car(V),T]);
           }
           return mcgrs(S,P|mat=1,slm=[TT,[[1,0]]]);
   }
   
 /*  /*
   F=0 : unify    F=0 : unify
   F=["add",S] :    F=["add",S] :
Line 7997  def mcgrs(G, R)
Line 9195  def mcgrs(G, R)
   F=["put",F,V] :    F=["put",F,V] :
   F=["get1",F,V] :    F=["get1",F,V] :
   F=["put1",F,V] :    F=["put1",F,V] :
     F=["max"] :
     F=["max",F.V] :
   F=["put1"] :    F=["put1"] :
   F=["val",F];    F=["val",F];
   F=["swap"];    F=["swap"];
Line 8041  def anal2sp(R,F)
Line 9241  def anal2sp(R,F)
                 return G;                  return G;
         }          }
         if(F[0]=="add") return append(R,F[1]);          if(F[0]=="add") return append(R,F[1]);
           if(F[0]=="max"){
                   if(length(F)==3) C=1;
                   else C=0;
                   M=-10^10;K=[-1];
                   for(I=0;R!=[];R=cdr(R),I++){
                           if(C>0&&car(R)[F[1]]!=F[2]) continue;
                           if(M<car(R)[0]){
                                   M=car(R)[0];K=[I,car(R)];
                           }
                   }
                   return K;
           }
         R=reverse(R);          R=reverse(R);
         if(F[0]=="sub"){          if(F[0]=="sub"){
                 for(S=F[1];S!=[];S=cdr(S))                  for(S=F[1];S!=[];S=cdr(S))
Line 8053  def anal2sp(R,F)
Line 9265  def anal2sp(R,F)
                 return G;                  return G;
         }          }
         if(F[0]=="+"){          if(F[0]=="+"){
                 for(G=[];R!=[];R=cdr(R))                  L=length(F);
                         G=cons([car(R)[0],car(R)[1]+F[1],car(R)[2]+F[2]],G);                  for(G=[];R!=[];R=cdr(R)){
                           for(S=[],I=L-1;I>0;I--) S=cons(car(R)[I]+F[I],S);
                           G=cons(cons(car(R)[0],S),G);
                   }
                 return G;                  return G;
         }          }
         if(F[0]=="*"){          if(F[0]=="*"){
                 for(G=[];R!=[];R=cdr(R))                  L=length(F);
                         G=cons([car(R)[0],car(R)[1]*F[1]+car(R)[2]*F[2]],G);                  for(G=[];R!=[];R=cdr(R)){
                           for(S=0,I=1;I<L;I++) S+=car(R)[I]*F[I];
                           G=cons([car(R)[0],S],G);
                   }
                 return G;                  return G;
         }          }
         if(F[0]=="mult"){          if(F[0]=="mult"){
Line 8117  def anal2sp(R,F)
Line 9335  def anal2sp(R,F)
  P=["get",L]   P=["get",L]
     L=n        for variable x_n      L=n        for variable x_n
     L=[m,n]    for residue [m,n]      L=[m,n]    for residue [m,n]
       L=[m,n,l]  for residue [m,n,l]
     L=[[m,n],[m',n']] for common spct      L=[[m,n],[m',n']] for common spct
    P=["eigen",I]   decomposition of A_I
  P=["get0",[m,n],[m',n']] for the sum of residues   P=["get0",[m,n],[m',n']] for the sum of residues
    P=["rest",[m,n]] restriction
  P=["swap",[m,n]] for symmetry   P=["swap",[m,n]] for symmetry
  P=["perm",[...]] for symmetry   P=["perm",[...]] for symmetry
  P=["deg"]   P=["deg"]
Line 8188  def mc2grs(G,P)
Line 9409  def mc2grs(G,P)
         }          }
         if(type(P)<2) return G;          if(type(P)<2) return G;
         F=0;          F=0;
         if(type(P)==7||(type(P)==4&&type(P[0])<4)) P=[P];          if(type(P)==7||(type(P)==4&&
                   (type(P[0])<4||(type(P[0])==4&&length(P[0])==2&&type(P[0][0])<4&&type(P[1])<4))
             )) P=[P];
         if((Dvi=getopt(dviout))!=1&&Dvi!=2&&Dvi!=-1) Dvi=0;          if((Dvi=getopt(dviout))!=1&&Dvi!=2&&Dvi!=-1) Dvi=0;
         Keep=(Dvi==2)?1:0;          Keep=(Dvi==2)?1:0;
         if(type(P)==4&&type(F=car(P))==7){          if(type(P)==4&&type(F=car(P))==7){
Line 8218  def mc2grs(G,P)
Line 9441  def mc2grs(G,P)
                         return R;                          return R;
                 }                  }
                 if(F=="show0"){                  if(F=="show0"){
                           if(type(Fig=getopt(fig))>0){
                                   PP=[[-1.24747,-5.86889],[1.24747,-5.86889],[3.52671,-4.8541],[5.19615,-3],
                                     [5.96713,-0.627171],[5.70634,1.8541],[4.45887,4.01478],[2.44042,5.48127],
                                     [0,6],[-2.44042,5.48127],[-4.45887,4.01478],[-5.70634,1.8541],
                                     [-5.96713,-0.627171],[-5.19615,-3],[-3.52671,-4.8541]];
                                   PL=[[1.8,-5.2],[5.7,-1.7],[3.2,5],[-3.6,4.7],[2.2,3],[-2.8,2.8],
                                     [-1.5,-1.4],[-3.2,-2.5],[0.76,-1.4],[-2,0.2]];
                                   PC=["black,dashed","green,dashed","red,dashed","blue,dashed",
                                           "black","cyan","green","blue","red","magenta"];
                                   N=["1","2","3","4","5","6","7","8","9","a","b","c","d","e","f"];
                                   LL=[[1,2,3],[4,5,6],[7,8,9],[10,11,12],[7,10,13],[4,11,14],[5,8,15],[1,12,15],
                                     [2,9,14],[3,6,13]];
                                   TB=str_tb("\\draw\n",TB);
                                   if(type(Fig)==4){
                                           if(type(car(Fig))==1){
                                                   PP=ptaffine(car(Fig)/12,PP);PL=ptaffine(car(Fig)/12,PL);
                                                   Fig=cdr(Fig);
                                           }
                                           if(Fig!=[]&&length(Fig)==10) PC=Fig;
                                   }
                                   for(R=mc2grs(G,"show0"|dviout=-1),I=0;R!="";I++){       /* ’¸“_ */
                                           J=str_chr(R,0,",");
                                           if(J>0){
                                                   S=str_cut(R,0,J-1);
                                                   R=str_cut(R,J+1,1000);
                                           }else{
                                                   S=R;R="";
                                           }
                                           T=(str_chr(S,0,"1")==0)?"":"[red]";
                                           str_tb(["node",T,"(",N[I],") at ",xypos(PP[I]),"{$",S,"$}\n"],TB);
                                   }
                                   for(S=PC,P=PL,I=0;I<4;I++){
                                           for(J=I+1;J<5;J++,S=cdr(S),P=cdr(P)){                   /* ü‚̔Ԇ */
                                                   SS=car(S);
                                                   if((K=str_chr(SS,0,","))>0) SS=sub_str(SS,0,K-1);
                                                   str_tb(["node[",SS,"] at ",xypos(car(P)),
                                                           "{$[",rtostr(I),rtostr(J),"]$}\n"],TB);
                                           }
                                   }
                                   str_tb(";\n",TB);
                                   for(I=0;I<10;I++){              /* ü */
                                           S=car(PC);P0=car(PC);L0=car(LL);PC=cdr(PC);LL=cdr(LL);
                                           C=[N[L0[0]-1],N[L0[1]-1],N[L0[2]-1]];
                                           str_tb(["\\draw[",S,"] (", C[0],")--(",C[1],") (",
                                                   C[0],")--(",C[2],") (",C[1],")--(",C[2],");\n"],TB);
                                   }
                                   R=str_tb(0,TB);
                                   if(TikZ==1&&Dvi!=-1) dviout(xyproc(R)|dviout=1,keep=Keep);
                                   return R;
                           }
                         for(S="",L=[];G!=[];G=cdr(G)){                          for(S="",L=[];G!=[];G=cdr(G)){
                                 for(TL=[],TG=cdr(car(G));TG!=[];TG=cdr(TG)) TL=cons(car(TG)[0],TL);                                  for(TL=[],TG=cdr(car(G));TG!=[];TG=cdr(TG)) TL=cons(car(TG)[0],TL);
                                 TL=msort(TL,[-1,0]);                                  TL=msort(TL,[-1,0]);
Line 8295  def mc2grs(G,P)
Line 9568  def mc2grs(G,P)
                                                 else S="A_{"+rtostr(T[0][0])+rtostr(T[0][1])+"}&"+S;                                                  else S="A_{"+rtostr(T[0][0])+rtostr(T[0][1])+"}&"+S;
                                         }                                          }
                                         L=ltotex(R|opt="GRS",pre=S);                                          L=ltotex(R|opt="GRS",pre=S);
                                           if(type(D=getopt(div))==1 || type(D)==4) L=divmattex(L,D);
                                         if(Dvi>0) dviout(L|eq=0,keep=Keep);                                          if(Dvi>0) dviout(L|eq=0,keep=Keep);
                                 }                                  }
                                 return L; /* get all spct */                                  return L; /* get all spct */
Line 8306  def mc2grs(G,P)
Line 9580  def mc2grs(G,P)
                                         if(I[0]>I[0]){S=I;I=J;J=S;};                                          if(I[0]>I[0]){S=I;I=J;J=S;};
                                         K=lsort(I,J,0);                                          K=lsort(I,J,0);
                                         if(length(K)==4){                                          if(length(K)==4){
                                                 S=sp2grs(G,["get0",[I,J]]);                                                  S=mc2grs(G,["get0",[I,J]]);
                                                 return anal2sp(S,[["*",1,1],0]);                                                  return anal2sp(S,[["*",1,1],0]);
                                         }                                          }
                                         I=lsort(K,lsort(I,J,2),1);                                          I=lsort(K,lsort(I,J,2),1);
                                         S=lsort([0,1,2,3,4],K,1);                                          S=lsort([0,1,2,3,4],K,1);
                                         D=sp2grs(G,"deg");                                          D=mc2grs(G,"deg");
                                         if(findin(4,S)<0) D=-D;                                          if(findin(4,S)<0) D=-D;
                                         J=sp2grs(G,["get0",[I,S]]);                                          J=mc2grs(G,["get0",[I,S]]);
                                         if(I[0]>S[0]) J=sp2grs(J,"swap");                                          if(I[0]>S[0]) J=sp2grs(J,"swap");
                                         return anal2sp(J,[["+",0,D],["*",-1,1]]);                                          return anal2sp(J,[["+",0,D],["*",-1,1]]);
                                 }                                  }
Line 8325  def mc2grs(G,P)
Line 9599  def mc2grs(G,P)
                                                 if(car(PG)[0]==T) return (F=="get")?car(PG):cdr(car(PG));                                                  if(car(PG)[0]==T) return (F=="get")?car(PG):cdr(car(PG));
                                         return [];      /* get common spct */                                          return [];      /* get common spct */
                                 }                                  }
                                   if(length(T)==3){
                                           T0=T;T=lsort([0,1,2,3,4],T,1);
                                           if(length(T)!=2) return [];
                                   }else T0=0;
                                 if(T[0]>T[1]) T=[T[1],T[0]];                                  if(T[0]>T[1]) T=[T[1],T[0]];
                                 for(FT=0,PG=G;PG!=[];PG=cdr(PG)){                                  for(FT=0,PG=G;PG!=[];PG=cdr(PG)){
                                         if(car(PG)[0][0]==T){                                          if(car(PG)[0][0]==T){
Line 8336  def mc2grs(G,P)
Line 9614  def mc2grs(G,P)
                                 }                                  }
                                 if(!FT) return [];                                  if(!FT) return [];
                                 L=anal2sp(cdr(car(PG)),[["get1",FT],0]);                                  L=anal2sp(cdr(car(PG)),[["get1",FT],0]);
                                   if(T0!=0){
                                           if((K=mc2grs(G,"deg"))!=0){
                                                   if(T[1]!=4) K=-K;
                                                   R=reverse(L);
                                                   for(L=[];R!=[];R=cdr(R)) L=cons([car(R)[0],car(R)[1]+K],L);
                                           }
                                           T=T0;
                                   }
                                 return (F=="get")?cons(T,L):L;                                  return (F=="get")?cons(T,L):L;
                         }                          }
                 }                  }
                   if(F=="rest"||F=="eigen"||F=="rest0"||F=="rest1"){
                           if(F!="eigen") G=mc2grs(G,"homog");
                           if(length(P)==1){
                                   for(R=[],I=0;I<4;I++){
                                           for(J=I+1;J<5;J++){
                                                   S=mc2grs(G,[F,[I,J]]);
                                                   if(S!=[]) R=cons(cons([I,J],S),R);
                                           }
                                   }
                                   R=reverse(R);
                                   if(Dvi){
                                           TB=str_tb(0,0);
                                           if(F=="rest0"||F=="rest1"){
                                                   for(T=R;;){
                                                           TT=car(T);
                                                           S=rtostr(car(TT)[0])+rtostr(car(TT)[1]);
                                                           str_tb(["[",S,"]","&: "],TB);
                                                           for(TR=[],TT=cdr(TT);TT!=[];TT=cdr(TT))
                                                                   TR=cons(car(TT)[1],TR);
                                                           for(TR=qsort(TR);TR!=[];TR=cdr(TR))
                                                                   str_tb([s2sp(car(TR)|short=1,std=-1),"\\ \\ "],TB);
                                                           if((T=cdr(T))==[]) break;
                                                           str_tb("\\\\\n",TB);
                                                   }
                                           }else{
                                                   TB=str_tb(0,0);
                                                   for(T=R;;){
                                                           TT=car(T);
                                                           S=rtostr(car(TT)[0])+rtostr(car(TT)[1]);
                                                           str_tb(["[",S,"]",":\\ "],TB);
                                                           for(TR=[],TT=cdr(TT);;){
                                                                   T0=car(TT);
                                                                   str_tb(["&",my_tex_form(car(T0)),"&&\\to\\ \n",
                                                                           ltotex(cdr(T0)|opt="GRS")],TB);
                                                                   if((TT=cdr(TT))==[]) break;
                                                                   str_tb("\\\\\n",TB);
                                                           }
                                                           if((T=cdr(T))==[]) break;
                                                           str_tb("\\allowdisplaybreaks\\\\\n",TB);
                                                   }
                                           }
                                           R=texbegin("align*",str_tb(0,TB));
                                           if(Dvi!=-1) dviout(R|keep=Keep);
                                   }
                                   return R;
                           }
                           I=P[1];
                           if(I[0]>I[1]) I=[I[1],I[0]];
                           L=lsort([0,1,2,3,4],I,1);
                           if(F=="rest"&&length(P)==3){
                                   J=P[2];if(J[0]>J[1]) J=[J[1],J[0]];
                                   L=lsort(L,J,1);
                                   if(length(L)!=1) return 0;
                                   return [mc2grs(G,["get0",I]),mc2grs(G,["get0",[I[0],J[0]],[I[1],J[1]]]),
                                           mc2grs(G,["get0",[I[0],J[1]],[I[1],J[0]]]),mc2grs(G,["get0",[I[0],I[1],L[0]]])];
                           }
                           L=[[L[0],L[1]],[L[0],L[2]],[L[1],L[2]]];
                           if(F!="eigen"){
                                   if(I==[0,4]) L=reverse(L);
                                   else{
                                           for(V=[],J=2;J>=0;J--){
                                                   if(L[J][0]==0) V=cons([L[J][1],J],V);
                                                   else{
                                                           for(K=4;K>=0;K--){
                                                                   if(findin(K,L[J])<0){
                                                                           V=cons([K,J],V);break;
                                                                   }
                                                           }
                                                   }
                                           }
                                           V=qsort(V);
                                           L=[L[V[0][1]],L[V[1][1]],L[V[2][1]]];
                                   }
                           }
                           for(LL=[],T=L;T!=[];T=cdr(T))
                                   LL=cons(mc2grs(G,["get0",[I,car(T)]]),LL);
                           LL=reverse(LL);
                           for(R=[],Q=mc2grs(G,["get0",I]);Q!=[];Q=cdr(Q)){
                                   for(T=[],J=2;J>=0;J--){
                                           V=anal2sp(LL[J],["get1",(I[0]<L[J][0])?1:2,car(Q)[1]]);
                                           if(F=="rest"){
                                                   if(I[0]==0){
                                                           if(I[1]!=4){
                                                                   if(L[J][1]!=4) V=anal2sp(V,["+",-car(Q)[1]]);
                                                           }else if (L[J][0]!=2) V=anal2sp(V,["+",-car(Q)[1]]);
                                                   }else if(L[J][0]!=0) V=anal2sp(V,["+",-car(Q)[1]]);
                                           }
                                           T=cons(V,T);
                                   }
                                   R=cons(cons(car(Q)[1],T),R);
                           }
                           if(F=="rest0"||F=="rest1"){
                                   for(L=[];R!=[];R=cdr(R)){
                                           TR=cdr(car(R));
                                           if(F=="rest1"&&chkspt(TR|opt="idx")==2) continue;
                                           L=cons([car(R)[0],s2sp(chkspt(TR|opt=6))],L);
                                   }
                                   R=reverse(L);
                           }
                           return R;
                   }
                 if(F=="deg"){                  if(F=="deg"){
                         for(S=I=0;I<3;I++){                          for(S=I=0;I<3;I++){
                                 for(J=I+1;J<4;J++){                                  for(J=I+1;J<4;J++){
Line 8349  def mc2grs(G,P)
Line 9736  def mc2grs(G,P)
                         }                          }
                         return S/L[0];                          return S/L[0];
                 }                  }
                 if(F=="spct"){                  if(F=="spct"||F=="spct1"){
                           K=(F=="spct")?5:6;
                         G=mc2grs(G,"get");                          G=mc2grs(G,"get");
                         M=newmat(5,5);                          M=newmat(5,K);
                         for(;G!=[];G=cdr(G)){                          for(;G!=[];G=cdr(G)){
                                 GT=car(G);I=GT[0][0];J=GT[0][1];                                  GT=car(G);I=GT[0][0];J=GT[0][1];
                                 for(S=0,L=[],GT=cdr(GT);GT!=[];GT=cdr(GT)){                                  for(S=0,L=[],GT=cdr(GT);GT!=[];GT=cdr(GT)){
Line 8368  def mc2grs(G,P)
Line 9756  def mc2grs(G,P)
                                         for(L=M[I][J];L!=[];L=cdr(L)) S+=car(L)^2;                                          for(L=M[I][J];L!=[];L=cdr(L)) S+=car(L)^2;
                                 }                                  }
                                 M[I][I]=S;                                  M[I][I]=S;
                                   if(K==6){
                                           for(S=[],J=4;J>=0;J--)
                                                   if(I!=J) S=cons(M[I][J],S);
                                           R=chkspt(S|opt=2);
                                           M[I][5]=((L=length(R))>1)?s2sp(R[L-2]|short=1):"";
                                   }
                         }                          }
                         if(Dvi){                          if(Dvi){
                                 S=[];                                  S=[];
                                 for(I=4;I>=0;I--){                                  for(I=4;I>=0;I--){
                                         L=[M[I][I]];                                          L=(K==6)?[M[I][5]]:[];
                                           L=cons(M[I][I],L);
                                         for(J=4;J>=0;J--){                                          for(J=4;J>=0;J--){
                                                 if(I==J) L=cons("",L);                                                  if(I==J) L=cons("",L);
                                                 else L=cons(s2sp([M[I][J]]),L);                                                  else L=cons(s2sp([M[I][J]]),L);
                                         }                                          }
                                         S=cons(L,S);                                          S=cons(L,S);
                                 }                                  }
                                 S=cons([x0,x1,x2,x3,x4,"idx"],S);                                  T=(K==6)?["reduction"]:[];
                                 M=ltotex(S|opt="tab",hline=[0,1,z],vline=[0,1,z-1,z],left=["","$x_0$","$x_1$","$x_2$","$x_3$","$x_4$"]);                                  S=cons(append([x0,x1,x2,x3,x4,"idx"],T),S);
                                   M=ltotex(S|opt="tab",hline=[0,1,z],
                                           vline=(K==6)?[0,1,z-2,z-1,z]:[0,1,z-1,z],
                                           left=["","$x_0$","$x_1$","$x_2$","$x_3$","$x_4$"]);
                                 if(Dvi>0) dviout(M|keep=Keep);                                  if(Dvi>0) dviout(M|keep=Keep);
                         }                          }
                         return M;                          return M;
Line 8643  def mcmgrs(G,P)
Line 10041  def mcmgrs(G,P)
         Keep=(Dvi==2)?1:0;          Keep=(Dvi==2)?1:0;
         if(type(P)==4 && type(F=car(P))==7){          if(type(P)==4 && type(F=car(P))==7){
                 if(F=="mult"){                  if(F=="mult"){
                         for(P=cdr(P);P!=[];P=cdr(P)) G=os_md.mc2grs(G,car(P)|option_list=getopt());                          for(P=cdr(P);P!=[];P=cdr(P)) G=mc2grs(G,car(P)|option_list=getopt());
                         return G;                          return G;
                 }                  }
                 if(F=="get"||F=="get0"){                  if(F=="get"||F=="get0"){
Line 8756  def mcmgrs(G,P)
Line 10154  def mcmgrs(G,P)
                                 L=cons(TL,L);                                  L=cons(TL,L);
                         }                          }
                         if(Dvi){                          if(Dvi){
                                 if(Dvi!=-1) dviout(S|eq=0);                                  if(Dvi!=-1) dviout(S|eq=0,keep=Keep);
                                 return S;                                  return S;
                         }                          }
                         return reverse(L);                          return reverse(L);
Line 9167  def str_str(S,T)
Line 10565  def str_str(S,T)
                 }else if(type(S)==4){                  }else if(type(S)==4){
                         for(; J<=LE; S=cdr(S),J++){                          for(; J<=LE; S=cdr(S),J++){
                                 if(car(S) != LP){                                  if(car(S) != LP){
                                         if(SJIS && (V=S[J])>128){                                          if(SJIS && (V=car(S))>128){
                                                 if(V<160 || (V>223 && V<240)) J++;                                                  if((V<160 || (V>223 && V<240))&&S!=[]) {
                                                           J++;S=cdr(S);
                                                   }
                                         }                                          }
                                         continue;                                          continue;
                                 }                                  }
Line 9671  def my_tex_form(S)
Line 11071  def my_tex_form(S)
                 }                  }
                 SS = cons(S[I], SS);                  SS = cons(S[I], SS);
         }          }
           SS=str_subst(SS,"\n\\\\\n\\end{pmatrix}","\n\\end{pmatrix}"|raw=1);
         SS=str_subst(SS,"\\\\\n\\end{pmatrix}","\n\\end{pmatrix}"|raw=1);          SS=str_subst(SS,"\\\\\n\\end{pmatrix}","\n\\end{pmatrix}"|raw=1);
         Subst=getopt(subst);          Subst=getopt(subst);
         Sub0=["{asin}","{acos}","{atan}"];          Sub0=["{asin}","{acos}","{atan}"];
Line 9824  def divmattex(S,T)
Line 11225  def divmattex(S,T)
         if(length(L0)>0) L=cons(reverse(L0),L);          if(length(L0)>0) L=cons(reverse(L0),L);
         L=lv2m(reverse(L));     /* get matrix */          L=lv2m(reverse(L));     /* get matrix */
         if(T==0) return L;          if(T==0) return L;
           if(type(T)==1) T=[T];
         Size=size(L);S0=Size[0];          Size=size(L);S0=Size[0];
         if(type(T[0])!=4){          if(type(T[0])!=4){
                 S1=Size[1];                  S1=Size[1];
Line 10084  def tocsv(L)
Line 11486  def tocsv(L)
                 if(type(LT)==5) LT=vtol(LT);                  if(type(LT)==5) LT=vtol(LT);
                 if(type(LT)<4) LT=[LT];                  if(type(LT)<4) LT=[LT];
                 for(N=0; LT!=[]; LT=cdr(LT),N++){                  for(N=0; LT!=[]; LT=cdr(LT),N++){
                         if(N) str_tb(", ",Tb);                          if(N) str_tb(",",Tb);
                         if((T=car(LT))==Null) continue;                          if((T=car(LT))==Null) continue;
                         if(type(T)==7){                          if(type(T)==7){
                                 K=str_len(T);                                  K=str_len(T);
Line 10095  def tocsv(L)
Line 11497  def tocsv(L)
                 }                  }
                 str_tb("\n",Tb);                  str_tb("\n",Tb);
         }          }
         return str_tb(0,Tb);          S=str_tb(0,Tb);
           if(type(EXE=getopt(exe))!=1&&EXE!=0&&type(EXE)!=7) return S;
           if(type(F)!=7){
                   fcat(-1,0);
                   F="risaout";
                   if(EXE>=2&&EXE<=9) F+=rtostr(EXE);
                   F=DIROUTD+F+".csv";
           }else F=S;
           if(EXE!=0 && access(F)) remove_file(F);
           fcat(F,S|exe=1);
           return 1;
 }  }
   
 def readcsv(F)  def readcsv(F)
Line 10108  def readcsv(F)
Line 11520  def readcsv(F)
                 else if(type(V)==1) V=[V];                  else if(type(V)==1) V=[V];
                 else V=[];                  else V=[];
         }          }
           Eq=getopt(eq);
         Sp=getopt(sp);          Sp=getopt(sp);
         if(type(T=getopt(col))!=1) T=0;          if(type(T=getopt(col))!=1) T=0;
         Null=getopt(null);          Null=getopt(null);
         if(type(Null)<0) Null="";          if(type(Null)<0) Null=(Eq==1)?0:"";
         while((S=get_line(ID))!=0){          while((S=get_line(ID))!=0){
                 S=strtoascii(S);                  S=strtoascii(S);
                 N=length(S);                  N=length(S);
Line 10139  def readcsv(F)
Line 11552  def readcsv(F)
                                 LT=cons(C,LT);continue;                                  LT=cons(C,LT);continue;
                         }                          }
                         LS=asciitostr(reverse(LT));                          LS=asciitostr(reverse(LT));
                         if(V==1||findin(++J,V)>=0) LS=(isdecimal(LS))?eval_str(LS):((LS=="")?Null:LS);                          if(V==1||findin(++J,V)>=0){
                                   if(Eq==1) LS=(LS=="")?Null:eval_str(LS);
                                   else LS=(isdecimal(LS))?eval_str(LS):((LS=="")?Null:LS);
                           }
                         if(!T || T==J) LL=cons(LS,LL);                          if(!T || T==J) LL=cons(LS,LL);
                         if(F==-2) while(++I<N && Sp!=1 && S[I]!=44);                          if(F==-2) while(++I<N && Sp!=1 && S[I]!=44);
                         F=0;LT=[];                          F=0;LT=[];
                 }                  }
                 if(I<=N && (Sp!=1 || length(LT)>0)){ /* lastline */                  if(I<=N && (Sp!=1 || length(LT)>0)){ /* lastline */
                         LS=asciitostr(reverse(LT));                          LS=asciitostr(reverse(LT));
                         if(findin(++J,V)>=0) LS=(isdecimal(LS))?eval_str(LS):((LS=="")?Null:LS);                          if(V==1||findin(++J,V)>=0){
                                   if(Eq==1) LS=(LS=="")?Null:eval_str(LS);
                                   else LS=(isdecimal(LS))?eval_str(LS):((LS=="")?Null:LS);
                           }
                         if(!T || T==J) LL=cons(LS,LL);                          if(!T || T==J) LL=cons(LS,LL);
                 }                  }
                 L=cons(reverse(LL),L);                  L=cons(reverse(LL),L);
         }          }
         close_file(ID);          close_file(ID);
         if(T) L=m2l(L|flat=1);          if(T) L=m2l(L|flat=1);
         return reverse(L);          L=reverse(L);
           return L;
 }  }
   
   def getline(ID)
   {
           if(isint(Maxlen=getopt(Max))>0) Maxlen=1024;
           if(type(CR=getopt(CR))!=4) CR=[13];
           if(type(LF=getopt(LF))!=4) LF=[10];
           S=[];
           for(I=0; I<1023; I++){
                   C=get_byte(ID);
                   if(C<0) return 0;
                   if(findin(C,CR)>=0) continue;
                   if(findin(C,LF)>=0) break;
                   S=cons(C,S);
           }
           return asciitostr(reverse(S));
   }
   
 def showbyshell(S)  def showbyshell(S)
 {  {
         Id = getbyshell(S);          Id = getbyshell(S);
Line 10175  def getbyshell(S)
Line 11611  def getbyshell(S)
                 Tmp=str_subst(DIROUT,["%HOME%","%ASIRROOT%"],[Home,get_rootdir()]);                  Tmp=str_subst(DIROUT,["%HOME%","%ASIRROOT%"],[Home,get_rootdir()]);
         Sep=isMs()?"\\":"/";          Sep=isMs()?"\\":"/";
         F=Tmp+Sep+"muldif.tmp";          F=Tmp+Sep+"muldif.tmp";
         if(type(S)<=1 && S>=0)  close_file(Id);          if(type(S)<=1 && S>=0)  close_file(S);
         remove_file(F);          remove_file(F);
         if(type(S)<=1) return -1;          if(type(S)<=1) return -1;
         shell(S+" > \""+F+"\"");          shell(S+" > \""+F+"\"");
Line 10622  def mtotex(M)
Line 12058  def mtotex(M)
   
 def sint(N,P)  def sint(N,P)
 {  {
     if( type(N)==1 ) {      if( type(N)==1 || N==0 ) {
                 NT=ntype(N);                  NT=ntype(N);
                 if((type(Opt=getopt(str))==1 || Opt==0) && Opt>=0 && P>=0){                  if((type(Opt=getopt(str))==1 || Opt==0) && Opt>=0 && P>=0){
                         if(Opt==2 || Opt==4 || Opt==0){                          if(Opt==2 || Opt==4 || Opt==0){
                                 if(N==0) return "0";                                  if(N==0) return (Opt>0)?"0":0;
                                 Pw=0;                                  Pw=0;
                                 if(NT==4){                                  if(NT==4){
                                         NN=abs(real(N));N1=abs(imag(N));                                          NN=abs(real(N));N1=abs(imag(N));
Line 10670  def sint(N,P)
Line 12106  def sint(N,P)
                                 N=-N;                                  N=-N;
                                 Neg="-";                                  Neg="-";
                         }else Neg="";                          }else Neg="";
                           N=rint(N*10^P)/10^P;
                         NN=floor(N);                          NN=floor(N);
                           NV=(N-NN+1)*10^P;
                         NS=rtostr(NN);                          NS=rtostr(NN);
                         if(P<=0) return Neg+NS;                          if(P<=0) return Neg+NS;
                         if(NN==0 && getopt(zero)==0) NS="";                          if(NN==0 && getopt(zero)==0) NS="";
                         return Neg+NS+"."+str_cut(rtostr(rint((N-NN+1)*10^P)),1,P);                          return Neg+NS+"."+str_cut(rtostr(NV),1,P);
                 }                  }
                 if(NT==4)                  if(NT==4)
                         return sint(real(N),P)+sint(imag(N),P)*@i;                          return sint(real(N),P)+sint(imag(N),P)*@i;
         X = rint( N*10^P );          X = rint( N*10^P );
         return ((X+1.0)-1.0)/10^P;          return deval(X/10^P);
         }          }
     if( (type(N)==2) || (type(N)==3) ){      if( (type(N)==2) || (type(N)==3) ){
                 NN = eval(N);                  NN = eval(N);
Line 10701  def frac2n(N)
Line 12139  def frac2n(N)
         if((T=type(N))<0) return N;          if((T=type(N))<0) return N;
         E=(getopt(big)==1)?eval(@e):0.1;          E=(getopt(big)==1)?eval(@e):0.1;
         if(T==1){          if(T==1){
                 if(ntype(N)==0) return (E+N)-E;                  if(ntype(N)==0) return (E*N)/E;
                 else if(ntype(N)!=4) return N;                  else if(ntype(N)!=4) return N;
                 else return (E*(1+@i)+N)-E*(1+@i);                  else return (E*(1+@i)*N)/(E*(1+@i));
         }          }
         if(T==3||T==2){          if(T==3||T==2){
                 N=red(N);                  N=red(N);
Line 10711  def frac2n(N)
Line 12149  def frac2n(N)
                 for(S=0,I=mydeg(Nm,V);I>=0;I--) S+=frac2n(mycoef(Nm,I,V))*V^I;                  for(S=0,I=mydeg(Nm,V);I>=0;I--) S+=frac2n(mycoef(Nm,I,V))*V^I;
                 return S/dn(N);                  return S/dn(N);
         }          }
         if(T<4) return (N+E)-E;          if(T<4) return (E*N)/E;
 #ifdef  USEMODULE  #ifdef  USEMODULE
         return mtransbys(os_md.frac2n,N,[]|option_list=getopt());          return mtransbys(os_md.frac2n,N,[]|option_list=getopt());
 #else  #else
Line 10830  def xyline(P,Q)
Line 12268  def xyline(P,Q)
   
 def xylines(P)  def xylines(P)
 {  {
 /*      mycat([P,getopt()]);    */  
         Lf=getopt(curve);          Lf=getopt(curve);
         if(type(Lf)!=1) Lf=0;          if(type(Lf)!=1) Lf=0;
         SS=getopt(opt);          SS=getopt(opt);
Line 10989  def saveproc(S,Out)
Line 12426  def saveproc(S,Out)
         }          }
 }  }
   
   def xygrid(X,Y)
   {
           for(RR=[],I=0,Z=X;I<2;I++){
                   U=Z[2];L=LL=[];M=Z[3];
                   if(Z[1]==1||Z[1]==-1){
                           if(type(M)==4) L=M;
                           else{
                                   if(U*(-dlog(1-1/20)/dlog(10))>=M){
                                           L=cons([1,2,1/10],L);
                                           LL=cons([1,2,1/2],LL);
                                   }else if(U*(-dlog(1-1/10)/dlog(10))>=M)
                                           L=cons([1,2,1/5],L);
                                   else if(U*(-dlog(1-1/4)/dlog(10))>=M)
                                           L=cons([1,2,1/2],L);
                                   if(U*(-dlog(1-1/50)/dlog(10))>=M){
                                           L=cons([2,5,1/10],L);
                                           LL=cons([2,5,1/2],LL);
                                   }else if(U*(-dlog(1-1/25)/dlog(10))>=M)
                                           L=cons([2,5,1/5],L);
                                   else if(U*(-dlog(1-1/10)/dlog(10))>=M)
                                           L=cons([2,5,1/2],L);
                                   if(U*(-dlog(1-1/100)/dlog(10))>=M){
                                           L=cons([5,10,1/10],L);
                                           LL=cons([5,10,1/2],LL);
                                   }
                                   else if(U*(-dlog(1-1/50)/dlog(10))>=M)
                                           L=cons([5,10,1/5],L);
                                   else if(U*(-dlog(1-1/20)/dlog(10))>=M)
                                           L=cons([5,10,1/2],L);
                                   L=cons(L,cons(LL,[[[1,10,1]]]));
                           }
                           R=scale(L|scale=U);
                           if(Z[1]==-1){
                                   for(LL=[];R!=[];R=cdr(R)){
                                           for(L=[],T=car(R);T!=[];T=cdr(T)) L=cons(U-car(T),L);
                                           LL=cons(reverse(L),LL);
                                   }
                                   R=reverse(LL);
                           }
                   }else if(Z[1]==0){
                           if(type(M)==4){
                                   R=scale(M|f=x,scale=U);
                           }else{
                                   V=0;
                                   if(U/10>=M) V=1/10;
                                   else if(U/5>=M) V=1/5;
                                   else if(U/2>=M) V=1/2;
                                   R=[];
                                   if(V>0){
                                           UU=U*V;
                                           for(R=[],J=UU;J<U;J+=UU) R=cons(J,R);
                                   }
                                   if(V==1/10) L=[U/2];
                                   else L=[];
                                   R=cons(R,cons(L,[[0,U]]));
                           }
                   }else if(type(Z[1])==4){
                           R=Z[1];
                           if(length(R)==0||type(R[0])!=4) R=[[],[],R];
                   }else return 0;
                   K=length(R);
                   S=newvect(K);
                   for(J=0;J<K;J++){
                           for(S[J]=[],JJ=0;JJ<=Z[0];JJ+=U){
                                   for(P=R[J];P!=[];P=cdr(P))
                                           if(car(P)+JJ<=Z[0]) S[J]=cons(car(P)+JJ,S[J]);
                           }
                   }
                   for(J=0;J<K;J++) S[J]=lsort(S[J],[],1);
                   for(U=[],J=K-1;J>0;J--){
                           U=lsort(S[J],U,0);S[J-1]=lsort(S[J-1],U,1);
                   }
                   RR=cons(vtol(S),RR);
                   Z=Y;
           }
           if((Raw=getopt(raw))==1) return RR;
           SS=[];
           if(type(Sf=getopt(shift))==7){
                   Sx=Sf[0];Sy=Sf[1];
           }else Sx=Sy=0;
           for(I=0;I<2;I++){
                   for(S0=[],L=RR[I];L!=[];L=cdr(L)){
                           for(S=[],T=car(L);T!=[];T=cdr(T)){
                                   if(S!=[]) S=cons(0,S);
                                   if(I==0){
                                           S=cons([X[0]+Sx,car(T)+Sy],S);
                                           S=cons([Sx,car(T)+Sy],S);
                                   }else{
                                           S=cons([car(T)+Sx,Y[0]+Sy],S);
                                           S=cons([car(T)+Sx,Sy],S);
                                   }
                           }
                           S0=cons(S,S0);
                   }
                   SS=cons(reverse(S0),SS);
           }
           SS=reverse(SS);
           if(Raw==2) return SS;
           if(length(Y)<5) T=[["",""]];
           else if(type(Y[4])==4) T=[Y[4]];
           else T=[Y[4],Y[4]];
           if(length(X[4])==4) T=cons([""],T);
           else if(type(X[4])==4) T=cons(X[4],T);
           else T=cons([X[4]],T);
           for(Sx=Sy=[],I=0;I<2;I++){
                   TT=T[I];
                   for(V=SS[I];V!=[];V=cdr(V)){
                           Op=car(TT);
                           if(length(TT)>1) TT=cdr(TT);
                           if(car(V)==[]) continue;
                           if(Op=="") S=xylines(car(V));
                           else S=xylines(car(V)|opt=Op);
                           if(I==0) Sx=cons(S,Sx);
                           else Sy=cons(S,Sy);
                   }
           }
           for(S="",Sx=reverse(Sx), Sy=reverse(Sy);Sx!=[]&&Sy!=[];){
                   if(Sx!=[]){
                           S+=car(Sx);Sx=cdr(Sx);
                   }
                   if(Sy!=[]){
                           S+=car(Sy);Sy=cdr(Sy);
                   }
           }
           return S;
   }
   
   
   def addIL(I,L)
   {
           if(I==0){
                   for(R=[];L!=[];L=cdr(L)) R=addIL(car(L),R);
                   return reverse(R);
           }
           if(type(In=getopt(in))==1){
                   if(In==-1){
                           J=JJ=I[1];I=I[0];
                           for(R=[];L!=[];L=cdr(L)){
                                   J=lmin([car(L)[0],JJ]);
                                   if(J>I) R=cons([I,J],R);
                                   I=lmax([car(L)[1],I]);
                           }
                           if(I<JJ) R=cons([I,JJ],R);
                           return reverse(R);
                   }else{
                           for(;L!=[];L=cdr(L)){
                                   if(car(L)[0]>I) return 0;
                                   if(car(L)[1]>=I){
                                           if(In==3) return car(L);
                                           if(In==1||(I!=car(L)[0]&&I!=car(L)[1])) return 1;
                                           return 2;
                                   }
                           }
                           return 0;
                   }
           }
           I0=car(I);I1=I[1];
           for(F=0,R=[];L!=[];L=cdr(L)){
                   if(I0>car(L)[1]){
                           R=cons(car(L),R);
                           continue;
                   }
                   if(I0<=car(L)[1]){
                           I0=lmin([I0,car(L)[0]]);
                           if(I1<car(L)[0]){
                                   R=cons([I0,I1],R);
                                   for( ;L!=[];L=cdr(L)) R=cons(car(L),R);
                                   F=1;
                                   break;
                           }
                           I1=lmax([I1,car(L)[1]]);
                   }
           }
           if(!F) R=cons([I0,I1],R);
           return reverse(R);
   }
   
   def xy2curve(F,N,Lx,Ly,Lz,A,B)
   {
           Raw=getopt(raw);
           if(type(Gap=getopt(gap))==4){
                   MG=Gap[1];Gap=car(Gap);
           }else MG=3;
           if(type(Gap)!=1 && Gap!=0) Gap=0.7;
           if(type(Dvi=getopt(dviout))<1) Dvi=0;
           OL=[["dviout",Dvi]];
           if(type(Opt=getopt(opt))<1) Opt=0;
           else OL=cons(["opt",Opt],OL);
           if(type(Sc=getopt(scale))!=1 && type(Sc)!=4) Sc=[1,1,1];
           else if(type(Sc)!=4) Sc=[Sc,Sc,Sc];
           else if(length(Sc)!=3) Sc=[Sc[0],Sc[1],Sc[1]];
           M=diagm(3,Sc);
           if(A!=0||B!=0){
                   if(type(A)==6) M=A;
                   else M=mrot([0,-B,-A]|deg=1)*M;
                   V=M*newvect(3,[x,y,z]);
                   Fx=compdf(V[0],[x,y,z],F);Fy=compdf(V[1],[x,y,z],F);Fz=compdf(V[2],[x,y,z],F);
           }else{
                   for(I=0;I<3;I++){
                           if(type(T=F[I])!=4) T=f2df(T);
                           if(type(T)==4) T=cons(car(T)*Sc[I],cdr(T));
                           else T*=Sc[I];
                           if(I==0) Fx=T;
                           else if(I==1) Fy=T;
                           else Fz=T;
                   }
           }
           if(Raw==5||!Gap)
                   return (Dvi||!Gap)? xygraph([Fy,Fz],N,Lx,Ly,Lz|option_list=OL):[Fx,Fy,Fz];
           R=xygraph([Fy,Fz],N,Lx,Ly,Lz|raw=2);
           R0=cdr(car(R));R1=R[1];
           for(LT=[];R0!=[];R0=cdr(R0),R1=cdr(R1))
                   if(car(R0)!=0) LT=cons([R1[0],R1[1]],LT);
           LT=reverse(LT);
           if(N<0){
                   Be=xylines(car(R)|curve=1,proc=3,close=-1);
                   LT=reverse(cdr(LT));
                   LT=reverse(cdr(LT));
           }
           else Be=xylines(car(R)|curve=1,proc=3);
           Be=cdr(cdr(Be));
           Be=lbezier(car(Be));
           if(Raw==4) return [Be,LT,Lx];
           X=ptcombz(Be,0,0);
           Var=(length(Lx)==3)?car(Lx):x;
           if(type(Eq=getopt(eq))!=1) Eq=0.01;
           if(TikZ==1){
                   Gap/=10;Eq/=10;
           }
           for(R=[],XT=X;XT!=[];XT=cdr(XT)){
                   V=car(XT);
                   U=LT[V[0][0]];
                   T=U[0]*V[1][0]+U[1]*(1-V[1][0]);
                   VV=myfdeval(Fx,[Var,T]);
                   U=LT[V[0][1]];
                   T=U[0]*V[1][1]+U[1]*(1-V[1][1]);
                   VV-=myfdeval(Fx,[Var,T]);
                   if(abs(VV)<Eq) continue;
                   I=(VV<0)?0:1;
                   R=cons([V[0][I],V[1][I],V[0][1-I],V[1][1-I]],R);
           }
           R=qsort(R);
           if(Raw==3) return [Be,R];
       Db=newvect(L=length(Be));
           for(I=0;I<L;I++) Db[I]=[];
           for(TR=R;TR!=[];TR=cdr(TR)){
           V1=ptbezier(Be,[I=car(TR)[0],P=car(TR)[1]])[1];
                   V2=ptbezier(Be,[car(TR)[2],car(TR)[3]])[1];
                   T=dsqrt(1-dvangle(V1,V2)^2);
                   if(T<1/MG) T=MG;
                   GP=Gap/T;
                   W=GP/dnorm(V1);
                   Db[I]=addIL([P-W,P+W],Db[I]);
                   if(P-W<0 && I>0) Db[I-1]=addIL([P-W+1,1],Db[I-1]);
                   if(P+W>1 && I+1<L) Db[I+1]=addIL([0,P+W-1],Db[I+1]);
           }
           Db=vtol(Db);
           for(Bf=[];Be!=[];Be=cdr(Be),Db=cdr(Db)){
                   if(car(Db)==[]) Bf=cons(car(Be),Bf);
                   else{
                           D=addIL([0,1],car(Db)|in=-1);
                           for(;D!=[];D=cdr(D))
                                   Bf=cons(tobezier(car(Be)|inv=car(D)),Bf);
                   }
           }
           Bf=reverse(Bf);
           if(Raw==2) return Bf;
           OL=[];
           if(Opt){
                   if(type(Opt)==4&&length(Opt)>1) OL=[["opt",Opt[0]],["cmd",Opt[1]]];
                   else OL=[["opt",Opt]];
           }else OL=[];
           S=xybezier(lbezier(Bf|inv=1)|option_list=OL);
           if(Raw==1||!Dvi) return S;
           return xyproc(S|dviout=Dvi);
   }
   
   def rungeKutta(F,N,Lx,Y,IY)
   {
           if((Pr=getopt(prec))==1){
                   One=eval(exp(0));
           }else{
                   One=deval(exp(0));Pr=0;
           }
           if(!isint(FL=getopt(mul))||!FL) FL=1;
           if(length(Lx)>2){
                   V=car(Lx);Lx=cdr(Lx);
           }else V=x;
           if(Pr==1) Lx=[eval(Lx[0]),eval(Lx[1])];
           else Lx=[deval(Lx[0]),deval(Lx[1])];
           if(type(Y)==4){
                   if((Sing=getopt(single))==1||type(F)!=4)
                           F=append(cdr(Y),[F]);
                   L=length(Y);
                   for(TF=[];F!=[];F=cdr(F))
                           TF=cons(f2df(car(F)),TF);
                   F=reverse(TF);
           }else{
                   L=1;
                   F=f2df(F);
           }
           if(getopt(val)==1) V1=1;
           else V1=0;
           if(FL>0) N*=FL;
           H=(Lx[1]-Lx[0])/N*One;H2=H/2;
           FV=findin(V,vars(F));
           K=newvect(4);
           if(L==1){
                   R=[[T=Lx[0],S=IY]];
                   if(!H) return R;
                   for(C=0;C<N;C++){
                           for(I=0;I<4;I++){
                                   if(I==0)      W=[[V,T],[Y,S]];
                                   else if(I==3) W=[[V,T+H],[Y,S+H*K[2]]];
                                   else          W=[[V,T+H2],[Y,S+H2*K[I-1]]];
                                   if(FV<0) W=cdr(W);
                                   K[I]=Pr?myfeval(F,W)*One:myfdeval(F,W);
                           }
                           S+=(K[0]+2*K[1]+2*K[2]+K[3])*H/6;T+=H;
                           if(FL>0&&!((C+1)%FL)) R=cons([deval(T),S],R);
                   }
           }else{
                   T=Lx[0];
                   R=[cons(T,V1?[car(IY)]:IY)];
                   S=ltov(IY);
                   if(!H) return R;
                   for(C=0;C<N;C++){
                           for(I=0;I<4;I++){
                                   if(I==0)      W=cons([V,T   ],lpair(Y,vtol(S)));
                                   else if(I==3) W=cons([V,T+H ],lpair(Y,vtol(S+H*K[2])));
                                   else          W=cons([V,T+H2],lpair(Y,vtol(S+H2*K[I-1])));
                                   if(FV<0) W=cdr(W);
                                   for(TK=[],TF=F;TF!=[];TF=cdr(TF)){
                                           TK=cons(Pr?myfeval(car(TF),W)*One:myfdeval(car(TF),W),TK);
                                   }
                                   K[I]=ltov(reverse(TK));
                           }
                           S+=(K[0]+2*K[1]+2*K[2]+K[3])*H/6;T+=H;
                           TS=vtol(S);
                           if(FL<0||(C+1)%FL) continue;
                           if(V1) TS=[car(TS)];
                           R=cons(cons(deval(T),TS),R);
                   }
           }
           L=(FL<0)?(V1?S[0]:S):reverse(R);
           return L;
   }
   
   def pwTaylor(F,N,Lx,Y,Ly,M)
   {
           if(!isint(FL=getopt(mul))||!FL) FL=1;
           if(getopt(val)==1) V1=1;
           else V1=0;
           if(isint(Er=getopt(er))&&Er>0){
                   Opt=delopt(getopt(),["er","mul"]);
                   L0=pwTaylor(F,N,Lx,Y,Ly,M|option_list=cons(["mul",FL*(Er+1)],Opt));
           }else Er=0;
           if(length(Lx)>2){
                   V=car(Lx);Lx=cdr(Lx);
           }else V=t;
           if(!isvar(T=getopt(var))) V=t;
           if((Pr=getopt(prec))==1){
                   One=eval(exp(0));
           }else{
                   One=deval(exp(0));Pr=0;
           }
           if(Pr==1) Lx=[eval(Lx[0]),eval(Lx[1])];
           else Lx=[deval(Lx[0]),deval(Lx[1])];
           if(type(Y)==4){
                   if(type(F)!=4)  F=append(cdr(Y),[F]);
           }else Y=[Y];
           if(type(Ly)!=4) Ly=[Ly];
           if(FL>0) N*=FL;
           H=(Lx[1]-Lx[0])/N*One;
           S=vtol(pTaylor(F,Y,M|time=V));
           LS=length(S);
           if(type(Vw=getopt(view))==4){
                   glib_window(car(Vw)[0], car(Vw)[1],car(Vw)[2],car(Vw)[3]);
                   if(length(Vw)>1 && (C=trcolor(Vw[1]))!=0) Opt=[["color",C]];
                   else Opt=[];
                   if(length(Vw)>2){
                           Mt=Vw[2];
                           if(LS==1){
                                   if(type(Mt)>1) Mt=0;
                           }else{
                                   if(type(Mt)!=6||size(Mt)!=[]) Mt=0;
                           }
                   }else Mt=0;
                   if(!Mt){
                           if(LS>1){
                                   Mt=newmat(2,LS);Mt[0][0]=Mt[1][1]=1;
                           }else Mt=1;
                           if(LS==1) glib_putpixel(Lx[0],Mt*Ly[0]|option_list=Opt);
                           else{
                                   YT=ptaffine(Mt,Ly);
                                   glib_putpixel(YT[0],YT[1]|option_list=Opt);
                           }
                   }
           }else Vw=0;
           R=[cons(T=Lx[0],Ly)];
           for(C=0,T+=H;C<N;T+=H,C++){
                   if(!C) SS=subst(S,V,H);
                   for(Dy=SS,TY=Y,TL=Ly;TY!=[];TY=cdr(TY),TL=cdr(TL))
                           Dy=subst(Dy,car(TY),One*car(TL));
   /*              Ly=subst(Dy,V,H); */
                   Ly=Dy;
                   if(FL<0||(C+1)%FL) continue;
                   if(Vw){
                           if(LS==1) glib_putpixel(deval(T),Mt*Ly[0]|option_list=Opt);
                           else{
                                   YT=ptaffine(Mt,Ly);
                                   glib_putpixel(YT[0],YT[1]|option_list=Opt);
                           }
                           continue;
                   }
                   TR=cons(deval(T),(V1)?[car(Ly)]:Ly);
                   R=cons(TR,R);
           }
           if(Vw) return 1;
           L=(FL<0)?((V1)?car(Ly):Ly):reverse(R);
           if(Er){
                   if(FL>0){
                           for(S=L,T=L0,D=[];S!=[];S=cdr(S),T=cdr(T)) D=cons(os_md.ladd(car(S),car(T),-1),D);
                           E=map(os_md.dnorm,reverse(D));F=map(os_md.nlog,E);
                   }else if(V1){
                           D=ladd(L,L0,-1);F=nlog(dnorm(D));
                   }else F=nlog(abs(L-L0));
                   return [L,F];
           }
           return L;
   }
   
 def xy2graph(F0,N,Lx,Ly,Lz,A,B)  def xy2graph(F0,N,Lx,Ly,Lz,A,B)
 {  {
         /* (x,y,z) -> ( -x sin A + y cos A, z cos B - x cos A sin B - y sin A sin B) */          /* (x,y,z) -> (z sin B + x cos A cos B + y sin A cos B,
               -x sin A + y cos A, z cos B - x cos A sin B - y sin A sin B) */
         if((Proc=getopt(proc))==1||Proc==2){          if((Proc=getopt(proc))==1||Proc==2){
                 OPT0=[["proc",3]];                  OPT0=[["proc",3]];
         }else{          }else{
Line 11425  def xy2graph(F0,N,Lx,Ly,Lz,A,B)
Line 13295  def xy2graph(F0,N,Lx,Ly,Lz,A,B)
         if(Dvi<0) return Lout;          if(Dvi<0) return Lout;
 }  }
   
   def orthpoly(N)
   {
           F=0;
           if(type(P=getopt(pol))==7){
                   for(L=["Le","Ge","Tc","2T","Ja","He","La","Se"];L!=[];L=cdr(L),F++)
                           if(str_str(P,car(L)|end=2)==0) break;
           }else P=0;
           if(type(D=N)==4) D=N[0];
           if(!isint(D)||D<0) return 0;
           if(F==0) return seriesHG([-D,D+1],[1],(1-x)/2,D);
           if(F==1) return red(seriesHG([-D,D+2*N[1]],[N[1]+1/2],(1-x)/2,D)*binom(D+2*N[1]-1,D));
           if(F==2) return seriesHG([-D,D],[1/2],(1-x)/2,D);
           if(F==3){
                   if(D==0) return 0;
                   return orthpoly([D-1,1]|pol="Ge");
           }
           if(F==4) return red(seriesHG([-D,D+N[1]],[N[2]],x,D));
           if(F==5){
                   for(S=I=1;I<=D;I+=2) S*=I;
                   if(iand(D,1)) return seriesHG([-(D-1)/2],[3/2],x^2/2,D-1)*x*S*(-1)^((D-1)/2);
                   else return seriesHG([-D/2],[1/2],x^2/2,D)*S*(-1)^(D/2);
           }
           if(F==6){
                   NN=(type(N)==4)?N[1]:0;
                   return red(seriesHG([-D],[NN+1],x,D)*binom(D+NN,D));
           }
           if(F==7){
                   NN=N[1];
                   for(S=1,I=1;I<=D;I++) S+=(-1)^I*binom(D,I)*binom(D+I,I)*sftpow(x,I)/sftpow(NN,I);
                   return S;
           }
           return 0;
   }
   
   def schurpoly(L)
   {
           N=length(L);
           for(R=[],I=1;L!=[];L=cdr(L),I++) R=cons(car(L)+N-I,R);
           L=reverse(R);
           if(type(X=getopt(var))!=4){
                   V=(type(X)>1)?X:"x";
                   for(X=[],I=0;I<N;I++) X=cons(makev([V,N-I]),X);
           }
           M=newmat(N,N);
           for(I=0;I<N;I++)
                   for(J=0;J<N;J++) M[I][J]=X[I]^L[J];
           P=det(M);
           for(I=0;I<N;I++)
                   for(J=I+1;J<N;J++) P=sdiv(P,X[I]-X[J]);
           return P;
   }
   
 def fouriers(A,B,X)  def fouriers(A,B,X)
 {  {
           if((Y=getopt(y))==0||type(Y)>0) Y=deval(Y);
           else Y=0;
           if((V=getopt(const))==0||type(V)>0){
                   V=myfeval(V,Y);
                   K=1;
           }else K=0;
         if(A!=[]&&type(car(A))>1){          if(A!=[]&&type(car(A))>1){
                 for(C=[],I=A[1];I>=0;I--) C=cons(myfeval(car(A),I),C);                  for(C=[],I=A[1];I>=K;I--) C=cons(myf2eval(car(A),I,Y),C);
                   if(K) C=cons(0,C);
                 A=C;                  A=C;
         }          }
           if(K){
                   if(A!=[]) A=cdr(A);
                   A=cons(V,A);
           }
         if(B!=[]&&type(car(B))>1){          if(B!=[]&&type(car(B))>1){
                 for(C=[],I=B[1];I>0;I--) C=cons(myfeval(car(B),I),C);                  for(C=[],I=B[1];I>0;I--) C=cons(myf2eval(car(B),I,Y),C);
                 B=C;                  B=C;
         }          }
         R=0;          L=length(B)+1;
           if(length(A)>=L) L=length(A)+1;
           if(type(Sum=getopt(sum))>0){
                   if(Sum==1) Sum=1-x;
                   else if(Sum==2) Sum=[(z__)/(3.1416*x),[z__,os_md.mysin,3.1416*x]];
                   else Sum=f2df(Sum);
                   C=[];
                   if(A!=[]){
                           C=cons(car(A),C);
                           A=cdr(A);
                   }
                   for(I=1;A!=[];A=cdr(A),I++) C=cons(car(A)*myf2eval(Sum,I/L,L),C);
                   A=reverse(C);
                   for(C=[],I=1;B!=[];B=cdr(B),I++) C=cons(car(B)*myf2eval(Sum,I/L,L),C);
                   B=reverse(C);
           }
         if(getopt(cpx)==1){          if(getopt(cpx)==1){
         if(type(X=eval(X))>1) return todf([os_md.fouriers,[["cpx",1]]],[[A],[B],[X]]);                  if(type(X=eval(X))>1) return todf([os_md.fouriers,[["cpx",1]]],[[A],[B],[X]]);
                 V=dexp(@i*X);                  V=dexp(@i*X);
                 for(C=A,P=1,I=0;C!=[];C=cdr(C),I++){                  for(C=A,P=1,I=0;C!=[];C=cdr(C),I++){
                         R+=car(C)*P;                          R+=S*car(C)*P;
                         P*=V;                          P*=V;
                 }                  }
                 V=dexp(-@i*X);                  V=dexp(-@i*X);
Line 11485  def mysin(Z)
Line 13433  def mysin(Z)
 def mytan(Z)  def mytan(Z)
 {  {
         if(type(Z=eval(Z))>1) return todf(os_md.mytan,[Z]);          if(type(Z=eval(Z))>1) return todf(os_md.mytan,[Z]);
         if((Im=imag(Z))==0) return dsin(Z);          if((Im=imag(Z))==0) return dtan(Z);
         V=myexp(2*Z*@i);          V=myexp(2*Z*@i);
         return @i*(1-V)/(1+V);          return @i*(1-V)/(1+V);
 }  }
Line 11493  def mytan(Z)
Line 13441  def mytan(Z)
 def mylog(Z)  def mylog(Z)
 {  {
         if(type(Z=eval(Z))>1) return todf(os_md.mylog,[Z]);          if(type(Z=eval(Z))>1) return todf(os_md.mylog,[Z]);
         if((Im=imag(Z))==0) return dlog(Z);          if(imag(Z)==0&&Z>=0) return dlog(Z);
         return dlog(dabs(Z))+@i*myarg(Z);          return dlog(dabs(Z))+@i*myarg(Z);
 }  }
   
   def nlog(X)
   {
           return mylog(X)/dlog(10);
   }
   
 def mypow(Z,R)  def mypow(Z,R)
 {  {
         if(type(Z=eval(Z))>1||type(R=eval(R))>1) return todf(os_md.mypow,[Z,R]);          if(type(Z=eval(Z))>1||type(R=eval(R))>1) return todf(os_md.mypow,[Z,R]);
Line 11956  def compdf(F,V,G)
Line 13909  def compdf(F,V,G)
 {  {
         FL=["abs","floor","rint","zeta","gamma","arg","real","imag","conj"];          FL=["abs","floor","rint","zeta","gamma","arg","real","imag","conj"];
         FS=[os_md.abs,floor,rint,os_md.zeta,os_md.gamma,os_md.myarg,real,imag,conj];          FS=[os_md.abs,floor,rint,os_md.zeta,os_md.gamma,os_md.myarg,real,imag,conj];
         if(type(V)==4){  
                 for(;V!=[];V=cdr(V),G=cdr(G)) F=compdf(F,car(V),car(G));  
                 return F;  
         }  
         if(type(F)==7){          if(type(F)==7){
                 if(str_str(F,"|")==0){                  if(str_str(F,"|")==0){
                         F="abs("+str_cut(F,1,str_len(F)-2)+")";                          F="abs("+str_cut(F,1,str_len(F)-2)+")";
Line 11980  def compdf(F,V,G)
Line 13929  def compdf(F,V,G)
         }          }
         if(type(F)!=4) F=f2df(F);          if(type(F)!=4) F=f2df(F);
         if(type(G)!=4) G=f2df(G);          if(type(G)!=4) G=f2df(G);
           if(V==G) return F;      /* subst(F(V),V,G) */
         VF=vars(F);VG=vars(G);          VF=vars(F);VG=vars(G);
           if(type(V)==4){
                   for(VT=[],VV=V;VV!=[];VV=cdr(VV)){
                           if(findin(car(VV),VF)>=0){
                                   X=makenewv(append(VF,VG));
                                   VF=cons(X,VF);
                                   F=mysubst(F,[car(VV),X]);
                                   VT=cons(X,VT);
                           }else VT=cons(car(VV),VT);
                   }
                   for(V=reverse(VT);V!=[];V=cdr(V),G=cdr(G)) F=compdf(F,car(V),car(G));
                   return F;
           }
         for(E=I=0;I<30;I++){          for(E=I=0;I<30;I++){
                 for(J=0;J<30;J++){                  for(J=0;J<30;J++){
                         X=makev(["z__",I,J]);                          X=makev(["z__",I,J]);
Line 11991  def compdf(F,V,G)
Line 13953  def compdf(F,V,G)
                 if(E) break;                  if(E) break;
         }          }
         if(!E) return 0;          if(!E) return 0;
         if(V==G) return F;      /* subst(F(V),V,G) */  
         if(type(G)<4) return mysubst(F,[V,G]);          if(type(G)<4) return mysubst(F,[V,G]);
         if(type(F)<4) F=[F]; /* return compdf([X,[X,0,F]],V,G); */          if(type(F)<4) F=[F]; /* return compdf([X,[X,0,F]],V,G); */
         F=mysubst(F,[V,X]);          F=mysubst(F,[V,X]);
Line 12255  def fcont(F,LX)
Line 14216  def fcont(F,LX)
         return reverse(L);          return reverse(L);
 }  }
   
   def xyplot(L,LX,LY)
   {
           LX=map(deval,LX);LY=map(deval,LY);
           Opt=getopt();Opt0=delopt(Opt,["dviout","proc"]);
           for(S="",L0=[],TL=L;TL!=[];TL=cdr(TL)){
                   TTL=map(deval,car(TL));
                   if(TTL[0]<LX[0]||TTL[0]>LX[1]||TTL[1]<LY[0]||TTL[1]>LY[1]){
                           S+=xylines(reverse(L0)|option_list=Opt0);
                           L0=[];
                   }else{
                           L0=cons(TTL,L0);
                   }
           }
           if(length(L0)>1) S+=xylines(reverse(L0)|option_list=Opt0);
           if(type(AX=getopt(ax))==4) S+=xygraph([0,0],0,LX,LX,LY|option_list=delopt(Opt0,"opt"));
           if(getopt(dviout)!=1) return S;
           xyproc(S|dviout=1);
   }
   
 def xygraph(F,N,LT,LX,LY)  def xygraph(F,N,LT,LX,LY)
 {  {
         if((Proc=getopt(proc))!=1&&Proc!=2&&Proc!=3) Proc=0;          if((Proc=getopt(proc))!=1&&Proc!=2&&Proc!=3) Proc=0;
Line 12439  def xygraph(F,N,LT,LX,LY)
Line 14419  def xygraph(F,N,LT,LX,LY)
                 }                  }
                 V=reverse(NV);                  V=reverse(NV);
         }          }
         if(getopt(raw)==1) return V;          if((Raw=getopt(raw))==1) return V;
           if(Raw==2) return [V,LT];
         OL=[["curve",1]];OLP=[];          OL=[["curve",1]];OLP=[];
         if(type(C=getopt(ratio))==1){          if(type(C=getopt(ratio))==1){
                 OL=cons(["ratio",C],OL);OLP=cons(["ratio",C],OLP);                  OL=cons(["ratio",C],OL);OLP=cons(["ratio",C],OLP);
Line 12566  def xygraph(F,N,LT,LX,LY)
Line 14547  def xygraph(F,N,LT,LX,LY)
                         if(length(Ax)>3){                          if(length(Ax)>3){
                                 D=Ax[3];                                  D=Ax[3];
                                 if(type(D)==1 && D>0){                                  if(type(D)==1 && D>0){
                                         I0=ceil((LY[0]-Ax[1])/D); I1=floor((LY[1]-Ax[0])/D);                                          I0=ceil((LY[0]-Ax[1])/D); I1=floor((LY[1]-Ax[1])/D);
                                         for(DD=[],I=I0; I<=I1; I++){                                          for(DD=[],I=I0; I<=I1; I++){
                                                 if(length(Ax)<5) DD=cons(I*D,DD);                                                  if(length(Ax)<5) DD=cons(I*D,DD);
                                                 else if(I!=0){                                                  else if(I!=0){
Line 12697  def polroots(L,V)
Line 14678  def polroots(L,V)
         Lim=Lim2=[];          Lim=Lim2=[];
         if(type(L)<4){          if(type(L)<4){
                 if(type(Lim=getopt(lim))==4){                  if(type(Lim=getopt(lim))==4){
                         if(type(Lim[0])!=4) Lim=[Lim];                          if(type(Lim[0])!=4){
                                   if(!isvar(Lim[0])) Lim=cons(V,[Lim]);
                                   Lim=[Lim];
                           }
                           if(!isvar(Lim[0][0])) Lim=[cons(V,Lim)];
                         Lim=delopt(Lim,V|inv=1);                          Lim=delopt(Lim,V|inv=1);
                         if(Lim!=[]){                          if(Lim!=[]){
                                 Lim=Lim[0];                                  Lim=Lim[0];
Line 12759  def polroots(L,V)
Line 14744  def polroots(L,V)
         if(SS==0&&INIT==1){          if(SS==0&&INIT==1){
                 SS=polroots(L,V|option_list=OL);                  SS=polroots(L,V|option_list=OL);
                 if(SS!=0) return SS;                  if(SS!=0) return SS;
                 for(C=0;SS==0&&C<4;C++){                  for(C=0;SS==0&&C<5;C++){
                         I=(C==0)?1:(iand(random(),0xff)-0x80);                          I=(C==0)?1:(iand(random(),0xff)-0x80);
                         for(LL=[],K=length(L)-1;K>=0;K--){                          for(LL=[],K=length(L)-1;K>=0;K--){
                                 for(Q=0,J=length(L)-1;J>=0;J--)                                  for(Q=0,J=length(L)-1;J>=0;J--)
Line 12785  def polroots(L,V)
Line 14770  def polroots(L,V)
         for(SS=[];R!=[];R=cdr(R)){          for(SS=[];R!=[];R=cdr(R)){
                 RS=(N==2)?[car(R)]:car(R);                  RS=(N==2)?[car(R)]:car(R);
                 for(I=0,L0=L[0];I<N-1;I++) L0=mysubst(L0,[V1[I],RS[I]]);                  for(I=0,L0=L[0];I<N-1;I++) L0=mysubst(L0,[V1[I],RS[I]]);
                   if(L0==0) return 0;
                 S0=polroots(L0,V[0]|option_list=OL);                  S0=polroots(L0,V[0]|option_list=OL);
                 if(type(S0)<2) return S0;                  if(type(S0)<2) return S0;
                 for(S=S0;S!=[];S=cdr(S)){                  for(S=S0;S!=[];S=cdr(S)){
Line 12977  def cutf(F,X,VV)
Line 14963  def cutf(F,X,VV)
                         if(car(V)!=[] && car(V)[0]<X) return myfeval(car(V)[1],Y);                          if(car(V)!=[] && car(V)[0]<X) return myfeval(car(V)[1],Y);
                         return myfeval(F,Y);                          return myfeval(F,Y);
                 }                  }
                 if(X>car(V)[0]) continue;                  if(car(V)==[]||X>car(V)[0]) continue;
                 if(X==car(V)[0]) return car(V)[1];                  if(X==car(V)[0]) return car(V)[1];
                 return myfeval(F,Y);                  return myfeval(F,Y);
         }          }
 }  }
   
 def fsum(F,L,X)  def fsum(F,L)
 {  {
           if(getopt(df)==1){
                   F=f2df(F);
           }else Sub=getopt(subst);
         if(type(L[0])==2){          if(type(L[0])==2){
                 X=L[0];                  X=L[0];
                 L=cdr(L);                  L=cdr(L);
Line 12992  def fsum(F,L,X)
Line 14981  def fsum(F,L,X)
         V=(length(L)>2)?L[2]:1;          V=(length(L)>2)?L[2]:1;
         for(R=0,I=L[0];;I+=V){          for(R=0,I=L[0];;I+=V){
                 if(V==0||(I-L[1])*V>0) return R;                  if(V==0||(I-L[1])*V>0) return R;
                 R+=os_md.myfeval(F,X?[X,I]:I);                  R+=(Sub==1)?subst(F,X?X:x,I):os_md.myfeval(F,X?[X,I]:I);
         }          }
 }  }
   
Line 13001  def periodicf(F,L,X)
Line 14990  def periodicf(F,L,X)
         if(type(L)==4) L=[eval(L[0]),eval(L[1])];          if(type(L)==4) L=[eval(L[0]),eval(L[1])];
         else L=eval(L);          else L=eval(L);
         if(isvar(X)){          if(isvar(X)){
                 Y=makenewv([X,V]);                  Y=makenewv([X,F]);
                 if(type(F)==5) return [Y,[Y,os_md.periodicf,[F],L,X]];                  Z=makenewv([X,Y,F]);
                 Z=makenewv([X,Y,V]);                  return [Z,[Z,os_md.periodicf,[mysubst(F,[x,Y])],(type(L)==4)?[L]:L,[[Y,X]]]];
                 return [Z,[Z,os_md.periodicf,[mysubst(F,[x,Y])],[L],[[Y,X]]]];  
         }          }
         X=eval(X);          if(type(X)==4){
         if(type(F)==5)                  V=X[0];
                 return myfeval(F[floor(X/L)%length(F)],X-floor(X/L)*L);                  X=X[1];
           }else V=x;
           if(type(F)==5){
                   X=eval(X);
                   return myfeval(F[floor(X/L)%length(F)],[V,X-floor(X/L)*L]);
           }
         if(type(L)==4){          if(type(L)==4){
                 if(type(X)==4){  
                         V=X[0];  
                         X=X[1];  
                 }else V=x;  
                 X-=floor((X-L[0])/(L[1]-L[0]))*(L[1]-L[0]);                  X-=floor((X-L[0])/(L[1]-L[0]))*(L[1]-L[0]);
                 return myfeval(F,[V,X]);                  return myfeval(F,[V,X]);
         }          }
Line 13410  def draw_bezier(ID,IDX,B)
Line 15399  def draw_bezier(ID,IDX,B)
         return 0;          return 0;
 }  }
   
   
   /*
   def redbezier(L)
   {
           V=newvect(4);ST=0;
           for(R=[],I=0,T=L;T=[];T=cdr(T){
                   if(type(car(T))<4){
                           F=0;
                           if(I==3)
                           if(car(T)==0){
                           }else if(car(T)==1){
                           }else if(car(T)==-1){
                                   if(I<3) V[I++]=ST;
                           }
                   }else if(I==3){
                           if(R==[] || car(R)!=1){
                                   R=cons(V[0],R);
                                   if(ST==0) ST=V[0];
                           }
                           for(J=1;J<3;J++) R=cons(V[J],R);
                           while((T=cdr(T))!=[]){
                                   R=cons(car(T),R);
                                   if(type(car(R))<4)
                           }
                   }else{
                           if(ST==0) ST=car(T);
                           V[I++]= car(T);
                   }
           }
   }
   */
   
 def lbezier(L)  def lbezier(L)
 {  {
         if((In=getopt(inv))==1||In==2||In==3){          if((In=getopt(inv))==1||In==2||In==3){
Line 13798  def xycirc(P,R)
Line 15819  def xycirc(P,R)
     return S+"}};\n";      return S+"}};\n";
 }  }
   
   def xypoch(W,H,R1,R2)
   {
           if(H>R1||2*H>R2){
                   errno(0);
                   return;
           }
           if(type(Ar=getopt(ar))!=1) Ar=TikZ?0.25:2.5;
           T1=dasin(H/R1);S1=R1*dcos(T1);
           T2=dasin(H/R2);S2=R2*dcos(T2);
           T3=dasin(2*H/R2);S3=R2*dcos(T3);
           S=xyline([R1,0],[W-R1,0]);
           S+=xyang(R1,[W,0],-@pi,@pi-T1);
           S+=xyline([S2,H],[W-S1,H]);
           S+=xyang(R2,[0,0],T2,2*@pi-T3);
           S+=xylines([[S3,-2*H],[W-H-R2,-2*H],[W-H-R2,2*H],[W-S3,2*H]]);
           S+=xyang(R2,[W,0],-@pi+T2,@pi-T3);
           S+=xyline([W-T2,-H],[W-T2,-H]);
           S+=xyang(R1,[0,0],0,2*@pi-T1);
           S+=xyline([W-S2,-H],[S1,-H]);
           if(Ar>0){
                   S+=xyang(Ar,[W/2,0],[0,0],8);
                   S+=xyang(Ar,[W/2,-2*H],[0,-2*H],8);
                   S+=xyang(Ar,[W/2-Ar,-H],[W,-H],8);
                   S+=xyang(Ar,[W/2-Ar,H],[W,H],8);
                   S+=xyang(Ar,[W-S3,2*H],[W-H-R2,2*H],8);
           }
           S+=xyput([R1,0,"$\\bullet$"]);
           S+=xyput([0,0,"$\\times$"]);
           S+=xyput([W,0,"$\\times$"]);
           if(TikZ) S=str_subst(S,";\n\\draw","\n");
           return S;
   }
   
 def ptaffine(M,L)  def ptaffine(M,L)
 {  {
         if(type(L)!=4&&type(L)!=5){          if(type(L)!=4&&type(L)!=5){
Line 14084  def ptcopy(L,V)
Line 16138  def ptcopy(L,V)
         }          }
 }  }
   
   def regress(L)
   {
           E=deval(exp(0));
           for(S0=T0=0,S=L;S!=[];S=cdr(S)){
                   S0+=car(S)[0]*E;T0+=car(S)[1]*E;
           }
           K=length(L);S0/=K;T0/=K;
           for(SS=TT=0,S=L;S!=[];S=cdr(S)){
                   SS+=(car(S)[0]-S0)^2*E;TT+=(car(S)[1]-T0)^2*E;
                   ST+=(car(S)[0]-S0)*(car(S)[1]-T0)*E;
           }
           if(!SS||!TT) return [];
           A=ST/SS;
           L=[A,A*S0-T0,ST/dsqrt(SS*TT),S0,dsqrt(SS/K),T0,dsqrt(TT/K)];
           if(isint(N=getopt(sint))){
                   R=reverse(L);
                   for(L=[];R!=[];R=cdr(R)) L=cons(sint(car(R),N|str=0),L);
           }
           return L;
   }
   
 def     average(L)  def     average(L)
 {  {
         L=os_md.m2l(L|flat=1);          if(getopt(opt)=="co"){
         M0=M1=car(L);                  S0=average(L[0]);V0=car(S0);
         for(I=SS=0, LT=L; LT!=[]; LT=cdr(LT), I++){                  S1=average(L[1]);V1=car(S1);
                 S+=(V=car(LT));                  L0=os_md.m2l(L[0]|flat=1);
                 SS+=V^2;                  L1=os_md.m2l(L[1]|flat=1);
                 if(V<M0)                M0=V;                  for(S=0;L0!=[];L0=cdr(L0),L1=cdr(L1))
                 else if(V>M1)   M1=V;                          S+=(car(L0)-V0)*(car(L1)-V1);
                   S/=S0[1]*S1[1]*S0[2];
                   S=[S,S0,S1];
           }else{
                   L=os_md.m2l(L|flat=1);
                   M0=M1=car(L);
                   for(I=SS=0, LT=L; LT!=[]; LT=cdr(LT), I++){
                           S+=(V=car(LT));
                           SS+=V^2;
                           if(V<M0)                M0=V;
                           else if(V>M1)   M1=V;
                   }
                   SS=dsqrt(SS/I-S^2/I^2);
                   S=[deval(S/I),SS,I,M0,M1];
         }          }
         SS=dsqrt(SS/I-S^2/I^2);  
         S=[deval(S/I),SS,I,M0,M1];  
         if(isint(N=getopt(sint))) S=sint(S,N);          if(isint(N=getopt(sint))) S=sint(S,N);
         return S;          return S;
 }  }
Line 14410  def ltotex(L)
Line 16496  def ltotex(L)
                 else Hline=subst(Hline,z,S);                  else Hline=subst(Hline,z,S);
                 for(VV=[],VT=Hline;VT!=[];VT=cdr(VT)){                  for(VV=[],VT=Hline;VT!=[];VT=cdr(VT)){
                         if(type(T=car(VT))==4 && T[1]>0){                          if(type(T=car(VT))==4 && T[1]>0){
                                 for(I=T[0];I<=CS;I+=T[1]) VV=cons(I,VV);                                  for(I=T[0];I<=S;I+=T[1]) VV=cons(I,VV);
                         }else VV=cons(T,VV);                          }else VV=cons(T,VV);
                 }                  }
                 Hline=qsort(VV);                  Hline=qsort(VV);
Line 14465  def ltotex(L)
Line 16551  def ltotex(L)
                 }                  }
                 str_tb("\\end{tabular}\n",Out);                  str_tb("\\end{tabular}\n",Out);
         }else if(Op==11){       /* graph */          }else if(Op==11){       /* graph */
                 Width=8; Hight=3; WRet=1/2; HMerg=0.2;                  if(type(Strip=getopt(strip))!=1) Strip=0;
                   if(type(MX=getopt(max))!=1)  MX=0;
                   if(type(ML=getopt(mult))!=1) ML=0;
                   if((REL=getopt(relative))!=1) REL=0;
                   CL=getopt(color);
                   OL=delopt(getopt(),["color","strip","mult"]);
                   if(ML==1&&type(CL)==4){
                           LL=L[1];L=L[0];K=length(L);S=T="";
                           if(!MX){
                                   MX=vector(length(L[0]));
                                   for(LT=L;LT!=[];LT=cdr(LT)){
                                           for(I=0,LTT=car(LT);LTT!=[];I++,LTT=cdr(LTT)){
                                                   if(REL==1) MX[I]+=car(LTT);
                                                   else if(MX[I]<car(LTT)) MX[I]=car(LTT);
                                           }
                                   }
                                   MX=lmax(MX);
                                   OL=cons(["max",MX],OL);
                           }
                           if(REL==1) MX=newvect(length(L[0]));
                           for(I=0;I<K;I++){
                                   for(R=[],J=length(L[I]);--J>=0;){
                                           if(REL==1){
                                                   R=cons([MX[J],V=MX[J]+L[I][J]],R);
                                                   MX[J]=V;
                                           }else R=cons([(!I)?0:L[I-1][J],L[I][J]],R);
                                   }
                                   OP=cons(["color",CL[I]],OL);
                                   S+=ltotex([R,LL]|option_list=cons(["value",0],cons(["strip",(!I)?1:2],OP)));
                                   T+=ltotex([R,LL]|option_list=cons(["strip",3],OP));
                           }
                           return(!Strip)?xyproc(S+T):(S+T);
                   }else if(!TikZ) CL=0;
                   if(type(Line=getopt(line))!=1){
                           if(type(Line)==4){
                                   if(type(Line[0])==1 && (type(Line[1])==7 || type(Line[1])==1)){
                                           Opt=Line[1]; Line=Line[0];
                                   }else if(ML==1){
                                           OL=delopt(OL,"line");
                                           LL=L[1];L=L[0];K=length(L);S="";
                                           if(!MX){
                                                   MX=newvect(length(L[0]));
                                                   for(LT=L;LT!=[];LT=cdr(LT)){
                                                           for(I=0,LTT=car(LT);LTT!=[];I++,LTT=cdr(LTT)){
                                                                   if(REL==1) MX[I]+=car(LTT);
                                                                   else if(MX[I]<car(LTT)) MX[I]=car(LTT);
                                                           }
                                                   }
                                                   MX=lmax(MX);
                                                   OL=cons(["max",MX],OL);
                                           }
                                           for(I=0;I<K;I++)
                                                   S+=ltotex([L[I],LL]|option_list
                                                           =cons(["line",Line[I]],cons(["strip",(!I)?1:2],OL)));
                                           return(!Strip)?xyproc(S):S;
                                   }
                           }else Line=0;
                   }else Opt="@{-}";
                   Width=8; Hight=3; WRet=1/2; HMerg=(getopt(horiz)==1)?0.3:0.2;
                 if(!TikZ){                  if(!TikZ){
                         Width*=10; Hight*=10; HMerg*=10;                          Width*=10; Hight*=10; HMerg*=10;
                 }                  }
                   VMerg=HMerg;
                   if(type(Shift=getopt(shift))!=1)
                           Shift=0;
                 if(type(V=getopt(size))==4){                  if(type(V=getopt(size))==4){
                         Width=V[0];Hight=V[1];                          Width=V[0];Hight=V[1];
                         if(Hight<0) Hight=-Hight*lmax(L[0]);  
                         if(length(V)>2) WRet=V[2];                          if(length(V)>2) WRet=V[2];
                         if(length(V)>3) HMerg=V[3];                          if(length(V)>3) VMerg=VMerg=V[3];
                           if(length(V)>4) HMerg=V[4];
                 }                  }
                 Val=getopt(value);                  Val=getopt(value);
                 if(!isint(Val)) Val=-1;                  if(!isint(Val)) Val=-1;
                 if(type(Shift=getopt(shift))!=1)  
                         Shift=0;  
                 if(type(Line=getopt(line))!=1){                  if(type(Line=getopt(line))!=1){
                         if(type(Line)==4 && type(Line[0])==1 && (type(Line[1])==7 || type(Line[1])==1)){                          if(type(Line)==4 && type(Line[0])==1 && (type(Line[1])==7 || type(Line[1])==1)){
                                 Opt=Line[1]; Line=Line[0];                                  Opt=Line[1]; Line=Line[0];
Line 14492  def ltotex(L)
Line 16637  def ltotex(L)
                                 if((S=car(LT))<=0) return 0;                                  if((S=car(LT))<=0) return 0;
                                 Sum+=S;                                  Sum+=S;
                         }                          }
                         for(R=[],LT=L;LT!=[];LT=cdr(LT))                          for(R=[],LT=L;LT!=[];LT=cdr(LT)) R=cons(car(LT)/Sum,R);
                                 R=cons(car(LT)/Sum,R);  
                         R=reverse(R);                          R=reverse(R);
                         Opt0=Opt*2/3;                          Opt0=Opt*2/3;
                         Out=str_tb(xyproc(1),0);                          Out=str_tb((Strip>0)?0:xyproc(1),0);
                         str_tb(xylines(ptpolygon(6,Opt)|close=1,curve=1),Out);                          if(type(CL)!=4) str_tb(xylines(ptpolygon(6,Opt)|close=1,curve=1),Out);
                         for(S=0,RT=R,LT=LL;RT!=[];RT=cdr(RT)){                          for(S=0,RT=R,LT=LL;RT!=[];RT=cdr(RT)){
                                 str_tb(xyline([0,0],[Opt*dsin(S*6.2832),Opt*dcos(S*6.2832)]),Out);                                  SS=S+RT[0];
                                 T=S+RT[0]/2;                                  if(type(CL)==4){
                                 S+=RT[0];                                          str_tb(xyang(Opt,[0,0],(0.25-SS)*6.2832,(0.25-S)*6.2832|ar=1,opt=car(CL)),Out);
                                           if(length(CL)>0) CL=cdr(CL);
                                   }else str_tb(xyline([0,0],[Opt*dsin(S*6.2832),Opt*dcos(S*6.2832)]),Out);
                                   T=(S+SS)/2;
                                   S=SS;
                                 if(LT!=[]){                                  if(LT!=[]){
                                         str_tb(xyput([Opt0*dsin(T*6.2832),Opt0*dcos(T*6.2832),SS]),Out);                                          str_tb(xyput([Opt0*dsin(T*6.2832),Opt0*dcos(T*6.2832),car(LT)]),Out);
                                         LT=cdr(LT);                                          LT=cdr(LT);
                                 }                                  }
                         }                          }
                         str_tb(xyproc(0),Out);                          if(!Strip) str_tb(xyproc(0),Out);
                         return str_tb(0,Out);                          return str_tb(0,Out);
                 }                  }
                 if(type(MX=getopt(max))!=1)  
                         MX=0;  
                 if(MX==0){                  if(MX==0){
                         for(MX=0,LT=L; LT!=[]; LT=cdr(LT))                          for(MX=0,LT=L; LT!=[]; LT=cdr(LT))
                                 if(car(LT)>MX) MX=car(LT);                                  if(car(LT)>MX) MX=car(LT);
Line 14520  def ltotex(L)
Line 16666  def ltotex(L)
                 S=length(L);                  S=length(L);
                 WStep=Width/S;                  WStep=Width/S;
                 WWStep=WStep*WRet;                  WWStep=WStep*WRet;
                 HStep=Hight/MX;                  HStep=(Hight<0)?-Hight:Hight/MX;
                 if(LL!=[]&&length(LL)==S-1) WS2=WStep/2;                  if(LL!=[]&&length(LL)==S-1) WS2=WStep/2;
                 else WS2=0;                  else WS2=0;
                 Out=str_tb(xyproc(1),0);                  Out=str_tb((Strip>0)?0:xyproc(1),0);
                 str_tb(xyline([0,0],[Width-WStep+WWStep,0]),Out);                  Hori=getopt(horiz);
                 if(TikZ) CL=getopt(color);                  if(Strip<2){
                 else CL=0;                          if(Hori==1)  str_tb(xyline([0,0],[0,Width-WStep+WWStep]),Out);
                           else str_tb(xyline([0,0],[Width-WStep+WWStep,0]),Out);
                   }
                 for(I=0,LT=L;LT!=[]; LT=cdr(LT),I++){                  for(I=0,LT=L;LT!=[]; LT=cdr(LT),I++){
                         XP=WStep*I; XPM=XP+WWStep/2; YP=(car(LT)-Shift)*HStep;                          XP=WStep*I; XPM=XP+WWStep/2;
                         if(Line!=0){                          if(type(LTT=car(LT))==4){
                                 if(I>0)                                  YP0=(car(LTT)-Shift)*HStep;YP=(LTT[1]-Shift)*HStep;
                                         str_tb(xyarrow([XPM-WStep,YPP],[XPM,YP]|opt=Opt),Out);                                  VL=LTT[1];
                                 if(Val!=0)                                  if(REL) VL-=LTT[0];
                                         str_tb(xyput([XPM,YP+HMerg,car(LT)]),Out);                          }else{
                                 if(Line==2)                                  YP0=0;YP=(LTT-Shift)*HStep;VL=LTT;
                                         str_tb(xyput([XPM,YP,"$\\bullet$"]),Out);                          }
                                 YPP=YP;                          if(Hori==1){
                         }else if(YP!=0 || Val==1){                                  if(Line!=0){
                                 if(CL) str_tb(xybox([[XP,0],[XP+WWStep,YP]]|color=CL),Out);                                          if(I>0)
                                 else str_tb(xybox([[XP,0],[XP+WWStep,YP]]),Out);                                                  str_tb(xyarrow([XPM,YP],[XPM-WStep,YPP]|opt=Opt),Out);
                                 if(Val!=0){                                          if(Val!=0)
                                          str_tb(xyput([XPM,(YP<0)?(YP-HMerg):(YP+HMerg),car(LT)]),Out);                                                  str_tb(xyput([YP+HMerg, XPM,car(LT)]),Out);
                                           if(Line==2)
                                                   str_tb(xyput([YP,XPM,"$\\bullet$"]),Out);
                                           YPP=YP;
                                   }else if(YP!=0 || Val==1){
                                           if(Strip!=3){
                                                   if(CL) str_tb(xybox([[YP,XP+WWStep], [YP0,XP]]|color=CL),Out);
                                                   else str_tb(xybox([[YP,XP+WWStep],[YP0,XP]]),Out);
                                           }
                                           if(Val!=0) str_tb(xyput([(YP<0||REL==1)?(YP-HMerg):(YP+HMerg),XPM,VL]),Out);
                                 }                                  }
                                   if(LL!=[]&&I<length(LL)&&Strip<2) str_tb(xyput([-VMerg,XPM+WS2,LL[I]]),Out);
                           }else{
                                   if(Line!=0){
                                           if(I>0)
                                                   str_tb(xyarrow([XPM-WStep,YPP],[XPM,YP]|opt=Opt),Out);
                                           if(Val!=0)
                                                   str_tb(xyput([XPM,YP+HMerg,car(LT)]),Out);
                                           if(Line==2)
                                                   str_tb(xyput([XPM,YP,"$\\bullet$"]),Out);
                                           YPP=YP;
                                   }else if(YP!=0 || Val==1){
                                           if(Strip!=3){
                                                   if(CL) str_tb(xybox([[XP,YP0],[XP+WWStep,YP]]|color=CL),Out);
                                                   else str_tb(xybox([[XP,YP0],[XP+WWStep,YP]]),Out);
                                           }
                                           if(Val!=0) str_tb(xyput([XPM,(YP<0||REL==1)?(YP-HMerg):(YP+HMerg),VL]),Out);
                                   }
                                   if(LL!=[]&&I<length(LL)&&Strip<2) str_tb(xyput([XPM+WS2,-VMerg,LL[I]]),Out);
                         }                          }
                         if(LL!=[]&&I<length(LL)) str_tb(xyput([XPM+WS2,-HMerg,LL[I]]),Out);  
                 }                  }
                 str_tb(xyproc(0),Out);                  if(!Strip)str_tb(xyproc(0),Out);
         }else if(Op==12){       /* coord */          }else if(Op==12){       /* coord */
                 Out=str_tb("(",0);                  Out=str_tb("(",0);
                 for(LT=L;;){                  for(LT=L;;){
Line 15655  def shiftop(M,S)
Line 17829  def shiftop(M,S)
         return [QQ,P,RS];          return [QQ,P,RS];
 }  }
   
   
   def shiftPfaff(A,B,G,X,M)
   {
           if(type(G)==4){
                   G0=G[1];G1=G[0];
           }
           if(type(G)==6){
                   G=map(red,G);
                   G0=llcm(G);G1=map(red,G0*G);
           }
           if(type(G)==3){
                   G=red(G);G0=dn(G);G1=nm(G);
           }
           if(type(M)==4){
                   M0=M[0];M1=M[1];
           }else{
                   M0=M;M1=0;
           }
           X=vweyl(X);
           D0=mydeg(G0,X[0]);D1=mydeg(G1,X[0]);
           if(M1>=0){
                   D=(D1-M1>D0)?D1-M1:D0;
                   G0=muldo(X[1]^D,G0,X);G1=muldo(X[1]^(D+M1),G1,X);
           }else{
                   D=(D0+M1>D1)?D0+M1:D1;
                   G0=muldo(X[1]^(D-M1),G0,X);G1==muldo(X[1]^D,G1,X);
           }
           G0=map(mc,G0,X,M0);G1=map(mc,G1,X,M0+M1);
           G0=appldo(G0,A,X|Pfaff=1);
           G1=sppldo(G1,B,X|Pfaff=1);
           return rmul(myinv(G0),G1);
   }
   
 def conf1sp(M)  def conf1sp(M)
 {  {
         if(type(M)==7) M=s2sp(M);          if(type(M)==7) M=s2sp(M);
Line 15744  def conf1sp(M)
Line 17951  def conf1sp(M)
         return P;          return P;
 }  }
   
   /* ((1)(1)) ((1))   111|11|21  [[ [2,[ [1,[1]],[1,[1]] ]], [1,[[1,[1]]]] ]]  */
   /* (11)(1),111  111|21,111 [[[2,[1,1]],[1,[1]]],[1,1,1]]  */
   def s2csp(S)
   {
           if(type(S)!=7){
                   U="";
                   if(type(N=getopt(n))>0){
                           for(D=0,S=reverse(S);S!=[];S=cdr(S),D++){
                                   if(D) U=","+U;
                                   T=str_subst(rtostr(car(S)),","," ");
                                   U=str_cut(T,1,str_len(T)-2)+U;
                           }
                           V=strtoascii(U);
                           for(R=[];V!=[];V=cdr(V)){
                                   if((CC=car(V))==91){    /* [ */
                                           if(length(V)>1 && V[1]==91) V=cdr(V);
                                           for(I=1;(CC=V[I])!=91&&CC!=93;I++);
                                           if(CC==91){
                                                   R=cons(40,R);   /* ( */
                                                   while(I--) V=cdr(V);
                                           }else{
                                                   V=cdr(V);
                                                   while(--I) R=cons(car(V),R);
                                           }
                                   }else if(CC==93){               /* ] */
                                           R=cons(41,R);
                                           if(length(V)>1 && V[1]==93) V=cdr(V);
                                   }else R=cons(CC,R);
                           }
                           return asciitostr(reverse(R));
                   }
                   for(;S!=[];S=cdr(S)){
                           if(U!="") U=U+",";
                           for(D=0,TU="",T=car(S);T!=[];D++){
                                   if(type(car(T))==4){
                                           R=lpair(T,0);
                                           T=R[0];R1=m2l(R[1]|flat=1);
                                   }else R1=[];
                                   if(D) TU="|"+TU;
                                   TU=s2sp([T])+TU;
                                   T=R1;
                           }
                                   U=U+TU;
                   }
                   return U;
           }
           S=strtoascii(S);
           if(type(N=getopt(n))>0){
                   S=ltov(S);
                   L=length(S);
                   R="";
                   for(I=J=N=0, V=[];J<L;J++){
                           if(S[J]==72) I=J;                       /* ( */
                           else if(S[J]>47&&S[J]<58) N=N*10+S[J]-48;
                           else{
                                   if(N>0){
                                           V=cons(N,V);
                                           N=0;
                                   }
                                   if(S[J]==41){   /* ) */
   
                                   }else if(S[J]==44){             /* , */
   
                                   }
                           }
                   }
           }
           for(P=TS=[],I=D=0; S!=[]; S=cdr(S)){
                   if((C=car(S))==44){                     /* , */
                           P=cons(D,P);D=0;
                   }else if(C==124){       /* | */
                           D++;C=44;
                   }
                   TS=cons(C,TS);
           }
           S=reverse(TS);
           P=reverse(cons(D,P));
           U=s2sp(asciitostr(S));
   
           for(R=[];P!=[];P=cdr(P),U=cdr(U)){
                   D=car(P);R0=car(U);
                   while(D--){
                           U=cdr(U);
                           for(U0=car(U),R2=[];U0!=[];U0=cdr(U0)){
                                   for(R1=[],N=car(U0);N>0;R0=cdr(R0)){
                                           R1=cons(car(R0),R1);
                                           if(type(car(R0))==4) N-=car(R0)[0];
                                           else N-=car(R0);
                                   }
                                   R2=cons([car(U0),reverse(R1)],R2);
                           }
                           R0=reverse(R2);
                   }
                   R=cons(R0,R);
           }
           return reverse(R);
   }
   
   
   def partspt(S,T)
   {
           if(length(S)>length(T)) return [];
           if(type(Op=getopt(opt))!=1) Op=0;
           else{
                   VS=ltov(S);
                   L=length(S)-1;
                   VT=ltov(qsort(T));
           }
           if(length(S)==length(T)){
                   if(S==T||qsort(S)==qsort(T)) R=S;
                   else return [];
           }else if(getopt(sort)==1){
                   S0=S1=[];
                   for(;S!=[]&&car(S)==car(T);S=cdr(S),T=cdr(T))
                           S0=cons(car(S),S0);
                   if(S!=[]&&car(S)<car(T)) return [];
                   S0=reverse(S0);
                   for(S=reverse(S),T=reverse(T);S!=[],car(S)==car(T);S=cdr(S),T=cdr(T))
                           S1=cons(car(S),S1);
                   if(car(S)!=[]&&car(S)<cat(T)) return [];
                   R=partspt(reverse(S),reverse(T));
                   if(S1!=[]){
                           for(R0=[];R!=[];R=cdr(R))
                                   R0=cons(append(car(R),S1),R0);
                           R=reverse(R0);
                   }
                   if(S0!=[]){
                           for(R0=[];R!=[];R=cdr(R))
                                   R0=cons(append(S0,car(R)),R0);
                           R=reverse(R0);
                   }
           }else{
             for(R=[];;){
                   for(I=J=P=0;I<L;I++){
                           P=VS[I];
                           X=100000;
                           while((P-=(Y=VT[J++]))>0){
                                   if(X<Y) break;
                                   X=Y;
                           }
                           if(X<Y||P<0) break;
                   }
                   if(!P&&X>=Y) R=cons(vtol(VT),R);
                   if(!vnext(VT)) break;
             }
           }
           if(Op){
                   for(W=[];R!=[];R=cdr(R)){
                           for(I=0,S=VS[0],K=U=[],TR=car(R);TR!=[];TR=cdr(TR)){
                                   K=cons(car(TR),K);
                                   if(!(S-=car(K))){
                                           U=cons([VS[I],reverse(K)],U);
                                           K=[];
                                           S=VS[++I];
                                           if(I==L){
                                                   U=cons([S,cdr(TR)],U);
                                                   break;
                                           }
                                   }
                           }
                           W=cons(reverse(U),W);
                   }
                   R=W;
                   if(iand(Op,1)){
                           for(R=[];W!=[];W=cdr(W))
                                   R=cons(reverse(qsort(car(W))),R);
                           R=lsort(R,[],1);
                   }
                   if(Op==3){
                           for(W=[];R!=[];R=cdr(R)){
                                   for(S=[],TR=car(R);TR!=[];TR=cdr(TR))
                                           S=append(S,car(TR)[1]);
                                   W=cons(S,W);
                           }
                           R=reverse(W);
                   }
           }
           return R;
   }
   
   #if 0
   def confspt(S,T)
   {
           R=[];
           LS=length(S);LT=length(T);
           if(LS<LT)  return R;
           if(LS==LT){
                   return(S==T)? return [[S,T]]:R;
           }
           R=[];
           for(ST=S,S0=T0=[],TT=T;ST!=[];ST=cdr(ST),TT=cdr(TT)){
                   if(car(ST)>car(TT)) return R;
                   if(car(ST)==car(TT){
                           S0=cons(car(ST));T0=cons(car(TT));
                           LS--;LT--;continue;
                   }
                   V=car(TT);D=LS-LT;
                   for(P=[ST],DD=D;DD>0;){
                           VD=V-car(car(ST));
                   }
           }
   }
   #endif
   
   def mcvm(N)
   {
     X=getopt(var);
     if((Z=getopt(z))!=1) Z=0;
     if(type(N)==4){
       if((K=length(N))==1&&isvar(X)) X=[X];
       if(type(X)!=4){
         for(X=[],I=0;I<K;I++) X=cons(asciitostr([97+I]),X);
         X=reverse(X);
       }
           if(getopt(e)==1){
             if(length(N)==4){
                   N=ltov(N);
                   if(N[1]<N[3]){
                           I=N[1];N[1]=N[3];N[3]=I;
                   }
                   if(N[2]<N[3]||N[2]>=N[1]+N[3]) return 0;
                   X=X[0];
                   for(R=[],I=1;I<N[3];I++) R=cons(makev([X[0],I]),R);
                   for(L=[],I=N[1];I<=N[2];I++) L=cons(makev([X[0],I]),L);
                   for(S=0,I=N[1];I<=N[2];I++){
                     V=makev([X[0],I]);
             S+=polbyroot(R,V)/polbyroot(lsort(L,V,1),V);
                     S=red(S);
                   }
                   return S;
         }
           }
           for(M=[],I=S=0;I<K;Z=0,I++){
                   M=cons(mcvm(N[I]|var=X[I],z=Z),M);
                   S+=N[I];
           }
           M=newbmat(K,K,reverse(M));
       NR=N;
           N=S;
     }else{
           if(type(X)==7) X=strtov(X);
           if(!isvar(X)) X=a;
       M=newmat(N,N);
       for(I=0;I<N;I++){
         V=makev([X,I+1]);
         for(J=0;J<=I;J++){
           R=polbyroot([1,J],V|var=X);
           if(Z==1) R*=V;
           M[I][J]=R;
         }
       }
     }
     if((Get=getopt(get))==1){
       for(R=[],I=0;I<N;I++){
         U=newmat(N,N);
         for(J=0;J<N;J++) U[J][J]=M[J][I];
         R=cons(rmul(rmul(myinv(M),U),M),R);
       }
       return reverse(R);
     }else if(Get==2||Get==3||Get==4){
           for(V=[],I=N;I>0;I--) V=cons(makev(["a0",I]),V);
       MI=myinv(M);
           V=ltov(V)*MI;
           for(R=[],I=0;I<N;I++){
         for(J=I+1;J<N;J++){
           K=newmat(N,N);
           K[I][I]=V[J];K[I][J]=-V[J];K[J][J]=V[I];K[J][I]=-V[I];
              R=cons(rmul(rmul(MI,K),M),R);
          }
           }
       R=reverse(R);
           if(Get==2||length(NR)!=2||Z==1) return R;
       for(V1=[],I=NR[0];I>0;I--) V1=cons(os_md.makev([X[0],I]),V1);
       for(V2=[],I=NR[1];I>0;I--) V2=cons(os_md.makev([X[1],I]),V2);
       R=subst(R,car(V1),0,car(V2),0);
       V1=subst(V1,car(V1),0);
       V2=subst(V2,car(V2),0);
       for(V=[],S=V1;S!=[];S=cdr(S)) for(T=V2;T!=[];T=cdr(T)) V=cons(car(T)-car(S),V);
       V=reverse(V);
       Mx=length(V);
       for(A0=[],I=J=NR[0]-1;J>=0;I+=--J) for(K=0;K<NR[1];K++,I++) A0=cons(R[I],A0);
       A0=reverse(A0);
       for(F0=[],T=1,I=Mx-1;I>=0;I--) F0=cons(1/(x-V[I]), F0);
       MV=confexp([F0,V]|sym=3);
       RR=newvect(Mx);
       for(K=0;K<Mx;K++) for(RR[K]=0,I=0;I<Mx;I++) RR[K]=map(red,RR[K]+MV[I][K]*A0[I]);
       RR0=mysubst(RR,[append(cdr(V1),cdr(V2)),vtol(newvect(Mx-2))]|lpair=1);
       RR0=vtol(RR0);
       return (Get==3)?[RR,RR0]:RR0;
     }
     return M;
   }
   
   def confexp(S)
   {
           if((Sym=getopt(sym))==1||Sym==2||Sym==3){
                   D=polbyroot(S[1],x);
                   for(R=[],T=S[0];T!=[];T=cdr(T)){
                           M=D*car(T);
                           if(type(M)>3) M=map(red,M);
                           else M=red(M);
                           R=cons(M,R);
                   }
                   R=reverse(R);
                   if(Sym==2) return R;
                   M=length(R);N=length(S[1]);
                   E=newmat(M,N);
                   for(I=0;I<M;I++){
                           for(J=0;J<N;J++) E[I][J]=mycoef(R[I],N-J-1,x);
                   }
                   if(Sym==3){
                           for(R=[],P=1,T=S[1];T!=[];T=cdr(T)) R=cons(P/=(x-car(T)),R);
                           R=confexp([reverse(R),S[1]]|sym=1);
                           return E*myinv(R);
                   }
                   return E;
           }
           if(type(S[0])==4){
                   for(E=[];S!=[];S=cdr(S)) E=cons(confexp(car(S),E));
                   return reverse(E);
           }
           V=x;E=[];
           for(P=0,Q=[],ST=S;ST!=[];ST=cdr(ST)){
                   Q=cons(car(ST)[0],Q);
                   P+=car(ST)[1]/(V-car(ST)[0]);
                   P=red(P);
           }
           P=red(P*polbyroot(Q,V));
           Q=cdr(reverse(Q));
           for(I=(length(W=Q));I>=0;I--){
                   C=mycoef(P,I,V);
                   P-=C*polbyroot(W,V);
                   W=cdr(W);
                   E=cons(red(C),E);
           }
           return reverse(E);
   }
   
 def pgen(L,VV)  def pgen(L,VV)
 {  {
         if(type(L[0])<4) L=[L];          if(type(L[0])<4) L=[L];
Line 15864  def newbmat(M,N,R)
Line 18409  def newbmat(M,N,R)
         S  = newvect(M);          S  = newvect(M);
         T  = newvect(N);          T  = newvect(N);
         IM = length(R);          IM = length(R);
           if(type(car(R))!=4 && M==N && M==IM){
                   for(RR=TR=[],I=0;I<M;I++){
                           for(TR=[R[I]],J=0;J<I;J++) TR=cons(0,TR);
                           RR=cons(TR,RR);
                   }
                   R=reverse(RR);
           }
         for(I = 0; I < IM; I++){          for(I = 0; I < IM; I++){
                 RI = R[I];                  RI = R[I];
                 JM = length(RI);                  JM = length(RI);
Line 16869  def integrate(P,X)
Line 19421  def integrate(P,X)
                         if(S!=RR) R=cons([[1,RR=S]],R);                          if(S!=RR) R=cons([[1,RR=S]],R);
                         for(V=FR=[];R!=[];R=cdr(R))                          for(V=FR=[];R!=[];R=cdr(R))
                                 if(car(R)!=FR) V=cons(FR=car(R),V);                                  if(car(R)!=FR) V=cons(FR=car(R),V);
                         Var=varargs(V|all=1)[1];                          Var=varargs(V|all=2);
                         for(S0=[x0,x1,x2,x3],S=[t,s,u,v,w];S0!=[]&&S!=[];){                          for(S0=[x0,x1,x2,x3],S=[t,s,u,v,w];S0!=[]&&S!=[];){
                                 if(findin(car(S0),Var)<0){                                  if(findin(car(S0),Var)<0){
                                         S0=cdr(S0); continue;                                          S0=cdr(S0); continue;
Line 18000  def linfrac01(X)
Line 20552  def linfrac01(X)
   
 def varargs(P)  def varargs(P)
 {  {
         if((All=getopt(all))!=1) All=0;          if((All=getopt(all))!=1&&All!=2) All=0;
         V=vars(P);          V=vars(P);
         for(Arg=FC=[];V!=[];V=cdr(V)){          for(Arg=FC=[];V!=[];V=cdr(V)){
                 if(vtype(CV=car(V))==0&&All==1){                  if(vtype(CV=car(V))==0&&All!=0){
                         Arg=lsort([CV],Arg,0);                          Arg=lsort([CV],Arg,0);
                 }                  }
                 if(vtype(CV)!=2) continue;                  if(vtype(CV)!=2) continue;
Line 18020  def varargs(P)
Line 20572  def varargs(P)
                         }                          }
                 }                  }
         }          }
         return [FC,Arg];          Arg=reverse(Arg);
           return (All==2)?Arg:[reverse(FC),Arg];
 }  }
   
 def pfargs(P,X)  def pfargs(P,X)
Line 18235  def ntable(F,II,D)
Line 20788  def ntable(F,II,D)
 {  {
         F=f2df(F|opt=-1);          F=f2df(F|opt=-1);
         Df=getopt(dif);          Df=getopt(dif);
         if(Df!=1) Df=0;  
         Str=getopt(str);          Str=getopt(str);
         L=[];T=II[1]-II[0];          if(Df!=1) Df=0;
           L=[];
           if(type(D)==4){
                   if(type(II[0])==4){
                           T1=II[0][1]-II[0][0];T2=II[1][1]-II[1][0];
                           for(L0=[],I=0;I<D[0];I++){
                                   for(R=[],J=0;J<D[1];J++)
                                           R=cons(myf2eval(F,II[0][0]+I*T1/D[0],II[1][0]+J*T2/D[1]),R);
                                   L=cons(reverse(R),L);L0=cons(II[0][0]+I*T1/D[0],L0);
                           }
                   }else{
                           for(T=II[1]-II[0],L0=[],I=0;I<D[0];I++){
                                   for(R=[],J=0;J<D[1];J++)
                                           R=cons(myfdeval(F,II[0]+I*T/D[0]+J*T/D[0]/D[1]),R);
                                   L=cons(reverse(R),L);L0=cons(II[0]+I*T/D[0],L0);
                           }
                   }
                   L=reverse(L);L0=reverse(L0);
                   if(type(Str)==4){
                           L0=mtransbys(os_md.sint,L0,[Str[0]]|str=1,zero=0);
                           L=mtransbys(os_md.sint,L,[Str[1]]|str=1,zero=0);
                           if(Df==1){
                                   for(DT=[],RT=L,I=0;RT!=[];){
                                           for(LT=[],TT=car(RT);TT!=[];TT=cdr(TT)){
                                                   VV=car(TT);
                                                   if((J=str_char(VV,0,"."))>=0){
                                                           if(J==0) VV=str_cut(VV,1,10000);
                                                           else VV=str_cut(VV,0,J-1)+str_cut(VV,J+1,10000);
                                                   }
                                                   V1=eval_str(VV);
                                                   if(I++) LT=cons(V1-V0,LT);
                                                   V0=V1;
                                           }
                                           DT=cons(LT,DT);
                                           if((RT=cdr(RT))==[]){
                                                   VE=rint(myfdeval(F,II[1])*10^Str[1]);
                                                   DT=cons([VE-V0],DT);
                                           }
                                   }
                                   for(I=0,D=[],TT=DT;TT!=[];TT=cdr(TT)){
                                           if(!I++) V=car(TT)[0];
                                           else{
                                                   T1=reverse(cons(V,car(TT)));
                                                   V=car(T1);
                                                   if(length(TT)>1) T1=cdr(T1);
                                                   D=cons(T1,D);
                                           }
                                   }
                                   for(DD=[],TT=D;TT!=[];TT=cdr(TT))
                                   DD=cons([os_md.lmin(car(TT)),os_md.lmax(car(TT))],DD);
                                   DD=reverse(DD);
                                   L=lsort(L,DD,"append");
                           }
                   }
                   L=lsort(L,L0,"cons");
                   if(type(Top=getopt(top))==4||getopt(TeX)==1){
                           if(type(Top)==4){
                                   K=length(L[0])-length(Top);
                                   if(K>0&&K<4){
                                           if(K>1){
                                                   Top=append(Top,["",""]);
                                                   K-=2;
                                           }
                                           if(K) Top=cons("",Top);
                                   }
                                   L=cons(Top,L);
                           }
                           if(type(H=getopt(hline))!=4) H=[0,1,z];
                           if(type(V=getopt(vline))!=4) V=[0,1,(DF)?z-2:z];
                           if(type(T=getopt(title))!=7) Out=ltotex(L|opt="tab",hline=H,vline=V);
                           else Out=ltotex(L|opt="tab",hline=H,vline=V,title=T);
                           if(Df) Out=str_subst(Out,"\\hline","\\cline{1-"+rtostr(length(L[0])-2)+"}");
                           return Out;
                   }
                   return L;
           }
         for(L=[],I=0;I<=D;I++){          for(L=[],I=0;I<=D;I++){
                 X=II[0]+I*T/D;                  X=II[0]+I*T/D;
                 L=cons([X,myfdeval(F,X)],L);                  L=cons([X,myfdeval(F,X)],L);
Line 18249  def ntable(F,II,D)
Line 20876  def ntable(F,II,D)
                 }                  }
                 L=reverse(LD);                  L=reverse(LD);
         }          }
         if(type(Str=getopt(str))==4){          if(type(Str)==4){
                 if(length(Str)==1) Str=[Str[0],Str[0]];                  if(length(Str)==1) Str=[Str[0],Str[0]];
                 if(Df==1 && length(Str)==2) Str=[Str[0],Str[1],Str[2]];                  if(Df==1 && length(Str)==2) Str=[Str[0],Str[1],Str[1]];
                 for(S=Str,Str=[];S!=[];S=cdr(S)){                  for(S=Str,Str=[];S!=[];S=cdr(S)){
                         if(type(car(S))!=4) Str=cons([car(S),3],Str);                          if(type(car(S))!=4) Str=cons([car(S),3],Str);
                         else Str=cons(car(S),Str);                          else Str=cons(car(S),Str);
Line 18339  def distpoint(L)
Line 20966  def distpoint(L)
   
 def keyin(S)  def keyin(S)
 {  {
         print(S,2);          mycat0(S,0);
         purge_stdin();          purge_stdin();
         S=get_line();          S=get_line();
         L=length(S=strtoascii(S));          L=length(S=strtoascii(S));
Line 18372  def init() {
Line 20999  def init() {
         }          }
         if(Id>=0){          if(Id>=0){
                 while((S=get_line(Id))!=0){                  while((S=get_line(Id))!=0){
                         if(type(P=str_str(S,LS))==4 && (P0=str_char(S,P[1]+5,"="))>0){                          if(type(P=str_str(S,LS))==4 && (P0=str_char(S,P[1]+4,"="))>0){
                                 if(P[0]<5){                                  if(P[0]<5){
                                         P0=str_chr(S,P0+1,"\"");                                          P0=str_chr(S,P0+1,"\"");
                                         if(P0>0){                                          if(P0>0){
Line 18396  def init() {
Line 21023  def init() {
                                         else if(P[0]==7)        TikZ=SV;                                          else if(P[0]==7)        TikZ=SV;
                                         else if(P[0]==8)        XYPrec=SV;                                          else if(P[0]==8)        XYPrec=SV;
                                         else if(P[0]==9)        XYcm=SV;                                          else if(P[0]==9)        XYcm=SV;
                                         else if(P[0]==10)       XYcm=Canvas;                                          else if(P[0]==10)       Canvas=SV;
                                 }                                  }
                         }                          }
                 }                  }

Legend:
Removed from v.1.8  
changed lines
  Added in v.1.58

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>