solution of system of quadratic eqns.



On 06/24/2012 01:35 PM, Barton Willis wrote:
>> Perhaps solve(), or algsys() can be made to recognize that the problem
>> requires a call to_poly_solve() and compute accordingly.
> Alternative:
>
> load(grobner);
>
> myalgsys(e,v) := (
>    e : map(lambda([s], rhs(s)-lhs(s)),e),
>    algsys(poly_reduced_grobner(e, block([listconstvars : false], listofvars(e))),v));
>
> Well maybe this should check for polynomial equations before calling  poly_reduced_grobner. This code  a bit redundant: algys
> doesn't know that the system is triangular.  I'm not sure that the author of the grobner package is interested in maintaining it.
> Maybe it is bug free :)
>
>> By the way, did you mean Maxima 5.27?
> Yes.
>
> --Barton
>
>
>
>
>
Yes, that procedure works to give a solution to the more general problem 
in Maxima 5.24. It's also worth noting that Aleksas' procedure also 
works for finding a solution to the general problem simply with solve() 
as well.

--

(%i1) load(grobner);

Loading maxima-grobner $Revision: 1.6 $ $Date: 2009/06/02 07:49:49 $
(%o1) /usr/share/maxima/5.24.0/share/contrib/Grobner/grobner.lisp
(%i2) myalgsys(e,v) := (
   e : map(lambda([s], rhs(s)-lhs(s)),e),
   algsys(poly_reduced_grobner(e, block([listconstvars : false], 
listofvars(e))),v));
(%o2) myalgsys(e, v) := (e : map(lambda([s], rhs(s) - lhs(s)), e),
algsys(poly_reduced_grobner(e, block([listconstvars : false], 
listofvars(e))),
v))
(%i3) e1: (x-a)^2 + (y-b)^2 = r1^2;
                                   2          2     2
(%o3)                      (y - b)  + (x - a)  = r1
(%i4) e2: (x-c)^2 + (y-d)^2 = r2^2;
                                   2          2     2
(%o4)                      (y - d)  + (x - c)  = r2
(%i5) myalgsys( [e1,e2], [x,y] );
                                 4        2      2 2              2
(%o5) [[x = - ((d - b) sqrt(- r2  + (2 r1  + 2 d  - 4 b d + 2 c  - 4 a c 
+ 2 b
       2    2     4       2              2              2      2 2    4
  + 2 a ) r2  - r1  + (2 d  - 4 b d + 2 c  - 4 a c + 2 b  + 2 a ) r1 - d
         3         2              2      2   2
  + 4 b d  + (- 2 c  + 4 a c - 6 b  - 2 a ) d
          2                3      2         4        3         2 2   2
  + (4 b c  - 8 a b c + 4 b  + 4 a  b) d - c  + 4 a c  + (- 2 b  - 6 a ) c
          2      3       4      2  2    4              2 2
  + (4 a b  + 4 a ) c - b  - 2 a  b  - a ) + (c - a) r2  + (a - c) r1
               2                        3      2     2    2 2    3
  + (- c - a) d  + (2 b c + 2 a b) d - c  + a c  + (a  - b ) c - a b - a )
      2              2              2      2
/(2 d  - 4 b d + 2 c  - 4 a c + 2 b  + 2 a ),
                       2               2    2            2 2    2
y = ((c - a) sqrt(- r2  + 2 r1 r2 - r1  + d  - 2 b d + c  - 2 a c + b  + a )
         2               2    2            2            2 2              2
  sqrt(r2  + 2 r1 r2 + r1  - d  + 2 b d - c  + 2 a c - b  - a ) + (b - d) r2
              2    3      2     2            2    2 2              3
  + (d - b) r1  + d  - b d  + (c  - 2 a c - b  + a ) d + b c  - 2 a b c + b
     2        2              2              2      2
  + a  b)/(2 d  - 4 b d + 2 c  - 4 a c + 2 b  + 2 a )],
                        4        2      2              2 2      2
[x = ((d - b) sqrt(- r2  + (2 r1  + 2 d  - 4 b d + 2 c  - 4 a c + 2 b  + 
2 a )
    2     4       2              2              2      2    2 4        3
  r2  - r1  + (2 d  - 4 b d + 2 c  - 4 a c + 2 b  + 2 a ) r1  - d  + 4 b d
          2              2      2   2         2                3 2
  + (- 2 c  + 4 a c - 6 b  - 2 a ) d  + (4 b c  - 8 a b c + 4 b  + 4 a  b) d
     4        3         2      2   2         2      3       4      2 2    4
  - c  + 4 a c  + (- 2 b  - 6 a ) c  + (4 a b  + 4 a ) c - b  - 2 a b  - a )
              2             2            2 3      2
  + (a - c) r2  + (c - a) r1  + (c + a) d  + (- 2 b c - 2 a b) d + c - a c
      2    2         2    3      2              2              2 2
  + (b  - a ) c + a b  + a )/(2 d  - 4 b d + 2 c  - 4 a c + 2 b  + 2 a ),
                         2               2    2 2            2    2
y = - ((c - a) sqrt(- r2  + 2 r1 r2 - r1  + d  - 2 b d + c  - 2 a c + b  
+ a )
         2               2    2            2            2 2              2
  sqrt(r2  + 2 r1 r2 + r1  - d  + 2 b d - c  + 2 a c - b  - a ) + (d - b) r2
              2    3      2       2            2    2 2              3
  + (b - d) r1  - d  + b d  + (- c  + 2 a c + b  - a ) d - b c  + 2 a b 
c - b
     2        2              2              2      2
  - a  b)/(2 d  - 4 b d + 2 c  - 4 a c + 2 b  + 2 a )]]
