[Date Index] [Thread Index] [Author Index]

Re: Re: Re: inconsistency with Inequality testing and Floor (long)

• To: mathgroup at smc.vnet.net
• Subject: [mg60155] Re: [mg60119] Re: [mg60079] Re: inconsistency with Inequality testing and Floor (long)
• From: Andrzej Kozlowski <akoz at mimuw.edu.pl>
• Date: Mon, 5 Sep 2005 02:27:48 -0400 (EDT)
• References: <200508251034.GAA10208@smc.vnet.net> <demmfd$rf1$1@smc.vnet.net> <200509010613.CAA08855@smc.vnet.net> <200509030606.CAA19076@smc.vnet.net>
• Sender: owner-wri-mathgroup at wolfram.com

On 3 Sep 2005, at 15:06, Andrzej Kozlowski wrote:

>
> On 1 Sep 2005, at 13:13, Richard J. Fateman wrote:
>
>
>>
>>
>> Andrzej Kozlowski wrote:
>>
>>
>>
>>
>>>
>>> You cannot expect inexact numbers, particularly "borderline  
>>> cases" as
>>> in this example, to obey the usual laws of arithmetic.
>>>
>>>
>>>
>>
>> These inconsistencies are not inherent in computer arithmetic,
>> even floating-point arithmetic. They represent choices made by the
>> designers of Mathematica.  In particular, choosing to allow two
>> distinct numbers to be numerically equal leads to problems. I
>> wonder if there is some compensating good reason for this choice
>> in Mathematica. I am not aware of any reason good enough to
>> justify this.
>>
>> Calling a bug a feature does not fix it.
>>
>> RJF
>>
>>
>>
>
>
> I can see a number of reasons for keeping numerical equality between
> approximate numbers in Mathematica -the main is that it makes many
> things simpler, in particular in  numerical polynomial algebra (for
> example in numerical Groebner basis). It is not strictly  necessary
> but it is convenient in interval based computations, which is what
> significance arithmetic essentially does. And by the way, while
> writing a review of Hans Stetter's "Numerical Polynomial Algebra", I
> have found that  Mathematica's significance arithmetic with its
> automatic precision tracking performs performs very well when
> combined with   backward error based techniques that Stetter uses; in
> fact quite a lot better than the other well known systems he was
> using. I could even send you some examples of Matheamtica's
> superiority but in view of the countless disputes on this topic that
> you have initiated and that never get anywhere,   I think it would be
> pointless. If you are interested you can read my review when it
> appears in Math Reviews as actually discusses Mathematica's numerics.
>
> Andrzej Kozlowski
>
>

It came to my mind that it might be a good idea to elucidate the  
above remarks somewhat. I have no intention to engage in any more  
polemics with RJF  but I think for the benefit of others it may be  
useful to try to explain what seem to me to be concrete benefits of  
the Mathematica approach to approximate computations. More precisely:  
I believe that this approach has shown itself to have clear  
advantages in at least one important area of computational  
mathematics known as "numerical polynomial algebra". I think this  
aspect of Mathematica has not been sufficiently advertised, probably  
because the subject is rather technical and difficult to explain so  
that the only person who can really do it properly is Daniel  
Lichtblau but natural modesty prevents him from putting as much  
stress on this as it deserves.  While I cannot match Daniel's  
knowledge and experience in this area I have the advantage of not  
needing to worry about engaging in self-praise when writing about  
this topic. On the other hand I do not want to try to write a very  
long essay o both because this is not the proper forum for this sort  
of thing and because trying to do so would expose the superficiality  
of my knowledge of this field. So I will be rather vague and if any  
one wants to know more I refer all questions to Daniel ;-)

