home *** CD-ROM | disk | FTP | other *** search
/ OS/2 Professional / OS2PRO194.ISO / os2 / progs / pari / pari_137 / src / mp.c < prev    next >
Text File  |  1992-05-20  |  37KB  |  1,985 lines
  1.  
  2. /*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~*/
  3. /*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~*/
  4. /*~                                    ~*/
  5. /*~               OPERATIONS DE BASE (NOYAU)            ~*/
  6. /*~                                    ~*/
  7. /*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~*/
  8. /*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~*/
  9.  
  10. # include "genpari.h"
  11.  
  12. GEN   mulsi(x,y)
  13.      long x;
  14.      GEN y;
  15. {
  16.   long s=signe(y),ly=lgef(y),i;
  17.   GEN z;
  18.   
  19.   if((!x)||(!s)) return gzero;
  20.   if(x<0) {s= -s;x= -x;}
  21.   z=cgeti(ly+1);hiremainder=0;
  22.   for(i=ly-1;i>=2;i--) z[i+1]=addmul(x,y[i]);
  23.   if(hiremainder) {z[2]=hiremainder;setlgef(z,ly+1);}
  24.   else {avma+=4;z[1]=z[0]-1;z++;setlgef(z,ly);}
  25.   setsigne(z,s);return z;
  26. }
  27.  
  28. GEN stoi(x)
  29.      long x;
  30. {
  31.   GEN y;
  32.   
  33.   if(!x) return gzero;
  34.   y=cgeti(3);
  35.   if(x>0) {y[1]=0x1000003;y[2]=x;}
  36.   else{y[1]=0xff000003;y[2]= -x;}
  37.   return y;
  38. }
  39.  
  40. GEN cgetg(x,y)
  41.      long x,y;
  42. {
  43.   unsigned long p1;
  44.   GEN z;
  45.   
  46.   p1=avma-(((unsigned short)x)<<2);
  47.   if(p1<bot) err(errpile);
  48.   avma=p1;z=(GEN)p1;z[0]=0x10000+x+(y<<24);
  49.   return z;
  50. }
  51.  
  52. GEN cgeti(x)
  53.      long x;
  54. {
  55.   unsigned long p1;
  56.   GEN z;
  57.   
  58.   p1=avma-(((unsigned short)x)<<2);
  59.   if(p1<bot) err(errpile);
  60.   avma=p1;z=(GEN)p1;z[0]=0x1010000+x;
  61.   return z;
  62. }
  63.  
  64. GEN cgetr(x)
  65.      long x;
  66. {
  67.   unsigned long p1;
  68.   GEN z;
  69.   
  70.   p1=avma-(((unsigned short)x)<<2);
  71.   if(p1<bot) err(errpile);
  72.   avma=p1;z=(GEN)p1;z[0]=0x2010000+x;
  73.   return z;
  74. }
  75.  
  76. GEN icopy(x)
  77.      GEN x;
  78. {
  79.   GEN y;
  80.   long lx=lgef(x),i;
  81.   
  82.   y=cgeti(lx);
  83.   for(i=1;i<lx;i++) y[i]=x[i];
  84.   return y;
  85. }
  86.  
  87. GEN rcopy(x)
  88.      GEN x;
  89. {
  90.   GEN y;
  91.   long lx=lg(x),i;
  92.   
  93.   y=cgetr(lx);
  94.   for(i=1;i<lx;i++) y[i]=x[i];
  95.   return y;
  96. }
  97.  
  98.  
  99. GEN negi(x)
  100.      GEN x;
  101. {
  102.   long s=signe(x);
  103.   GEN y;
  104.   
  105.   if(!s) return gzero;
  106.   y=icopy(x);setsigne(y,-s);
  107.   return y;
  108. }
  109.  
  110. GEN negr(x)
  111.      GEN x;
  112. {
  113.   GEN y;
  114.   
  115.   y=rcopy(x);setsigne(y,-signe(x));
  116.   return y;
  117. }
  118.  
  119.  
  120. GEN absi(x)
  121.      GEN x;
  122. {
  123.   GEN y;
  124.   long s=signe(x);
  125.   
  126.   if(!s) return gzero;
  127.   y=icopy(x);setsigne(y,1);return y;
  128. }
  129.  
  130. GEN absr(x)
  131.      GEN x;
  132. {
  133.   GEN y;
  134.   long s=signe(x);
  135.   
  136.   y=rcopy(x);
  137.   if(s) setsigne(y,1);
  138.   return y;
  139. }
  140.  
  141. int expi(x)
  142.      GEN x;
  143. {
  144.   long lx=x[1]&0xffff;
  145.   
  146.   return lx==2 ? -8388608 : ((lx-2)<<5)-bfffo(x[2])-1;
  147. }
  148.  
  149. long itos(x)
  150.      GEN x;
  151. {
  152.   long s=signe(x),p2;
  153.   unsigned long p1;
  154.   
  155.   if(!s) return 0;
  156.   if(lgef(x)>3) err(affer2);
  157.   p1=x[2];if(p1>=0x80000000) err(affer2);
  158.   p2=(s>0)?p1:(-((long)p1));return p2;
  159. }
  160.  
  161. void mpaff(x,y)
  162.      GEN x,y;
  163. {
  164.   long tx=typ(x),ty=typ(y);
  165.   if(tx==1)
  166.     {
  167.       if(ty==1) affii(x,y);else affir(x,y);
  168.     }
  169.   else
  170.     {
  171.       if(ty==1) affri(x,y);else affrr(x,y);
  172.     }
  173. }
  174.  
  175. void affsi(s,x)
  176.      long s;
  177.      GEN x;
  178. {
  179.   long lx;
  180.   
  181.   if(!s) {x[1]=2;return;}
  182.   lx=lg(x);if(lx<3) err(affer1);
  183.   if(s>0) {x[1]=0x1000003;x[2]=s;}
  184.   else {x[1]=0xff000003;x[2]= -s;}
  185. }
  186.  
  187. void affsr(s,x)
  188.      long s;
  189.      GEN x;
  190. {
  191.   long l,i,d;
  192.   
  193.   if(!s)
  194.     {
  195.       l=((x[0]&0xffff)-2)<<5;x[1]=0x800000-l;x[2]=0;
  196.     }
  197.   else
  198.     {
  199.       d=1;if(s<0) {d= -1;s= -s;}
  200.       l=bfffo(s);setexpo(x,31-l);setsigne(x,d);
  201.       x[2]=(s<<l);for(i=3;i<lg(x);i++) x[i]=0;
  202.     }
  203. }
  204.  
  205. void affii(x,y)
  206.      GEN x,y;
  207. {
  208.   long lx=lgef(x),i;
  209.   
  210.   if(x==y) return;
  211.   if(lg(y)<lx) err(affer3);
  212.   for(i=1;i<lx;i++) y[i]=x[i];
  213. }
  214.  
  215. void affir(x,y)
  216.      GEN x,y;
  217. {
  218.   long lx=lgef(x),ly=lg(y),s,i,l,k;
  219.   
  220.   s=signe(x);
  221.   if(!s)
  222.     {
  223.       y[1]=0x800000-((ly-2)<<5);y[2]=0;
  224.     }
  225.   else
  226.     {
  227.       setsigne(y,s);l=((x[1]&0xffff)-2)<<5;s=bfffo(x[2]);
  228.       if(s)
  229.     {
  230.       setexpo(y,l-s-1);
  231.       if(lx<=ly)
  232.         {
  233.           for(i=lx;i<ly;i++) y[i]=0;k=0;
  234.           for(i=lx-1;i>=2;i--) {y[i]=shiftl(x[i],s)+k;k=hiremainder;}
  235.         }
  236.       else
  237.         {
  238.           shiftl(x[ly],s);k=hiremainder;
  239.           for(i=ly-1;i>=2;i--) {y[i]=shiftl(x[i],s)+k;k=hiremainder;}
  240.         }
  241.     }
  242.       else
  243.     {
  244.       setexpo(y,l-1);
  245.       if(lx<=ly)
  246.         {
  247.           for(i=lx;i<ly;i++) y[i]=0;
  248.           for(i=lx-1;i>=2;i--) y[i]=x[i];
  249.         }
  250.       else for(i=ly-1;i>=2;i--) y[i]=x[i];
  251.     }
  252.     }
  253. }
  254.  
  255. void affrr(x,y)
  256.      GEN x,y;
  257. {
  258.   long lx=lg(x),ly=lg(y),i;
  259.   
  260.   if(x==y) return;
  261.   y[1]=x[1];
  262.   if(!(x[1]&0xff000000)) y[2]=0;
  263.   else
  264.     {
  265.       if(lx<=ly)
  266.     {
  267.       for(i=2;i<lx;i++) y[i]=x[i];
  268.       for(i=lx;i<ly;i++) y[i]=0;
  269.     }
  270.       else for(i=2;i<ly;i++) y[i]=x[i];
  271.     }
  272. }
  273.  
  274. GEN shifts(x,y)
  275.      long x,y;
  276. {
  277.   long t[3];
  278.   
  279.   if(!x) return gzero;
  280.   t[0]=0x1010003;
  281.   if(x>0) {t[1]=0x1000003;t[2]=x;}
  282.   else {t[1]=0xff000003;t[2]= -x;}
  283.   return shifti(t,y);
  284. }
  285.  
  286. GEN shifti(x,n)
  287.      GEN x;
  288.      long n;
  289. {
  290.   long lx=lgef(x),i,s=signe(x),d,m,p1,p2,k;
  291.   GEN y;
  292.   
  293.   if(!s) return gzero;
  294.   if(!n) return icopy(x);
  295.   if(n>0)
  296.     {
  297.       d=n>>5;m=n&31;
  298.       if(m)
  299.     {
  300.       p1=shiftl(x[2],m);p2=hiremainder;k=0;
  301.       if(p2)
  302.         {
  303.           y=cgeti(lx+d+1);for(i=lx+1;i<=lx+d;i++) y[i]=0;
  304.           for(i=lx;i>=4;i--) {y[i]=shiftl(x[i-1],m)+k;k=hiremainder;}
  305.           y[3]=p1+k;y[2]=p2;
  306.         }
  307.       else
  308.         {
  309.           y=cgeti(lx+d);for(i=lx;i<lx+d;i++) y[i]=0;
  310.           for(i=lx-1;i>=3;i--) {y[i]=shiftl(x[i],m)+k;k=hiremainder;}
  311.           y[2]=p1+k;
  312.         }
  313.     }
  314.       else
  315.     {
  316.       y=cgeti(lx+d);for(i=lx;i<lx+d;i++) y[i]=0;
  317.       for(i=lx-1;i>=2;i--) y[i]=x[i];
  318.     }
  319.     }
  320.   else
  321.     {
  322.       n= -n;d=n>>5;m=n&31;if(lx<d+3) return gzero;
  323.       if(!m)
  324.     {
  325.       y=cgeti(lx-d);for(i=2;i<lx-d;i++) y[i]=x[i];
  326.     }
  327.       else 
  328.     {
  329.       m=32-m;d++;
  330.       p1=shiftl(x[2],m);
  331.       if(hiremainder)
  332.         {
  333.           y=cgeti(lx-d+1);y[2]=hiremainder;
  334.           for(i=3;i<=lx-d;i++)
  335.         {
  336.           p2=shiftl(x[i],m);y[i]=p1+hiremainder;p1=p2;
  337.         }
  338.         }
  339.       else
  340.         {
  341.           if(lx==d+2) return gzero;
  342.           y=cgeti(lx-d);
  343.           for(i=3;i<=lx-d;i++) 
  344.         {
  345.           p2=shiftl(x[i],m);y[i-1]=p1+hiremainder;p1=p2;
  346.         }
  347.         }
  348.     }
  349.     }   
  350.   y[1]=y[0];setsigne(y,s);return y;
  351. }
  352.  
  353. GEN shiftr(x,n)
  354.      GEN x;
  355.      long n;
  356. {
  357.   long l;
  358.   GEN y;
  359.   
  360.   y=rcopy(x);l=expo(x)+n;
  361.   if(l>0x7fffff||l<-0x800000) err(shier2);
  362.   setexpo(y,l);return y;
  363. }
  364.  
  365. GEN mptrunc(x)
  366.      GEN x;
  367. {
  368.   long e,i,s,lx=lg(x),p1,p2,m;
  369.   unsigned long d;
  370.   GEN y;
  371.   
  372.   if(typ(x)==1) return icopy(x);
  373.   s=signe(x);if(!s) return gzero;
  374.   e=expo(x);if(e<0) return gzero;
  375.   d=e>>5;m=e&31;if(d>=lx-2) err(truer2);
  376.   y=cgeti(d+3);y[1]=y[0];setsigne(y,s);
  377.   if(m==31) for(i=2;i<=d+2;i++) y[i]=x[i];
  378.   else
  379.     {
  380.       m++;p1=0;
  381.       for(i=2;i<=d+2;i++)
  382.     {
  383.       p2=shiftl(x[i],m);y[i]=hiremainder+p1;p1=p2;
  384.     }
  385.     }
  386.   return y;
  387. }
  388.  
  389. GEN mpent(x)
  390.      GEN x;
  391. {
  392.   long e,i,lx=lg(x),m,f,p1,p2;
  393.   unsigned long d;
  394.   GEN y,z;
  395.   
  396.   if(typ(x)==1) return icopy(x);
  397.   if(signe(x)>=0) return mptrunc(x);
  398.   e=expo(x);if(e<0) {y=cgeti(3);y[2]=1;y[1]=0xff000003;return y;}
  399.   d=e>>5;m=e&31;if(d>=lx-2) err(truer2);
  400.   y=cgeti(d+3);y[1]=0xff000003+d;
  401.   if(m==31) 
  402.     {
  403.       for(i=2;i<=d+2;i++) y[i]=x[i];
  404.       while((i<lx)&&(!x[i])) i++;
  405.       f=(i<lx);
  406.     }    
  407.   else
  408.     {
  409.       m++;p1=0;
  410.       for(i=2;i<=d+2;i++)
  411.     {
  412.       p2=shiftl(x[i],m);y[i]=hiremainder+p1;p1=p2;
  413.     }
  414.       if(p1) f=1;
  415.       else
  416.     {
  417.       while((i<lx)&&(!x[i])) i++;
  418.       f=(i<lx);
  419.     }
  420.     }
  421.   if(f)
  422.     {
  423.       for(i=d+2;(i>=2)&&(y[i]==0xffffffff);i--) y[i]=0;
  424.       if(i>=2) y[i]++;
  425.       else
  426.     {
  427.       z=y;y=cgeti(1);*y=(*z)+1;y[1]=z[1]+1;
  428.     }
  429.     }
  430.   return y;
  431. }
  432.  
  433. int mpcmp(x,y)
  434.      GEN x,y;
  435. {
  436.   if(typ(x)==1) return (typ(y)==1) ? cmpii(x,y) : cmpir(x,y);
  437.   return (typ(y)==1) ? -cmpir(y,x) : cmprr(x,y);
  438. }
  439.  
  440. int cmpsi(x,y)
  441.      long x;
  442.      GEN y;
  443. {
  444.   ulong p;
  445.   
  446.   if(!x) return -signe(y);
  447.   if(x>0)
  448.     {
  449.       if(signe(y)<=0) return 1;
  450.       if(lgef(y)>3) return -1;
  451.       p=y[2];if(p==x) return 0;
  452.       return (p<(ulong)x) ? 1 : -1;
  453.     }
  454.   else
  455.     {
  456.       if(signe(y)>=0) return -1;
  457.       if(lgef(y)>3) return 1;
  458.       p=y[2];if(p== -x) return 0;
  459.       return (p<(ulong)(-x)) ? -1 : 1;
  460.     }
  461. }
  462.  
  463. int cmpsr(x,y)
  464.      long x;
  465.      GEN y;      
  466. {
  467.   int p;
  468.   long av;
  469.   GEN z;
  470.   
  471.   if(!x) return -signe(y);
  472.   av=avma;z=cgetr(3);affsr(x,z);
  473.   p=cmprr(z,y);avma=av;return p;
  474. }    
  475.  
  476. int cmpii(x,y)
  477.      GEN x,y;
  478. {
  479.   long sx=signe(x),sy=signe(y),lx,ly,i;
  480.   
  481.   if(sx<sy) return -1;
  482.   if(sx>sy) return 1;
  483.   if(!sx) return 0;
  484.   lx=lgef(x);ly=lgef(y);
  485.   if(lx>ly) return sx;
  486.   if(lx<ly) return -sx;
  487.   for(i=2;(i<lx)&&(x[i]==y[i]);i++);
  488.   if(i==lx) return 0;
  489.   return ((ulong)x[i]>(ulong)y[i]) ? sx : -sx;
  490. }
  491.  
  492. int cmpir(x,y)
  493.      GEN x,y;
  494. {
  495.   long av=avma;
  496.   int p;
  497.   GEN z;
  498.   
  499.   if(!signe(x)) return -signe(y);
  500.   z=cgetr(lg(y));affir(x,z);
  501.   p=cmprr(z,y);avma=av;return p;
  502. }
  503.  
  504. int cmprr(x,y)
  505.      GEN x,y;
  506. {
  507.   long sx=signe(x),sy=signe(y),ex,ey,lx,ly,lz,i;
  508.   
  509.   if(sx<sy) return -1;
  510.   if(sx>sy) return 1;
  511.   if(!sx) return 0;
  512.   ex=expo(x);ey=expo(y);
  513.   if(ex>ey) return sx;
  514.   if(ex<ey) return -sx;
  515.   lx=lg(x);ly=lg(y);lz=(lx<ly)?lx:ly;
  516.   for(i=2;(i<lz)&&(x[i]==y[i]);i++);
  517.   if(i<lz) return ((ulong)x[i]>(ulong)y[i]) ? sx : -sx;
  518.   if(lx>=ly)
  519.     {
  520.       while((i<lx)&(!x[i])) i++;
  521.       return (i==lx) ? 0 : sx;
  522.     }
  523.   else
  524.     {
  525.       while((i<ly)&(!y[i])) i++;
  526.       return (i==ly) ? 0 : -sx;
  527.     }
  528. }    
  529.  
  530. GEN mpadd(x,y)
  531.      GEN x,y;
  532. {
  533.   if(typ(x)==1) return (typ(y)==1) ? addii(x,y) : addir(x,y);
  534.   return (typ(y)==1) ? addir(y,x) : addrr(x,y);
  535. }
  536.  
  537. GEN addss(x,y)
  538.      long x,y;
  539. {
  540.   long t[3];
  541.   
  542.   if(!x) return stoi(y);
  543.   t[0]=0x1010003;
  544.   if(x>0) {t[1]=0x1000003;t[2]=x;} else {t[1]=0xff000003;t[2]= -x;}
  545.   return addsi(y,t);
  546. }
  547.  
  548.  
  549. GEN addsi(x,y)
  550.      long x;
  551.      GEN y;
  552. {
  553.   long sx,sy,ly,p,i;
  554.   GEN z;
  555.   
  556.   if(!x) return icopy(y);
  557.   sy=signe(y);if(!sy) return stoi(x);
  558.   if(x<0) {sx= -1;x= -x;} else sx=1;
  559.   ly=lgef(y);
  560.   if(sx==sy)
  561.     {
  562.       p=addll(x,y[ly-1]);
  563.       if(overflow)
  564.     {
  565.       z=cgeti(ly+1);z[ly]=p;
  566.       for(i=ly-1;(i>2)&&(y[i-1]==0xffffffff);i--) z[i]=0;
  567.       if(i>2)
  568.         {
  569.           z[i]=y[i-1]+1;i--;while(i>=3) {z[i]=y[i-1];i--;}
  570.           z[2]=z[1]=z[0]-1;z++;avma+=4;
  571.         }
  572.       else {z[2]=1;z[1]=z[0];}
  573.     }
  574.       else
  575.     {
  576.       z=cgeti(ly);z[ly-1]=p;for(i=1;i<ly-1;i++) z[i]=y[i];
  577.     }
  578.       setsigne(z,sx);
  579.     }
  580.   else
  581.     {
  582.       if(ly==3)
  583.     {
  584.       if((ulong)y[2]>(ulong)x)
  585.         {
  586.           z=cgeti(3);z[1]=(sy<<24)+3;z[2]=y[2]-x;return z;
  587.         }
  588.       if(y[2]==x) return gzero;
  589.       z=cgeti(3);z[1]=((-sy)<<24)+3;z[2]=x-y[2];return z;
  590.     }
  591.       p=subll(y[ly-1],x);
  592.       if(overflow)
  593.     {
  594.       z=cgeti(ly);z[ly-1]=p;
  595.       for(i=ly-2;!(y[i]);i--) z[i]=0xffffffff;
  596.       z[i]=y[i]-1;
  597.       if((i>2)||z[i]) {i--;for(;i>=1;i--) z[i]=y[i];}
  598.       else
  599.         {
  600.           z[2]=z[1]=z[0]-1;z++;avma+=4;setsigne(z,sy);
  601.         }
  602.     }
  603.       else
  604.     {
  605.       z=cgeti(ly);z[ly-1]=p;for(i=1;i<ly-1;i++) z[i]=y[i];
  606.     }    
  607.     }
  608.   return z;
  609. }
  610.  
  611. GEN addii(x,y)
  612.      GEN x,y;
  613. {
  614.   long sx=signe(x),sy=signe(y),sz,lx=lgef(x),ly=lgef(y),i,j,p;
  615.   GEN z;
  616.   
  617.   if(!sx) return icopy(y);
  618.   if(!sy) return icopy(x);
  619.   if(lx<ly) {z=x;x=y;y=z;sz=sx;sx=sy;sy=sz;sz=lx;lx=ly;ly=sz;}
  620.   if(sx==sy)
  621.     {
  622.       z=cgeti(lx+1);overflow=0;
  623.       for(i=ly-1,j=lx-1;i>=2;i--,j--) z[j+1]=addllx(x[j],y[i]);
  624.       if(overflow)
  625.     {
  626.       for(;(j>=2)&&(x[j]==0xffffffff);j--) z[j+1]=0;
  627.       if(j>=2)
  628.         {
  629.           z[j+1]=x[j]+1;j--;
  630.           for(;j>=2;j--) z[j+1]=x[j];
  631.           z[1]=z[0]-1;z[2]=x[1];z++;avma+=4;
  632.         }
  633.       else {z[2]=1;z[1]=x[1]+1;}
  634.     }
  635.       else
  636.     {
  637.       for(;j>=2;j--) z[j+1]=x[j];
  638.       z[1]=z[0]-1;z[2]=x[1];z++;avma+=4;
  639.     }
  640.     }
  641.   else
  642.     {
  643.       if(lx==ly)
  644.     {
  645.       setsigne(y,1);setsigne(x,1);p=cmpii(x,y);
  646.       setsigne(y,sy);setsigne(x,sx);if(!p) return gzero;
  647.       if(p<0) {z=x;x=y;y=z;sz=sx;sx=sy;sy=sz;}
  648.     }
  649.       z=cgeti(lx);overflow=0;
  650.       for(i=ly-1,j=lx-1;i>=2;i--,j--) z[j]=subllx(x[j],y[i]);
  651.       if(overflow)
  652.     {
  653.       for(;!(x[j]);j--) z[j]=0xffffffff;
  654.       z[j]=x[j]-1;j--;
  655.       for(;j>=2;j--) z[j]=x[j];
  656.     }
  657.       else
  658.     {
  659.       for(;j>=2;j--) z[j]=x[j];
  660.     }
  661.       if(z[2]) z[1]=x[1];
  662.       else
  663.     {
  664.       for(j=3;(j<lx)&&(!z[j]);j++);
  665.       i=j-2;z[i+1]=z[i]=z[0]-i;z+=i;avma+=(i<<2);setsigne(z,sx);
  666.     }
  667.     }
  668.   return z;
  669. }      
  670.  
  671. GEN addsr(x,y)
  672.      long x;
  673.      GEN y;
  674. {
  675.   long p[3];
  676.   
  677.   if(!x) return rcopy(y);
  678.   p[0]=0x1010003;
  679.   if(x>0) {p[1]=0x1000003;p[2]=x;} else {p[1]=0xff000003;p[2]= -x;}
  680.   return addir(p,y);
  681. }
  682.  
  683. GEN addir(x,y)
  684.      GEN x,y;
  685. {
  686.   long l,e,ly,av,i,l1;
  687.   GEN z;
  688.   
  689.   if(!signe(x)) return rcopy(y);
  690.   if(!signe(y))
  691.     {
  692.       l=lgef(x)-(expo(y)>>5);if((l<3)||(l>32767)) err(adder3);
  693.       z=cgetr(l);affir(x,z);return z;
  694.     }
  695.   else
  696.     {
  697.       e=expo(y)-expi(x);ly=lg(y);
  698.       if(e>0)
  699.     {
  700.       l=ly-(e>>5);if(l<=2) return rcopy(y);
  701.     }
  702.       else
  703.     { 
  704.       l=ly+((-e)>>5)+1;if(l>32767) err(adder3);
  705.     }
  706.       av=avma;z=cgetr(l);affir(x,z);l1=av-avma;l=l1>>2;
  707.       z=addrr(z,y);
  708.       for(i=lg(z)-1;i>=0;i--) z[i+l]=z[i];z+=l;avma+=l1;
  709.     }
  710.   return z;
  711. }
  712.  
  713. GEN addrr(x,y)
  714.      GEN x,y;
  715. {
  716.   long sx=signe(x),sy=signe(y),lx=lg(x),ly=lg(y),lz,ex=expo(x),ey=expo(y),sz;
  717.   long av0=avma,e,l,i,d,m,flag,lp1,lp2,av,k,j,cex,f2;
  718.   GEN z,p1,p2;
  719.   
  720.   if(!sy)
  721.     {
  722.       if(!sx) {e=(ex>=ey)?ex:ey;z=cgetr(3);z[2]=0;z[1]=e+0x800000;return z;}
  723.       e=ex-ey;
  724.       if(e<=0) {z=cgetr(3);z[2]=0;z[1]=ey+0x800000;return z;}
  725.       l=(e>>5)+3;if(l>lx) l=lx;z=cgetr(l);
  726.       for(i=1;i<l;i++) z[i]=x[i];return z;
  727.     }
  728.   e=ey-ex;
  729.   if(!sx)
  730.     {
  731.       if(e<=0) {z=cgetr(3);z[2]=0;z[1]=ex+0x800000;return z;}
  732.       l=(e>>5)+3;if(l>ly) l=ly;z=cgetr(l);
  733.       for(i=1;i<l;i++) z[i]=y[i];return z;
  734.     }
  735.   if(e)
  736.     {
  737.       if(e<0) {z=x;x=y;y=z;lz=lx;lx=ly;ly=lz;ey=ex;e= -e;sz=sx;sx=sy;sy=sz;}
  738.       d=e>>5;m=e&31;
  739.       if(d>=ly-2) return rcopy(y);
  740.       l=d+lx;
  741.       if(l>=ly)
  742.     {
  743.       flag=1;p1=cgetr(ly);lp1=ly;lp2=ly-d;
  744.     }
  745.       else
  746.     {
  747.       flag=0;p1=cgetr(l+1);lp2=lx+1;lp1=l+1;
  748.     }
  749.       av=avma;
  750.       if(m)
  751.     {
  752.       p2=cgetr(lp2);m=32-m;
  753.       if(flag) {shiftl(x[lp2-1],m);k=hiremainder;}
  754.       else k=0;
  755.       for(i=lp2-1;i>=3;i--) 
  756.         {
  757.           p2[i]=shiftl(x[i-1],m)+k;k=hiremainder;
  758.         }
  759.       p2[2]=k;
  760.     }
  761.       else p2=x;
  762.     }
  763.   else
  764.     {
  765.       l=(lx>ly)?ly:lx;p1=cgetr(l);av=avma;lp2=lp1=l;flag=2;p2=x;m=0;
  766.     }
  767.   if(sx==sy)
  768.     {
  769.       overflow=0;
  770.       if(m+flag) for(i=lp1-1,j=lp2-1;j>=2;i--,j--) p1[i]=addllx(p2[j],y[i]);
  771.       else 
  772.     {
  773.       p1[lp1-1]=y[lp1-1];
  774.       for(i=lp1-2,j=lp2-2;j>=2;i--,j--) p1[i]=addllx(p2[j],y[i]);
  775.     }
  776.       if(overflow)
  777.     {
  778.       for(;(i>=2)&&(y[i]==0xffffffff);i--) p1[i]=0;
  779.       if(i>=2) {cex=0;p1[i]=y[i]+1;while(i>=3) {i--;p1[i]=y[i];}}
  780.       else 
  781.         {
  782.           cex=1;k=0x80000000;if(ey==0x7fffff) err(adder4);
  783.           for(i=2;i<lp1;i++) {p1[i]=shiftlr(p1[i],1)+k;k=hiremainder;}
  784.         }
  785.     }
  786.       else {cex=0;for(;i>=2;i--) p1[i]=y[i];}
  787.       p1[1]=(sx<<24)+ey+cex+0x800000;
  788.       avma=av;return p1;
  789.     }
  790.   else 
  791.     {
  792.       if(!e) 
  793.     {
  794.       for(i=2;(i<l)&&(p2[i]==y[i]);i++);
  795.       if(i==l)
  796.         {
  797.           e=ex-((l-2)<<5)+0x800000;if(e<0) err(adder5);
  798.           if(e>=0x1000000) err(adder4);
  799.           avma=av0;z=cgetr(3);z[2]=0;z[1]=e;return z;
  800.         }
  801.       else
  802.         {
  803.           f2=(((ulong)y[i])>((ulong)p2[i]))?1:0;
  804.         }
  805.     }
  806.       else f2=1;
  807.       if(f2)
  808.     {
  809.       overflow=0;
  810.       if(m+flag) for(i=lp1-1,j=lp2-1;j>=2;i--,j--) p1[i]=subllx(y[i],p2[j]);
  811.       else 
  812.         {
  813.           p1[lp1-1]=y[lp1-1];
  814.           for(i=lp1-2,j=lp2-2;j>=2;i--,j--) p1[i]=subllx(y[i],p2[j]);
  815.         }
  816.       if(overflow)
  817.         {
  818.           for(;(i>=2)&&(!y[i]);i--) p1[i]=0xffffffff;
  819.           p1[i]=y[i]-1;while(i>=3) {i--;p1[i]=y[i];}
  820.         }
  821.       else for(;i>=2;i--) p1[i]=y[i];
  822.     }
  823.       else
  824.     {
  825.       overflow=0;
  826.       if(m+flag) for(i=lp1-1;i>=2;i--) p1[i]=subllx(p2[i],y[i]);
  827.       else 
  828.         {
  829.           p1[lp1-1]=subllx(0,y[lp1-1]);
  830.           for(i=lp1-2;i>=2;i--) p1[i]=subllx(p2[i],y[i]);
  831.         }
  832.     }
  833.       for(i=2;!p1[i];i++);j=i-2;avma=av+(j<<2);p1[j]=p1[0]-j;p1+=j;
  834.       m=bfffo(p1[2]);e=ey-(j<<5)-m+0x800000;
  835.       if(e<0) err(adder5);
  836.       p1[1]=f2 ? (sy<<24)+e : (sx<<24)+e;
  837.       if(m)
  838.     {
  839.       k=0;for(i=lp1-1-j;i>=2;i--) {p1[i]=shiftl(p1[i],m)+k;k=hiremainder;}
  840.     }
  841.       return p1;
  842.     }
  843. }
  844.  
  845. void addssz(x,y,z)
  846.      long x,y;
  847.      GEN z;
  848. {
  849.   long av=avma;
  850.   GEN p1;
  851.   
  852.   if(typ(z)==1) gops2ssz(addss,x,y,z);
  853.   else
  854.     {
  855.       p1=cgetr(lg(z));affsr(x,p1);p1=addrs(p1,y);
  856.       affrr(p1,z);avma=av;
  857.     }
  858. }
  859.  
  860. GEN mpsub(x,y)
  861.      GEN x,y;
  862. {
  863.   if(typ(x)==1) return (typ(y)==1) ? subii(x,y) : subir(x,y);
  864.   return (typ(y)==1) ? subri(x,y) : subrr(x,y);
  865. }
  866.  
  867. GEN subii(x,y)
  868.      GEN x,y;
  869. {
  870.   long s=signe(y);
  871.   GEN z;
  872.   
  873.   if(x==y) return gzero;
  874.   setsigne(y,-s);z=addii(x,y);setsigne(y,s);
  875.   return z;
  876. }
  877.  
  878. GEN subrr(x,y)
  879.      GEN x,y;
  880. {
  881.   long s=signe(y);
  882.   GEN z;
  883.   
  884.   if(x==y)
  885.     {
  886.       z=cgetr(3);z[2]=0;z[1]=0x800000-(lg(x)<<5);return z;
  887.     }
  888.   setsigne(y,-s);z=addrr(x,y);setsigne(y,s);return z;
  889. }
  890.  
  891. GEN subsi(x,y)
  892.      long x;
  893.      GEN y;
  894. {
  895.   long s=signe(y);
  896.   GEN z;
  897.   
  898.   setsigne(y,-s);z=addsi(x,y);setsigne(y,s);return z;
  899. }
  900.  
  901. GEN subsr(x,y)
  902.      long x;
  903.      GEN y;
  904. {
  905.   long s=signe(y);
  906.   GEN z;
  907.   
  908.   setsigne(y,-s);z=addsr(x,y);setsigne(y,s);return z;
  909. }
  910.  
  911. GEN subss(x,y)
  912.      long x,y;
  913. {
  914.   return addss(-y,x);
  915. }
  916.  
  917.  
  918. GEN subir(x,y)
  919.      GEN x,y;
  920. {
  921.   long s=signe(y);
  922.   GEN z;
  923.   
  924.   setsigne(y,-s);z=addir(x,y);setsigne(y,s);return z;
  925. }
  926.  
  927. GEN subri(x,y)
  928.      GEN x,y;
  929. {
  930.   long s=signe(y);
  931.   GEN z;
  932.   
  933.   setsigne(y,-s);z=addir(y,x);setsigne(y,s);return z;
  934. }
  935.  
  936. void subssz(x,y,z)
  937.      long x,y;
  938.      GEN z;
  939. {
  940.   long av=avma;
  941.   GEN p1;
  942.   
  943.   if(typ(z)==1) gops2ssz(addss,x,-y,z);
  944.   else
  945.     {
  946.       p1=cgetr(lg(z));affsr(x,p1);p1=addrs(p1,-y);
  947.       affrr(p1,z);avma=av;
  948.     }
  949. }
  950.  
  951. GEN mpmul(x,y)
  952.      GEN x,y;
  953. {
  954.   if(typ(x)==1) return (typ(y)==1) ? mulii(x,y) : mulir(x,y);
  955.   return (typ(y)==1) ? mulir(y,x) : mulrr(x,y);
  956. }
  957.  
  958. GEN mulss(x,y)
  959.      long x,y;
  960. {
  961.   long s,p1;
  962.   GEN z;
  963.   
  964.   if((!x)||(!y)) return gzero;
  965.   s=1;if(x<0) {s= -1;x= -x;} if(y<0) {s= -s;y= -y;}
  966.   p1=mulll(x,y);
  967.   if(hiremainder) {z=cgeti(4);z[2]=hiremainder;z[3]=p1;}
  968.   else {z=cgeti(3);z[2]=p1;}
  969.   z[1]=z[0];setsigne(z,s);return z;
  970. }
  971.  
  972. GEN mulsr(x,y)
  973.      long x;
  974.      GEN y;
  975. {
  976.   long lx,i,k,s,p1,p2,e;
  977.   GEN z;
  978.   
  979.   if(!x) return gzero;
  980.   s=signe(y);if(x<0) {s= -s;x= -x;}
  981.   if(!s)
  982.     {
  983.       p1=bfffo(x);e=y[1]+31-p1;if(e>=0x1000000) err(muler2);
  984.       z=cgetr(3);z[1]=e;z[2]=0;
  985.     }
  986.   else
  987.     {
  988.       if(x==1) {z=rcopy(y);setsigne(z,s);return z;}
  989.       lx=lg(y);z=cgetr(lx);setsigne(z,s);
  990.       p2=mulll(x,y[lx-1]);
  991.       for(i=lx-2;i>=2;i--) z[i+1]=addmul(x,y[i]);
  992.       z[2]=hiremainder;p1=bfffo(hiremainder);
  993.       if(p1)
  994.     {
  995.       shiftl(p2,p1);k=hiremainder;
  996.       for(i=lx-1;i>=2;i--)
  997.         {
  998.           z[i]=shiftl(z[i],p1)+k;k=hiremainder;
  999.         }
  1000.     }
  1001.       e=32-p1+expo(y);if(e>=0x800000) err(muler2);
  1002.       setexpo(z,e);
  1003.     }  
  1004.   return z;
  1005. }
  1006.  
  1007. GEN mulii(x,y)
  1008.      GEN x,y;
  1009. {
  1010.   long i,j,lx=lgef(x),ly=lgef(y),sx,sy,lz,p1,p2;
  1011.   GEN z;
  1012.   
  1013.   sx=signe(x);if(!sx) return gzero;
  1014.   sy=signe(y);if(!sy) return gzero;
  1015.   if(sy<0) sx= -sx;
  1016.   if(lx>ly) {z=x;x=y;y=z;lz=lx;lx=ly;ly=lz;}
  1017.   lz=lx+ly-2;if(lz>=0x10000) err(muler1);
  1018.   z=cgeti(lz);z[1]=z[0];setsigne(z,sx);
  1019.   p1=x[lx-1];hiremainder=0;
  1020.   for(i=ly-1;i>=2;i--) z[lx+i-2]=addmul(p1,y[i]);
  1021.   z[lx-1]=hiremainder;
  1022.   for(j=lx-2;j>=2;j--)
  1023.     {
  1024.       p1=x[j];hiremainder=0;
  1025.       for(i=ly-1;i>=2;i--)
  1026.     {
  1027.       p2=addmul(p1,y[i]);
  1028.       z[i+j-1]=addll(p2,z[i+j-1]);hiremainder+=overflow;
  1029.     }
  1030.       z[j]=hiremainder;
  1031.     }
  1032.   if(!(z[2]))
  1033.     {
  1034.       z[2]=z[1]-1;z[1]=z[0]-1;z++;avma+=4;
  1035.     }
  1036.   return z;
  1037. }
  1038.  
  1039. GEN mulrr(x,y)
  1040.      GEN x,y;
  1041. {
  1042.   long i,j,lx=lg(x),ly=lg(y),sx=signe(x),sy=signe(y),ex=expo(x),ey=expo(y);
  1043.   long e,flag,garde,p1,p2,lz;
  1044.   GEN z;
  1045.   
  1046.   e=ex+ey+0x800000;if(e>=0xffffff) err(muler4);
  1047.   if(e<0) err(muler5);
  1048.   if((!sx)||(!sy)) {z=cgetr(3);z[2]=0;z[1]=e;return z;}
  1049.   if(sy<0) sx= -sx;
  1050.   if(lx>ly) {lz=ly;z=x;x=y;y=z;flag=1;} else {lz=lx;flag=(lx!=ly);}
  1051.   z=cgetr(lz);z[1]=(sx<<24)+e;
  1052.   if(flag) mulll(x[2],y[lz]);else hiremainder=0;
  1053.   if(lz==3)
  1054.     {
  1055.       garde=flag ? addmul(x[2],y[2]) : mulll(x[2],y[2]);
  1056.       if((long)hiremainder<0) {z[2]=hiremainder;z[1]++;}
  1057.       else {z[2]=(garde<0)?(hiremainder<<1)+1:(hiremainder<<1);}
  1058.       return z;
  1059.     }
  1060.   p1=x[lz-1];garde=hiremainder;
  1061.   if(p1)
  1062.     {
  1063.       mulll(p1,y[3]);p2=addmul(p1,y[2]);
  1064.       garde=addll(p2,garde);z[lz-1]=overflow+hiremainder;
  1065.     }
  1066.   else z[lz-1]=0;
  1067.   for(j=lz-2;j>=3;j--)
  1068.     {
  1069.       p1=x[j];
  1070.       if(p1)
  1071.     {
  1072.       mulll(p1,y[lz+2-j]);
  1073.       p2=addmul(p1,y[lz+1-j]);
  1074.       garde=addll(p2,garde);hiremainder+=overflow;
  1075.       for(i=lz-j;i>=2;i--)
  1076.         {
  1077.           p2=addmul(p1,y[i]);
  1078.           z[i+j-1]=addll(p2,z[i+j-1]);hiremainder+=overflow;
  1079.         }
  1080.       z[j]=hiremainder;
  1081.     }
  1082.       else z[j]=0;
  1083.     }
  1084.   p1=x[2];p2=mulll(p1,y[lz-1]);
  1085.   garde=addll(p2,garde);hiremainder+=overflow;
  1086.   for(i=lz-2;i>=2;i--)
  1087.     {
  1088.       p2=addmul(p1,y[i]);
  1089.       z[i+1]=addll(p2,z[i+1]);hiremainder+=overflow;
  1090.     }
  1091.   z[2]=hiremainder;
  1092.   if((long)hiremainder>0)
  1093.     {
  1094.       overflow=(garde<0)?1:0;
  1095.       for(i=lz-1;i>=2;i--) {p1=z[i];z[i]=addllx(p1,p1);}
  1096.     }
  1097.   else z[1]++;
  1098.   return z;
  1099. }
  1100.  
  1101. GEN mulir(x,y)
  1102.      GEN x,y;
  1103. {
  1104.   long sx=signe(x),sy,av,lz,ey,e,garde,p1,p2,i,j;
  1105.   GEN z,temp;
  1106.   
  1107.   if(!sx) return gzero;
  1108.   sy=signe(y);ey=expo(y);
  1109.   if(!sy)
  1110.     {
  1111.       z=cgetr(3);z[2]=0;e=expi(x)+ey+0x800000;if(e>=0x1000000) err(muler6);
  1112.       z[1]=e;return z;
  1113.     }
  1114.   lz=lg(y);if(sy<0) sx= -sx;
  1115.   z=cgetr(lz);setsigne(z,sx);av=avma;
  1116.   temp=cgetr(lz+1);affir(x,temp);x=y;y=temp;
  1117.   e=expo(y)+ey+0x800000;if(e>=0xffffff) err(muler4);
  1118.   if(e<0) err(muler5);
  1119.   z[1]=(sx<<24)+e;
  1120.   mulll(x[2],y[lz]);
  1121.   if(lz==3)
  1122.     {
  1123.       garde=addmul(x[2],y[2]);
  1124.       if((long)hiremainder<0) {z[2]=hiremainder;z[1]++;}
  1125.       else {z[2]=(garde<0)?(hiremainder<<1)+1:(hiremainder<<1);}
  1126.       avma=av;return z;
  1127.     }
  1128.   garde=hiremainder;
  1129.   p1=x[lz-1];mulll(p1,y[3]);p2=addmul(p1,y[2]);
  1130.   garde=addll(p2,garde);z[lz-1]=overflow+hiremainder;
  1131.   for(j=lz-2;j>=3;j--)
  1132.     {
  1133.       p1=x[j];mulll(p1,y[lz+2-j]);
  1134.       p2=addmul(p1,y[lz+1-j]);
  1135.       garde=addll(p2,garde);hiremainder+=overflow;
  1136.       for(i=lz-j;i>=2;i--)
  1137.     {
  1138.       p2=addmul(p1,y[i]);
  1139.       z[i+j-1]=addll(p2,z[i+j-1]);hiremainder+=overflow;
  1140.     }
  1141.       z[j]=hiremainder;
  1142.     }
  1143.   p1=x[2];p2=mulll(p1,y[lz-1]);
  1144.   garde=addll(p2,garde);hiremainder+=overflow;
  1145.   for(i=lz-2;i>=2;i--)
  1146.     {
  1147.       p2=addmul(p1,y[i]);
  1148.       z[i+1]=addll(p2,z[i+1]);hiremainder+=overflow;
  1149.     }
  1150.   z[2]=hiremainder;
  1151.   if((long)hiremainder>0)
  1152.     {
  1153.       overflow=(garde<0)?1:0;
  1154.       for(i=lz-1;i>=2;i--) {p1=z[i];z[i]=addllx(p1,p1);}
  1155.     }
  1156.   else z[1]++;
  1157.   avma=av;return z;
  1158. }
  1159.  
  1160. void mulssz(x,y,z)
  1161.      long x,y;
  1162.      GEN z;
  1163. {
  1164.   long av=avma;
  1165.   GEN p1;
  1166.   
  1167.   if(typ(z)==1) gops2ssz(mulss,x,y,z);
  1168.   else
  1169.     {
  1170.       p1=cgetr(lg(z));affsr(x,p1);p1=mulsr(y,p1);
  1171.       mpaff(p1,z);avma=av;
  1172.     }
  1173. }
  1174.  
  1175. GEN convi(x)
  1176.      GEN x;
  1177. {
  1178.   long lx,av=avma,lz;
  1179.   GEN z,p1,p2;  
  1180.   
  1181.   if(!signe(x))
  1182.     {
  1183.       z=cgeti(3);z[1]= -1;z[2]=0;avma=av;return z+3;
  1184.     }
  1185.   p1=absi(x);lx=lgef(p1);lz=((lx-2)*15)/14+3;z=cgeti(lz);z[1]= -1;
  1186.   for(p2=z+2;signe(p1);p2++) *p2=divisii(p1,1000000000,p1);
  1187.   avma=av;return p2;
  1188. }
  1189.  
  1190. GEN confrac(x)
  1191.      GEN x;
  1192. {
  1193.   long lx=lg(x),ex= -expo(x)-1,ex1,av=avma,ly,ey;
  1194.   long lr,nbdec,k,i,j;
  1195.   GEN y,res;
  1196.   
  1197.   ey=((lx-2)<<5)+ex;ly=(ey+63)>>5;y=cgeti(ly);ex1=ex>>5; /* 95 dans mp.s faux? */
  1198.   for(i=0;i<ex1;i++) y[i]=0;
  1199.   ex&=31;
  1200.   if(!ex) for(j=2;j<lx;j++) y[i++]=x[j];
  1201.   else
  1202.     {
  1203.       k=0;
  1204.       for(j=2;j<lx;j++) {y[i++]=shiftlr(x[j],ex)+k;k=hiremainder;}
  1205.       y[ly-2]=k;
  1206.     }
  1207.   y[ly-1]=0;
  1208.   nbdec=ey*0.30103+1;lr=(nbdec+17)/9;res=cgeti(lr);
  1209.   *res=nbdec;
  1210.   for(j=1;j<lr;j++)
  1211.     {
  1212.       hiremainder=0;
  1213.       for(i=ly-1;i>=0;i--) y[i]=addmul(y[i],1000000000);
  1214.       res[j]=hiremainder;
  1215.     }
  1216.   avma=av;return res;
  1217. }
  1218.  
  1219. void mulsii(x,y,z)
  1220.      long x;
  1221.      GEN y,z;
  1222. {
  1223.   long av=avma;
  1224.   GEN p1;
  1225.   
  1226.   p1=mulsi(x,y);affii(p1,z);avma=av;
  1227. }
  1228.  
  1229. void addsii(x,y,z)
  1230.      long x;
  1231.      GEN y,z;
  1232. {
  1233.   long av=avma;
  1234.   GEN p1;
  1235.   
  1236.   p1=addsi(x,y);affii(p1,z);avma=av;
  1237. }
  1238.  
  1239. long divisii(x,y,z)
  1240.      long y;
  1241.      GEN x,z;
  1242. {
  1243.   long av=avma,k;
  1244.   GEN p1;
  1245.   
  1246.   p1=divis(x,y);affii(p1,z);avma=av;
  1247.   k=hiremainder;return k;
  1248. }
  1249.  
  1250.  
  1251. /*
  1252. static int tv[37]={0, 0, 25, 1, 22, 26, 31, 2, 15, 23, 29, 27, 10, 0, 12, 3, 6, 16, 0, 24, 21, 30, 14, 28, 9, 11, 5, 0, 20, 13, 8, 4, 19, 7, 18, 17, 0};
  1253.  
  1254. #define vals(x) (x?tv[(x^(x-1))%37]:-1)
  1255. */
  1256.  
  1257.  
  1258. long vals(x)
  1259.      long x;
  1260. {
  1261.   unsigned short int y,z;
  1262.   int s;
  1263.  
  1264.   if(!x) return -1;
  1265.   y=x;if(!y) {s=16;y=((ulong)x)>>16;} else s=0;
  1266.   z=y&255;if(!z) {s+=8;z=y>>8;}
  1267.   y=z&15;if(!y) {s+=4;y=z>>4;}
  1268.   z=y&3;if(!z) {s+=2;z=y>>2;}
  1269.   return (z&1) ? s : s+1;
  1270. }
  1271.  
  1272. long vali(x)
  1273.      GEN x;
  1274. {
  1275.   long i,lx=lgef(x);
  1276.   
  1277.   if(!signe(x)) return -1;
  1278.   for(i=lx-1;(i>=2)&&(!x[i]);i--);
  1279.   return ((lx-1-i)<<5)+vals(x[i]);
  1280. }
  1281.  
  1282. GEN mpdiv(x,y)
  1283.      GEN x,y;
  1284. {
  1285.   if(typ(x)==1) return (typ(y)==1) ? divii(x,y) : divir(x,y);
  1286.   return (typ(y)==1) ? divri(x,y) : divrr(x,y);
  1287. }
  1288.  
  1289. GEN divss(x,y)
  1290.      long x,y;
  1291. {
  1292.   long p1;
  1293.   
  1294.   if(!y) err(diver1);
  1295.   hiremainder=0;p1=divll((ulong)abs(x),(ulong)abs(y));
  1296.   if(y<0) {hiremainder= -((long)hiremainder);p1= -p1;}
  1297.   if(x<0) p1= -p1;
  1298.   return stoi(p1);
  1299. }
  1300.  
  1301. GEN modss(x,y)
  1302.      long x,y;
  1303. {
  1304.   long y1;
  1305.   
  1306.   if(!y) err(moder1);
  1307.   hiremainder=0;divll(abs(x),y1=abs(y));
  1308.   if(!hiremainder) return gzero;
  1309.   return (((long)hiremainder)<0) ? stoi(y1-hiremainder) : stoi(hiremainder);
  1310. }
  1311.  
  1312. GEN resss(x,y)
  1313.      long x,y;
  1314. {
  1315.   if(!y) err(reser1);
  1316.   hiremainder=0;divll(abs(x),abs(y));
  1317.   return (y<0) ? stoi(-((long)hiremainder)) : stoi(hiremainder);
  1318. }
  1319.  
  1320. GEN dvmdss(x,y,z)
  1321.      long x,y;
  1322.      GEN *z;
  1323. {
  1324.   GEN p1;
  1325.  
  1326.   p1=divss(x,y);*z=stoi(hiremainder);
  1327.   return p1;
  1328. }
  1329.  
  1330. void dvmdssz(x,y,z,t)
  1331.      long x,y;
  1332.      GEN z,t;
  1333. {
  1334.   long av=avma;
  1335.   GEN p1;
  1336.  
  1337.   p1=divss(x,y);affsi(hiremainder,t);
  1338.   mpaff(p1,z);avma=av;
  1339. }
  1340.  
  1341. GEN dvmdsi(x,y,z)
  1342.      long x;
  1343.      GEN y,*z;
  1344. {
  1345.   GEN p1;
  1346.   p1=divsi(x,y);*z=stoi(hiremainder);
  1347.   return p1;
  1348. }
  1349.  
  1350. void dvmdsiz(x,y,z,t)
  1351.      long x;
  1352.      GEN t,y,z;
  1353. {
  1354.   long av=avma;
  1355.   GEN p1;
  1356.   
  1357.   p1=divsi(x,y);affsi(hiremainder,t);
  1358.   mpaff(p1,z);avma=av;
  1359. }
  1360.  
  1361. GEN dvmdis(x,y,z)
  1362.      long y;
  1363.      GEN x,*z;
  1364. {
  1365.   GEN p1;
  1366.   p1=divis(x,y);*z=stoi(hiremainder);
  1367.   return p1;
  1368. }
  1369.  
  1370. void dvmdisz(x,y,z,t)
  1371.      long y;
  1372.      GEN x,t,z;
  1373. {
  1374.   long av=avma;
  1375.   GEN p1;
  1376.   
  1377.   p1=divis(x,y);affsi(hiremainder,t);
  1378.   mpaff(p1,z);avma=av;
  1379. }
  1380.  
  1381. void dvmdiiz(x,y,z,t)
  1382.      GEN x,y,z,t;
  1383. {
  1384.   long av=avma;
  1385.   GEN p1,p2;
  1386.  
  1387.   p1=dvmdii(x,y,&p2);mpaff(p1,z);mpaff(p2,t);
  1388.   avma=av;
  1389. }
  1390.  
  1391. GEN divsi(x,y)
  1392.      long x;
  1393.      GEN y;
  1394. {
  1395.   long s=signe(y),ly=lgef(y),p1;
  1396.  
  1397.   if(!s) err(diver2);
  1398.   if((!x)||(ly>3)||(y[2]<0)) {hiremainder=x;return gzero;}
  1399.   hiremainder=0;p1=divll(abs(x),y[2]);
  1400.   if(signe(y)<0) {hiremainder= -((long)hiremainder);p1= -p1;}
  1401.   if(x<0) p1= -p1;
  1402.   return stoi(p1);
  1403. }
  1404.  
  1405. GEN ressi(x,y)
  1406.      long x;
  1407.      GEN y;
  1408. {
  1409.   divsi(x,y);return stoi(hiremainder);
  1410. }
  1411.  
  1412. GEN modsi(x,y)
  1413.      long x;
  1414.      GEN y;
  1415. {
  1416.   long s;
  1417.   GEN p1;
  1418.   
  1419.   divsi(x,y);
  1420.   if(!hiremainder) return gzero;
  1421.   if(x>0) return stoi(hiremainder);
  1422.   else
  1423.     {
  1424.       s=signe(y);setsigne(y,1);p1=addsi(hiremainder,y);
  1425.       setsigne(y,s);return p1;
  1426.     }
  1427. }
  1428.  
  1429. GEN divis(y,x)
  1430.      long x;
  1431.      GEN y;
  1432. {
  1433.   long s=signe(y),ly=lgef(y),i,d;
  1434.   GEN z;
  1435.   
  1436.   if(!x) err(diver4);
  1437.   if(!s) {hiremainder=0;return gzero;}
  1438.   if(x<0) {s= -s;x= -x;} 
  1439.   if((ulong)x>(ulong)y[2])
  1440.     {
  1441.       if(ly==3) {hiremainder=itos(y);return gzero;}
  1442.       else {z=cgeti(ly-1);d=1;hiremainder=y[2];}
  1443.     }
  1444.   else {z=cgeti(ly);d=0;hiremainder=0;}
  1445.   for(i=d+2;i<ly;i++) z[i-d]=divll(y[i],x);
  1446.   z[1]=z[0];setsigne(z,s);if(s<0) hiremainder= -((long)hiremainder);
  1447.   return z;
  1448. }
  1449.  
  1450. GEN modis(x,y)
  1451.      long y;
  1452.      GEN x;
  1453. {
  1454.   divis(x,y);if(!hiremainder) return gzero;
  1455.   return (signe(x)>0) ? stoi(hiremainder) : stoi(abs(y)+hiremainder);
  1456. }
  1457.  
  1458. GEN resis(x,y)
  1459.      long y;
  1460.      GEN x;
  1461. {
  1462.   divis(x,y);return stoi(hiremainder);
  1463. }
  1464.      
  1465. void divisz(x,y,z)
  1466.      long y;
  1467.      GEN x,z;
  1468. {
  1469.   long av=avma;
  1470.   GEN p1;
  1471.   
  1472.   if(typ(z)==1) gops2gsz(divis,x,y,z);
  1473.   else
  1474.     {
  1475.       p1=cgetr(lg(z));affir(x,p1);p1=divrs(p1,y);
  1476.       affrr(p1,z);avma=av;
  1477.     }
  1478. }
  1479.  
  1480. void divsiz(x,y,z)
  1481.      long x;
  1482.      GEN y,z;
  1483. {
  1484.   long av=avma,lz;
  1485.   GEN p1,p2;
  1486.   
  1487.   if(typ(z)==1) gops2sgz(divsi,x,y,z);
  1488.   else
  1489.     {
  1490.       lz=lg(z);p1=cgetr(lz);p2=cgetr(lz);affsr(x,p1);affir(y,p2);
  1491.       p1=divrr(p1,p2);affrr(p1,z);avma=av;
  1492.     }
  1493. }
  1494.  
  1495. void divssz(x,y,z)
  1496.      long x,y;
  1497.      GEN z;
  1498. {
  1499.   long av=avma;
  1500.   GEN p1;
  1501.   
  1502.   if(typ(z)==1) gops2ssz(divss,x,y,z);
  1503.   else
  1504.     {
  1505.       p1=cgetr(lg(z));affsr(x,p1);p1=divrs(p1,y);
  1506.       affrr(p1,z);avma=av;
  1507.     }
  1508. }
  1509.  
  1510. GEN divir(x,y)
  1511.      GEN x,y;
  1512. {
  1513.   GEN xr,z;
  1514.   long av,ly;
  1515.   
  1516.   if(!signe(y)) err(diver5);
  1517.   if(!signe(x)) return gzero;
  1518.   ly=lg(y);z=cgetr(ly);av=avma;affir(x,xr=cgetr(ly+1));
  1519.   xr=divrr(xr,y);affrr(xr,z);avma=av;return z;
  1520. }
  1521.  
  1522. GEN divri(x,y)
  1523.      GEN x,y;
  1524.  
  1525. {
  1526.   GEN yr,z;
  1527.   long av,lx,ex,s=signe(y);
  1528.  
  1529.   if(!s) err(diver8);
  1530.   if(!signe(x))
  1531.   {
  1532.     ex=expo(x)-expi(y)+0x800000;
  1533.     if(ex<0) err(diver12);
  1534.     z=cgetr(3);z[1]=ex;z[2]=0;return z;
  1535.   }
  1536.   if((lg(y)==3)&&(y[2]>0)) return (s>0) ? divrs(x,y[2]) : divrs(x,-y[2]);
  1537.   lx=lg(x);z=cgetr(lx);av=avma;affir(y,yr=cgetr(lx+1));
  1538.   yr=divrr(x,yr);affrr(yr,z);avma=av;return z;
  1539. }
  1540.  
  1541. void diviiz(x,y,z)
  1542.      GEN x,y,z;
  1543. {
  1544.   long av=avma,lz;
  1545.   GEN p1,p2;
  1546.   
  1547.   if(typ(z)==1) {p1=divii(x,y);affii(p1,z);avma=av;}
  1548.   else
  1549.     {
  1550.       lz=lg(z);p1=cgetr(lz);p2=cgetr(lz);affir(x,p1);affir(y,p2);
  1551.       p1=divrr(p1,p2);affrr(p1,z);avma=av;
  1552.     }
  1553. }
  1554.      
  1555. void divrrz(x,y,z)
  1556.      GEN x,y,z;
  1557. {
  1558.   long av=avma;
  1559.   GEN p1;
  1560.  
  1561.   p1=divrr(x,y);mpaff(p1,z);avma=av;
  1562. }
  1563.  
  1564. void mpdivz(x,y,z)
  1565.      GEN x,y,z;
  1566. {
  1567.   long av=avma,lz;
  1568.   GEN p1,p2;
  1569.  
  1570.   if(typ(z)==1)
  1571.     {
  1572.       if(typ(x)==2||typ(y)==2) err(divzer1);
  1573.       p1=divii(x,y);affii(p1,z);avma=av;
  1574.     }
  1575.   else
  1576.     {
  1577.       if(typ(x)==1)
  1578.     {
  1579.       if(typ(y)==2) {p1=divir(x,y);mpaff(p1,z);avma=av;}
  1580.       else
  1581.         {
  1582.           lz=lg(z);p1=cgetr(lz);p2=cgetr(lz);affir(x,p1);affir(y,p2);
  1583.           p1=divrr(p1,p2);affrr(p1,z);avma=av;
  1584.         }
  1585.     }
  1586.       else
  1587.     {
  1588.       if(typ(y)==2) {p1=divrr(x,y);affrr(p1,z);avma=av;}
  1589.       else {p1=divri(x,y);affrr(p1,z);avma=av;}
  1590.     }
  1591.     }
  1592. }
  1593.  
  1594. GEN divsr(x,y)
  1595.      long x;
  1596.      GEN y;
  1597. {
  1598.   long av,ly;
  1599.   GEN p1,z;
  1600.  
  1601.   if(!signe(y)) err(diver3);
  1602.   if(!x) return gzero;
  1603.   ly=lg(y);z=cgetr(ly);av=avma;p1=cgetr(ly+1);affsr(x,p1);p1=divrr(p1,y);
  1604.   affrr(p1,z);avma=av;return z;
  1605. }
  1606.  
  1607. GEN modii(x,y)
  1608.      GEN x,y;
  1609. {
  1610.   long av=avma,tetpil;
  1611.   GEN p1;
  1612.  
  1613.   p1=dvmdii(x,y,-1);
  1614.   if(signe(p1)>=0) return p1;
  1615.   tetpil=avma;p1=(signe(y)>0) ? addii(p1,y) : subii(p1,y);
  1616.   return gerepile(av,tetpil,p1);
  1617. }
  1618.  
  1619. void modiiz(x,y,z)
  1620.      GEN x,y,z;
  1621. {
  1622.   long av=avma;
  1623.   GEN p1;
  1624.  
  1625.   p1=modii(x,y);affii(p1,z);avma=av;
  1626. }
  1627.  
  1628. void resiiz(x,y,z)
  1629.      GEN x,y,z;
  1630. {
  1631.   long av=avma;
  1632.   GEN p1;
  1633.  
  1634.   p1=resii(x,y);affii(p1,z);avma=av;
  1635. }
  1636.  
  1637. GEN divrs(x,y)
  1638.      long y;
  1639.      GEN x;
  1640. {
  1641.   long i,k,lx,ex,garde,sh,s=signe(x);
  1642.   GEN z;
  1643.  
  1644.   if(!y) err(diver6);
  1645.   if(!s)
  1646.     {
  1647.       z=cgetr(3);z[2]=0;z[1]=x[1]-31+bfffo(y);
  1648.       if(z[1]<0) err(diver7);return z;
  1649.     }
  1650.   if(y<0) {s= -s;y= -y;}
  1651.   if(y==1) {z=rcopy(x);setsigne(z,s);return z;}
  1652.   z=cgetr(lx=lg(x));setsigne(z,s);hiremainder=0;
  1653.   for(i=2;i<lx;i++) z[i]=divll(x[i],y);
  1654.   garde=divll(0,y);sh=bfffo(z[2]);ex=expo(x)-sh;if((-ex)>0x800000) err(diver7);
  1655.   setexpo(z,ex);shiftl(garde,sh);k=hiremainder;
  1656.   for(i=lx-1;i>=2;i--) {z[i]=shiftl(z[i],sh)+k;k=hiremainder;}
  1657.   return z;
  1658. }
  1659.  
  1660. int mpdivis(x,y,z)
  1661.      GEN x,y,z;
  1662. {
  1663.   long av=avma;
  1664.   GEN p1,p2;
  1665.  
  1666.   p1=dvmdii(x,y,&p2);
  1667.   if(signe(p2)) {avma=av;return 0;}
  1668.   affii(p1,z);avma=av;return 1;
  1669. }
  1670.  
  1671. int divise(x,y)
  1672.      GEN x,y;
  1673. {
  1674.   long av=avma;
  1675.   GEN p1;
  1676.  
  1677.   p1=dvmdii(x,y,-1);avma=av;
  1678.   return signe(p1) ? 0 : 1;
  1679. }
  1680.  
  1681.  
  1682. GEN dvmdii(x,y,z)
  1683.      GEN x,y,*z;
  1684. {
  1685.   long av,av2,lx,ly,lz,i,j,dec,sh,k,k1,sx=signe(x),sy=signe(y);
  1686.   long saux,k3,k4,av1,flk4;
  1687.   ulong si,qp;
  1688.   GEN p1,p2,p3,p4;
  1689.   
  1690.   if(!sy) err(dvmer1);
  1691.   if(!sx)
  1692.     {
  1693.       if(((long)z==0xffffffff)||((long)z==0)) return gzero;
  1694.       *z=gzero;return gzero;
  1695.     }
  1696.   lx=lgef(x);ly=lgef(y);lz=lx-ly;
  1697.   if(lz<0)
  1698.     {
  1699.       if((long)z==0xffffffff) return icopy(x);
  1700.       if(z==0) return gzero;
  1701.       *z=icopy(x);return gzero;
  1702.     }
  1703.   av=avma;if(sx<0) sy= -sy;
  1704.   if(ly==3)
  1705.     {
  1706.       si=y[2];
  1707.       if(si>(ulong)x[2]) {p1=cgeti(lx-1);hiremainder=x[2];dec=1;}
  1708.       else {p1=cgeti(lx);hiremainder=0;dec=0;}
  1709.       for(i=2+dec;i<lx;i++) p1[i-dec]=divll(x[i],si);
  1710.       if((long)z==0xffffffff)
  1711.     {
  1712.       avma=av;if(!hiremainder) return gzero;
  1713.       p2=cgeti(3);p2[1]=(sx<<24)+3;p2[2]=hiremainder;return p2;
  1714.     }
  1715.       if(lx!=(dec+2)) {p1[1]=p1[0];setsigne(p1,sy);} else {avma=av;p1=gzero;}
  1716.       if(z==0) return p1;
  1717.       if(!hiremainder) *z=gzero;
  1718.       else {p2=cgeti(3);p2[1]=(sx<<24)+3;p2[2]=hiremainder;*z=p2;}
  1719.       return p1;
  1720.     }
  1721.   else
  1722.     {
  1723.       p1=cgeti(lx);
  1724.       sh=bfffo(y[2]);
  1725.       if(sh)
  1726.     {
  1727.       p2=cgeti(ly);k=shiftl(y[2],sh);
  1728.       for(i=3;i<ly;i++) 
  1729.         {
  1730.           k1=shiftl(y[i],sh);p2[i-1]=k+hiremainder;k=k1;
  1731.         }
  1732.       p2[ly-1]=k;k=0;
  1733.       for(i=2;i<lx;i++)
  1734.         {
  1735.           k1=shiftl(x[i],sh);p1[i-1]=k+hiremainder;k=k1;
  1736.         }
  1737.       p1[lx-1]=k;
  1738.     }
  1739.       else {p1[1]=0;for(j=2;j<lx;j++) p1[j]=x[j];p2=y;}
  1740.       si=p2[2];saux=p2[3];
  1741.       for(i=1;i<=lz+1;i++)
  1742.     {
  1743.       if(p1[i]==si) 
  1744.         {
  1745.           qp=0xffffffff;k=addll(si,p1[i+1]);
  1746.         }
  1747.       else
  1748.         {
  1749.           hiremainder=p1[i];qp=divll(p1[i+1],si);
  1750.           overflow=0;k=hiremainder;
  1751.         }
  1752.       if(!overflow)
  1753.         {
  1754.           k1=mulll(qp,saux);k3=subll(k1,p1[i+2]);k+=overflow;
  1755.           flk4=((ulong)hiremainder>(ulong)k);k4=subll(hiremainder,k);
  1756.           while(flk4) {qp--;k3=subll(k3,saux);k4-=overflow;flk4=((ulong)k4>(ulong)si);k4=subll(k4,si);}
  1757.         }
  1758.       hiremainder=0;
  1759.       for(j=ly-1;j>=2;j--)
  1760.         {
  1761.           k1=addmul(qp,p2[j]);
  1762.           p1[i+j-1]=subll(p1[i+j-1],k1);hiremainder+=overflow;
  1763.         }
  1764.       if((ulong)p1[i]<(ulong)hiremainder)
  1765.         {
  1766.           overflow=0;qp--;
  1767.           for(j=ly-1;j>=2;j--) p1[i+j-1]=addllx(p1[i+j-1],p2[j]);
  1768.         }
  1769.       p1[i]=qp;
  1770.     }
  1771.       av1=avma;
  1772.       if((long)z!=0xffffffff)
  1773.     {
  1774.       if(p1[1])
  1775.         {
  1776.           p3=cgeti(lz+3);for(j=2;j<=lz+2;j++) p3[j]=p1[j-1];
  1777.         }
  1778.       else
  1779.         {
  1780.           p3=cgeti(lz+2);
  1781.           if(!lz) sy=0;else for(j=2;j<=lz+1;j++) p3[j]=p1[j];
  1782.         }
  1783.       if(lg(p3)<3) p3[1]=2;else {p3[1]=p3[0];setsigne(p3,sy);}
  1784.     }
  1785.       if(z!=0)
  1786.     {
  1787.       for(j=lz+2;(j<lx)&&(!p1[j]);j++);
  1788.       if(j==lx) p4=gzero;
  1789.       else
  1790.         {
  1791.           p4=cgeti(lx-j+2);p4[1]=p4[0];
  1792.           if(!sh) for(i=2;j<lx;j++,i++) p4[i]=p1[j];
  1793.           else
  1794.         {
  1795.           hiremainder=0;k1=shiftlr(p1[j++],sh);k=hiremainder;
  1796.           if(k1) {p4[2]=k1;dec=1;} else {p4[1]=p4[0]-1;p4++;avma+=4;p4[1]=p4[0];dec=0;}
  1797.           for(i=2+dec;j<lx;j++,i++)
  1798.             {
  1799.               p4[i]=shiftlr(p1[j],sh)+k;k=hiremainder;
  1800.             }
  1801.         }
  1802.           setsigne(p4,sx);
  1803.         }
  1804.     }
  1805.       if((long)z==0xffffffff) return gerepile(av,av1,p4);
  1806.       if((long)z==0) return gerepile(av,av1,p3);
  1807.       av2=avma;dec=lpile(av,av1,0)>>2;*z=adecaler(p4,av1,av2)?p4+dec:p4;
  1808.       return adecaler(p3,av1,av2)?p3+dec:p3;
  1809.     }
  1810. }
  1811.  
  1812. GEN divrr(x,y)
  1813.      GEN x,y;
  1814. {
  1815.   long sx=signe(x),sy=signe(y),lx,ly,lz,ex,ex1,i,z0;
  1816.   long ldif,y0,y1,si,saux,qp,k,k1,k3,k4,j,flk4;
  1817.   GEN z;
  1818.   
  1819.   if(!sy) err(diver9);
  1820.   ex=expo(x)-expo(y)+0x800000;
  1821.   if(ex<=0) err(diver10);
  1822.   if(ex>=0x1000000) err(diver11);
  1823.   if(!sx)
  1824.     {
  1825.       z=cgetr(3);z[1]=ex;z[2]=0;return z;
  1826.     }
  1827.   lx=lg(x);ly=lg(y);lz=(lx<=ly)?lx:ly;
  1828.   z=cgetr(lz);if(sy<0) sx= -sx;
  1829.   ex1=(sx<<24)+ex;
  1830.   if(ly==3)
  1831.     {
  1832.       i=x[2];si=(lx>3)?x[3]:0;
  1833.       if((ulong)i<(ulong)y[2])
  1834.     {
  1835.       hiremainder=i;z[2]=divll(si,y[2]);
  1836.       z[1]=ex1-1;return z;
  1837.     }
  1838.       else
  1839.     {
  1840.       hiremainder=((ulong)i)>>1;
  1841.       z[2]=(i&1)?divll((((ulong)si)>>1)|(0x80000000),y[2]):divll(((ulong)si)>>1,y[2]);
  1842.       z[1]=ex1;return z;
  1843.     }
  1844.     }
  1845.   z0= *z;*z=0;
  1846.   for(i=2;i<=lz-1;i++) z[i-1]=x[i];
  1847.   z[lz-1]=(lx>lz) ? x[lz] : 0;
  1848.   ldif=ly-lz;if(!ldif) {y0=y[lz];y[lz]=0;}
  1849.   if(ldif<=1) {y1=y[lz+1];y[lz+1]=0;}
  1850.   si=y[2];saux=y[3];
  1851.   for(i=0;i<lz-1;i++)
  1852.     {
  1853.       if(z[i]==si) 
  1854.     {
  1855.       qp=0xffffffff;k=addll(si,z[i+1]);
  1856.     }
  1857.       else
  1858.     {
  1859.       hiremainder=z[i];qp=divll(z[i+1],si);
  1860.       overflow=0;k=hiremainder;
  1861.     }
  1862.       if(!overflow)
  1863.     {
  1864.       k1=mulll(qp,saux);k3=subll(k1,z[i+2]);k+=overflow;
  1865.       flk4=((ulong)hiremainder>(ulong)k);k4=subll(hiremainder,k);
  1866.       while(flk4) {qp--;k3=subll(k3,saux);k4-=overflow;flk4=((ulong)k4>(ulong)si);k4=subll(k4,si);}
  1867.     }
  1868.       mulll(qp,y[lz+1-i]);
  1869.       for(j=lz-i;j>=2;j--)
  1870.     {
  1871.       k1=addmul(qp,y[j]);
  1872.       z[i+j-1]=subll(z[i+j-1],k1);hiremainder+=overflow;
  1873.     }
  1874.       if((ulong)z[i]<(ulong)hiremainder)
  1875.     {
  1876.       overflow=0;qp--;
  1877.       for(j=lz-i;j>=2;j--) z[i+j-1]=addllx(z[i+j-1],y[j]);
  1878.     }
  1879.       z[i]=qp;
  1880.     }
  1881.   if(!ldif) y[lz]=y0;if(ldif<=1) y[lz+1]=y1;
  1882.   for(j=lz-1;j>=2;j--) z[j]=z[j-1];
  1883.   if(*z)
  1884.     {
  1885.       k=0x80000000;
  1886.       for(j=2;j<lz;j++) {z[j]=shiftlr(z[j],1)+k;k=hiremainder;}
  1887.     }
  1888.   else ex1--;
  1889.   z[1]=ex1;*z=z0;return z;
  1890. }
  1891.  
  1892. double rtodbl(x)
  1893.      GEN x;
  1894. {
  1895.   double y,ma,co=32.0;
  1896.   long ex,s=signe(x);
  1897.   
  1898.   if((!s)||((ex=expo(x))< -1023)) return 0.0;
  1899.   if(ex>=0x3ff) err(rtodber);
  1900.   ma=(ulong)x[2]+(ulong)x[3]/exp2(co);
  1901.   y=exp2(ex+1-co);
  1902.   return (s<0) ? -y*ma : y*ma;
  1903. }
  1904.  
  1905. GEN dbltor(x)
  1906.      double x;
  1907. {
  1908.   GEN z;
  1909.   double co=32.0;
  1910.   ulong y;
  1911.   long ex;
  1912.   int n;
  1913.   
  1914.   if(x==0) {z=cgetr(3);z[2]=0;z[1]=0x800000-308;return z;}
  1915.   z=cgetr(4);if(x<0) {setsigne(z,-1);x= -x;} else setsigne(z,1);
  1916.   ex=floor(log2(x));setexpo(z,ex);x=x*exp2(-ex+co-1);y=x;x-=y;z[2]=y;
  1917.   y=exp2(co)*x;z[3]=y;
  1918.  
  1919.   if(!(z[2]&(1<<31)))
  1920.     {
  1921.       n=bfffo(z[2]);z[2]=shiftl(z[2],n);z[3]=shiftl(z[3],n);
  1922.       z[2]+=hiremainder;setexpo(z,expo(z)-n);
  1923.     }
  1924.   return z;
  1925. }
  1926.  
  1927. double  gtodouble(x)
  1928.      
  1929.      GEN     x;
  1930.      
  1931. {
  1932.   GEN     x1;
  1933.   long  t=typ(x);
  1934.   static long reel4[4]={0x2010004,0,0,0};
  1935.   
  1936.   if (t==2) x1=x;
  1937.   else gaffect(x,x1=(GEN)reel4);
  1938.   return rtodbl(x1);
  1939. }
  1940.   
  1941. GEN gerepile(l,p,q)
  1942.     GEN l,p,q;
  1943.  
  1944. {
  1945.   long av,declg,tl;
  1946.   GEN ll,pp,l1,l2,l3;
  1947.  
  1948.   declg=(long)l-(long)p;if(declg<=0) return q;
  1949.   for(ll=l,pp=p;pp>(GEN)avma;) *--ll= *--pp;
  1950.   av=(long)ll;
  1951.   while((ll<l)||((ll==l)&&(long)q))
  1952.   {
  1953.     l2=ll+lontyp[tl=typ(ll)];
  1954.     if(tl==10) {l3=ll+lgef(ll);ll+=lg(ll);if(l3>ll) l3=l2;}
  1955.     else {ll+=lg(ll);l3=ll;} 
  1956.     for(;l2<l3;l2++) 
  1957.       {
  1958.     l1=(GEN)(*l2);
  1959.     if((l1<l)&&(l1>=(GEN)avma))
  1960.       {
  1961.         if(l1<p) *l2+=declg;
  1962.         else
  1963.           if(ll<=l) err(gerper);
  1964.       }
  1965.       }
  1966.   }
  1967.   if((!((long)q))||((q<p)&&(q>=(GEN)avma)))
  1968.   {
  1969.     avma=av;return q+(declg>>2);
  1970.   }
  1971.   else {avma=av;return q;}
  1972. }
  1973.  
  1974. void cgiv(x)
  1975.      GEN x;
  1976. {
  1977.   long p;
  1978.  
  1979.   if((p=pere(x))==255) return;
  1980.   if((x!=(GEN)avma)||(p>1)) {setpere(x,p-1);return;}
  1981.   do x+=lg(x);while(!pere(x));
  1982.   avma=(long)x;
  1983.   return;
  1984. }
  1985.