Groups | Search | Server Info | Keyboard shortcuts | Login | Register [http] [https] [nntp] [nntps]

Groups > comp.lang.forth > #13188 > unrolled thread

intersection of circles

Back to article view | Back to comp.lang.forth

Contents

  intersection of circles Krishna Myneni <krishna.myneni@ccreweb.org> - 2012-06-22 17:31 -0700
    Re: intersection of circles mhx@iae.nl (Marcel Hendrix) - 2012-06-23 12:26 +0200
      Re: intersection of circles Krishna Myneni <krishna.myneni@ccreweb.org> - 2012-06-23 04:24 -0700
    Re: intersection of circles Krishna Myneni <krishna.myneni@ccreweb.org> - 2012-06-23 04:28 -0700
    Re: intersection of circles Bernd Paysan <bernd.paysan@gmx.de> - 2012-06-23 22:12 +0200
      Re: intersection of circles Krishna Myneni <krishna.myneni@ccreweb.org> - 2012-06-23 16:06 -0700
        Re: intersection of circles Bernd Paysan <bernd.paysan@gmx.de> - 2012-06-24 02:11 +0200

#13188 — intersection of circles

Surprisingly, the following problem isn't easily handled by a couple
of scientific computing packages:

If solution(s) exist, find the intersection points of two circles in
the x-y plane. This requires solving a system of two qudratic
equations in two variables:

1)  (x - a)^2 + (y - b)^2 = r1^2

2)  (x - c)^2 + (y - d)^2 = r2^2

where the points (a, b) and (c, d) are the centers of two circles of
radii r1 and r2.

For example, if I take a = b = 0, c = d = 0.5, and r1 = r2 = 1, the
two intersection points may be found using Wolfram Alpha. However,
using the solve() function in Maxima, which is supposed to handle
systems of both linear and nonlinear equations, I obtain a result
indicating there is no solution (see below). I believe Matlab requires
some sort of optimization toolbox to solve a system of nonlinear
equations, so I don't know how direct is the procedure in Matlab. I
expect Mathematica can handle it, given that Wolfram Alpha finds the
solution.

The problem should be solvable with a quadratic equation solver (see
FSL algorithm 63). This problem may make a nice Forth exercise. Before
attempting the solution, the program should first check to see if
solutions exist -- the two circles may not intersect, or may intersect
at a single point, or intersect at two points.

Krishna

--
Example in Maxima:

$ maxima

Maxima 5.24.0 http://maxima.sourceforge.net
using Lisp GNU Common Lisp (GCL) GCL 2.6.7 (a.k.a. GCL)
Distributed under the GNU Public License. See the file COPYING.
Dedicated to the memory of William Schelter.
The function bug_report() provides bug reporting information.
(%i1) e1: x^2 + y^2 = 1;
                                       2     2
(%o1)                             y  + x  = 1
(%i2) e2: (x - 0.5)^2 + (y - 0.5)^2 = 1;
                                        2               2
(%o2)                     (y - 0.5)  + (x - 0.5)  = 1
(%i3) solve( [e1, e2], [x,y] );

rat: replaced -0.5 by -1/2 = -0.5

rat: replaced -0.5 by -1/2 = -0.5
(%o3)                                 []

The last line indicates that Maxima was not able to find a solution.
It does however find a solution for the following,

(%i4) e2: (x - 0.5)^2 + y^2 = 1;
                                    2             2
(%o4)                         y  + (x - 0.5)  = 1
(%i5) solve( [e1, e2], [x, y] );

rat: replaced -0.5 by -1/2 = -0.5
                         1        sqrt(15)        1      sqrt(15)
(%o5)          [[x = -, y = - --------], [x = -, y = --------]]
                          4           4            4         4

Formatting may be messed up in the output lines ...

[toc] | [next] | [standalone]

#13193

Krishna Myneni <krishna.myneni@ccreweb.org> writes Re: intersection of circles

> Surprisingly, the following problem isn't easily handled by a couple
> of scientific computing packages:

> If solution(s) exist, find the intersection points of two circles in
> the x-y plane. This requires solving a system of two quadratic
> equations in two variables:

> 1)  (x - a)^2 + (y - b)^2 = r1^2

> 2)  (x - c)^2 + (y - d)^2 = r2^2

I guess 'requires' should be 'could be done by'?

> where the points (a, b) and (c, d) are the centers of two circles of
> radii r1 and r2.

[..]

-marcel

-- ------------------

0e FVALUE a    0e FVALUE b    0e FVALUE r1  
0e FVALUE c    0e FVALUE d    0e FVALUE r2  

0e FVALUE alpha  
0e FVALUE beta   

