R J Cano on Sat, 22 Apr 2017 18:27:22 +0200


[Date Prev] [Date Next] [Thread Prev] [Thread Next] [Date Index] [Thread Index]

Re: 4 Gigabytes of RAM are not enough!


 Hi!,

 I had a wrong result... some things cannot be grouped inside a single
"while" loop as I shown previously.

 But I fixed it... a little bit slower but correct.

 Cheers,

 Remy

 P.S.: Bonus feature, but beware with overflows when using w=1 option.

Notice if 2<=m<n, how the first m!-1 terms given by Ticket1b(m) are
the same corresponding values for Ticket1b(n), but multiplied/divided
by B^(n-m);
/* (PARI) R. J. Cano, Apr 22 2017; Note: For a best performance, use this script with GP2C. */

Ticket1b(n,{B=10},{w=0})={ /* Computes a "final representation" (by default in decimal) for the first n!-1 terms of A217626, storing them as a "t_VECSMALL" object */

  n=max(2,n);
  my(x0:vecsmall,x:vecsmall,y:vecsmall,i:small,j:small,k:small,u:small,t:small,Q:vecsmall,p:small,bb:small,m:small);
  x=Vecsmall(vector(n,i,i-1));
  m=n!-1;
  Q=Vecsmall(0,m);
  k=0;
  bb=B;

  while(k++<=#Q,

    x0=x;

   if(w,
     
      /* Implements Schmuel Zak's algorithm -- Start */
      u=1;
      r=k;
      while(u++<k,if(r%u,break,r\=u)); /* These 3 initial lines are an adaptation from the Charles Greathouse IV 's A055881 program. */
      j=u;
      j\=2;
      i=0;
      while(i++<=j,
       t=x[i];
        x[i]=x[1+u-i];
        x[1+u-i]=t
      )
      /* Implements Schmuel Zak's algorithm -- End */
      
   ,
   
      /* Implements Narayana Pandit's lex-order algorithm -- Start */
      i=#x-1;
      while((x[i]>=x[i+1])&&i,i--);
      j=#x;
      while((x[i]>=x[j])&&(j>i),j--);
      t=x[j];
      x[j]=x[i];
      x[i]=t;
      u=#x;
      u-=i;
      u\=2;
      u++;
      while(u--,
        t=x[i+u];
        x[i+u]=x[1+#x-u];
        x[1+#x-u]=t;
      )
      /* Implements Narayana Pandit's lex-order algorithm -- End */
      
    );

    y=x;
    
    r=0;
    while(r++<#y, y[r  ]-= x0[r] );
    r=0;
    while(r++<#y, y[r+1]+=  y[r] );
    p=1;
    r=0;
    while(r++<#y, Q[k  ]+=  y[#y-r]*p; p*=bb )

  );

  return(Q);

}