Ruud H.G. van Tol on Wed, 29 Nov 2023 22:16:01 +0100


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

Re: log_int_rat




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?
Using logint(2^n,3/2) would be much slower even if it worked, because that would
force you to compute 2^n.

Really, the only defect of your function is that you do not take the accuracy
into account so

? a054414(2^1000)
   *** floor: precision too low in truncr (precision loss in truncation).

but this is easily fixed:

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

? a054414(2^1000)
%8 = 29032646699514618648207963388262277036465203935797085941659905401800597040332207967488444244921882872990640968909754332584729580377856680178651318617855116668095773049095433042501905586769687829752260479426676858529835781762375806913040422676074996080564646620120844538289072080537948895898995900461439
   ***   last result computed in 0 ms.

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]


Thanks for the broad answers, Karim's as well; I will certainly use parts of it.

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));

Quotes from the page:
<>
a(n) = 1 + floor(n/(1-log(2)/log(3)))
Also, numbers n such that the first digit in ternary expansion on 2^n is 2.
</>

-- Ruud