-- Shift a,b to origin, then rotate c-a, d-b to x-axis  
: transform-> ( -- )
	d b F- FSQR  c a F- FSQR F+ FSQRT TO alpha
	d b F-       c a F- FATAN2        TO beta ; 

: solution' ( F: -- x y )
	alpha F0= ABORT" SOLUTION' :: The origins of the circles should be distinct"
	r1 FSQR  r2 FSQR F- alpha FSQR F+   alpha F2* F/ ( F: -- x )
	r1 FSQR  FOVER FSQR F- FDUP F0< ABORT" SOLUTION' :: the circles do not touch" FSQRT ; 

: transform<- ( F: x y -- x1 y1 x2 y2 )
	F2DUP FSWAP FATAN2       FLOCAL gamma
	FSQR FSWAP FSQR F+ FSQRT FLOCAL rr 
	beta gamma F+ FSINCOS rr F* a F+  FSWAP rr F* b F+
	beta gamma F- FSINCOS rr F* a F+  FSWAP rr F* b F+ ; 

: 2CIRCLES ( a b r1  c d r2 -- x1 y1 x2 y2 )
	TO r2 TO d TO c  TO r1 TO b TO a
	transform-> solution' transform<- ;

        A few results:

	FORTH> 0e 0e 1e  0.5e 0.5e 1e    2CIRCLES f2swap fswap f. f. fswap f. f.  -0.411438 0.911438 0.911438 -0.411438  ok
	FORTH> 0e 0e 1e  0.5e 1e   1e    2CIRCLES f2swap fswap f. f. fswap f. f.  -0.491620 0.870810 0.991620  0.129190  ok
	FORTH> 0e 0e 1e  0.5e 0e   1e    2CIRCLES f2swap fswap f. f. fswap f. f.   0.250000 0.968246 0.250000 -0.968246  ok
	FORTH> 0e 0e 1e  0.5e 0e   0.25e 2CIRCLES f2swap fswap f. f. fswap f. f.
	**** RETURN STACK DUMP ****
		SOLUTION'                                                  (+$000000BE)
		2CIRCLES                                                   (+$000000C9)
	Error -2
	SOLUTION' :: the circles do not touch ?
	FORTH> 0e 0e   0.5e 1e    0e   0.5e 2CIRCLES f2swap fswap f. f. fswap f. f.   0.500000  0.000000  0.500000  0.000000  ok
	FORTH> 0e 0e   1e   0.5e  0e   1e   2CIRCLES f2swap fswap f. f. fswap f. f.   0.250000  0.968246  0.250000 -0.968246  ok
	FORTH> 0e 0e   1e  -0.5e  0e   1e   2CIRCLES f2swap fswap f. f. fswap f. f.  -0.250000 -0.968246 -0.250000  0.968246  ok
	FORTH> 0e 0e   1e  -0.5e  0.5e 1e   2CIRCLES f2swap fswap f. f. fswap f. f.  -0.911438 -0.411438  0.411438  0.911438  ok
	FORTH> 0e 0e   1e  -0.5e -0.5e 1e   2CIRCLES f2swap fswap f. f. fswap f. f.   0.411438 -0.911438 -0.911438  0.411438  ok
	FORTH> 0e 0e   1e   0.5e -0.5e 1e   2CIRCLES f2swap fswap f. f. fswap f. f.   0.911438  0.411438 -0.411438 -0.911438  ok
	FORTH> 0e 0.1e 1e   0.5e -0.5e 1e   2CIRCLES f2swap fswap f. f. fswap f. f.   0.957223  0.389352 -0.457223 -0.789352  ok

[toc] | [prev] | [next] | [standalone]

#13195

On Jun 23, 5:26 am, m...@iae.nl (Marcel Hendrix) wrote:
> Krishna Myneni <krishna.myn...@ccreweb.org> writes Re: intersection of circles
>
> > Surprisingly, the following problem isn't easily handled by a couple
> > of scientific computing packages:
> > If solution(s) exist, find the intersection points of two circles in
> > the x-y plane. This requires solving a system of two quadratic
> > equations in two variables:
> > 1)  (x - a)^2 + (y - b)^2 = r1^2
> > 2)  (x - c)^2 + (y - d)^2 = r2^2
>
> I guess 'requires' should be 'could be done by'?
>

Nice. I hadn't considered doing the rotation to further simplify the
equations. I think the "requires" in the statement is valid -- all
solutions must obey equations 1) and 2). There may be several methods,
however, of finding the solutions of 1) and 2). Your elegant procedure
consists of first doing the translation transformation:

x' = x - a, y' = y - b

