[<prev] [next>] [<thread-prev] [day] [month] [year] [list]
Message-ID: <20260428082645.GZ3520958@port70.net>
Date: Tue, 28 Apr 2026 10:26:45 +0200
From: Szabolcs Nagy <nsz@...t70.net>
To: Damian McGuckin <damianm@....com.au>
Cc: MUSL <musl@...ts.openwall.com>
Subject: Re: Implementation of pow(x, y)
* Damian McGuckin <damianm@....com.au> [2026-04-28 14:44:05 +1000]:
> Is there any documentation of the algorithm in the latest pow(x,y)?
no, only in the source. i can collect the important bits:
pow(x,y) = exp(y*log(x))
special cases:
(x < 0x1p-1022 or inf or nan) or
(|y| < 0x1p-65 or |y| >= 0x1p63 or nan)
for details see first if block in pow()
log(x):
given x = 2^k z where z is in [OFF,2*OFF)
log(x) = k ln2 + log(c) + log(z/c)
log(z/c) = poly(z*invc - 1)
log(c) = logc + logc_tail
7bits of z is used to lookup invc, logc, logc_tail
and z*invc-1 is in (-1/128,1/128) and computed exactly
log(x) is computed as hi + lo to 0x1p-15 ulp precision
for details see pow_data.c and log_inline() in pow.c
y*log(x):
ehi + elo = y*(hi + lo)
exp(hi+lo):
given hi+lo = ln2/N*k + r, with int k and r in [-ln2/2N, ln2/2N]
exp(hi+lo) = 2^(k/N) * exp(r), with exp(r) in [2^(-1/2N),2^(1/2N)]
k = round(N/ln2 * hi)
r = x - ln2/N*k + lo
2^(k/N) = scale * (1 + tail)
exp(hi+lo) ~= scale + scale*(tail + exp(r) - 1)
tail+exp(r)-1 = poly(r)
for details see exp_inline() in pow.c
> Are there any plans for rootn(x, n) and pown(x, n) routines or even powr(x,
> y);
i dont have plans.
Powered by blists - more mailing lists
Confused about mailing lists and their use?
Read about mailing lists on Wikipedia
and check out these
guidelines on proper formatting of your messages.
