James Amundson wrote:
> ---------- Forwarded message ----------
> Date: Tue, 22 Feb 2005 10:59:33 -0500 (EST)
> From: George Gesslein II
> To: maxima@www.math.utexas.edu
> Subject: Two bugs in one
>
> Hello,
>
> The following is problem that I had with version 5.9.0 of Maxima
> and still with V5.9.1:
>
> (%i1) x^(2^.5);
> 1.414213562373095
> (%o1) x
> (%i2) ratsimp(%);
>
> RAT replaced 1.414213562373095 by 8119//5741 = 1.414213551646055
> 8119/5741
> (%o2) x
> (%i3)
>
> First, it shouldn't have converted 2^.5 to float,
Floating point numbers are "contagious" over integers, by design.
Most people seem to want that. If you prefer exact rationals,
I suggest you use x^(2^(1/2)). There is a long defense of the
strategy in Macsyma in the documentation, which you may not agree
with. I disagree with some of it myself, but not this part.
There is a flag "keepfloat" which affects how rat() works with
respect to floating point numbers. Some operations that would
ordinarily be routine for rational functions, like removing common
polynomial factors by GCD, are complicated by keeping floating point
numbers around, so most people do not use rat() or ratsimp() with
floats. Or they use it for polynomials only.
The usual use of floats in Maxima is to defer their
use as long as possible to avoid introducing approximations.
Sometimes this is not a good idea and is much more computationally
expensive than using floats earlier, so you have to be judicious
and it certainly
> shouldn't convert the float to 8119/5741.
Well, it differs from the sqrt(2) by about 10^(-8), so it is good
to about single precision.
>
> I can explain how to better validate floats to fractions,
> by using an epsilon value. Here is my C code:
>
> if (fabs(k3 - k4) > (small_epsilon * fabs(k3))) {
> return false; /* fraction is too inaccurate */
> }
>
> k3 = the original float;
> k4 = the new float obtained by dividing out the fraction;
> small_epsilon = 0.000000000000005; a good value for doubles.
>
>
> Regards,
> George
>
> George John Gesslein II
> http://www.mathomatic.com
> "The light of souls is a warm, little glow."
ANSI Common Lisp has two functions, rational and rationalize which
can be used to convert floats to rationals. Maxima may use its
own, based on its non-ANSI past.
The relative error between two values r1, r2 is
abs((r1-r2)/r1) if r1 is not zero.
Your program asks if the relative error is less than epsilon.
This does not tell you how to find a good rational approximation
to a floating point number. I think you will find such an algorithm
in the Common Lisp literature.
rationalize finds a compact representation that is pretty good to the
given precision. rational finds a possibly clumsier ratio that attempts
to recover "all" the bits. This is always possible if start with a ratio,
float it, and then convert back. NOte that you have big integers in CL.
since a float internally in usual IEEE standard arithmetic is
fraction X 2^exponent , by computing an appropriate power of 2,
positive or negative, the float can be made into a rational number.
A typical common lisp includes the following constants:
single-float-epsilon value: 5.960465e-8
double-float-epsilon value: 1.1102230246251568d-16
Thanks for your interest.
RJF