followed by a rotation transformation,

x'' =  x'*cos(phi) + y'*sin(phi)
y'' = -x'*sin(phi) + y'*cos(phi)

where phi = atan( (d-b)/(c-a) )

to equations (1) and (2), to obtain a simpler set of equations to
solve. I need to go through the algebra to see what the equations look
like in the double-primed coordinates.

Krishna


>
> 0e FVALUE a    0e FVALUE b    0e FVALUE r1
> 0e FVALUE c    0e FVALUE d    0e FVALUE r2
>
> 0e FVALUE alpha
> 0e FVALUE beta
>
> -- Shift a,b to origin, then rotate c-a, d-b to x-axis
> : transform-> ( -- )
>         d b F- FSQR  c a F- FSQR F+ FSQRT TO alpha
>         d b F-       c a F- FATAN2        TO beta ;
>
> : solution' ( F: -- x y )
>         alpha F0= ABORT" SOLUTION' :: The origins of the circles should be distinct"
>         r1 FSQR  r2 FSQR F- alpha FSQR F+   alpha F2* F/ ( F: -- x )
>         r1 FSQR  FOVER FSQR F- FDUP F0< ABORT" SOLUTION' :: the circles do not touch" FSQRT ;
>
> : transform<- ( F: x y -- x1 y1 x2 y2 )
>         F2DUP FSWAP FATAN2       FLOCAL gamma
>         FSQR FSWAP FSQR F+ FSQRT FLOCAL rr
>         beta gamma F+ FSINCOS rr F* a F+  FSWAP rr F* b F+
>         beta gamma F- FSINCOS rr F* a F+  FSWAP rr F* b F+ ;
>
> : 2CIRCLES ( a b r1  c d r2 -- x1 y1 x2 y2 )
>         TO r2 TO d TO c  TO r1 TO b TO a
>         transform-> solution' transform<- ;
>
>         A few results:
>
>         FORTH> 0e 0e 1e  0.5e 0.5e 1e    2CIRCLES f2swap fswap f. f. fswap f. f.  -0.411438 0.911438 0.911438 -0.411438  ok
>         FORTH> 0e 0e 1e  0.5e 1e   1e    2CIRCLES f2swap fswap f. f. fswap f. f.  -0.491620 0.870810 0.991620  0.129190  ok
>         FORTH> 0e 0e 1e  0.5e 0e   1e    2CIRCLES f2swap fswap f. f. fswap f. f.   0.250000 0.968246 0.250000 -0.968246  ok
>         FORTH> 0e 0e 1e  0.5e 0e   0.25e 2CIRCLES f2swap fswap f. f. fswap f. f.
>         **** RETURN STACK DUMP ****
>                 SOLUTION'                                                  (+$000000BE)
>                 2CIRCLES                                                   (+$000000C9)
>         Error -2
>         SOLUTION' :: the circles do not touch ?
...

[toc] | [prev] | [next] | [standalone]

#13196

On Jun 22, 7:31 pm, Krishna Myneni <krishna.myn...@ccreweb.org> wrote:
...
> For example, if I take a = b = 0, c = d = 0.5, and r1 = r2 = 1, the
> two intersection points may be found using Wolfram Alpha. However,
> using the solve() function in Maxima, which is supposed to handle
> systems of both linear and nonlinear equations, I obtain a result
> indicating there is no solution (see below).
...
> --
> Example in Maxima:
>
> $ maxima
>
> Maxima 5.24.0http://maxima.sourceforge.net
> using Lisp GNU Common Lisp (GCL) GCL 2.6.7 (a.k.a. GCL)
> Distributed under the GNU Public License. See the file COPYING.
> Dedicated to the memory of William Schelter.
> The function bug_report() provides bug reporting information.
> (%i1) e1: x^2 + y^2 = 1;
>                                    2    2
> (%o1)                             y  + x  = 1
> (%i2) e2: (x - 0.5)^2 + (y - 0.5)^2 = 1;
>                                    2            2
> (%o2)                     (y - 0.5)  + (x - 0.5)  = 1
> (%i3) solve( [e1, e2], [x,y] );
>
> rat: replaced -0.5 by -1/2 = -0.5
>
> rat: replaced -0.5 by -1/2 = -0.5
> (%o3)                                 []
>
> The last line indicates that Maxima was not able to find a solution.