(%i6)
--

Also, without loading grobner and defining myalgsys,

--
(%i1)  e1: (x-a)^2 + (y-b)^2 = r1^2;
                                   2          2     2
(%o1)                      (y - b)  + (x - a)  = r1
(%i2) e2: (x-c)^2 + (y-d)^2 = r2^2;
                                   2          2     2
(%o2)                      (y - d)  + (x - c)  = r2
(%i3) solve( [e1,e2-e1], [x,y] );
                                 4        2      2 2              2
(%o3) [[x = - ((d - b) sqrt(- r2  + (2 r1  + 2 d  - 4 b d + 2 c  - 4 a c 
+ 2 b
       2    2     4       2              2              2      2 2    4
  + 2 a ) r2  - r1  + (2 d  - 4 b d + 2 c  - 4 a c + 2 b  + 2 a ) r1 - d
         3         2              2      2   2
  + 4 b d  + (- 2 c  + 4 a c - 6 b  - 2 a ) d
          2                3      2         4        3         2 2   2
  + (4 b c  - 8 a b c + 4 b  + 4 a  b) d - c  + 4 a c  + (- 2 b  - 6 a ) c
          2      3       4      2  2    4              2 2
  + (4 a b  + 4 a ) c - b  - 2 a  b  - a ) + (c - a) r2  + (a - c) r1
               2                        3      2     2    2 2    3
  + (- c - a) d  + (2 b c + 2 a b) d - c  + a c  + (a  - b ) c - a b - a )
      2              2              2      2
/(2 d  - 4 b d + 2 c  - 4 a c + 2 b  + 2 a ),
                       2               2    2            2 2    2
y = ((c - a) sqrt(- r2  + 2 r1 r2 - r1  + d  - 2 b d + c  - 2 a c + b  + a )
         2               2    2            2            2 2              2
  sqrt(r2  + 2 r1 r2 + r1  - d  + 2 b d - c  + 2 a c - b  - a ) + (b - d) r2
              2    3      2     2            2    2 2              3
  + (d - b) r1  + d  - b d  + (c  - 2 a c - b  + a ) d + b c  - 2 a b c + b
     2        2              2              2      2
  + a  b)/(2 d  - 4 b d + 2 c  - 4 a c + 2 b  + 2 a )],
                        4        2      2              2 2      2
[x = ((d - b) sqrt(- r2  + (2 r1  + 2 d  - 4 b d + 2 c  - 4 a c + 2 b  + 
2 a )
    2     4       2              2              2      2    2 4        3
  r2  - r1  + (2 d  - 4 b d + 2 c  - 4 a c + 2 b  + 2 a ) r1  - d  + 4 b d
          2              2      2   2         2                3 2
  + (- 2 c  + 4 a c - 6 b  - 2 a ) d  + (4 b c  - 8 a b c + 4 b  + 4 a  b) d
     4        3         2      2   2         2      3       4      2 2    4
  - c  + 4 a c  + (- 2 b  - 6 a ) c  + (4 a b  + 4 a ) c - b  - 2 a b  - a )
              2             2            2 3      2
  + (a - c) r2  + (c - a) r1  + (c + a) d  + (- 2 b c - 2 a b) d + c - a c
      2    2         2    3      2              2              2 2
  + (b  - a ) c + a b  + a )/(2 d  - 4 b d + 2 c  - 4 a c + 2 b  + 2 a ),
                         2               2    2 2            2    2
y = - ((c - a) sqrt(- r2  + 2 r1 r2 - r1  + d  - 2 b d + c  - 2 a c + b  
+ a )
         2               2    2            2            2 2              2
  sqrt(r2  + 2 r1 r2 + r1  - d  + 2 b d - c  + 2 a c - b  - a ) + (d - b) r2
              2    3      2       2            2    2 2              3
  + (b - d) r1  - d  + b d  + (- c  + 2 a c + b  - a ) d - b c  + 2 a b 
c - b
     2        2              2              2      2
  - a  b)/(2 d  - 4 b d + 2 c  - 4 a c + 2 b  + 2 a )]]
(%i4)
--

Krishna