set n 123456789......(lotsa digits).....0123456789
expr {sqrt($n)} -> 1.1111111061111111e+49
converting this to wide or int, gives me utter garbage
(something modulo 2**64 or 2**32)
In my case, I'd even be satisfied with something like:
111111110611111110000000000000000000...0000
thus, a large integer in "some" vicinity of
sqrt($n), but I don't even get that :-(
Of course,
11111111061111110993611110581861108108154842000984
(obtained with "bc") would be the ultimate goal.
% expr entier(1.1111111061111111e+49)
11111111061111110659303881096431740515280041803776
about a few minutes after I wrote my previous posting,
I idly tried: info commands ::tcl::mathfunc::*
and noticed the "entier" (the name didn't convey anything
at all to me, and it is not yet mentioned in the mathfunc
page), so I just tried:
% expr {entier(1.0e49)}
9999999999999999464902769475481793196872414789632
I would have hoped for sthg like: 1[string repeat 0 49]
I know, "decimal vs binary" and all that... :-(
While it is exactly what I asked for, I now notice that
my question was wrong in the first place.
PS:
meanwhile I changed the algorithm, and now it needs
an exact-to-the-last-integer-digit sqrt (and it is
now three times faster than before) and my "almost exact"
square-root is now done like this:
namespace eval ::tcl::mathfunc {
proc isqrt {n} {
string map {\\\n {}} [exec bc -q << [format sqrt(%lld)\n $n]]
}
}
(if this bc just weren't so stubborn about breaking lines...
oh, there is also "dc", but how do I efficiently read in
a "P"-style output from "dc" ?)
I actually don't need anything *after* the decimal point.
I downloaded it anyway... I'll see if it can still help me.
Thanks.
We already have square-root-of-bignum in the library, but it's
not exported. Primarily for reasons of controlling "code bloat",
the bignum stuff hasn't introduced any new math functions.
And sqrt() needs to preserve its old semantics, otherwise
code that expects the result of sqrt(int) to be a floating point
number will break rather badly.
What I'm planning to do - in the fairly near future - is come
up with an extension that can be [package require]d that adds
the missing functionality:
- direct base64 and base256 encoding/decoding of bignums
- greatest common divisor, least common multiple,
extended Euclid algorithm
- integer square (squaring is faster than multiplying)
and square root
- modular inverse
- modular exponentiation
- Legendre-Jacobi symbol
- primality testing
- random prime generation
- Possibly, optimizations for reduction by fixed moduli
(Barrett, Montgomery, diminished-radix)
Is there anything obvious missing from this list that would
be needed for your application?
If, in the meantime, you need a better square root urgently,
the double-precision square root followed by a couple of rounds
of Newton-Raphson ought to work. (Ugly, I know, but I haven't
anything better to offer just at the moment.)
--
73 de ke9tv/2, Kevin
> What I'm planning to do - in the fairly near future - is come
> up with an extension that can be [package require]d that adds
> the missing functionality:
> - direct base64 and base256 encoding/decoding of bignums
base256 would be "binary encoding", or something else ?
At the moment I was rather surprised, that specifying
a number as a decimal string appeared to be faster than
using [scan $hexstring %llx]. Are bignums internally
decimal-based, or is this just another "not-yet-optimised" ?
In my opinion, those (below) I marked with a (*) are
simple or common enough, to not have to require
package-requiring, those with (**) are deeper into
number-theory or computationally more complex, so
having to package.require them is fine:
> - greatest common divisor, least common multiple,
> extended Euclid algorithm
cool. (*)
I'm not sure about the "extended Euclid algorithm":
If it is not a means to do gcd&lcm, I'd say: (**)
> - integer square (squaring is faster than multiplying)
> and square root
Extra function for squaring? Couldn't that be just an internal
optimization of "$bignum**2" ?
Integer-squareroot: cool (*) (suggestion for name: isqrt)
> - modular inverse
> - modular exponentiation
In my eyes this is really advanced stuff ... (**)
> - Legendre-Jacobi symbol
Isn't "==" just that?
> - primality testing
(**)
Since it may take a loooooong time for certain arguments,
it is perhaps better to do this outside of expr. Then it
could also be event-loop-compatible...
> - random prime generation
(**)
> - Possibly, optimizations for reduction by fixed moduli
> (Barrett, Montgomery, diminished-radix)
(**) or (*) depending on the "code-bloat" this would really cause.
> Is there anything obvious missing from this list that would
> be needed for your application?
For my "application", all I actually need (and what isn't already
there) is:
isqrt,
faster conversion (if internally appropriate) from/to
(big-endian)binary(base=256),bin(base=2),hex strings.
> If, in the meantime, you need a better square root urgently,
> the double-precision square root followed by a couple of rounds
> of Newton-Raphson ought to work. (Ugly, I know, but I haven't
> anything better to offer just at the moment.)
I currently use external program "dc" for square-rooting.
on my maschine it takes about 4 milliseconds per call.
calling out to "dc" and having it return the original
number takes about 2 milliseconds, so I'm looking forward
to how long a builtin isqrt would take. I'd be very surprised,
if Newton-Raphson implemented in pure Tcl could compete with
calling out to bc...
Btw. having dc return a hex-number and then [scan ... %llx] takes
6 millisecs, and having dc return base256-binary fails due
to non-binarysafe exec.
> string map {\\\n {}} [exec bc -q << [format sqrt(%lld)\n $n]]
Meanwhile changed to:
string map {\\\n {}} [exec dc -e ${n}vp]]
... [scan [string map ... [exec dc -e 16o${n}vp]] %llx] ...
turned out to be slower... why is conversion from
decimal faster than from hexadecimal?
> oh, there is also "dc", but how do I efficiently read in
> a "P"-style output from "dc" ?)
This is void anyway, due to deficiencies of [exec] :-(
You might want to look at adapting the paper-and-pencil
method of taking square roots. The algorithm adapts
really nicely to binary numbers and only involves shifts
subtractions and compares. Check out http://wiki.tcl.tk/SquareRoot.
Keith
I now use the snippet, someone recently posted, using
Newton's method, and it is indeed about 25 times faster
than calling out to "dc" :-)
I'll further hold my breath till isqrt gets checked in on
CVS, and blame Kevin Kenny, in case I suffocate ;-)
Not wanting to take responsibility for asphyxiating
Andreas, I've posted a better stopgap (immune to floating
point overflows) at the bottom of http://wiki.tcl.tk/9352.
I hope to have isqrt in by the end of 2005, but with the
holidays coming up, I make no promises.
Thanks for caring :-)
> I hope to have isqrt in by the end of 2005, but with the
> holidays coming up, I make no promises.
I'm somewhat sitting on needles, wanting to know
an estimate of the speed of the final isqrt.
How much faster will it be, e.g. compared to
that recently added version in wiki?
e.g. for 30, 100 and 1000 digit integers (if I may be so bold)
Surely. (Note that I had to fix a performance bug in libtommath
before doing this...) Times are in microseconds on a 1.6 GHz
Pentium-M with a --enable-symbols build. The test constants were
[string repeat 5 $numberOfDigits].
Number of Wiki Built-in
digits version isqrt
30 59 47
100 243 78
300 751 312
1000 4627 2322
Other n-digit constants gave different results (the initial estimate
for the square root is developed differently in the C code) but there
were consistent speedups in the built-in function.
Note that libtommath uses Newton-Raphson internally for its square
root. Prior to my tweak, it was actually slower than the Tcl version,
owing to a poor initial estimate of the square root. It could be
improved still further by replacing the Newton-Raphson square root
with a Zimmermann one (http://www.inria.fr/rrrt/rr-3805.html).
One of these years, I might even get around to it.