From: murison@cfa203.HARVARD.EDU (Marc A. Murison, RG)
Newsgroups: sci.math.symbolic
Subject: Branching frustrations (was: Integral Tables)
Date: 8 Oct 92 15:11:57 GMT
Organization: Smithsonian Astrophysical Observatory, Cambridge, MA,  USA

Here's an example of a branching problem that I stubbed
my toe on recently with Maple V (release 1.1 for DOS).  
Given the expression

> foo := (r+2)*L - (r+1)*d + sqrt(L)*sqrt(d+r^2*(L-d));

                                    1/2       2      2   1/2
    foo := (r + 2) L - (r + 1) d + L    (d + r  L - r  d)

solve for r:

> bar := [solve( foo, r )];

                  4 L - d
      bar := [-1, -------]
                     d

Check by substituting the solutions back into foo, and 
at first glance the solution r=-1 seems wrong, while 
the second one is obviously right:

> simplify( subs(r=bar[1],foo) );
                                      2 L
> simplify( subs(r=bar[2],foo) );
                                       0
But what's happening is that Maple is considering both 
branches of the sqrt.  r = -1 is a solution of the 
equation when the sign of the sqrt term is negative 
(even though I explicitely wrote it positive!).  If you 
get rid of the sqrt by squaring,

> foo := (op(1,foo)+op(2,foo))^2 = op(3,foo)^2;

                                   2           2      2
     foo := ((r + 2) L - (r + 1) d)  = L (d + r  L - r  d)

which is, I suppose, what Maple must be doing, then the 
(same) solutions

> bar := map(simplify,[solve( foo, r )]);
                 4 L - d
     bar := [-1, -------]
                    d
are both valid:

> simplify( subs(r=bar[1],foo) );

                 2    2
                L  = L

> factor( subs(r=bar[2],foo) );

        2            2    2            2
       L  (4 L - 3 d)    L  (4 L - 3 d)
       --------------- = ---------------
               2                 2
              d                 d

When the equation was squared, the fact that I 
originally wanted a specific branch was lost to Maple.  
I suppose you can live with this kind of thing you 
happen to be doing calculations step by painful step 
interactively (i.e., not automating things via 
procedures).  But my original expression

                            1/2       2      2   1/2
   (r + 2) L - (r + 1) d + L    (d + r  L - r  d)    = 0

and many others like it are the intermediate results of 
some Maple procedures I've written to do some lengthy 
calculations.  And of course the procedures merrily 
munch right along and the wrong turns promptly get 
buried by further calculations and are then tough to 
track down. It's very discouraging to get "wrong" 
answers like r = -1, especially since Maple does so 
many other (very difficult!) things wonderfully.  The 
moral of the story is that you must be very very 
careful.  

In hindsight this is pretty obvious.  But I'm finding 
that when you're exploring unknown territory in a 
research project, trying to use a CAS as a tool rather 
than as a toy, the particular pitfalls that await you 
are annoyingly hard to detect until after you figure 
out you got nailed.

Marc A. Murison                | #include <fstream.h>
Smithsonian Astrophysical Obs. | void main( void ) { 
Cambridge, MA 02138            | ifstream in("obquote.txt",ios::in); 
murison@cfa.harvard.edu        | char c; while(in.get(c)) cout.put(c); }
Newsgroups: sci.math.symbolic