Re: log_int_rat

["Ruud H.G. van Tol" <rvtol@isolution.nl>]

On 2024-01-22 04:26, Ruud H.G. van Tol wrote:
> On 2023-11-29 22:15, Ruud H.G. van Tol wrote:
> > On 2023-11-29 21:17, Bill Allombert wrote:
> > > On Wed, Nov 29, 2023 at 07:44:11PM +0100, Ruud H.G. van Tol wrote:
> > > > ? a054414(n) = 1 + n + floor( n * log(2) / log(3/2) );
> > > >
> > > > ? [ a054414(n) |n<-[0..20]]
> > > > % [1, 3, 6, 9, 11, 14, 17, 19, 22, 25, 28, 30, 33, 36, 38, 41, 44, 
> > > > 47, 49,
> > > > 52, 55]
> > > >
> > > > Is there a cleaner way, similar to logint, to do that 
> > > > floor-expression?
> > > [...]
> > > If being slow is not an issue, I suppose you can do:
> > >   a054414(n)=my(a=2^n,b=1);while(a>=b,a*=2;b*=3;n++);n
> > > and
> > > table(n)=my(a=1/2,b=1);vector(n+1,i,a*=2;while(a>=b,a*=2;b*=3);valuation(a,2)); 
> > >
> > >
> > > ? table(20)
> > > %12 = [1,3,6,9,11,14,17,19,22,25,28,30,33,36,38,41,44,47,49,52,55]
> >
> > [...]
> > The https://oeis.org/A054414 page doesn't have a PARI-prog yet, so I 
> > will at least propose
> > alist(n=60) = my(a=1/2, b=1); vector(n+1, i, a*=2; 
> > while(a>=b,a*=2;b*=3); valuation(a, 2));
> > [...]
>
> Now proposed it as:
>
> alist(N=61) = {
>   my(a=1/2, b=1, r=-1);
>   vector(N, i, a*=4; b*=3; r+=2; if(a>b, a*=2; b*=3; r++); r);
> }

On the non-discrete code branch:

  a(n) = localbitprec(logint(n+2,2)); 1 + n + floor(n*log(2)/log(3/2));

with test:

  {default(parisizemax,2^30); my(v=alist(2^20)); for(n=1,#v, 
(v[n]!=a(n-1)) && print(n);break);}

I received this report:

> that code wrong result a(83690) = 226760 (in a 32-bit pari)
> with no hint that it didn't know right from wrong.


? my(v=alist(83700)); v[83691]
cpu time = 485 ms, real time = 490 ms.
% 226759

-- Ruud