A challenging 3 equations and 3 unknowns
pnachtwey
circle using three ultrasonic sensors. The thread is here:
http://www.plctalk.net/qanda/showthread.php?p=316629#post316629
My wxMaxima doesn't seem to be getting the job done.
eq1: (x1-x0)^2+(y1-y0)^2-(r0+r1)^2;
eq2: (x2-x0)^2+(y2-y0)^2-(r0+r2)^2;
eq3: (x3-x0)^2+(y3-y0)^2-(r0+r3)^2;
eq4: solve([eq1,eq2,eq3],[x0,y0,r0]);
grind(%);
where the center of the inner circle is x0,y0 and the radius is r0.
x1,y1 and x2,y2 and x3,y3 are the locations of the ultrasonic sensors.
r1, r2, r3 are the distances between the sensors and the inner circle.
It looks easy but my wxMaxima has been working on the solution for
about a half hour and so far there is no solution.
Perhaps there is a better way to approach this but I can't see how.
Thanks
Peter Nachtwey
Roman Pearce
> The is post on this forum about finding the location and size of a
> circle using three ultrasonic sensors. The thread is here:http://www.plctalk.net/qanda/showthread.php?p=316629#post316629
>
> My wxMaxima doesn't seem to be getting the job done.
>
> eq1: (x1-x0)^2+(y1-y0)^2-(r0+r1)^2;
> eq2: (x2-x0)^2+(y2-y0)^2-(r0+r2)^2;
> eq3: (x3-x0)^2+(y3-y0)^2-(r0+r3)^2;
> eq4: solve([eq1,eq2,eq3],[x0,y0,r0]);
Maple 12 can quickly get a triangular form (F4+FGLM) but the system
doesn't split. I obtain the following Groebner basis in lex order
with y0 > x0 > r0:
sys := {(x1-x0)^2+(y1-y0)^2-(r0+r1)^2,(x2-x0)^2+(y2-y0)^2-(r0+r2)^2,
(x3-x0)^2+(y3-y0)^2-(r0+r3)^2}:
vars := {x0,y0,r0};
G := Groebner[Basis](sys,'tord',order=plex,variables=vars):
[-6*x3^2*x1^2*x2^2+2*y1^2*x1^2*y2^2-2*x1*x2^3*y3^2-2*x1^3*y3^2*x2+4*x3^2*x1^2*r2^2-2*r1^2*x1^2*y2^2-2*r1^2*y3^2*x1^2+2*y1^2*y3^2*x1^2-2*y3^2*y2^2*r3^2-2*x1^3*y2^2*x2-2*y2^4*x1*x3+2*y3^2*y2^2*x3^2+2*y2^2*y3^2*x2^
2-2*x3^2*y2^2*r3^2-2*y1^2*x2^3*x3-2*r1^4*y2*y3-2*y3^2*x2^2*r3^2-2*x1^3*y2^2*x3+
(-8*y1*r3^2*y2+8*x3*y1^2*x2-8*y1^2*r3*r2-8*r3*x1^2*r2+8*x3^2*y2*y1-8*y1*y3*r2^2+8*y1*x2^2*y3-8*y2*x1*y3*x2-8*y3*r1^2*y2+8*x1^2*y2*y3
+8*x1*y2^2*x3-8*x1*r2^2*x3-8*x2*r1^2*x3+8*y3^2*x1*x2-8*r3^2*x1*x2-4*y3^2*x2^2-4*x1^2*y2^2-4*y3^2*x1^2-4*x3^2*y2^2+4*x2^2*r1^2+4*r1^2*y2^2+4*r3^2*y2^2+4*r3^2*x1^2+4*r3^2*x2^2+4*r2^2*x1^2-8*x3*y2*y1*x2-8*x3*y2*y3*
x1+8*x3*y2*y3*x2-8*x3*y1*x1*y2+8*x3*y1*y3*x1-8*x3*y1*y3*x2+8*x1*y2*y1*x2-8*y1*x2*y3*x1-4*x3^2*y1^2-4*y1^2*x2^2+4*y1^2*r2^2+8*r1*x1*x3*r2+4*r1^2*x3^2+4*y3^2*r1^2+4*y3^2*r2^2+8*r2*x1*x3*r3+8*r3*x1*x2*r2-8*r2*x2*x3
*r3+4*x3^2*r2^2-8*r3*y2*y3*r2+8*y1*r2*y2*r3+8*y1*r3*y3*r2+4*y1^2*r3^2-8*r1*r2*x3^2-8*r1*y2^2*r3-8*y3^2*r1*r2-8*r1*x2^2*r3-8*y1*r1*y2*r2+8*r2*y2*y3*r1+8*r3*y2*y3*r1+8*r1*x2*x3*r3+8*r1*r2*x3*x2-8*r1*x1*x2*r2+8*r1*
x1*x2*r3-8*r1*x1*x3*r3+8*y1*y3*r1*r2-8*y1*y3*r1*r3+8*y1*r1*y2*r3)
*r0^2-2*x3^2*x2^2*r2^2-2*y1^4*x3*x2+2*x3^2*x2^2*y3^2-2*x1^4*y2*y3-2*y1^3*x2^2*y3+2*y2^2*y1^2*x2^2+2*x3^2*x1^2*y3^2+2*r2^2*x3^2*y2*y1-2*x3^2*x2^2*
r3^2-2*r1^2*y2^2*r3^2+2*y2^2*y3^2*x1^2-2*y3^2*x2^2*r2^2+2*x1^2*x3^2*y2^2+2*x3^2*x2^2*y2^2+2*y1^2*x3^2*y2^2+2*x1^2*y2^2*x2^2-2*y1^2*x2^2*r2^2+2*y1^2*x2^2*x1^2-2*x3^2*y2^3*y3-2*y1*y2^3*x3^2+2*x1^2*y3^2*x2^2-2*y1^3
*x3^2*y2+4*y1*x1^2*x3*y3*x2+2*y1^2*y3^2*x2^2-2*y1*y2^3*x1^2-2*y2^3*x1^2*y3-2*y1^3*y2*x2^2-2*x3^2*x1^2*r3^2+y1^4*x3^2+y1^4*x2^2+y2^4*x1^2+y2^4*x3^2-4*y1^2*x1^2*y2*y3-2*y1^2*x1*y2^2*x3-2*y1^2*y3^2*x1*x2+4*y1^3*x3*
y2*x2+4*y1^3*x3*y3*x2+2*y1^2*y2*x2^2*y3+2*y1*y2^2*x1^2*y3+4*y1*y2^3*x1*x3+2*y1*y2*y3^2*x2^2+2*y1*y2*y3^2*x1^2-2*y1^2*y2^2*x3*x2-2*y1^2*y2^2*x1*x2-4*y2^2*y1*x2^2*y3-2*y2^2*y3^2*x1*x2+4*y2^3*x3*y3*x1+4*y1^2*y2*x1*
y3*x2+4*y1^2*x3*y2*y3*x1-8*y1^2*x3*y2*y3*x2+4*y1*y2^2*x1*y3*x2-8*y1*y2*y3^2*x1*x2-8*y1*y2^2*x3*y3*x1+4*y1*y2^2*x3*y3*x2+2*y1^2*y3*y2^3+4*y1*x2*x3*y3^2*y2-2*r1^2*x3^2*y1^2+y3^2*x2^4-2*x2^3*x3^3-2*x1^3*x2^3+y3^2*
y2^4+r1^4*y2^2+x3^2*x2^4+y3^4*y2^2+x1^4*y3^2+y3^4*x2^2+r3^4*x1^2-2*x1^3*x3^3+r2^4*x1^2+x3^4*x1^2+x3^4*y2^2+y3^4*x1^2+r3^4*y2^2+4*y3*r2^2*y2^2*y1-2*r1^2*y1^2*x2^2-2*r2^2*x1^2*y2^2-2*r2^2*x3^2*y2^2-2*r3^2*y2^3*y1+
4*r3^2*y2^2*y1^2+2*r1^2*y2^3*y1-2*y3*r1^2*y2^3+4*y3^2*r1^2*y2^2+2*y1^3*r2^2*y2-2*r3^2*y2*y1^3-2*y3^2*r1^2*y1^2-2*r1^2*y2^2*y1^2-2*y3*r2^2*y1^3+4*y3^2*r2^2*y1^2-2*y3^2*r2^2*y2^2-2*y1^2*r2^2*y2^2+4*r1^2*x3*y1^2*x2
+2*r1^2*x3^2*y2*y1-4*r1^2*x3*y2*y1*x2+4*r2^2*x1*y2^2*x3-4*r2^2*x3*y1*x1*y2+4*y3*r1^2*y1^2*y2+2*r1^2*x2^2*y1*y2+2*r2^2*x1^2*y2*y1-2*r3^2*y1*x2^2*y2-2*r3^2*y1*x1^2*y2-2*r1^2*y2^2*y1*y3-2*y1^2*r2^2*y2*y3-2*r2^2*x1*
x3*y1^2-2*r1^2*x2*x3*y2^2-2*r3^2*y2^2*x1*x2+8*r3^2*y2*x1*x2*y1-2*r3^2*y1^2*x1*x2-2*y3^2*r1^2*y1*y2-2*y3^2*r2^2*y1*y2+x1^4*y2^2+y3^2*r2^4+x3^4*x2^2+x2^4*x1^2+r3^4*x2^2+y3^2*r1^4-2*y3^3*y2^3+2*x1^2*y2^2*x2*x3-2*x1
^2*y2*x2^2*y3+4*x1^3*x3*y2*y3+4*x1^3*y2*y3*x2+2*r1^2*x1*y2^2*x2+4*r1^2*x1^2*y2*y3+2*r1^2*x1*y2^2*x3-2*r1^2*y2*x2^2*y3+2*r1^2*y3^2*x1*x2-4*r1^2*x3*y2*y3*x1+8*r1^2*x3*y2*y3*x2-8*x1^2*x3*y2*y3*x2-4*r1^2*y2*x1*y3*x2
-2*y3^2*x1^2*r3^2-2*y3^3*y2*x1^2-2*y3^3*y2*r1^2-2*y3^3*y2*x2^2+2*y3^3*y2*r2^2+2*x1^3*x2*x3^2-2*x1^3*x2*r3^2+2*x1^3*x2*r2^2-2*x1^3*x3*y3^2+2*x1^3*x3*r3^2+2*x1^3*x3*x2^2-2*x1^3*x3*r2^2+4*x2^2*x1^2*r3^2-2*x2^2*x1^2
*r2^2-2*x2^3*x3*y3^2+2*x2^3*x3*r3^2+4*x2^2*r1^2*x3^2-2*x2^2*r1^2*r3^2-2*x1*y2^2*x3^3-2*x3^2*y2*y3*x1^2-2*x3^2*y2*y3*r1^2-2*x3^2*y2*y3*x2^2+2*x3^2*y2*y3*r2^2+2*r1^2*y2*y3*r2^2+2*x1^2*r2^2*y2*y3+2*x2^2*y3*r3^2*y2+
2*x1*y2^2*x2*x3^2+2*x1^2*y2*y3*r3^2-2*x1*y2^2*x3*y3^2+2*x1*y2^2*x3*r3^2-4*x1*y2^2*x3*x2^2+4*x1*y3^3*x2*y2+2*x1*y3^2*x2*r2^2+2*x1*r1^2*x3*y3^2-2*x1*r1^2*x3*r3^2-2*x1*r1^2*x3*x2^2+2*x1*r1^2*x3*r2^2-2*x1*r1^2*x2*x3
^2+2*x1*r1^2*x2*r3^2-2*x1*r1^2*x2*r2^2+4*x1*y3*x2*x3^2*y2-4*x1*y3*x2*r3^2*y2+2*x1*x2^2*x3*y3^2-2*x1*x2^2*x3*r3^2+4*x1*x2^2*x3*r2^2-2*x1*r2^2*x3*y3^2+2*x1*r2^2*x3*r3^2-2*x2*r1^2*x3*y3^2+2*x2*r1^2*x3*r3^2+2*x2*r1^
2*x3*r2^2-2*x2*x3*y2^2*y3^2+2*x2*x3*y2^2*r3^2+2*x2*x1^2*x3*y3^2-2*x2*x1^2*x3*r3^2-2*x2*x1^2*x3*r2^2+2*x2*r2^2*x3*y3^2-2*x2*r2^2*x3*r3^2-4*x3^2*x1*x2*y3^2+4*x3^2*x1*x2*r3^2-2*x3^2*x1*x2*r2^2+4*y3^2*x1*x2*r3^2+2*
r3^2*x1*x2*r2^2+2*y3*r1^2*y2*r3^2-2*y3*r2^2*y2*r3^2+4*r1^2*x1^2*x3*x2+4*x1*y2*x3*y3*x2^2-4*x1*y2*x3*y3*r2^2+2*r3^2*y2^3*y3+2*x1*r1^2*x3^3+2*x1*r1^2*x2^3+2*x1*x2^2*x3^3-2*x1*x2^4*x3-2*x1*r2^2*x3^3-2*x1*r2^4*x3-2*
x2*r1^2*x3^3-2*x2^3*r1^2*x3-2*x2*x3^3*y2^2+2*x2*x1^2*x3^3+2*x2^3*x1^2*x3+2*x2*r2^2*x3^3-2*x3^4*x1*x2+2*x3^2*x1*x2^3-2*y3^4*x1*x2-2*r3^4*x1*x2-2*r3^2*x1*x2^3-2*r3^2*x1^2*r2^2-2*y3^2*r1^2*r2^2-2*x1^4*x3*x2-2*r1^4*
x3*x2-2*r1^2*x1^2*x3^2-2*r1^2*x1^2*x2^2-2*r1^2*r2^2*x3^2+x1^4*x3^2+x1^4*x2^2+r1^4*x3^2+r1^4*x2^2+r2^4*x3^2-2*y1^4*y3*y2-2*y1^3*x3^2*y3+2*y1^3*r3^2*y3+2*y1^3*y3*y2^2+2*y1^3*y3^2*y2-6*y1^2*y3^2*y2^2+2*y1^2*x3^2*x2
^2+2*y1^2*x3^2*x1^2+y1^4*y3^2+y1^4*y2^2-2*y1^3*y3^3-2*y1^3*y2^3+y1^2*x2^4+y1^2*y2^4+y1^2*r2^4+y1^2*x3^4+y1^2*y3^4+y1^2*r3^4+2*y1^2*x2*x3*r3^2-2*y1^2*x2*x3*y3^2+2*y1^2*x1*x2*r2^2+2*y1^2*x3^2*y2*y3+2*y1^2*x1*x2^2*
x3-4*y1^2*x2*x1^2*x3+2*y1^2*x2*r2^2*x3+2*y1^2*x3^2*x1*x2-2*y1^2*x1*x3*y3^2+2*y1^2*x1*x3*r3^2-2*y1^2*y3*r3^2*y2-2*y1*x1^2*x2^2*y3-4*y1*x3^2*y2*y3^2-2*y1*r1^2*y2*r2^2+4*y1*x3^2*y2*r3^2-2*y1*x1^2*x2^2*y2+2*y1*x2^2*
r1^2*y3+4*y1*x1*x2^3*y3+2*y1*r1^2*y2*r3^2-2*y1*x3^2*x1^2*y3+4*y1*x2^3*x3*y3+4*y1*x2^2*r2^2*y3-2*y1*x3^2*y2*x1^2-2*y1*x3^2*y2*x2^2+4*y1*x2*x3^3*y2-2*y1*x1^2*r2^2*y3-2*y1*x3^2*x2^2*y3+4*y1*y3^2*y2*r3^2+2*y1*x3^2*
y2^2*y3+2*y1*x2^2*y3*r3^2+4*y1*x1*x3^3*y2+2*y1*x3^2*y3*r1^2-2*y1*x3^2*y3*r2^2+2*y1*r3^2*y3*x1^2-2*y1*r3^2*y3*r1^2-2*y1*r3^2*y3*y2^2+2*y1*r3^2*y3*r2^2+2*y1*r2^2*r3^2*y2+2*y1*r2^2*y3*r1^2+4*y1*x1*y3^3*x2+2*y1^2*x3
^2*y3^2-2*y1^2*x3^2*r3^2-2*y1^2*x1*x3^3-2*y1^2*x1*x2^3+2*y1^2*y3^3*y2-2*y1^2*x2*x3^3-2*y1^2*y3^2*r3^2-2*y1^2*r3^2*r2^2-2*y1*x2^4*y3+2*y1*y3^2*y2^3-2*y1*y3*y2^4-2*y1*y3^3*x1^2-2*y1*x3^4*y2-2*y1*y3^4*y2-2*y1*y3^3*
x2^2+2*y1*y3^3*y2^2+2*y1*y3^3*r1^2-2*y1*y3^3*r2^2-2*y1*r3^4*y2-2*y1*r2^4*y3+4*y1*x1*x2^2*x3*y2-4*y1*x1*r2^2*y3*x2+4*y1*x2*x1^2*x3*y2-8*y1*x3^2*y2*x1*x2-4*y1*r1^2*x3*y3*x2-4*y1*x2*x3*r3^2*y2-4*y1*x2*x3*y3*r2^2+4*
y1*x1*x3*y3^2*y2-4*y1*x1*x3*r3^2*y2-8*y1*x1*x3*y3*x2^2+8*y1*x1*x3*y3*r2^2+4*y1*x1*y3*x2*x3^2-4*y1*x1*y3*x2*r3^2+(4*x1*r1^2*x2*r3-4*x1*r1^2*x2*r2-4*x1*r1^2*x3*r3+4*y1^2*r3^3+4*y1^2*r2^3+4*r2^3*x3^2+4*r1^3*x2^2+4*
r1^3*x3^2+4*x1*r1^2*x3*r2+4*y3^2*r2^3+4*x1*y3^2*x2*r2+4*y3^2*r1^3+4*r2^3*x1^2+4*r3^3*x2^2+4*r3^3*x1^2+4*r1^2*y2*y3*r2+4*r3^3*y2^2+4*x3^2*y2*y3*r2-4*x3^2*y2^2*r3+4*r1^3*y2^2-4*x3^2*y2*y3*r1+4*r1*y3^2*x1*x2+4*r1*
x1*y2^2*x3+8*r1*x1^2*y2*y3+4*r1*y2^2*x1*x2-4*r1*y2^2*x3*x2-4*r1*y2*x2^2*y3-4*y1*r2*x1^2*y3+4*y1*r1^2*y2*r3+4*y1^2*x2*x3*r3-4*y1*r1^2*y2*r2+4*y1*r2*x1^2*y2-4*y3^2*y2^2*r3+4*x1*y2^2*x3*r3+4*x1^2*y2*y3*r3+4*x1^2*y2
*y3*r2-4*x1*y2^2*x2*r3+4*x2^2*y3*r3*y2-4*y3^2*x2^2*r3+4*r3^2*y2*y3*r1-4*y3*r2^2*y2*r3+4*y3*r1^2*y2*r3+4*r3^2*x1*x2*r2+8*y3^2*x1*x2*r3-4*y3^2*x2^2*r2+8*x3^2*x1*x2*r3-4*x3^2*x1*x2*r2+8*x1*y2^2*x3*r2+8*y1^2*r1*x3*
x2-4*y1*y3^2*y2*r2-4*r2*x2*x3*r3^2+4*r2*x2*x3*y3^2+4*r3*x1*x2*r2^2+4*r1*x1*x3*r2^2-4*r1*x1*x3*x2^2-4*r1*x1*x3*r3^2-4*x2^2*x1^2*r2+4*r1*x1*x3*y3^2-4*r1*x1*x2*r2^2+4*r1*x1*x2*r3^2-4*r1*x1*x2*x3^2+4*r1*x2*x3*r3^2-4
*r1*x2*x3*y3^2+4*r2*x1*x3*r3^2-4*r2*x1*x3*y3^2-4*x2*r2^2*x3*r3-4*x2*x1^2*x3*r3-4*x3^2*y2^2*r2+4*x1^3*x3*r3-4*x2*x1^2*x3*r2+4*x2*x3*y2^2*r3+4*x2*r1^2*x3*r3+4*x2*r1^2*x3*r2-4*x1^3*x3*r2+4*x1*r2^2*x3*r3-4*x1*x2^2*
x3*r3+8*x1*x2^2*x3*r2+4*r2^2*r1*x3*x2+8*x1^2*r1*x3*x2+4*r2^2*y2*y3*r1-4*x1^3*x2*r3+4*x1^3*x2*r2+4*y3^3*y2*r2-4*y3^3*y2*r1-4*y3^2*x1^2*r3-4*x1^2*y2^2*r2-4*x3^2*x1^2*r3+8*x3^2*x1^2*r2-8*r1^3*y2*y3-4*r3^2*y2*y3*r2-
\
4*r1^2*y2^2*r3-4*x1^2*r1*x2^2+8*y1*y2^2*y3*r2+4*y1*x2^2*y3*r3-4*x1^2*r1*x3^2-4*y3^2*r2^2*r1+8*y1*x2^2*y3*r2+4*y3*y2^3*r3-4*y3^2*r1^2*r2+8*y1*y3^2*y2*r3-4*y1*y2*x2^2*r3+8*y1*x3^2*y2*r3-4*r2^2*x1^2*r3+4*y1*x3^2*y2
*r2+4*y1*r1*y2*x2^2+4*y1*r1*x3^2*y2-4*r3^2*x1^2*r2-4*r3^2*y2^2*r1-8*r3^3*x1*x2+4*y1*r1*x2^2*y3-8*x1*r2^3*x3-4*y2^3*y3*r1-4*y1*y2^2*y3*r1-4*y1^2*y3*r3*y2-4*x2^2*r1^2*r3+4*x2^3*x3*r3-4*y1^2*y2*y3*r2+8*y1^2*y2*y3*
r1+4*y1^2*x1*x3*r3-4*y1^2*x1*x3*r2+8*x2^2*x1^2*r3+4*y1^2*x2*x3*r2-4*y1^2*x1*x2*r3+4*y1^2*x1*x2*r2-4*y1*y3^2*r1*y2-4*r1*x2^2*r3^2+8*x2^2*r1*x3^2-4*y1*x3^2*y3*r2+4*y1*x3^2*y3*r1-4*x2^3*r1*x3-8*r1^3*x3*x2-4*r1^2*r2
*x3^2-4*y1^2*y2^2*r1+4*y1*y3*r2^2*r3-4*y1*y3*r1^2*r3+4*y1*y3*r1^2*r2-4*y1*y3*y2^2*r3+4*y1*r3^2*y2*r2+4*y1*r2^2*y3*r1+4*y1*r2^2*r3*y2+4*y1^3*y2*r2-4*y1^3*r3*y2+4*y1^3*y3*r3-4*y1^3*y3*r2-4*r3*x1*x2^3-4*y1*r2^2*r1*
y2+4*y1*r3^2*y3*r2-4*y1*r3^2*y3*r1+4*r1*x1*x3^3+4*r1*x1*x2^3-4*r1*x2*x3^3+4*y1*r3^2*r1*y2-4*r2*x1*x3^3-4*y1^2*y3^2*r3+8*y1^2*y3^2*r2-4*y1^2*x3^2*r3+8*y1^2*y2^2*r3-4*y1^2*x2^2*r2-4*y1^2*x2^2*r1-4*y1^2*y2^2*r2-4*
y1^2*r1*x3^2+4*y1*r3*x1^2*y3-4*y1*r3*x1^2*y2-4*y1*y3^3*r2+4*y1*y3^3*r1-4*y1*y2^3*r3+4*y1*y2^3*r1-4*y1^2*y3^2*r1-4*y1^2*r2^2*r3-4*y1^2*r3^2*r2-8*y1*r2^3*y3-8*y1*r3^3*y2-8*r1*y2*x1*y3*x2-8*y1*x2*x3*r3*y2-8*y1*y2*
x1*x3*r2+16*y1*x1*x3*y3*r2-8*y1*x1*x3*r3*y2-8*y1*x1*y3*x2*r3-8*y1*x1*x2*y3*r2-8*y1*x2*x3*y3*r2+16*y1*x1*x2*r3*y2-8*y1*r1*x3*y3*x2-8*x1*y2*x3*y3*r2-8*x1*y3*x2*r3*y2+16*r1*x3*y2*y3*x2-8*r1*x3*y2*y3*x1-8*y1*r1*x3*
y2*x2-4*y3^2*y2^2*r2-4*x3^2*x2^2*r3+8*y3^2*y2^2*r1-4*r2^2*r1*x3^2-4*x3^2*x2^2*r2-4*r1*y3^2*x1^2-4*r1*x1^2*y2^2+4*r2*x2*x3^3)
*r0, -x1^2*y2-y1^2*y2+r1^2*y2-x3^2*y1+x3^2*y2-y3^2*y1+y3^2*y2+r3^2*y1-
r3^2*y2+y1*x2^2+
y1*y2^2-y1*r2^2+y3*x1^2+y3*y1^2-y3*r1^2-y3*x2^2-y3*y2^2+y3*r2^2+
(-2*x3*y2+2*x3*y1+2*x1*y2-2*y1*x2-2*y3*x1+2*y3*x2)
*x0+(2*r1*y2+2*r3*y1-2*r3*y2-2*y1*r2-2*y3*r1+2*y3*r2)*r0, x3*y2^2-
x1*y2^2+y1^2*x2-x3*y1^2-x1*x2^2
+x1*r2^2-x2*r1^2+x2*x1^2-x1*r3^2+r1^2*x3+x2^2*x3-r2^2*x3-
x3^2*x2+r3^2*x2+x1*x3^2-x1^2*x3+x1*y3^2-y3^2*x2+
(-2*x3*y2+2*x3*y1+2*x1*y2-2*y1*x2-2*y3*x1+2*y3*x2)*y0+
(-2*x1*r3+2*x2*r3-2*x3*r2+2*x1*r2+2*r1*x3-2*r1*x2)*
r0]
Daniel Lichtblau
You probably do not want to do this in a purely symbolic manner (can
certainly be done but you don't want to look at the result, not even
as wallpaper).
First observe that you can make the problem linear in two variables
just by subtracting one equation from the other two (say). So solve
for the two variables in terms of the third. Now you can substitute
into any one of the polynomials to get a quadratic in that third. It
is not at all difficult to then solve symbolically for that third
variable, but the solution, even after some massaging, will still be
fairly unwieldy.
An alternative that keeps size quite well contained is to do each part
(linear and quadratic) after first substituting specific values for
the positional and distance parameters. I show a version of this in
Mathematica below.
(* Here are the parameter values I use. *)
xrep = Thread[{x[1],x[2],x[3]} -> {1,3,4}];
yrep = Thread[{y[1],y[2],y[3]} -> {2,1,3}];
rrep = Thread[{r[1],r[2],r[3]} -> {.2,.1,.3}];
reps = Join[xrep,yrep,rrep];
(*Here are the three polynomials. *)
polys = Table[(x[j]-x[0])^2 + (y[j]-y[0])^2 - (r[0]+r[j])^2, {j,3}];
(*Subtract last from first two. *)
newpolys = Most[polys]-Last[polys];
(* Substitute parameter values and solve for center coordinates. *)
soln1 = First[Solve[(newpolys/.reps)==0, {x[0],y[0]}]];
(* Find a positive solution for the radius. I simply use a local root
finder,
and give a positive initial value in the (naive) hope that it will
converge to
a positive solution. No check is made for this. *)
solnr = FindRoot[Numerator[Together[polys[[1]] /. soln1] /. reps]==0,
{r[0],1}];
result = {x[0],y[0],r[0]} /. soln1 /.solnr
Out[46]= {2.498, 2.34723, 1.33772}
Of course there is the (perhaps awkward?) possibility that there are
two viable solutions (viable meaning having a positive radius). For
example, if we use for our radial displacements
rrep2 = Thread[{r[1],r[2],r[3]} -> {1.2,.1,2.3}];
then we have two such solutions. If this is an issue in rpactice, then
it might be best just to numerically solve the original system,
skipping all the intermediate steps. In Mathematica this could be done
with NSolve as below.
reps2 = Join[xrep,yrep,rrep2];
In[52]:= NSolve[polys/.reps2]
Out[52]= {{r[0] -> 3.15225, x[0] -> 2.258, y[0] -> -2.16648},
{r[0] -> 0.95346, x[0] -> 2.258, y[0] -> 0.252194}}
Daniel Lichtblau
Wolfram Research
pnachtwey
Thanks guys, I did try eliminating using the eliminate function but
even that was slow and messy. I agree that this application requires
an iterative approach and that is practically impossible to do on a
PLC. I was hoping that the solution would be only a little more
difficult than finding the center of a circle given three points but I
can see it isn't.
Peter Nachtwey
Robert H. Lewis
> and size of a
> circle using three ultrasonic sensors. The thread is
> here:
> http://www.plctalk.net/qanda/showthread.php?p=316629#p
> ost316629
>
> My wxMaxima doesn't seem to be getting the job done.
>
> eq1: (x1-x0)^2+(y1-y0)^2-(r0+r1)^2;
> eq2: (x2-x0)^2+(y2-y0)^2-(r0+r2)^2;
> eq3: (x3-x0)^2+(y3-y0)^2-(r0+r3)^2;
> eq4: solve([eq1,eq2,eq3],[x0,y0,r0]);
>
This sort of thing is very easy to do symbolically with resultants. I only have a minute now, but you might look at
http://fordham.academia.edu/RobertLewis/Papers/82784/Apollonius-Problems-in-Biochemistry
Robert H. Lewis
Fordham University
Daniel Lichtblau
It's not hard to do the symbolic elimination to obtain, say, the
radius. I show a way to do it below (and Robert Lewis already posted
to the effect that one can use e.g. Dixon resultants for the same
purpose). But the result is quite large and possibly not so useful.
On the bright side, it is quadratic in r[0] (we already knew this
because we know there are two solutions). So it should be efficient to
solve for the two candidate radii, once given actual values for the
parameters. It is just a question of whether this symbolic approach
will beat purely numeric methods that begin with numeric equations and
solve each time parameters change.
In[8]:= Timing[rpoly = First[GroebnerBasis[polys, r[0], {x[0],y[0]},
CoefficientDomain->RationalFunctions,
MonomialOrder->EliminationOrder]];]
Out[8]= {0.192971, Null}
In[9]:= InputForm[rpoly]
Out[9]//InputForm=
r[2]^4*x[1]^2 - 2*r[2]^2*r[3]^2*x[1]^2 + r[3]^4*x[1]^2 -
2*r[1]^2*r[2]^2*x[1]*x[2] + 2*r[1]^2*r[3]^2*x[1]*x[2] +
2*r[2]^2*r[3]^2*x[1]*x[2] - 2*r[3]^4*x[1]*x[2] + 2*r[2]^2*x[1]^3*x[2]
-
2*r[3]^2*x[1]^3*x[2] + r[1]^4*x[2]^2 - 2*r[1]^2*r[3]^2*x[2]^2 +
r[3]^4*x[2]^2 - 2*r[1]^2*x[1]^2*x[2]^2 - 2*r[2]^2*x[1]^2*x[2]^2 +
4*r[3]^2*x[1]^2*x[2]^2 + x[1]^4*x[2]^2 + 2*r[1]^2*x[1]*x[2]^3 -
2*r[3]^2*x[1]*x[2]^3 - 2*x[1]^3*x[2]^3 + x[1]^2*x[2]^4 +
2*r[1]^2*r[2]^2*x[1]*x[3] - 2*r[2]^4*x[1]*x[3] - 2*r[1]^2*r[3]^2*x[1]
*x[3] +
2*r[2]^2*r[3]^2*x[1]*x[3] - 2*r[2]^2*x[1]^3*x[3] + 2*r[3]^2*x[1]^3*x
[3] -
2*r[1]^4*x[2]*x[3] + 2*r[1]^2*r[2]^2*x[2]*x[3] + 2*r[1]^2*r[3]^2*x[2]
*x[3] -
2*r[2]^2*r[3]^2*x[2]*x[3] + 4*r[1]^2*x[1]^2*x[2]*x[3] -
2*r[2]^2*x[1]^2*x[2]*x[3] - 2*r[3]^2*x[1]^2*x[2]*x[3] - 2*x[1]^4*x[2]
*x[3] -
2*r[1]^2*x[1]*x[2]^2*x[3] + 4*r[2]^2*x[1]*x[2]^2*x[3] -
2*r[3]^2*x[1]*x[2]^2*x[3] + 2*x[1]^3*x[2]^2*x[3] - 2*r[1]^2*x[2]^3*x
[3] +
2*r[3]^2*x[2]^3*x[3] + 2*x[1]^2*x[2]^3*x[3] - 2*x[1]*x[2]^4*x[3] +
r[1]^4*x[3]^2 - 2*r[1]^2*r[2]^2*x[3]^2 + r[2]^4*x[3]^2 -
2*r[1]^2*x[1]^2*x[3]^2 + 4*r[2]^2*x[1]^2*x[3]^2 - 2*r[3]^2*x[1]^2*x[3]
^2 +
x[1]^4*x[3]^2 - 2*r[1]^2*x[1]*x[2]*x[3]^2 - 2*r[2]^2*x[1]*x[2]*x[3]^2
+
4*r[3]^2*x[1]*x[2]*x[3]^2 + 2*x[1]^3*x[2]*x[3]^2 + 4*r[1]^2*x[2]^2*x
[3]^2 -
2*r[2]^2*x[2]^2*x[3]^2 - 2*r[3]^2*x[2]^2*x[3]^2 - 6*x[1]^2*x[2]^2*x[3]
^2 +
2*x[1]*x[2]^3*x[3]^2 + x[2]^4*x[3]^2 + 2*r[1]^2*x[1]*x[3]^3 -
2*r[2]^2*x[1]*x[3]^3 - 2*x[1]^3*x[3]^3 - 2*r[1]^2*x[2]*x[3]^3 +
2*r[2]^2*x[2]*x[3]^3 + 2*x[1]^2*x[2]*x[3]^3 + 2*x[1]*x[2]^2*x[3]^3 -
2*x[2]^3*x[3]^3 + x[1]^2*x[3]^4 - 2*x[1]*x[2]*x[3]^4 + x[2]^2*x[3]^4
+
r[2]^4*y[1]^2 - 2*r[2]^2*r[3]^2*y[1]^2 + r[3]^4*y[1]^2 +
2*r[2]^2*x[1]*x[2]*y[1]^2 - 2*r[3]^2*x[1]*x[2]*y[1]^2 -
2*r[1]^2*x[2]^2*y[1]^2 - 2*r[2]^2*x[2]^2*y[1]^2 + 2*x[1]^2*x[2]^2*y[1]
^2 -
2*x[1]*x[2]^3*y[1]^2 + x[2]^4*y[1]^2 - 2*r[2]^2*x[1]*x[3]*y[1]^2 +
2*r[3]^2*x[1]*x[3]*y[1]^2 + 4*r[1]^2*x[2]*x[3]*y[1]^2 +
2*r[2]^2*x[2]*x[3]*y[1]^2 + 2*r[3]^2*x[2]*x[3]*y[1]^2 -
4*x[1]^2*x[2]*x[3]*y[1]^2 + 2*x[1]*x[2]^2*x[3]*y[1]^2 -
2*x[2]^3*x[3]*y[1]^2 - 2*r[1]^2*x[3]^2*y[1]^2 - 2*r[3]^2*x[3]^2*y[1]
^2 +
2*x[1]^2*x[3]^2*y[1]^2 + 2*x[1]*x[2]*x[3]^2*y[1]^2 +
2*x[2]^2*x[3]^2*y[1]^2 - 2*x[1]*x[3]^3*y[1]^2 - 2*x[2]*x[3]^3*y[1]^2
+
x[3]^4*y[1]^2 + x[2]^2*y[1]^4 - 2*x[2]*x[3]*y[1]^4 + x[3]^2*y[1]^4 -
2*r[1]^2*r[2]^2*y[1]*y[2] + 2*r[1]^2*r[3]^2*y[1]*y[2] +
2*r[2]^2*r[3]^2*y[1]*y[2] - 2*r[3]^4*y[1]*y[2] + 2*r[2]^2*x[1]^2*y[1]
*y[2] -
2*r[3]^2*x[1]^2*y[1]*y[2] + 8*r[3]^2*x[1]*x[2]*y[1]*y[2] +
2*r[1]^2*x[2]^2*y[1]*y[2] - 2*r[3]^2*x[2]^2*y[1]*y[2] -
2*x[1]^2*x[2]^2*y[1]*y[2] - 4*r[2]^2*x[1]*x[3]*y[1]*y[2] -
4*r[3]^2*x[1]*x[3]*y[1]*y[2] - 4*r[1]^2*x[2]*x[3]*y[1]*y[2] -
4*r[3]^2*x[2]*x[3]*y[1]*y[2] + 4*x[1]^2*x[2]*x[3]*y[1]*y[2] +
4*x[1]*x[2]^2*x[3]*y[1]*y[2] + 2*r[1]^2*x[3]^2*y[1]*y[2] +
2*r[2]^2*x[3]^2*y[1]*y[2] + 4*r[3]^2*x[3]^2*y[1]*y[2] -
2*x[1]^2*x[3]^2*y[1]*y[2] - 8*x[1]*x[2]*x[3]^2*y[1]*y[2] -
2*x[2]^2*x[3]^2*y[1]*y[2] + 4*x[1]*x[3]^3*y[1]*y[2] +
4*x[2]*x[3]^3*y[1]*y[2] - 2*x[3]^4*y[1]*y[2] + 2*r[2]^2*y[1]^3*y[2]
-
2*r[3]^2*y[1]^3*y[2] - 2*x[2]^2*y[1]^3*y[2] + 4*x[2]*x[3]*y[1]^3*y[2]
-
2*x[3]^2*y[1]^3*y[2] + r[1]^4*y[2]^2 - 2*r[1]^2*r[3]^2*y[2]^2 +
r[3]^4*y[2]^2 - 2*r[1]^2*x[1]^2*y[2]^2 - 2*r[2]^2*x[1]^2*y[2]^2 +
x[1]^4*y[2]^2 + 2*r[1]^2*x[1]*x[2]*y[2]^2 - 2*r[3]^2*x[1]*x[2]*y[2]^2
-
2*x[1]^3*x[2]*y[2]^2 + 2*x[1]^2*x[2]^2*y[2]^2 + 2*r[1]^2*x[1]*x[3]*y
[2]^2 +
4*r[2]^2*x[1]*x[3]*y[2]^2 + 2*r[3]^2*x[1]*x[3]*y[2]^2 -
2*x[1]^3*x[3]*y[2]^2 - 2*r[1]^2*x[2]*x[3]*y[2]^2 +
2*r[3]^2*x[2]*x[3]*y[2]^2 + 2*x[1]^2*x[2]*x[3]*y[2]^2 -
4*x[1]*x[2]^2*x[3]*y[2]^2 - 2*r[2]^2*x[3]^2*y[2]^2 -
2*r[3]^2*x[3]^2*y[2]^2 + 2*x[1]^2*x[3]^2*y[2]^2 +
2*x[1]*x[2]*x[3]^2*y[2]^2 + 2*x[2]^2*x[3]^2*y[2]^2 - 2*x[1]*x[3]^3*y
[2]^2 -
2*x[2]*x[3]^3*y[2]^2 + x[3]^4*y[2]^2 - 2*r[1]^2*y[1]^2*y[2]^2 -
2*r[2]^2*y[1]^2*y[2]^2 + 4*r[3]^2*y[1]^2*y[2]^2 + 2*x[1]^2*y[1]^2*y[2]
^2 -
2*x[1]*x[2]*y[1]^2*y[2]^2 + 2*x[2]^2*y[1]^2*y[2]^2 -
2*x[1]*x[3]*y[1]^2*y[2]^2 - 2*x[2]*x[3]*y[1]^2*y[2]^2 +
2*x[3]^2*y[1]^2*y[2]^2 + y[1]^4*y[2]^2 + 2*r[1]^2*y[1]*y[2]^3 -
2*r[3]^2*y[1]*y[2]^3 - 2*x[1]^2*y[1]*y[2]^3 + 4*x[1]*x[3]*y[1]*y[2]^3
-
2*x[3]^2*y[1]*y[2]^3 - 2*y[1]^3*y[2]^3 + x[1]^2*y[2]^4 -
2*x[1]*x[3]*y[2]^4 + x[3]^2*y[2]^4 + y[1]^2*y[2]^4 +
2*r[1]^2*r[2]^2*y[1]*y[3] - 2*r[2]^4*y[1]*y[3] - 2*r[1]^2*r[3]^2*y[1]
*y[3] +
2*r[2]^2*r[3]^2*y[1]*y[3] - 2*r[2]^2*x[1]^2*y[1]*y[3] +
2*r[3]^2*x[1]^2*y[1]*y[3] - 4*r[2]^2*x[1]*x[2]*y[1]*y[3] -
4*r[3]^2*x[1]*x[2]*y[1]*y[3] + 2*r[1]^2*x[2]^2*y[1]*y[3] +
4*r[2]^2*x[2]^2*y[1]*y[3] + 2*r[3]^2*x[2]^2*y[1]*y[3] -
2*x[1]^2*x[2]^2*y[1]*y[3] + 4*x[1]*x[2]^3*y[1]*y[3] - 2*x[2]^4*y[1]*y
[3] +
8*r[2]^2*x[1]*x[3]*y[1]*y[3] - 4*r[1]^2*x[2]*x[3]*y[1]*y[3] -
4*r[2]^2*x[2]*x[3]*y[1]*y[3] + 4*x[1]^2*x[2]*x[3]*y[1]*y[3] -
8*x[1]*x[2]^2*x[3]*y[1]*y[3] + 4*x[2]^3*x[3]*y[1]*y[3] +
2*r[1]^2*x[3]^2*y[1]*y[3] - 2*r[2]^2*x[3]^2*y[1]*y[3] -
2*x[1]^2*x[3]^2*y[1]*y[3] + 4*x[1]*x[2]*x[3]^2*y[1]*y[3] -
2*x[2]^2*x[3]^2*y[1]*y[3] - 2*r[2]^2*y[1]^3*y[3] + 2*r[3]^2*y[1]^3*y
[3] -
2*x[2]^2*y[1]^3*y[3] + 4*x[2]*x[3]*y[1]^3*y[3] - 2*x[3]^2*y[1]^3*y[3]
-
2*r[1]^4*y[2]*y[3] + 2*r[1]^2*r[2]^2*y[2]*y[3] + 2*r[1]^2*r[3]^2*y[2]
*y[3] -
2*r[2]^2*r[3]^2*y[2]*y[3] + 4*r[1]^2*x[1]^2*y[2]*y[3] +
2*r[2]^2*x[1]^2*y[2]*y[3] + 2*r[3]^2*x[1]^2*y[2]*y[3] - 2*x[1]^4*y[2]
*y[3] -
4*r[1]^2*x[1]*x[2]*y[2]*y[3] - 4*r[3]^2*x[1]*x[2]*y[2]*y[3] +
4*x[1]^3*x[2]*y[2]*y[3] - 2*r[1]^2*x[2]^2*y[2]*y[3] +
2*r[3]^2*x[2]^2*y[2]*y[3] - 2*x[1]^2*x[2]^2*y[2]*y[3] -
4*r[1]^2*x[1]*x[3]*y[2]*y[3] - 4*r[2]^2*x[1]*x[3]*y[2]*y[3] +
4*x[1]^3*x[3]*y[2]*y[3] + 8*r[1]^2*x[2]*x[3]*y[2]*y[3] -
8*x[1]^2*x[2]*x[3]*y[2]*y[3] + 4*x[1]*x[2]^2*x[3]*y[2]*y[3] -
2*r[1]^2*x[3]^2*y[2]*y[3] + 2*r[2]^2*x[3]^2*y[2]*y[3] -
2*x[1]^2*x[3]^2*y[2]*y[3] + 4*x[1]*x[2]*x[3]^2*y[2]*y[3] -
2*x[2]^2*x[3]^2*y[2]*y[3] + 4*r[1]^2*y[1]^2*y[2]*y[3] -
2*r[2]^2*y[1]^2*y[2]*y[3] - 2*r[3]^2*y[1]^2*y[2]*y[3] -
4*x[1]^2*y[1]^2*y[2]*y[3] + 4*x[1]*x[2]*y[1]^2*y[2]*y[3] +
2*x[2]^2*y[1]^2*y[2]*y[3] + 4*x[1]*x[3]*y[1]^2*y[2]*y[3] -
8*x[2]*x[3]*y[1]^2*y[2]*y[3] + 2*x[3]^2*y[1]^2*y[2]*y[3] -
2*y[1]^4*y[2]*y[3] - 2*r[1]^2*y[1]*y[2]^2*y[3] + 4*r[2]^2*y[1]*y[2]
^2*y[3] -
2*r[3]^2*y[1]*y[2]^2*y[3] + 2*x[1]^2*y[1]*y[2]^2*y[3] +
4*x[1]*x[2]*y[1]*y[2]^2*y[3] - 4*x[2]^2*y[1]*y[2]^2*y[3] -
8*x[1]*x[3]*y[1]*y[2]^2*y[3] + 4*x[2]*x[3]*y[1]*y[2]^2*y[3] +
2*x[3]^2*y[1]*y[2]^2*y[3] + 2*y[1]^3*y[2]^2*y[3] - 2*r[1]^2*y[2]^3*y
[3] +
2*r[3]^2*y[2]^3*y[3] - 2*x[1]^2*y[2]^3*y[3] + 4*x[1]*x[3]*y[2]^3*y[3]
-
2*x[3]^2*y[2]^3*y[3] + 2*y[1]^2*y[2]^3*y[3] - 2*y[1]*y[2]^4*y[3] +
r[1]^4*y[3]^2 - 2*r[1]^2*r[2]^2*y[3]^2 + r[2]^4*y[3]^2 -
2*r[1]^2*x[1]^2*y[3]^2 - 2*r[3]^2*x[1]^2*y[3]^2 + x[1]^4*y[3]^2 +
2*r[1]^2*x[1]*x[2]*y[3]^2 + 2*r[2]^2*x[1]*x[2]*y[3]^2 +
4*r[3]^2*x[1]*x[2]*y[3]^2 - 2*x[1]^3*x[2]*y[3]^2 - 2*r[2]^2*x[2]^2*y
[3]^2 -
2*r[3]^2*x[2]^2*y[3]^2 + 2*x[1]^2*x[2]^2*y[3]^2 - 2*x[1]*x[2]^3*y[3]
^2 +
x[2]^4*y[3]^2 + 2*r[1]^2*x[1]*x[3]*y[3]^2 - 2*r[2]^2*x[1]*x[3]*y[3]^2
-
2*x[1]^3*x[3]*y[3]^2 - 2*r[1]^2*x[2]*x[3]*y[3]^2 +
2*r[2]^2*x[2]*x[3]*y[3]^2 + 2*x[1]^2*x[2]*x[3]*y[3]^2 +
2*x[1]*x[2]^2*x[3]*y[3]^2 - 2*x[2]^3*x[3]*y[3]^2 + 2*x[1]^2*x[3]^2*y
[3]^2 -
4*x[1]*x[2]*x[3]^2*y[3]^2 + 2*x[2]^2*x[3]^2*y[3]^2 -
2*r[1]^2*y[1]^2*y[3]^2 + 4*r[2]^2*y[1]^2*y[3]^2 - 2*r[3]^2*y[1]^2*y[3]
^2 +
2*x[1]^2*y[1]^2*y[3]^2 - 2*x[1]*x[2]*y[1]^2*y[3]^2 +
2*x[2]^2*y[1]^2*y[3]^2 - 2*x[1]*x[3]*y[1]^2*y[3]^2 -
2*x[2]*x[3]*y[1]^2*y[3]^2 + 2*x[3]^2*y[1]^2*y[3]^2 + y[1]^4*y[3]^2 -
2*r[1]^2*y[1]*y[2]*y[3]^2 - 2*r[2]^2*y[1]*y[2]*y[3]^2 +
4*r[3]^2*y[1]*y[2]*y[3]^2 + 2*x[1]^2*y[1]*y[2]*y[3]^2 -
8*x[1]*x[2]*y[1]*y[2]*y[3]^2 + 2*x[2]^2*y[1]*y[2]*y[3]^2 +
4*x[1]*x[3]*y[1]*y[2]*y[3]^2 + 4*x[2]*x[3]*y[1]*y[2]*y[3]^2 -
4*x[3]^2*y[1]*y[2]*y[3]^2 + 2*y[1]^3*y[2]*y[3]^2 + 4*r[1]^2*y[2]^2*y
[3]^2 -
2*r[2]^2*y[2]^2*y[3]^2 - 2*r[3]^2*y[2]^2*y[3]^2 + 2*x[1]^2*y[2]^2*y[3]
^2 -
2*x[1]*x[2]*y[2]^2*y[3]^2 + 2*x[2]^2*y[2]^2*y[3]^2 -
2*x[1]*x[3]*y[2]^2*y[3]^2 - 2*x[2]*x[3]*y[2]^2*y[3]^2 +
2*x[3]^2*y[2]^2*y[3]^2 - 6*y[1]^2*y[2]^2*y[3]^2 + 2*y[1]*y[2]^3*y[3]
^2 +
y[2]^4*y[3]^2 + 2*r[1]^2*y[1]*y[3]^3 - 2*r[2]^2*y[1]*y[3]^3 -
2*x[1]^2*y[1]*y[3]^3 + 4*x[1]*x[2]*y[1]*y[3]^3 - 2*x[2]^2*y[1]*y[3]^3
-
2*y[1]^3*y[3]^3 - 2*r[1]^2*y[2]*y[3]^3 + 2*r[2]^2*y[2]*y[3]^3 -
2*x[1]^2*y[2]*y[3]^3 + 4*x[1]*x[2]*y[2]*y[3]^3 - 2*x[2]^2*y[2]*y[3]^3
+
2*y[1]^2*y[2]*y[3]^3 + 2*y[1]*y[2]^2*y[3]^3 - 2*y[2]^3*y[3]^3 +
x[1]^2*y[3]^4 - 2*x[1]*x[2]*y[3]^4 + x[2]^2*y[3]^4 + y[1]^2*y[3]^4 -
2*y[1]*y[2]*y[3]^4 + y[2]^2*y[3]^4 +
r[0]^2*(4*r[2]^2*x[1]^2 - 8*r[2]*r[3]*x[1]^2 + 4*r[3]^2*x[1]^2 -
8*r[1]*r[2]*x[1]*x[2] + 8*r[1]*r[3]*x[1]*x[2] + 8*r[2]*r[3]*x[1]*x
[2] -
8*r[3]^2*x[1]*x[2] + 4*r[1]^2*x[2]^2 - 8*r[1]*r[3]*x[2]^2 +
4*r[3]^2*x[2]^2 + 8*r[1]*r[2]*x[1]*x[3] - 8*r[2]^2*x[1]*x[3] -
8*r[1]*r[3]*x[1]*x[3] + 8*r[2]*r[3]*x[1]*x[3] - 8*r[1]^2*x[2]*x[3]
+
8*r[1]*r[2]*x[2]*x[3] + 8*r[1]*r[3]*x[2]*x[3] - 8*r[2]*r[3]*x[2]*x
[3] +
4*r[1]^2*x[3]^2 - 8*r[1]*r[2]*x[3]^2 + 4*r[2]^2*x[3]^2 + 4*r[2]^2*y
[1]^2 -
8*r[2]*r[3]*y[1]^2 + 4*r[3]^2*y[1]^2 - 4*x[2]^2*y[1]^2 +
8*x[2]*x[3]*y[1]^2 - 4*x[3]^2*y[1]^2 - 8*r[1]*r[2]*y[1]*y[2] +
8*r[1]*r[3]*y[1]*y[2] + 8*r[2]*r[3]*y[1]*y[2] - 8*r[3]^2*y[1]*y[2]
+
8*x[1]*x[2]*y[1]*y[2] - 8*x[1]*x[3]*y[1]*y[2] - 8*x[2]*x[3]*y[1]*y
[2] +
8*x[3]^2*y[1]*y[2] + 4*r[1]^2*y[2]^2 - 8*r[1]*r[3]*y[2]^2 +
4*r[3]^2*y[2]^2 - 4*x[1]^2*y[2]^2 + 8*x[1]*x[3]*y[2]^2 - 4*x[3]^2*y
[2]^2 +
8*r[1]*r[2]*y[1]*y[3] - 8*r[2]^2*y[1]*y[3] - 8*r[1]*r[3]*y[1]*y[3]
+
8*r[2]*r[3]*y[1]*y[3] - 8*x[1]*x[2]*y[1]*y[3] + 8*x[2]^2*y[1]*y[3]
+
8*x[1]*x[3]*y[1]*y[3] - 8*x[2]*x[3]*y[1]*y[3] - 8*r[1]^2*y[2]*y[3]
+
8*r[1]*r[2]*y[2]*y[3] + 8*r[1]*r[3]*y[2]*y[3] - 8*r[2]*r[3]*y[2]*y
[3] +
8*x[1]^2*y[2]*y[3] - 8*x[1]*x[2]*y[2]*y[3] - 8*x[1]*x[3]*y[2]*y[3]
+
8*x[2]*x[3]*y[2]*y[3] + 4*r[1]^2*y[3]^2 - 8*r[1]*r[2]*y[3]^2 +
4*r[2]^2*y[3]^2 - 4*x[1]^2*y[3]^2 + 8*x[1]*x[2]*y[3]^2 -
4*x[2]^2*y[3]^2) + r[0]*(4*r[2]^3*x[1]^2 - 4*r[2]^2*r[3]*x[1]^2 -
4*r[2]*r[3]^2*x[1]^2 + 4*r[3]^3*x[1]^2 - 4*r[1]^2*r[2]*x[1]*x[2] -
4*r[1]*r[2]^2*x[1]*x[2] + 4*r[1]^2*r[3]*x[1]*x[2] +
4*r[2]^2*r[3]*x[1]*x[2] + 4*r[1]*r[3]^2*x[1]*x[2] +
4*r[2]*r[3]^2*x[1]*x[2] - 8*r[3]^3*x[1]*x[2] + 4*r[2]*x[1]^3*x[2]
-
4*r[3]*x[1]^3*x[2] + 4*r[1]^3*x[2]^2 - 4*r[1]^2*r[3]*x[2]^2 -
4*r[1]*r[3]^2*x[2]^2 + 4*r[3]^3*x[2]^2 - 4*r[1]*x[1]^2*x[2]^2 -
4*r[2]*x[1]^2*x[2]^2 + 8*r[3]*x[1]^2*x[2]^2 + 4*r[1]*x[1]*x[2]^3 -
4*r[3]*x[1]*x[2]^3 + 4*r[1]^2*r[2]*x[1]*x[3] + 4*r[1]*r[2]^2*x[1]*x
[3] -
8*r[2]^3*x[1]*x[3] - 4*r[1]^2*r[3]*x[1]*x[3] + 4*r[2]^2*r[3]*x[1]*x
[3] -
4*r[1]*r[3]^2*x[1]*x[3] + 4*r[2]*r[3]^2*x[1]*x[3] - 4*r[2]*x[1]^3*x
[3] +
4*r[3]*x[1]^3*x[3] - 8*r[1]^3*x[2]*x[3] + 4*r[1]^2*r[2]*x[2]*x[3]
+
4*r[1]*r[2]^2*x[2]*x[3] + 4*r[1]^2*r[3]*x[2]*x[3] -
4*r[2]^2*r[3]*x[2]*x[3] + 4*r[1]*r[3]^2*x[2]*x[3] -
4*r[2]*r[3]^2*x[2]*x[3] + 8*r[1]*x[1]^2*x[2]*x[3] -
4*r[2]*x[1]^2*x[2]*x[3] - 4*r[3]*x[1]^2*x[2]*x[3] -
4*r[1]*x[1]*x[2]^2*x[3] + 8*r[2]*x[1]*x[2]^2*x[3] -
4*r[3]*x[1]*x[2]^2*x[3] - 4*r[1]*x[2]^3*x[3] + 4*r[3]*x[2]^3*x[3]
+
4*r[1]^3*x[3]^2 - 4*r[1]^2*r[2]*x[3]^2 - 4*r[1]*r[2]^2*x[3]^2 +
4*r[2]^3*x[3]^2 - 4*r[1]*x[1]^2*x[3]^2 + 8*r[2]*x[1]^2*x[3]^2 -
4*r[3]*x[1]^2*x[3]^2 - 4*r[1]*x[1]*x[2]*x[3]^2 - 4*r[2]*x[1]*x[2]*x
[3]^2 +
8*r[3]*x[1]*x[2]*x[3]^2 + 8*r[1]*x[2]^2*x[3]^2 - 4*r[2]*x[2]^2*x[3]
^2 -
4*r[3]*x[2]^2*x[3]^2 + 4*r[1]*x[1]*x[3]^3 - 4*r[2]*x[1]*x[3]^3 -
4*r[1]*x[2]*x[3]^3 + 4*r[2]*x[2]*x[3]^3 + 4*r[2]^3*y[1]^2 -
4*r[2]^2*r[3]*y[1]^2 - 4*r[2]*r[3]^2*y[1]^2 + 4*r[3]^3*y[1]^2 +
4*r[2]*x[1]*x[2]*y[1]^2 - 4*r[3]*x[1]*x[2]*y[1]^2 - 4*r[1]*x[2]^2*y
[1]^2 -
4*r[2]*x[2]^2*y[1]^2 - 4*r[2]*x[1]*x[3]*y[1]^2 + 4*r[3]*x[1]*x[3]*y
[1]^2 +
8*r[1]*x[2]*x[3]*y[1]^2 + 4*r[2]*x[2]*x[3]*y[1]^2 +
4*r[3]*x[2]*x[3]*y[1]^2 - 4*r[1]*x[3]^2*y[1]^2 - 4*r[3]*x[3]^2*y[1]
^2 -
4*r[1]^2*r[2]*y[1]*y[2] - 4*r[1]*r[2]^2*y[1]*y[2] +
4*r[1]^2*r[3]*y[1]*y[2] + 4*r[2]^2*r[3]*y[1]*y[2] +
4*r[1]*r[3]^2*y[1]*y[2] + 4*r[2]*r[3]^2*y[1]*y[2] - 8*r[3]^3*y[1]*y
[2] +
4*r[2]*x[1]^2*y[1]*y[2] - 4*r[3]*x[1]^2*y[1]*y[2] +
16*r[3]*x[1]*x[2]*y[1]*y[2] + 4*r[1]*x[2]^2*y[1]*y[2] -
4*r[3]*x[2]^2*y[1]*y[2] - 8*r[2]*x[1]*x[3]*y[1]*y[2] -
8*r[3]*x[1]*x[3]*y[1]*y[2] - 8*r[1]*x[2]*x[3]*y[1]*y[2] -
8*r[3]*x[2]*x[3]*y[1]*y[2] + 4*r[1]*x[3]^2*y[1]*y[2] +
4*r[2]*x[3]^2*y[1]*y[2] + 8*r[3]*x[3]^2*y[1]*y[2] + 4*r[2]*y[1]^3*y
[2] -
4*r[3]*y[1]^3*y[2] + 4*r[1]^3*y[2]^2 - 4*r[1]^2*r[3]*y[2]^2 -
4*r[1]*r[3]^2*y[2]^2 + 4*r[3]^3*y[2]^2 - 4*r[1]*x[1]^2*y[2]^2 -
4*r[2]*x[1]^2*y[2]^2 + 4*r[1]*x[1]*x[2]*y[2]^2 - 4*r[3]*x[1]*x[2]*y
[2]^2 +
4*r[1]*x[1]*x[3]*y[2]^2 + 8*r[2]*x[1]*x[3]*y[2]^2 +
4*r[3]*x[1]*x[3]*y[2]^2 - 4*r[1]*x[2]*x[3]*y[2]^2 +
4*r[3]*x[2]*x[3]*y[2]^2 - 4*r[2]*x[3]^2*y[2]^2 - 4*r[3]*x[3]^2*y[2]
^2 -
4*r[1]*y[1]^2*y[2]^2 - 4*r[2]*y[1]^2*y[2]^2 + 8*r[3]*y[1]^2*y[2]^2
+
4*r[1]*y[1]*y[2]^3 - 4*r[3]*y[1]*y[2]^3 + 4*r[1]^2*r[2]*y[1]*y[3]
+
4*r[1]*r[2]^2*y[1]*y[3] - 8*r[2]^3*y[1]*y[3] - 4*r[1]^2*r[3]*y[1]*y
[3] +
4*r[2]^2*r[3]*y[1]*y[3] - 4*r[1]*r[3]^2*y[1]*y[3] +
4*r[2]*r[3]^2*y[1]*y[3] - 4*r[2]*x[1]^2*y[1]*y[3] +
4*r[3]*x[1]^2*y[1]*y[3] - 8*r[2]*x[1]*x[2]*y[1]*y[3] -
8*r[3]*x[1]*x[2]*y[1]*y[3] + 4*r[1]*x[2]^2*y[1]*y[3] +
8*r[2]*x[2]^2*y[1]*y[3] + 4*r[3]*x[2]^2*y[1]*y[3] +
16*r[2]*x[1]*x[3]*y[1]*y[3] - 8*r[1]*x[2]*x[3]*y[1]*y[3] -
8*r[2]*x[2]*x[3]*y[1]*y[3] + 4*r[1]*x[3]^2*y[1]*y[3] -
4*r[2]*x[3]^2*y[1]*y[3] - 4*r[2]*y[1]^3*y[3] + 4*r[3]*y[1]^3*y[3]
-
8*r[1]^3*y[2]*y[3] + 4*r[1]^2*r[2]*y[2]*y[3] + 4*r[1]*r[2]^2*y[2]*y
[3] +
4*r[1]^2*r[3]*y[2]*y[3] - 4*r[2]^2*r[3]*y[2]*y[3] +
4*r[1]*r[3]^2*y[2]*y[3] - 4*r[2]*r[3]^2*y[2]*y[3] +
8*r[1]*x[1]^2*y[2]*y[3] + 4*r[2]*x[1]^2*y[2]*y[3] +
4*r[3]*x[1]^2*y[2]*y[3] - 8*r[1]*x[1]*x[2]*y[2]*y[3] -
8*r[3]*x[1]*x[2]*y[2]*y[3] - 4*r[1]*x[2]^2*y[2]*y[3] +
4*r[3]*x[2]^2*y[2]*y[3] - 8*r[1]*x[1]*x[3]*y[2]*y[3] -
8*r[2]*x[1]*x[3]*y[2]*y[3] + 16*r[1]*x[2]*x[3]*y[2]*y[3] -
4*r[1]*x[3]^2*y[2]*y[3] + 4*r[2]*x[3]^2*y[2]*y[3] +
8*r[1]*y[1]^2*y[2]*y[3] - 4*r[2]*y[1]^2*y[2]*y[3] -
4*r[3]*y[1]^2*y[2]*y[3] - 4*r[1]*y[1]*y[2]^2*y[3] +
8*r[2]*y[1]*y[2]^2*y[3] - 4*r[3]*y[1]*y[2]^2*y[3] - 4*r[1]*y[2]^3*y
[3] +
4*r[3]*y[2]^3*y[3] + 4*r[1]^3*y[3]^2 - 4*r[1]^2*r[2]*y[3]^2 -
4*r[1]*r[2]^2*y[3]^2 + 4*r[2]^3*y[3]^2 - 4*r[1]*x[1]^2*y[3]^2 -
4*r[3]*x[1]^2*y[3]^2 + 4*r[1]*x[1]*x[2]*y[3]^2 + 4*r[2]*x[1]*x[2]*y
[3]^2 +
8*r[3]*x[1]*x[2]*y[3]^2 - 4*r[2]*x[2]^2*y[3]^2 - 4*r[3]*x[2]^2*y[3]
^2 +
4*r[1]*x[1]*x[3]*y[3]^2 - 4*r[2]*x[1]*x[3]*y[3]^2 -
4*r[1]*x[2]*x[3]*y[3]^2 + 4*r[2]*x[2]*x[3]*y[3]^2 - 4*r[1]*y[1]^2*y
[3]^2 +
8*r[2]*y[1]^2*y[3]^2 - 4*r[3]*y[1]^2*y[3]^2 - 4*r[1]*y[1]*y[2]*y[3]
^2 -
4*r[2]*y[1]*y[2]*y[3]^2 + 8*r[3]*y[1]*y[2]*y[3]^2 + 8*r[1]*y[2]^2*y
[3]^2 -
4*r[2]*y[2]^2*y[3]^2 - 4*r[3]*y[2]^2*y[3]^2 + 4*r[1]*y[1]*y[3]^3 -
4*r[2]*y[1]*y[3]^3 - 4*r[1]*y[2]*y[3]^3 + 4*r[2]*y[2]*y[3]^3)
Daniel Lichtblau
Wolfram Research
Axel Vogt
What do you get for x0=y0 = 0? May be I got trapped, but I was
thinking of quadratic normal forms (though I do not have code
for simultaneous situations):
Each is equivalent to x[1]^2+x[2]^2-x[3]^2 over the Reals, so one
may look at x[1]^2+x[2]^2-(x[3]-r[i])^2, r[1]=0 and the other are
free.
May be I got it wrong ... but what for that special case?
Robert H. Lewis
I computed the 3 resultants. Each takes about 0.019 seconds. This probably duplicates Dan Lichtblau's result, but here is the resultant for x0 in nested form:
((4y2^2 + (-8y1)y2 + 4x2^2 + (-8x1)x2 + 4y1^2 + 4x1^2)r3^2 + (((-8y2 + 8y1)y3 +
(-8x2 + 8x1)x3 + (8y1)y2 + (8x1)x2 - 8y1^2 - 8x1^2)r2 + ((8y2 - 8y1)y3 + (8x2
- 8x1)x3 - 8y2^2 + (8y1)y2 - 8x2^2 + (8x1)x2)r1)r3 + (4y3^2 + (-8y1)y3 + 4x3^2 +
(-8x1)x3 + 4y1^2 + 4x1^2)r2^2 + ((-8y3^2 + (8y2 + 8y1)y3 - 8x3^2 + (8x2 + 8x1)
x3 + (-8y1)y2 + (-8x1)x2)r1)r2 + (4y3^2 + (-8y2)y3 + 4x3^2 + (-8x2)x3 + 4y2^2
+ 4x2^2)r1^2 + (-4x2^2 + (8x1)x2 - 4x1^2)y3^2 + (((8x2 - 8x1)y2 + (-8y1)x2 + (8x1)
y1)x3 + ((-8x1)x2 + 8x1^2)y2 + (8y1)x2^2 + ((-8x1)y1)x2)y3 + (-4y2^2 + (8y1)y2
- 4y1^2)x3^2 + ((8x1)y2^2 + ((-8y1)x2 + (-8x1)y1)y2 + (8y1^2)x2)x3 + (-4x1^2)y2^2
+ (((8x1)y1)x2)y2 + (-4y1^2)x2^2)x0^2 + (((-4x2 + 4x1)r2 + (4x2 - 4x1)r1)r3^3 +
((4x3 + 4x2 - 8x1)r2^2 + ((-8x3 + 4x2 + 4x1)r1)r2 + (4x3 - 8x2 + 4x1)r1^2 +
((4x2 - 4x1)y2 + (-4y1)x2 + (4x1)y1)y3 + (-4y2^2 + (8y1)y2 - 4y1^2)x3 + (-4x2)
y2^2 + ((4y1)x2 + (4x1)y1)y2 - 4x2^3 + (4x1)x2^2 + (4x1^2)x2 + (-4x1)y1^2 - 4x1^3)
r3^2 + ((-4x3 + 4x1)r2^3 + ((4x3 - 8x2 + 4x1)r1)r2^2 + ((4x3 + 4x2 - 8x1)r1^2 +
(4x2 - 4x1)y3^2 + ((16x1)y2 + (-8y1)x2 + (-8x1)y1)y3 + (4x2 - 4x1)x3^2 + (4y2^2
+ (-8y1)y2 + 4x2^2 + 4y1^2 - 4x1^2)x3 + (-4x1)y2^2 + ((-8x1)y1)y2 + (-4x1)x2^2 +
(4y1^2 - 4x1^2)x2 + (8x1)y1^2 + 8x1^3)r2 + (-4x3 + 4x2)r1^3 + ((-4x2 + 4x1)y3^2 +
((-8x2 - 8x1)y2 + (16y1)x2)y3 + (-4x2 + 4x1)x3^2 + (4y2^2 + (-8y1)y2 - 4x2^2 +
4y1^2 + 4x1^2)x3 + (8x2 + 4x1)y2^2 + ((-8y1)x2)y2 + 8x2^3 + (-4x1)x2^2 + (-4y1^2
- 4x1^2)x2)r1)r3 + ((4x3 - 4x1)r1)r2^3 + ((-8x3 + 4x2 + 4x1)r1^2 + (-4x3 - 4x2)
y3^2 + ((4y2 + 4y1)x3 + (-4x1)y2 + (8y1)x2 + (4x1)y1)y3 - 4x3^3 + (4x1)x3^2 +
((-4y1)y2 + 4x1^2)x3 + ((4x1)y1)y2 + (-4y1^2)x2 + (-4x1)y1^2 - 4x1^3)r2^2 + ((4x3
- 4x2)r1^3 + ((8x3 + 4x2 + 4x1)y3^2 + ((-8y2 - 8y1)x3 + (-8x1)y2 + (-8y1)x2)
y3 + 8x3^3 + (-4x2 - 4x1)x3^2 + (-4y2^2 + (16y1)y2 - 4x2^2 - 4y1^2 - 4x1^2)x3 + (4x1)
y2^2 + (4x1)x2^2 + (4y1^2 + 4x1^2)x2)r1)r2 + ((-4x3 - 4x1)y3^2 + ((4y2 + 4y1)x3
+ (4x2 + 8x1)y2 + (-4y1)x2)y3 - 4x3^3 + (4x2)x3^2 + ((-4y1)y2 + 4x2^2)x3 + (-4x2
- 4x1)y2^2 + ((4y1)x2)y2 - 4x2^3)r1^2 + ((-4x2 + 4x1)y2 + (4y1)x2 + (-4x1)y1)y3^3
+ ((4y2^2 + (-8y1)y2 + 4y1^2)x3 + (4x2 - 8x1)y2^2 + ((4y1)x2 + (4x1)y1)y2 +
4x2^3 + (-4x1)x2^2 + (-8y1^2 - 4x1^2)x2 + (4x1)y1^2 + 4x1^3)y3^2 + (((-4x2 + 4x1)
y2 + (4y1)x2 + (-4x1)y1)x3^2 + (-4y2^3 + (4y1)y2^2 + (-4x2^2 + 4y1^2 + 4x1^2)y2
+ (4y1)x2^2 - 4y1^3 + (-4x1^2)y1)x3 + (4x1)y2^3 + ((-8y1)x2 + (4x1)y1)y2^2 + ((4x1)
x2^2 + (4y1^2 + 4x1^2)x2 + (-8x1)y1^2 - 8x1^3)y2 + (-8y1)x2^3 + ((4x1)y1)x2^2 + (4y1^3
+ (4x1^2)y1)x2)y3 + (4y2^2 + (-8y1)y2 + 4y1^2)x3^3 + ((-4x1)y2^2 + ((4y1)x2 +
(4x1)y1)y2 + (-4y1^2)x2)x3^2 + ((4y1)y2^3 + (-8y1^2 - 4x1^2)y2^2 + ((4y1)x2^2 + 4y1^3
+ (4x1^2)y1)y2 + (-4y1^2)x2^2)x3 + ((-4x1)y1)y2^3 + ((4y1^2)x2 + (4x1)y1^2 + 4x1^3)
y2^2 + (((-4x1)y1)x2^2 + (-4y1^3 + (-4x1^2)y1)x2)y2 + (4y1^2)x2^3)x0 + (r2^2 + (-2r1)
r2 + r1^2 - y2^2 + (2y1)y2 - y1^2)r3^4 + (-2r2^3 + (2r1)r2^2 + (2r1^2 + 2y2^2 + (-4y1)
y2 + 2x2^2 + 2y1^2 - 2x1^2)r2 - 2r1^3 + (2y2^2 + (-4y1)y2 - 2x2^2 + 2y1^2 + 2x1^2)
r1)r3^3 + (r2^4 + (2r1)r2^3 + (-6r1^2 - 2y3^2 + (2y2 + 2y1)y3 - 2x3^2 - 2y2^2 + (2y1)
y2 - 2x2^2 - 2y1^2 + 4x1^2)r2^2 + (2r1^3 + (4y3^2 + (-4y2 - 4y1)y3 + 4x3^2 - 2y2^2 +
(8y1)y2 - 2x2^2 - 2y1^2 - 2x1^2)r1)r2 + r1^4 + (-2y3^2 + (2y2 + 2y1)y3 - 2x3^2 - 2y2^2
+ (2y1)y2 + 4x2^2 - 2y1^2 - 2x1^2)r1^2 + (2y2^2 + (-4y1)y2 + 2y1^2)y3^2 + (-2y2^3 +
(2y1)y2^2 + (-2x2^2 + 2y1^2 + 2x1^2)y2 + (2y1)x2^2 - 2y1^3 + (-2x1^2)y1)y3 + (2y2^2
+ (-4y1)y2 + 2y1^2)x3^2 + y2^4 + (-2y1)y2^3 + (2x2^2 + 2y1^2)y2^2 + ((-2y1)x2^2 - 2y1^3
+ (-2x1^2)y1)y2 + x2^4 + (-2x1^2)x2^2 + y1^4 + (2x1^2)y1^2 + x1^4)r3^2 + ((-2r1)r2^4
+ (2r1^2 + 2y3^2 + (-4y1)y3 + 2x3^2 + 2y1^2 - 2x1^2)r2^3 + (2r1^3 + (-2y3^2 + (-4y2
+ 8y1)y3 - 2x3^2 + 4y2^2 + (-4y1)y2 + 4x2^2 - 2y1^2 - 2x1^2)r1)r2^2 + (-2r1^4 + (-2y3^2
+ (8y2 - 4y1)y3 - 2x3^2 - 2y2^2 + (-4y1)y2 - 2x2^2 + 4y1^2 + 4x1^2)r1^2 + (-2y2^2 +
(4y1)y2 - 2x2^2 - 2y1^2 + 2x1^2)y3^2 + ((4y1)y2^2 + (-8y1^2 - 8x1^2)y2 + (4y1)x2^2 +
4y1^3 + (4x1^2)y1)y3 + (-2y2^2 + (4y1)y2 - 2x2^2 - 2y1^2 + 2x1^2)x3^2 + (-2y1^2 +
2x1^2)y2^2 + (4y1^3 + (4x1^2)y1)y2 + (-2y1^2 + 2x1^2)x2^2 - 2y1^4 + (-4x1^2)y1^2 - 2x1^4)
r2 + (2y3^2 + (-4y2)y3 + 2x3^2 + 2y2^2 - 2x2^2)r1^3 + ((-2y2^2 + (4y1)y2 + 2x2^2
- 2y1^2 - 2x1^2)y3^2 + (4y2^3 + (-8y1)y2^2 + (4x2^2 + 4y1^2 + 4x1^2)y2 + (-8y1)x2^2)
y3 + (-2y2^2 + (4y1)y2 + 2x2^2 - 2y1^2 - 2x1^2)x3^2 - 2y2^4 + (4y1)y2^3 + (-4x2^2 - 2y1^2
- 2x1^2)y2^2 + ((4y1)x2^2)y2 - 2x2^4 + (2y1^2 + 2x1^2)x2^2)r1)r3 + (r1^2 - y3^2 + (2y1)
y3 - y1^2)r2^4 + (-2r1^3 + (2y3^2 + (-4y1)y3 - 2x3^2 + 2y1^2 + 2x1^2)r1)r2^3 + (r1^4
+ (-2y3^2 + (2y2 + 2y1)y3 + 4x3^2 - 2y2^2 + (2y1)y2 - 2x2^2 - 2y1^2 - 2x1^2)r1^2 +
y3^4 + (-2y2 - 2y1)y3^3 + (2x3^2 + 2y2^2 + (2y1)y2 + 2x2^2 + 2y1^2)y3^2 + ((-2y2
- 2y1)x3^2 + (-4y1)y2^2 + (2y1^2 + 2x1^2)y2 + (-4y1)x2^2 - 2y1^3 + (-2x1^2)y1)y3 +
x3^4 + ((2y1)y2 - 2x1^2)x3^2 + (2y1^2)y2^2 + (-2y1^3 + (-2x1^2)y1)y2 + (2y1^2)x2^2 +
y1^4 + (2x1^2)y1^2 + x1^4)r2^2 + ((2y3^2 + (-4y2)y3 - 2x3^2 + 2y2^2 + 2x2^2)r1^3 +
(-2y3^4 + (4y2 + 4y1)y3^3 + (-4x3^2 - 2y2^2 + (-8y1)y2 - 2x2^2 - 2y1^2 - 2x1^2)y3^2 +
((4y2 + 4y1)x3^2 + (4y1)y2^2 + (4y1^2 + 4x1^2)y2 + (4y1)x2^2)y3 - 2x3^4 + (2y2^2 +
(-8y1)y2 + 2x2^2 + 2y1^2 + 2x1^2)x3^2 + (-2y1^2 - 2x1^2)y2^2 + (-2y1^2 - 2x1^2)x2^2)r1)
r2 + (-y3^2 + (2y2)y3 - y2^2)r1^4 + (y3^4 + (-2y2 - 2y1)y3^3 + (2x3^2 + 2y2^2 + (2y1)
y2 + 2y1^2 + 2x1^2)y3^2 + ((-2y2 - 2y1)x3^2 - 2y2^3 + (2y1)y2^2 + (-2x2^2 - 4y1^2 - 4x1^2)
y2 + (2y1)x2^2)y3 + x3^4 + ((2y1)y2 - 2x2^2)x3^2 + y2^4 + (-2y1)y2^3 + (2x2^2 +
2y1^2 + 2x1^2)y2^2 + ((-2y1)x2^2)y2 + x2^4)r1^2 + (-y2^2 + (2y1)y2 - y1^2)y3^4 + (2y2^3
+ (-2y1)y2^2 + (2x2^2 - 2y1^2 - 2x1^2)y2 + (-2y1)x2^2 + 2y1^3 + (2x1^2)y1)y3^3 + (
(-2y2^2 + (4y1)y2 - 2y1^2)x3^2 - y2^4 + (-2y1)y2^3 + (-2x2^2 + 6y1^2 + 4x1^2)y2^2 +
((-2y1)x2^2 - 2y1^3 + (-2x1^2)y1)y2 - x2^4 + (4y1^2 + 2x1^2)x2^2 - y1^4 + (-2x1^2)y1^2
- x1^4)y3^2 + ((2y2^3 + (-2y1)y2^2 + (2x2^2 - 2y1^2 - 2x1^2)y2 + (-2y1)x2^2 + 2y1^3 +
(2x1^2)y1)x3^2 + (2y1)y2^4 + (-2y1^2 - 2x1^2)y2^3 + ((4y1)x2^2 - 2y1^3 + (-2x1^2)y1)y2^2
+ ((-2y1^2 - 2x1^2)x2^2 + 2y1^4 + (4x1^2)y1^2 + 2x1^4)y2 + (2y1)x2^4 + (-2y1^3 + (-2x1^2)
y1)x2^2)y3 + (-y2^2 + (2y1)y2 - y1^2)x3^4 + ((-2y1)y2^3 + (4y1^2 + 2x1^2)y2^2 + ((-2y1)
x2^2 - 2y1^3 + (-2x1^2)y1)y2 + (2y1^2)x2^2)x3^2 + (-y1^2)y2^4 + (2y1^3 + (2x1^2)y1)y2^3
+ ((-2y1^2)x2^2 - y1^4 + (-2x1^2)y1^2 - x1^4)y2^2 + ((2y1^3 + (2x1^2)y1)x2^2)y2 + (-y1^2)
x2^4
It's just quadratic, so this ought to be quite useful.
Don't you need a 3D version of this? With four spheres maybe? That would be a bit more challenging, but still easy.
Robert H. Lewis
Fordham University
pnachtwey
> > > The is post on this forum about finding the location
> > > and size of a
> > > circle using three ultrasonic sensors. The thread is
> > > here:
>
> >http://www.plctalk.net/qanda/showthread.php?p=316629#post316629
>
> > > My wxMaxima doesn't seem to be getting the job
> > done.
>
> > > eq1: (x1-x0)^2+(y1-y0)^2-(r0+r1)^2;
> > > eq2: (x2-x0)^2+(y2-y0)^2-(r0+r2)^2;
> > > eq3: (x3-x0)^2+(y3-y0)^2-(r0+r3)^2;
> > > eq4: solve([eq1,eq2,eq3],[x0,y0,r0]);
>
> > This sort of thing is very easy to do symbolically
> > with resultants. I only have a minute now, but you
> > might look at
>
wxWidgets or Java Swing or similar.
If you read the first post I posted a link to the problem on www.plcs.net.
This for for someone else and I doubt they need a 3d version. I have
been e-mail the complete solution. At least I think it is. It is so
long it is hard to make heads or tails out of it. In any case the OP
on the plc forum said this must run on a PLC so the symbolic solution
is way too complex. I agree with Daniel Lichtblau that a practical
solution will require iteration.
Thanks for everyone's time and effort but now we have proven that a
symbolic or exact solution is not practical on the intended hardware.
Peter Nachtwey
pnachtwey
I was assuming the x0=y0=0 would be the inner circle being exactly in
the center. In this case r1=r2=r3. The three ultrasonic sensors are
spaced 120 degrees apart if you see the link in the original post.
Peter Nachtwey
Robert H. Lewis
> location
> > > > and size of a
> > > > circle using three ultrasonic sensors. The thread is
> > > > here:
> >
> >
> >http://www.plctalk.net/qanda/showthread.php?p=316629#post316629
> >
> > > > My wxMaxima doesn't seem to be getting the job
> > > done.
> >
> > > > eq1: (x1-x0)^2+(y1-y0)^2-(r0+r1)^2;
> > > > eq2: (x2-x0)^2+(y2-y0)^2-(r0+r2)^2;
> > > > eq3: (x3-x0)^2+(y3-y0)^2-(r0+r3)^2;
> > > > eq4: solve([eq1,eq2,eq3],[x0,y0,r0]);
> >
> > > This sort of thing is very easy to do
> symbolically
> > > with resultants. I only have a minute now, but
> you
> > > might look at
> >
> >
> >http://fordham.academia.edu/RobertLewis/Papers/82784/
> Apollonius-Probl...
> >
> > > Robert H. Lewis
> > > Fordham University
> >
> > I computed the 3 resultants. Each takes about
> 0.019 seconds. This probably duplicates Dan
> Lichtblau's result, but here is the resultant for x0
> > + ((-2y1^2)x2^2 - y1^4 + (-2x1^2)y1^2 - x1^4)y2^2 +
> ((2y1^3 + (2x1^2)y1)x2^2)y2 + (-y1^2)
> > x2^4
> >
> > It's just quadratic, so this ought to be quite
> useful.
> Someone should write a universal graphical interface
> to Fermat using
> wxWidgets or Java Swing or similar.
Please do!
> If you read the first post I posted a link to the
> problem on www.plcs.net. This for for someone else and I doubt they need a 3d
> version. I have been e-mail the complete solution. At least I think
> it is. It is so long it is hard to make heads or tails out of it. In
> any case the OP on the plc forum said this must run on a PLC so the
> symbolic solution is way too complex. I agree with Daniel Lichtblau
> that a practical solution will require iteration.
>
> Peter Nachtwey
I just read the post you quote from the other forum. What is a PLC? I gather this is some kind
of very small processor?
If, as the original post indicates, the three points are on an equilateral triangle, that could probably simplify the formulas.
pnachtwey
> > version. I have been e-mail the complete solution. At least I think
> > it is. It is so long it is hard to make heads or tails out of it. In
> > any case the OP on the plc forum said this must run on a PLC so the
> > symbolic solution is way too complex. I agree with Daniel Lichtblau
> > that a practical solution will require iteration.
>
> > Peter Nachtwey
>
> I just read the post you quote from the other forum. What is a PLC? I gather this is some kind
> of very small processor?
>
is used in industrial control systems. Most programming is done in
ladder logic. Here is a text from a respected professor in controls.
Note this has little to do with symbolic processing. I recommend that
you just skim it to get the idea.
http://claymore.engineer.gvsu.edu/~jackh/books/plcs/pdf/plcbook5_1.pdf
http://sites.google.com/site/automatedmanufacturingsystems/
PLCs are powerful but not as powerful as PC because PLCs must be able
to run without a fan in very hot conditions so processing power is
limited by how much heat that can be dissipated. Low power micro
controllers like ARM and Atom processors are used. They have firmware
that simulates programming using relays and while they can do math.
These are the kinds of computers that electricians and industrial
engineers use to keep the factories working.
I monitor the sites because my company makes motion controllers that
are used in industrial control. My questions usually deal with motion
profile generation.
> If, as the original post indicates, the three points are on an equilateral triangle, that could probably simplify the formulas.
I am not sure how. The equations seem pretty simple as they are.
However, that is one of the reasons I posted the question here because
I know there are better mathematicians here than I am that would tell
me if the three equations and three unknowns could be simplified. I
did try using the 'assume()' function to specify that all the
distances are positive. I saw this problem like you do. The inner
circle fits between three other circles with radii of the distances
from the ultra sonic transducers. Your comment about the 3D solution
finding the location of a sphere in a tetrahedron is interesting but I
don't think we need to go there.
Back in the 80s we would solve problems like this iteratively. We
didn't know about LFBGS or Levenberg-Marquardt back then and symbolic
solutions were a dream. The log would move end wise between the three
sensors and a circular cross section would be approximated for every
foot of travel. Computers, HP1000s or PDP-11, would when find the best
way to cut the log.
Peter Nachtwey
Robert H. Lewis
> <rle...@fordham.edu> wrote:
> > > > > > The is post on this forum about finding the
> > > location
> > > > > > and size of a
> > > > > > circle using three ultrasonic sensors. The thread is
> > > > > > here:
> > >
> >http://www.plctalk.net/qanda/showthread.php?p=316629#post316629
> >
> > > > > > My wxMaxima doesn't seem to be getting the job
> > > > > done.
> >
> > > > > > eq1: (x1-x0)^2+(y1-y0)^2-(r0+r1)^2;
> > > > > > eq2: (x2-x0)^2+(y2-y0)^2-(r0+r2)^2;
> > > > > > eq3: (x3-x0)^2+(y3-y0)^2-(r0+r3)^2;
> > > > > > eq4: solve([eq1,eq2,eq3],[x0,y0,r0]);
> >
> >
....
> > I just read the post you quote from the other forum. What is a PLC? I gather this is some kind
> > of very small processor?
> >
> A PLC, programmable logic controller, is an
> industrial computer that
> is used in industrial control systems. Most
> programming is done in ladder logic. Here is a text from a respected
> professor in controls. Note this has little to do with symbolic processing.
> I recommend that you just skim it to get the idea.
> http://claymore.engineer.gvsu.edu/~jackh/books/plcs/pdf/plcbook5_1.pdf
> http://sites.google.com/site/automatedmanufacturingsystems/
.....
> Peter Nachtwey
I wish I had time to get into this subject, but I don't. I'll just say one more thing.
The resultant I posted for x0 has 593 terms. Evidently plugging in values for x1, x2, ... , r3 is too complex for your PLC system. By a translation, one can assume that one point is the origin, say x1=0=y1. That reduces the polynomial to 215 terms. Optimizing this for numerical computing should be a routine task.
Looks like you folks ought to hire more mathematicians! There are certainly lots of unemployed ones these days.
Robert H. Lewis
Fordham University
Dave Rusin
pnachtwey <pnac...@gmail.com> wrote:
>>>>>> eq1: (x1-x0)^2+(y1-y0)^2-(r0+r1)^2;
>>>>>> eq2: (x2-x0)^2+(y2-y0)^2-(r0+r2)^2;
>>>>>> eq3: (x3-x0)^2+(y3-y0)^2-(r0+r3)^2;
>>>>>> eq4: solve([eq1,eq2,eq3],[x0,y0,r0]);
>> If, as the original post indicates, the three points are on an
>> equilateral triangle, that could probably simplify the formulas.
>I am not sure how. The equations seem pretty simple as they are.
Oh, it helps a lot because the equations simplify. If you have the
luxury of drawing the coordinate system after the diagram is presented
to you, then you can do it with the origin at point #1, the first
axis passing through point #2 and with a unit of measure equal to
half the distance between points 1 and 2. Then the equations are
simpler because x1=y1=y2=0, x2=2, x3=1, y3^2=3 . In this case a
computer-algebra system works much faster to determine the relationship
between the four radii. Moreover, the symmetry of the diagram
assures us that the polynomial relating r0 to the other three radii
will be symmetric in r1, r2, r3, and hence expressible in terms of
S = r1 + r2 + r3
T = r1*r2 + r2*r3 + r3*r1
U = r1*r2*r3
I get the defining relation to be simply 0 =
(S^4-4*S^2-4*S^2*T+6*U*S+T^2+8*T+16) +(2*S^3-7*S*T-4*S+9*U)*X +(S^2-3*T-3)*X^2
where X = 2*r0 . (Of course the actual presentation of the quadratic
relation for r0 is a matter of taste -- there are other sets of
generators for the symmetric polynomials besides {S,T,U} .)
So I think your problem is very easy to solve now:
1. Choose units so the equilateral triangle has sides of length 2
2. Compute S,T,U as above
3. Solve for r0 = X/2 as above.
As has been noted earlier in this thread, once r0 is known, you can solve
for x0, y0 by solving the two linear equations eq2-eq1 = eq3-eq1 = 0.
In the coordinate system I described above, the solution is
x0 = ( (4+r1^2-r2^2) + (r1-r2) *X )/ 4
y0 = ( (4+r1^2+r2^2-2*r3^2) + (r1+r2-2*r3)*X )/(4*sqrt(3))
but you can also easily solve the linear system in a
previously-defined coordinate system, too.
dave
Robert H. Lewis
> eq2: (x2-x0)^2+(y2-y0)^2-(r0+r2)^2;
> eq3: (x3-x0)^2+(y3-y0)^2-(r0+r3)^2;
>
> radius is r0. x1,y1 and x2,y2 and x3,y3 are the locations of the
> ultrasonic sensors. r1, r2, r3 are the distances between the sensors and
> the inner circle. It looks easy but my wxMaxima has been working on the
> solution for about a half hour and so far there is no solution.
I couldn't resist generalizing this to 3 dimensions. For 4 spheres in 3-space the equations are
eq1 := (x1-x0)^2+(y1-y0)^2+(z1-z0)^2-(r0+r1)^2;
eq2 := (x2-x0)^2+(y2-y0)^2+(z2-z0)^2-(r0+r2)^2;
eq3 := (x3-x0)^2+(y3-y0)^2+(z3-z0)^2-(r0+r3)^2;
eq4 := (x4-x0)^2+(y4-y0)^2+(z4-z0)^2-(r0+r4)^2;
We want to solve for x0, y0, z0, r0. One could reasonably assume that one center is the origin and plug in, say, x1=y1=z1=0. I did NOT do that.
The resultant for x0 has 18366 terms. This takes 10 meg RAM and 0.7 seconds on a year-old Imac. I'm using the Dixon-EDF method as explained in the reference below, running on Fermat.
If I plug in x1=y1=z1=0 at the start and run it, it takes 0.08 seconds. The resultant has 4214 terms.
About 7 years ago we had a discussion here about similar problems, Apollonius-type sphere and ellipsoid problems. A lot can happen in 7 years. I wonder how Groebner basis and other techniques do now on this particular problem?
Robert H. Lewis
Fordham University
http://fordham.academia.edu/RobertLewis/Papers/82784/Apollonius-Problems-in-Biochemistry
http://fordham.academia.edu/RobertLewis/Papers/80768/Algorithmic-Search-for-Flexibility-using-Resultants-of-Polynomial-Systems
Dave Rusin
Robert H. Lewis <rle...@fordham.edu> wrote:
>I couldn't resist generalizing this to 3 dimensions. For 4 spheres in 3-space the equations are
>
>eq1 := (x1-x0)^2+(y1-y0)^2+(z1-z0)^2-(r0+r1)^2;
>eq2 := (x2-x0)^2+(y2-y0)^2+(z2-z0)^2-(r0+r2)^2;
>eq3 := (x3-x0)^2+(y3-y0)^2+(z3-z0)^2-(r0+r3)^2;
>eq4 := (x4-x0)^2+(y4-y0)^2+(z4-z0)^2-(r0+r4)^2;
Since the original poster wanted us to consider the case in which the
points were equidistant, I looked at that special case for dimensions
3 through 10: normalizing so that the distance between any two points
is always 1, the relationship between r0 and the other radii r1, r2, ...
can be expressed as
0 = (8*(d+1)*s2-4*d*s1^2+2*(d+1))*r0^2
+ (-12*(d+1)*s3+4*s1+(12*d+4)*s1*s2-4*d*s1^3)*r0
+ (-d+(4*d+4)*s4-4*s2+(2-2*d)*s2^2-4*(d+1)*s3*s1+(4*s2*d+2)*s1^2-d*s1^4)
where s1, s2, s3, s4 are the first four (only!) of the elementary
symmetric functions in the r_i :
s1 = r1+r2+r3+...
s2 = r1*r2 + r1*r3 + ...
s3 = r1*r2*r3 + r1*r2*r4 + ...
s4 = r1*r2*r3*r4 + r1*r2*r3*r5 + ...
(In the planar case d=2, there are only three radii and so s4 = 0 .)
(Just as a sanity check: when d=1 this relation states that for
positive values of the r_i, we must have r0 = (1 - r1 - r2)/2
unless (r1-r2)^2 = 1, which is exactly what the geometry forces.)
dave
Robert H. Lewis
> <16050249.8831.1238689766096.JavaMail.jakarta@nitrogen
> .mathforum.org>,
> Robert H. Lewis <rlewis at fordham <dot> edu> wrote:
>
> >I couldn't resist generalizing this to 3 dimensions.
> For 4 spheres in 3-space the equations are
> >
> >eq1 := (x1-x0)^2+(y1-y0)^2+(z1-z0)^2-(r0+r1)^2;
> >eq2 := (x2-x0)^2+(y2-y0)^2+(z2-z0)^2-(r0+r2)^2;
> >eq3 := (x3-x0)^2+(y3-y0)^2+(z3-z0)^2-(r0+r3)^2;
> >eq4 := (x4-x0)^2+(y4-y0)^2+(z4-z0)^2-(r0+r4)^2;
>
> Since the original poster wanted us to consider the
> case in which the points were equidistant, I looked at that special
> case for dimensions 3 through 10: normalizing so that the distance
That VERY special case rarely holds in applications. I am asking here about the general case with the equations above. I would appreciate knowing how Groebner bases and various computer algebra systems handle those equations.
pnachtwey
> In article <b8cbcb94-b30b-4041-950c-212b14f0b...@j18g2000prm.googlegroups.com>,
>
> pnachtwey <pnacht...@gmail.com> wrote:
> >>>>>> eq1: (x1-x0)^2+(y1-y0)^2-(r0+r1)^2;
> >>>>>> eq2: (x2-x0)^2+(y2-y0)^2-(r0+r2)^2;
> >>>>>> eq3: (x3-x0)^2+(y3-y0)^2-(r0+r3)^2;
> >>>>>> eq4: solve([eq1,eq2,eq3],[x0,y0,r0]);
> >> If, as the original post indicates, the three points are on an
> >> equilateral triangle, that could probably simplify the formulas.
> >I am not sure how. The equations seem pretty simple as they are.
>
> Oh, it helps a lot because the equations simplify. If you have the
> luxury of drawing the coordinate system after the diagram is presented
> to you, then you can do it with the origin at point #1, the first
> axis passing through point #2 and with a unit of measure equal to
> half the distance between points 1 and 2. Then the equations are
> simpler because x1=y1=y2=0, x2=2, x3=1, y3^2=3 . In this case a
> computer-algebra system works much faster to determine the relationship
> between the four radii.
I tried that. You are right. It is much faster. Now my wxMaxima
returns with no answer much quicker :) I have had time to try this
on my Mathcad to see if it is better.
> Moreover, the symmetry of the diagram
> assures us that the polynomial relating r0 to the other three radii
> will be symmetric in r1, r2, r3, and hence expressible in terms
I don't see that. The inner circle does not need to be at the center
of the three sensors. The inner circle can be off center.
Peter Nachtwey
pnachtwey
Where does one find information on the Dixon EDF? I did a search and
your name is on every article. You seem to be 'the man' in on this
topic. I saw the code on your website but it has few comments so it is
hard follow.
Peter Nachtwey
Roman Pearce
> > >eq1 := (x1-x0)^2+(y1-y0)^2+(z1-z0)^2-(r0+r1)^2;
> > >eq2 := (x2-x0)^2+(y2-y0)^2+(z2-z0)^2-(r0+r2)^2;
> > >eq3 := (x3-x0)^2+(y3-y0)^2+(z3-z0)^2-(r0+r3)^2;
> > >eq4 := (x4-x0)^2+(y4-y0)^2+(z4-z0)^2-(r0+r4)^2;
>
> I would appreciate knowing how Groebner bases and various computer
> algebra systems handle those equations.
For solving {x0,y0,z0,r0} in terms of all the other variables Maple's
Groebner bases basically suck. It takes 45 sec to run F4 and 62
seconds to run FGLM using rational function arithmetic (Maple code).
Evaluation and interpolation would obviously be much better, and the
slimgb algorithm would probably also be much better.
Robert H. Lewis
> > I couldn't resist generalizing this to 3
> dimensions. For 4 spheres in 3-space the equations are
> >
> > eq1 := (x1-x0)^2+(y1-y0)^2+(z1-z0)^2-(r0+r1)^2;
> > eq2 := (x2-x0)^2+(y2-y0)^2+(z2-z0)^2-(r0+r2)^2;
> > eq3 := (x3-x0)^2+(y3-y0)^2+(z3-z0)^2-(r0+r3)^2;
> > eq4 := (x4-x0)^2+(y4-y0)^2+(z4-z0)^2-(r0+r4)^2;
......
> > About 7 years ago we had a discussion here about
> similar problems, Apollonius-type sphere and
> ellipsoid problems. A lot can happen in 7 years. I
> wonder how Groebner basis and other techniques do now
> on this particular problem?
> >
> > Robert H. Lewis
> > Fordham University
> >
> >
> http://fordham.academia.edu/RobertLewis/Papers/82784/A
> pollonius-Probl...http://fordham.academia.edu/RobertLe
> wis/Papers/80768/Algorithmic-Sear...
>
> Where does one find information on the Dixon EDF? I
> did a search and your name is on every article. You seem to be 'the
> man' in on this topic. I saw the code on your website but it has few
> comments so it is hard follow.
>
> Peter Nachtwey
It's been published at
Heuristics to Accelerate the Dixon Resultant, Mathematics and Computers in Simulation 77, Issue 4, p. 400 -- 407, April 2008.
It's also in the paper I mentioned in an earlier post here, Algorithmic Search for Flexibility Using Resultants of Polynomial Systems (with E. Coutsias)