Polynomial Algebra is primarily concerned with the solution of  
systems of polynomial equations. It is one of the oldest areas of  
mathematics going back back to the ancient Babylonians. Linear  
Algebra is a special case where all polynomials have degree =1.
I think few people who know about modern symbolic algebra would deny  
that Polynomial Algebra is one of the most important areas, even if  
one does not include Linear Algebra as part of it. Until recently  
there has been a very big difference in the way Linear Algebra  
(degree 1 case) and Polynomial Algebra (degree >=2 case) has been  
implemented in CAS systems. In the case of linear algebra there are  
usually quite different algorithms for solving the corresponding  
problems in a symbolic (or exact) and in numerical setting. The  
reason for this is primarily performance: when working with  
approximate (floating point) numbers specially designed "approximate"  
algorithms perform much better than the "exact ones". The situation  
for polynomial systems of degree >=2 has been quite different. Until  
relatively recently there existed essentially no (multivariate)  
Numerical Polynomial Algebra in the same sense as Numerical Linear  
Algebra  mentioned above. There are several reasons for this, the  
main one being the dominance of the field of multivariate polynomial  
algebra by the famous Groebner Basis algorithm due to Bruno  
Buchberger. This is one of the most important algorithms in all of  
computer algebra but it has one striking weakness: it is numerically  
unstable. In other words, when used with floating point numbers it  
may produce correct answers for certain values but wildly wrong ones  
for others. I don't want to go into the reason for this but they are  
caused by the need to choose a "term order". In any case, what this  
meant in practise was that faced with a problem in numerical  
polynomial algebra requiring a Groebner basis almost all CAS would do  
the following: they would either rationalise the coefficients and  
solve the exact Groebner bassi problem using Buchberger"s algorithm  
or just go ahead and apply it directly to (fixed precision) floating  
point data. The first approach involves a big performance hit (in  
fact many problems can"t be solved int his way at all) while the  
second risks returning a completely meaningless answer. This has been  
the situation and as far as I know still is the situation with all  
the standard CAS systems except one. Before I reveal the name of the  
exception let me state  that it is not the old and now defunct CAS  
system that RJF was involved with. As we are not supposed to mention  
the names of other CAS here let me just refer to it as the FF system  
( you can think of your own reasons for this choice of name).
  I would never even consider suggesting that RJF is anything but  
altruistic and objective in his periodic assaults on Mathematica but   
but I also note that he has sent me in the past personal messages  
with remarks like "The FF System can do this while Mathematica  
can't". For some reason I completely forgot to reply "Matheamtica has  
a Namerical Groebner Basis and has had it for nearly a decade, how  
about FF ?".
Indeed, Daniel Lichtblau succeeded many years ago in implementing a  
numerical-symbolic Groebner basis, which is used by NSolve and  
perhaps some other polynomial algebra functions. It works pretty well  
as I have had the chance to check myself while writing a review of a  
new book on "Numerical Polynomial Algebra" by Hans Stetter. It means  
that NSolve can deal with problems that other systems cannot without  
expert  human intervention. And the key to Daniel's implementation is  
the Mathematica model of big num arithmetic, based on what is known  
as "significance arithmetic". The key role is played by the fact that  
thanks to significance aritmetic Mathematica can perform automatic  
tracking of Precision/Accuracy. RJF has in the past written  
dismissively of this as of something that can be of use only to  
"sloppy people". In a sense he is of course right. Someone skilled in  
polynomial algebra and numerical analysis (which I suspect is still a  
set containing only a few individuals) does not need a numerical  
Groebner basis, as Hans Stetter impressively demonstrates in his  
book. But not all of us have the time, the sill an the patience to  
overcome our sloppiness at dealing with rather esoteric numerical  
issues.
This may be a good point to say a few words about this issue of  
equality between approximate numbers and why it is useful to have it  
in Mathematica. The basic fact is that (over simplifying somewhat)  
algebra works with equalities while analysis inequalities. To see  
what i mean compare the two equivalent (over the real numbers)  
statements: a=b and  b-epsilon<a<b+epsilon for every epsilon >0. The  
first belongs to algebra, the second to analysis. From the point of  
view of analysis equality is too "unstable", even the slightest  
perturbation of the parameters will break it. Thus analysis are  
generally experts in dealing with inequalities, algebraist in solving  
equations. Numerical analysis is really a subspecies of analysis and  
thus it prefers inequalities to equalities. Suppose however we are  
interested in "numerical algebra'. From the point of view of a  
computer algebra has a big advantage over analysis: it's algorithms  
being of discrete nature are much easier to implement as computer  
programs. But if we are going to do algebra we need to deal with  
equality.  So what can be the meaning of an equality a==b where a and  
be are inexact numbers, for example how should we interpret an  
equation f[x]=0 where f is a polynomial whose coefficients are  
inexact numbers? There are two basic approaches. The more common one  
is to interpret an equation as always involving exact numbers only.  
Instead we talk of an approximate root of such an equation by which  
we mean an exact root of another exact equation whose coefficients  
are very close to the coefficients of our original equation. We can  
then use ordinary algebra to find such a root, but when  veryfying if  
it is a 'close enough" root or trying to refine it we use methods  
which essentially belong to analysis. If x1 is our "approximate root"  
we do not talk of the equality f[x1]==0 but we talk of the "validity"  
of x1 as a root of the original equation. x1 is a valid root if it is  
a root of a polynomial equation whose coefficients are sufficiently  
close to those of the original equation. This involves dealing with  
inequalities.
Mathematica does not take this approach. Instead it represents  
approximate numbers as intervals with centre at some fixed value  
whose lenght reflects the precision or accuracy of our data. These  
intervals can then be directly substituted into various algebraic  
formulae; they can be added , multiplied etc. In other words we can  
do "algebra' with them. Since we want to do algebra with such  
quantities we need also the notion of equality. Of course the rules  
of our "algebra" as were as of our equality will not be quite the  
same as those we are used to and we will end up with some unintuitive  
results (wand this is made somewhat worse by the fact that in order  
to improve performance significance arithmetic is only an  
approximation to interval arithmetic). So there is a price that has  
to be paid. But the gains such as automatic tracing of precision and  
ultimately numerical Groebner basis are, in my opinion, worth this  
price.
Finally, I have to say that recently a great deal of progress has  
been made in the area of numerical polynomial algebra. In particular  
Kondratieff, a student of Stetter, has described (in his PhD thesis)  
a fixed precision implementation of numerical Groebner basis. This  
approach may even be in some ways superior than Mathematica's but as  
much as I have understood of it seems considerably more complicated  
and even less intuitive than significance arithmetic.

Andrzej Kozlowski

• References:
  • Re: inconsistency with Inequality testing and Floor
    • From: "Richard J. Fateman" <fateman@eecs.berkeley.edu>
  • Re: Re: inconsistency with Inequality testing and Floor
    • From: Andrzej Kozlowski <andrzej@akikoz.net>