Tell me more ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.
up vote 4 down vote favorite
2

I am currently having problems with some floating points.

I have a function, which gives as an intermediate result (for example)

-(10000. - 10000. a) E^(-16.4157 kp^2 x^0.277)/((-1. + a) a x^0.252525 xj^4.47047)

It gets passed through into another function, which will fill in the value a = 1. Then Mathematica reacts with a 1/Sqrt[0] error. However, as you can see, the factor (10000. - 10000. a) should cancel with the factor (-1.+ a) in the denominator to give 10000, the result is thus well-defined.

Do you have any idea how to let Mathematica cancel these factors out? I have tried Simplify and FullSimplify with Assumptions -> a!=1 but it doesn't work. I cannot change much to the function itself, as it is an intermediate result (and it is just an example for this post; other intermediate results also occur, sometimes with the same problem, sometimes they work fine).

share|improve this question
3  
Floating point numbers are different than rational numbers. Since floating point numbers don't really follow the rules of algebra, you'll find that Simplify and the like will not do much at all wit them. - Run the calculations with increased precision using SetPrecision. (still essentially the same problem) I think there are two options: - Use the function Rationalize to turn the floating point numbers into rational numbers and then use FullSimplify. (I'm not a big fan of this because it ignores the fact that the numbers of a limited precision) – Searke Mar 1 '12 at 19:47
@Searke "Floats don't really follow the rules of algebra"?! – David Mar 1 '12 at 19:56
4  
@David That is right. Addition/multiplication is not even associative. – Szabolcs Mar 1 '12 at 19:57
1  
@David for example, (1 + $MaxMachineNumber) == $MaxMachineNumber evaluates to True – acl Mar 1 '12 at 20:02
@David try 1.*10^16 + 1. - 1.*10^16 which returns different answers depending on where 1. is, however -10.^16 + (10.^16 + 1.) and -10.^16 + (1. + 10.^16) return the same incorrect answer. This indicates that addition is not associative, but commutative. Also, given g[a_, b_, c_] := a^2 - 2 a b + b^2 - c^2, g[10.^10, 10.^10, 10.] == 0 implying that mma is doing speculative processing as the first three terms simplify to $(a - b)^2$ and should cancel. – rcollyer Mar 1 '12 at 20:43
show 3 more comments

3 Answers

up vote 5 down vote accepted

The problem is that symbolic computation doesn't mix well with floating point. If you use Rationalize first to convert to an exact representation, Simplify (or the faster Cancel) will be able to do it's job.

share|improve this answer
up vote 5 down vote

Use Limit, Mathematica's implementation of $\lim$:

$\displaystyle\lim_{x\to a}f(x) = $ Limit[f[x], x -> a]

Applied to your equation, this results in

Limit[
    -(((10000. - 10000. a) E^(-16.4157 kp^2 x^0.277)/((-1. + a) a x^0.252525 xj^4.47047))),
    a -> 1
]

(10000. E^(-16.4157 kp^2 x^0.277))/(x^0.252525 xj^4.47047)

The error you are talking about is not caused by floating point arithmetic, it's how Mathematica (or any other language for that matter) handles computations in general. If you substitute $a=1$ then you've got a division by zero, and once this occurs you're out of the domain of, well, mathematics; therefore, you can't cancel the two factors anymore. (Actually, your function is not defined at $a=1$, it has a removable singularity there.)

share|improve this answer
The function is well defined at a = 1. The singularity is artificial, and has more to do with the order of evaluation than a being a property of the function, unlike in $\sin(x)/x$. Since mma evaluates the terms in isolation from one another, i.e. $(a - 1)^{-1}$ is evaluated separately from $1/a$, etc., this is likely to occur. – rcollyer Mar 1 '12 at 20:08
1  
$\frac xx$ stands for the product of x and its multiplicative inverse, which does not exist for $x=0$. It's not well-defined or even defined at that point, although the limit exists of course. – David Mar 1 '12 at 20:10
On the surface, I'd disagree with you, as it is plain that $\frac xx$ simplifies to $1$ prior to the limit being taken. But, upon further review, it appears that removable singularities may in fact correspond to how we evaluate things, as opposed to being intrinsic. After all, $\frac\sin{x} = 1 - x^2/3! + x^4/5!\ \cdots$, implying that the singularity has more to due with the evaluation process than actually being present. – rcollyer Mar 1 '12 at 20:34
It's the same thing for $\sin(x)/x$, not defined at $x=0$. It's of course very convenient to always assume that the removable singularity is removed (by defining the additional point using the limit). Mathematically what you do when canceling fractions is choosing another element of the equivalence class of two numbers, but if the denominator is zero you don't have a rational number in the first place, meaning no arithmetic is defined. (Also see Wikipedia, Constructing the rationals) – David Mar 1 '12 at 20:41
2  
This, of course, reminds me why I hated my analysis course. – rcollyer Mar 1 '12 at 20:44
up vote 5 down vote

Well...

[Warning: This act is performed by trained professionals. Do not try it at home.]

In[24]:= expr = -(10000. - 10000. a) E^(-16.4157 kp^2 x^0.277)/
  ((-1. + a) a x^0.252525 xj^4.47047);

In[25]:= Cancel[expr]
Out[25]//InputForm= 10000./(a*E^(16.4157*kp^2*x^0.277)*x^0.252525*xj^4.47047)

Okay, I cheated. I'm using the development version of Mathematica.

I will add some emphasis to remarks others have made. Doing computer algebra with approximate numbers is almost never a safe bet. Even to the extent that it might be supported, things can go wrong. And do. On a remarkably frequent basis.

share|improve this answer
+1, for cheek, alone. :) – rcollyer Mar 1 '12 at 23:03

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.