I followed up with a question on the Maxima mailing list about why its
solve() function doesn't work with the above example, and the response
was that solve() has a bug, believed to not exist in earlier versions
of Maxima (although that hasn't been verified, I think).

KM

[toc] | [prev] | [next] | [standalone]

#13202

Krishna Myneni wrote:

> Surprisingly, the following problem isn't easily handled by a couple
> of scientific computing packages:
> 
> If solution(s) exist, find the intersection points of two circles in
> the x-y plane. This requires solving a system of two qudratic
> equations in two variables:
> 
> 1)  (x - a)^2 + (y - b)^2 = r1^2
> 
> 2)  (x - c)^2 + (y - d)^2 = r2^2
> 
> where the points (a, b) and (c, d) are the centers of two circles of
> radii r1 and r2.

I had that problem a bit over a month ago, when I solved the Triceps 2 
geometry.  There, intersection of two circles is one of the steps to be 
done.  Fortunately, I didn't know that this was a "hard" problem, which 
required a quadratic solver, and just decided to use the law of cosines.

You know:  The distance of the two points r3²=(a-c)²+(b-d)².  You know 
r1 and r2.  These are the three sides of a triangle, which could be on 
either side (for Triceps 2, only one side is possible, as gravity goes 
downwards).

cos alpha=(r1²+r2²-r3²)/(2*r1*r2)

You need to add the angle of the rectangular triangle where r3 is the 
hypothenuse, and the a and b sides are parallel to x and y, but this is 
left as an exersice to the reader.

-- 
Bernd Paysan
"If you want it done right, you have to do it yourself"
http://bernd-paysan.de/

[toc] | [prev] | [next] | [standalone]

#13205

On Jun 23, 3:12 pm, Bernd Paysan <bernd.pay...@gmx.de> wrote:
> Krishna Myneni wrote:
> > Surprisingly, the following problem isn't easily handled by a couple
> > of scientific computing packages:
>
> > If solution(s) exist, find the intersection points of two circles in
> > the x-y plane. This requires solving a system of two qudratic
> > equations in two variables:
>
> > 1)  (x - a)^2 + (y - b)^2 = r1^2
>
> > 2)  (x - c)^2 + (y - d)^2 = r2^2
>
> > where the points (a, b) and (c, d) are the centers of two circles of
> > radii r1 and r2.
>
> I had that problem a bit over a month ago, when I solved the Triceps 2
> geometry.  There, intersection of two circles is one of the steps to be
> done.  Fortunately, I didn't know that this was a "hard" problem, which
> required a quadratic solver, and just decided to use the law of cosines.
>

It's always true that a geometric problem has an algebraic
counterpart, and the converse is often true. In such cases, sometimes
the problem is easier to solve by applying geometry rather than
solving the algebraic counterpart. Even the algebraic version is not
"hard", just tedious. What surprised me was that some high level
computation packages were either unable to solve the algebraic
problem, or did not offer even a rudimentary solver for a system of
the simplest types of nonlinear equations (without resort to an add-on
toolkit). One of the advantages of a computer algebra system is that
we can easily extend the equations to three dimensions (e.g. the
intersection of three spheres). My visualization is not quite good
enough to handle such a problem geometrically.

> You know:  The distance of the two points r3²=(a-c)²+(b-d)².  You know
> r1 and r2.  These are the three sides of a triangle, which could be on
> either side (for Triceps 2, only one side is possible, as gravity goes
> downwards).
>
> cos alpha=(r1²+r2²-r3²)/(2*r1*r2)
>
> You need to add the angle of the rectangular triangle where r3 is the
> hypothenuse, and the a and b sides are parallel to x and y, but this is
> left as an exersice to the reader.
>

Will have to work this out. Thanks.

Krishna

> --
> Bernd Paysan
> "If you want it done right, you have to do it yourself"http://bernd-paysan.de/

[toc] | [prev] | [next] | [standalone]

#13207

Krishna Myneni wrote:
> One of the advantages of a computer algebra system is that
> we can easily extend the equations to three dimensions (e.g. the
> intersection of three spheres). My visualization is not quite good
> enough to handle such a problem geometrically.

Triceps 2 actually computes intersection of a circle with a sphere.  
Fortunately, I have the hardware here, so visualizing it is to some 
extent just looking at it :-).

So three spheres, that's two spheres intersecting into a circle (same 
triangle as before, but its degree of freedom is a plane, i.e. it can 
rotate around the (a,b)-(c,d) axis), and then a circle intersecting a 
sphere, giving again two possible solutions.  To intersect circle and 
sphere, you need to cut the sphere with the circle's plane, so it's a 
circle&circle problem again (see above).

Don't ask me for intersecting four 4D hyperspheres.

-- 
Bernd Paysan
"If you want it done right, you have to do it yourself"
http://bernd-paysan.de/

[toc] | [prev] | [standalone]

Back to top | Article view | comp.lang.forth

csiph-web