Line circle intersection: Difference between revisions

From Rosetta Code
Content added Content deleted
m (→‎{{header|Haskell}}: One possible consolidation of the main – OP may prefer to revert.)
m (added a blank line)
Line 22: Line 22:
*See [https://math.stackexchange.com/questions/228841/how-do-i-calculate-the-intersections-of-a-straight-line-and-a-circle Math Stack Exchange] for development of the formulae needed.
*See [https://math.stackexchange.com/questions/228841/how-do-i-calculate-the-intersections-of-a-straight-line-and-a-circle Math Stack Exchange] for development of the formulae needed.
*See [https://mathworld.wolfram.com/Circle-LineIntersection.html Wolfram] for the formulae needed.
*See [https://mathworld.wolfram.com/Circle-LineIntersection.html Wolfram] for the formulae needed.



=={{header|Go}}==
=={{header|Go}}==

Revision as of 21:41, 16 March 2020

Line circle intersection is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.

In plane geometry, a line (or segment) may intersect a circle at 0, 1 or 2 points.

Task

Implement a method (function, procedure etc.) in your language which takes as parameters:

  • the starting point of a line;
  • the point where the line ends;
  • the center point of a circle;
  • the circle's radius; and
  • whether the line is a segment or extends to infinity beyond the above points.


The method should return the intersection points (if any) of the circle and the line.

Illustrate your method with some examples (or use the Go examples below).

References


Go

<lang go>package main

import (

   "fmt"
   "math"

)

const eps = 1e-14 // say

type point struct{ x, y float64 }

func (p point) String() string {

   // hack to get rid of negative zero
   // compiler treats 0 and -0 as being same 
   if p.x == 0 {
       p.x = 0
   }
   if p.y == 0 {
       p.y = 0
   }
   return fmt.Sprintf("(%g, %g)", p.x, p.y)

}

func sq(x float64) float64 { return x * x }

// Returns the intersection points (if any) of a circle, center 'cp' with radius 'r', // and either an infinite line containing the points 'p1' and 'p2' // or a segment drawn between those points. func intersects(p1, p2, cp point, r float64, segment bool) []point {

   var res []point
   x0, y0 := cp.x, cp.y
   x1, y1 := p1.x, p1.y
   x2, y2 := p2.x, p2.y
   A := y2 - y1
   B := x1 - x2
   C := x2*y1 - x1*y2
   a := sq(A) + sq(B)
   var b, c float64
   var bnz = true
   if math.Abs(B) >= eps { // if B isn't zero or close to it
       b = 2 * (A*C + A*B*y0 - sq(B)*x0)
       c = sq(C) + 2*B*C*y0 - sq(B)*(sq(r)-sq(x0)-sq(y0))
   } else {
       b = 2 * (B*C + A*B*x0 - sq(A)*y0)
       c = sq(C) + 2*A*C*x0 - sq(A)*(sq(r)-sq(x0)-sq(y0))
       bnz = false
   }
   d := sq(b) - 4*a*c // discriminant
   if d < 0 {
       // line & circle don't intersect
       return res
   }
   // checks whether a point is within a segment
   within := func(x, y float64) bool {
       d1 := math.Sqrt(sq(x2-x1) + sq(y2-y1)) // distance between end-points
       d2 := math.Sqrt(sq(x-x1) + sq(y-y1))   // distance from point to one end
       d3 := math.Sqrt(sq(x2-x) + sq(y2-y))   // distance from point to other end
       delta := d1 - d2 - d3
       return math.Abs(delta) < eps // true if delta is less than a small tolerance
   }
   var x, y float64
   fx := func() float64 { return -(A*x + C) / B }
   fy := func() float64 { return -(B*y + C) / A }
   rxy := func() {
       if !segment || within(x, y) {
           res = append(res, point{x, y})
       }
   }
   if d == 0 {
       // line is tangent to circle, so just one intersect at most
       if bnz {
           x = -b / (2 * a)
           y = fx()
           rxy()
       } else {
           y = -b / (2 * a)
           x = fy()
           rxy()
       }
   } else {
       // two intersects at most
       d = math.Sqrt(d)
       if bnz {
           x = (-b + d) / (2 * a)
           y = fx()
           rxy()
           x = (-b - d) / (2 * a)
           y = fx()
           rxy()
       } else {
           y = (-b + d) / (2 * a)
           x = fy()
           rxy()
           y = (-b - d) / (2 * a)
           x = fy()
           rxy()
       }
   }
   return res

}

func main() {

   cp := point{3, -5}
   r := 3.0
   fmt.Println("The intersection points (if any) between:")
   fmt.Println("\n  A circle, center (3, -5) with radius 3, and:")
   fmt.Println("\n    a line containing the points (-10, 11) and (10, -9) is/are:")
   fmt.Println("     ", intersects(point{-10, 11}, point{10, -9}, cp, r, false))
   fmt.Println("\n    a segment starting at (-10, 11) and ending at (-11, 12) is/are")
   fmt.Println("     ", intersects(point{-10, 11}, point{-11, 12}, cp, r, true))
   fmt.Println("\n    a horizontal line containing the points (3, -2) and (7, -2) is/are:")
   fmt.Println("     ", intersects(point{3, -2}, point{7, -2}, cp, r, false))
   cp = point{0, 0}
   r = 4.0
   fmt.Println("\n  A circle, center (0, 0) with radius 4, and:")
   fmt.Println("\n    a vertical line containing the points (0, -3) and (0, 6) is/are:")
   fmt.Println("     ", intersects(point{0, -3}, point{0, 6}, cp, r, false))
   fmt.Println("\n    a vertical segment starting at (0, -3) and ending at (0, 6) is/are:")
   fmt.Println("     ", intersects(point{0, -3}, point{0, 6}, cp, r, true))
   cp = point{4, 2}
   r = 5.0
   fmt.Println("\n  A circle, center (4, 2) with radius 5, and:")
   fmt.Println("\n    a line containing the points (6, 3) and (10, 7) is/are:")
   fmt.Println("     ", intersects(point{6, 3}, point{10, 7}, cp, r, false))
   fmt.Println("\n    a segment starting at (7, 4) and ending at (11, 8) is/are:")
   fmt.Println("     ", intersects(point{7, 4}, point{11, 8}, cp, r, true))

}</lang>

Output:
The intersection points (if any) between:

  A circle, center (3, -5) with radius 3, and:

    a line containing the points (-10, 11) and (10, -9) is/are:
      [(6, -5) (3, -2)]

    a segment starting at (-10, 11) and ending at (-11, 12) is/are
      []

    a horizontal line containing the points (3, -2) and (7, -2) is/are:
      [(3, -2)]

  A circle, center (0, 0) with radius 4, and:

    a vertical line containing the points (0, -3) and (0, 6) is/are:
      [(0, 4) (0, -4)]

    a vertical segment starting at (0, -3) and ending at (0, 6) is/are:
      [(0, 4)]

  A circle, center (4, 2) with radius 5, and:

    a line containing the points (6, 3) and (10, 7) is/are:
      [(8, 5) (1, -2)]

    a segment starting at (7, 4) and ending at (11, 8) is/are:
      [(8, 5)]

Haskell

<lang haskell>import Data.Tuple.Curry

main :: IO () main =

 mapM_ putStrLn $
 concatMap
   (("" :) . uncurryN task)
   [ ((-10, 11), (10, -9), ((3, -5), 3))
   , ((-10, 11), (-11, 12), ((3, -5), 3))
   , ((3, -2), (7, -2), ((3, -5), 3))
   , ((3, -2), (7, -2), ((0, 0), 4))
   , ((0, -3), (0, 6), ((0, 0), 4))
   , ((6, 3), (10, 7), ((4, 2), 5))
   , ((7, 4), (11, 18), ((4, 2), 5))
   ]

task :: (Double, Double)

    -> (Double, Double)
    -> ((Double, Double), Double)
    -> [String]

task pt1 pt2 circle@(pt3@(a3, b3), r) = [line, segment]

 where
   xs = map fun $ lineCircleIntersection pt1 pt2 circle
   ys = map fun $ segmentCircleIntersection pt1 pt2 circle
   to x = (fromIntegral . round $ 100 * x) / 100
   fun (x, y) = (to x, to y)
   yo = show . fun
   start = "Intersection: Circle " ++ yo pt3 ++ " " ++ show (to r) ++ " and "
   end = yo pt1 ++ " " ++ yo pt2 ++ ": "
   line = start ++ "Line " ++ end ++ show xs
   segment = start ++ "Segment " ++ end ++ show ys

segmentCircleIntersection

 :: (Double, Double)
 -> (Double, Double)
 -> ((Double, Double), Double)
 -> [(Double, Double)]

segmentCircleIntersection pt1 pt2 circle =

 filter (go p1 p2) $ lineCircleIntersection pt1 pt2 circle
 where
   [p1, p2]
     | pt1 < pt2 = [pt1, pt2]
     | otherwise = [pt2, pt1]
   go (x, y) (u, v) (i, j)
     | x == u = y <= j && j <= v
     | otherwise = x <= i && i <= u

lineCircleIntersection

 :: (Double, Double)
 -> (Double, Double)
 -> ((Double, Double), Double)
 -> [(Double, Double)]

lineCircleIntersection (a1, b1) (a2, b2) ((a3, b3), r) = go delta

 where
   (x1, x2) = (a1 - a3, a2 - a3)
   (y1, y2) = (b1 - b3, b2 - b3)
   (dx, dy) = (x2 - x1, y2 - y1)
   drdr = dx * dx + dy * dy
   d = x1 * y2 - x2 * y1
   delta = r * r * drdr - d * d
   sqrtDelta = sqrt delta
   (sgnDy, absDy) = (sgn dy, abs dy)
   u1 = (d * dy + sgnDy * dx * sqrtDelta) / drdr
   u2 = (d * dy - sgnDy * dx * sqrtDelta) / drdr
   v1 = (-d * dx + absDy * sqrtDelta) / drdr
   v2 = (-d * dx - absDy * sqrtDelta) / drdr
   go x
     | 0 > x = []
     | 0 == x = [(u1 + a3, v1 + b3)]
     | otherwise = [(u1 + a3, v1 + b3), (u2 + a3, v2 + b3)]

sgn :: Double -> Double sgn x

 | 0 > x = -1
 | otherwise = 1</lang>
Output:
Intersection: Circle (3.0,-5.0) 3.0 and Line (-10.0,11.0) (10.0,-9.0): [(3.0,-2.0),(6.0,-5.0)]
Intersection: Circle (3.0,-5.0) 3.0 and Segment (-10.0,11.0) (10.0,-9.0): [(3.0,-2.0),(6.0,-5.0)]

Intersection: Circle (3.0,-5.0) 3.0 and Line (-10.0,11.0) (-11.0,12.0): [(3.0,-2.0),(6.0,-5.0)]
Intersection: Circle (3.0,-5.0) 3.0 and Segment (-10.0,11.0) (-11.0,12.0): []

Intersection: Circle (3.0,-5.0) 3.0 and Line (3.0,-2.0) (7.0,-2.0): [(3.0,-2.0)]
Intersection: Circle (3.0,-5.0) 3.0 and Segment (3.0,-2.0) (7.0,-2.0): [(3.0,-2.0)]

Intersection: Circle (0.0,0.0) 4.0 and Line (3.0,-2.0) (7.0,-2.0): [(3.46,-2.0),(-3.46,-2.0)]
Intersection: Circle (0.0,0.0) 4.0 and Segment (3.0,-2.0) (7.0,-2.0): [(3.46,-2.0)]

Intersection: Circle (0.0,0.0) 4.0 and Line (0.0,-3.0) (0.0,6.0): [(0.0,4.0),(0.0,-4.0)]
Intersection: Circle (0.0,0.0) 4.0 and Segment (0.0,-3.0) (0.0,6.0): [(0.0,4.0)]

Intersection: Circle (4.0,2.0) 5.0 and Line (6.0,3.0) (10.0,7.0): [(8.0,5.0),(1.0,-2.0)]
Intersection: Circle (4.0,2.0) 5.0 and Segment (6.0,3.0) (10.0,7.0): [(8.0,5.0)]

Intersection: Circle (4.0,2.0) 5.0 and Line (7.0,4.0) (11.0,18.0): [(7.46,5.61),(5.03,-2.89)]
Intersection: Circle (4.0,2.0) 5.0 and Segment (7.0,4.0) (11.0,18.0): [(7.46,5.61)]

Perl

<lang perl>use strict; use warnings; use feature 'say'; use List::Util 'sum';

sub find_intersection {

  my($P1, $P2, $center, $radius) = @_;
  my @d = ($$P2[0] -     $$P1[0], $$P2[1] -     $$P1[1]);
  my @f = ($$P1[0] - $$center[0], $$P1[1] - $$center[1]);
  my $a = sum map { $_**2 } @d;
  my $b = 2 * ($f[0]*$d[0] + $f[1]*$d[1]);
  my $c = sum(map { $_**2 } @f) - $radius**2;
  my $D =  $b**2 - 4*$a*$c;
  return unless $D >= 0;
  my ($t1, $t2) = ( (-$b - sqrt $D) / (2*$a), (-$b + sqrt $D) / (2*$a) );
  return unless $t1 >= 0 and $t1 <= 1 or $t2 >= 0 and $t2  <= 1;
  my ($dx, $dy) = ($$P2[0] - $$P1[0], $$P2[1] - $$P1[1]);
  return ([$dx*$t1 + $$P1[0], $dy*$t1 + $$P1[1]],
          [$dx*$t2 + $$P1[0], $dy*$t2 + $$P1[1]])

}

my @data = (

  [ [-10, 11], [ 10, -9], [3, -5], 3 ],
  [ [-10, 11], [-11, 12], [3, -5], 3 ],
  [ [  3, -2], [  7, -2], [3, -5], 3 ],
  [ [  3, -2], [  7, -2], [0,  0], 4 ],
  [ [  0, -3], [  0,  6], [0,  0], 4 ],
  [ [  6,  3], [ 10,  7], [4,  2], 5 ],
  [ [  7,  4], [ 11, 18], [4,  2], 5 ],

);

sub rnd { map { sprintf('%.2f', $_) =~ s/\.00//r } @_ }

for my $d (@data) {

  my @solution = find_intersection @$d[0] , @$d[1] , @$d[2], @$d[3];
  say 'For input: ' . join ', ', (map { '('. join(',', @$_) .')' } @$d[0,1,2]), @$d[3];
  say 'Solutions: ' . (@solution > 1 ? join ', ', map { '('. join(',', rnd @$_) .')' } @solution : 'None');
  say ;

}</lang>

Output:
For input: (-10,11), (10,-9), (3,-5), 3
Solutions: (3,-2), (6,-5)

For input: (-10,11), (-11,12), (3,-5), 3
Solutions: None

For input: (3,-2), (7,-2), (3,-5), 3
Solutions: (3,-2), (3,-2)

For input: (3,-2), (7,-2), (0,0), 4
Solutions: (-3.46,-2), (3.46,-2)

For input: (0,-3), (0,6), (0,0), 4
Solutions: (0,-4), (0,4)

For input: (6,3), (10,7), (4,2), 5
Solutions: (1,-2), (8,5)

For input: (7,4), (11,18), (4,2), 5
Solutions: (5.03,-2.89), (7.46,5.61)

Phix

Translation of: Go
Translation of: zkl

<lang Phix>constant epsilon = 1e-14 -- say atom cx, cy, r, x1, y1, x2, y2

function sq(atom x) return x*x end function

function within(atom x, y)

   --
   -- checks whether a point is within a segment
   -- ie:    <-------d1------->
   --        <--d2--><---d3--->    -- within, d2+d3 ~= d1
   --   x1,y1^      ^x,y      ^x2,y2
   -- vs:
   --  <-d2->
   --  <-----------d3--------->    -- not "", d2+d3 > d1
   --  ^x,y - and obviously ditto when x,y is (say) out here^
   --
   -- (obviously only works when x,y is on the same line as x1,y1 to x2,y2)
   --
   atom d1 := sqrt(sq(x2-x1) + sq(y2-y1)), -- distance between end-points
        d2 := sqrt(sq(x -x1) + sq(y -y1)), -- distance from point to one end
        d3 := sqrt(sq(x2-x ) + sq(y2-y )), -- distance from point to other end
        delta := (d2 + d3) - d1
   return abs(delta) < epsilon -- true if delta is less than a small tolerance

end function

function pf(atom x,y)

   return sprintf("(%g,%g)",{x,y})

end function

function intersects(bool bSegment) -- -- Returns the intersection points (if any) of a circle, center (cx,cy) with radius r, -- and line containing the points (x1,y1) and (x2,y2) being either infinite or limited -- to the segment drawn between those points. --

   sequence res = {}
   atom A = y2 - y1,       sqA = sq(A),
        B = x1 - x2,       sqB = sq(B),
        C = x2*y1 - x1*y2, sqC = sq(C),
        sqr = r*r-cx*cx-cy*cy,
        a := sqA + sqB,
        b, c
   bool bDivA = false
   if abs(B)<epsilon then  -- B is zero or close to it
       b = 2 * (B*C + A*B*cx - sqA*cy)
       c = sqC + 2*A*C*cx - sqA*sqr
       bDivA = true -- (and later divide by A instead!)
   else
       b = 2 * (A*C + A*B*cy - sqB*cx)
       c = sqC + 2*B*C*cy - sqB*sqr
   end if
   atom d := b*b - 4*a*c -- discriminant
   if d>=0 then -- (-ve means line & circle do not intersect)
       d = sqrt(d)
       atom ux,uy, vx,vy
       if bDivA then       
           {uy,vy} = sq_div(sq_sub({+d,-d},b),2*a)
           {ux,vx} = sq_div(sq_sub(sq_mul(-B,{uy,vy}),C),A)
       else
           {ux,vx} = sq_div(sq_sub({+d,-d},b),2*a)
           {uy,vy} = sq_div(sq_sub(sq_mul(-A,{ux,vx}),C),B)
       end if
       if not bSegment or within(ux,uy) then
           res = append(res,pf(ux,uy))
       end if
       if d!=0 and (not bSegment or within(vx,vy)) then
           res = append(res,pf(vx,vy))
       end if              
   end if
   return res

end function

-- cx cy r x1 y1 x2 y2 bSegment constant tests = {{3,-5,3,{{-10,11, 10,-9,false},

                          {-10,11,-11,12,true},
                          {  3,-2,  7,-2,false}}},
                 {0, 0,4,{{  0,-3,  0, 6,false},
                          {  0,-3,  0, 6,true}}},
                 {4, 2,5,{{  6, 3, 10, 7,false},
                          {  7, 4, 11, 8,true}}}}

for t=1 to length(tests) do

   {cx, cy, r, sequence lines} = tests[t]
   string circle = sprintf("Circle at %s radius %d",{pf(cx,cy),r})
   for l=1 to length(lines) do
       {x1, y1, x2, y2, bool bSegment} = lines[l]
       sequence res = intersects(bSegment)
       string ls = iff(bSegment?"segment":"   line"),
              at = iff(length(res)?"intersect at "&join(res," and ")
                                  :"do not intersect")
       printf(1,"%s and %s %s to %s %s.\n",{circle,ls,pf(x1,y1),pf(x2,y2),at})
       circle = repeat(' ',length(circle))
   end for

end for</lang>

Output:
Circle at (3,-5) radius 3 and    line (-10,11) to (10,-9) intersect at (6,-5) and (3,-2).
                          and segment (-10,11) to (-11,12) do not intersect.
                          and    line (3,-2) to (7,-2) intersect at (3,-2).
Circle at (0,0) radius 4 and    line (0,-3) to (0,6) intersect at (0,4) and (0,-4).
                         and segment (0,-3) to (0,6) intersect at (0,4).
Circle at (4,2) radius 5 and    line (6,3) to (10,7) intersect at (8,5) and (1,-2).
                         and segment (7,4) to (11,8) intersect at (8,5).

Raku

(formerly Perl 6) Extend solution space to 3D. Reference: this SO question and answers <lang perl6>#!/usr/bin/env perl6

sub LineCircularOBJintersection(@P1, @P2, @Centre, \Radius) {

  my @d = @P2 »-« @P1 ;           # d
  my @f = @P1 »-« @Centre ;       # c
  my \a = [+] @d»²;               # d dot d
  my \b = 2 * ([+] @f »*« @d);    # 2 * f dot d
  my \c = ([+] @f»²) - Radius²;   # f dot f - r²
  my \Δ =  b²-(4*a*c);            # discriminant
  if (Δ < 0) {
     return [];
  } else {
     my (\t1,\t2) = (-b - Δ.sqrt)/(2*a), (-b + Δ.sqrt)/(2*a);
     if 0 ≤ t1|t2 ≤ 1 {
        return @P1 »+« ( @P2 »-« @P1 ) »*» t1, @P1 »+« ( @P2 »-« @P1 ) »*» t2
     } else {
        return []
     }
  }

}

my \DATA = [

  [ <-10 11>, < 10 -9>, <3 -5>, 3 ],
  [ <-10 11>, <-11 12>, <3 -5>, 3 ],
  [ <  3 -2>, <  7 -2>, <3 -5>, 3 ],
  [ <  3 -2>, <  7 -2>, <0  0>, 4 ],
  [ <  0 -3>, <  0  6>, <0  0>, 4 ],
  [ <  6  3>, < 10  7>, <4  2>, 5 ],
  [ <  7  4>, < 11 18>, <4  2>, 5 ],
  [ <5  2 −2.26 >, <0.77 2 4>, <1 4 0>, 4 ]

];

for DATA {

  my @solution = LineCircularOBJintersection $_[0] , $_[1] , $_[2], $_[3];
  say "For data set: ", $_;
  say "Solution(s) is/are: ", @solution.Bool ?? @solution !! "None";

}</lang>

Output:
For data set: [(-10 11) (10 -9) (3 -5) 3]
Solution(s) is/are: [(3 -2) (6 -5)]
For data set: [(-10 11) (-11 12) (3 -5) 3]
Solution(s) is/are: None
For data set: [(3 -2) (7 -2) (3 -5) 3]
Solution(s) is/are: [(3 -2) (3 -2)]
For data set: [(3 -2) (7 -2) (0 0) 4]
Solution(s) is/are: [(-3.4641016151377544 -2) (3.4641016151377544 -2)]
For data set: [(0 -3) (0 6) (0 0) 4]
Solution(s) is/are: [(0 -4) (0 4)]
For data set: [(6 3) (10 7) (4 2) 5]
Solution(s) is/are: [(1 -2) (8 5)]
For data set: [(7 4) (11 18) (4 2) 5]
Solution(s) is/are: [(5.030680985703315 -2.892616550038399) (7.459885052032535 5.60959768211387)]
For data set: [(5 2 −2.26) (0.77 2 4) (1 4 0) 4]
Solution(s) is/are: [(4.2615520237084015 2 -1.1671668246843006) (1.13386504516801 2 3.461514141193441)]

zkl

Translation of: Go

<lang zkl>const EPS=1e-14; // a close-ness to zero

   // p1,p2 are (x,y), circle is ( (x,y),r )

fcn intersectLineCircle(p1,p2, circle, segment=False) // assume line {

  cx,cy := circle[0].apply("toFloat");
  r     := circle[1].toFloat();
  x1,y1 := p1.apply("toFloat");   x2,y2 := p2.apply("toFloat");
  A,B,C,a := (y2 - y1),  (x1 - x2),  (x2*y1 - x1*y2),  (A*A + B*B);
  b,c,bnz := 0,0,True;
  if(B.closeTo(0,EPS)){ // B is zero or close to it
     b   = 2.0 * (B*C + A*B*cx - A*A*cy);
     c   = C*C + 2.0*A*C*cx - A*A*(r*r - cx*cx - cy*cy);
     bnz = False
  }else{
     b   = 2.0*( A*C + A*B*cy - B*B*cx );
     c   = C*C + 2.0*B*C*cy - B*B*( r*r - cx*cx - cy*cy );
  }
  d := b*b - 4.0*a*c;		// discriminant
  if(d<0.0){			// no real solution?  zero --> one solution
     if (d>-0.005) d=0.0;	// close enough to zero
     else return(T);		// no intersection
  }
  d=d.sqrt();
  reg ux,uy, vx,vy;
  if(bnz){
     ux,vx = (-b + d) / (2*a),  (-b - d) / (2*a);
     uy,vy = -(A*ux + C) / B,   -(A*vx + C) / B;
  }else{
     uy,vy = (-b + d) / (2*a), (-b - d) / (2*a);
     ux,vx = -(B*uy + C) / A,  -(B*vy + C) / A;
  }
  if(segment){
     within:='wrap(x,y){  // is (x,y) within segment p1 p2?

d1:=( (x2 - x1).pow(2) + (y2 - y1).pow(2) ).sqrt(); d2:=( (x - x1).pow(2) + (y - y1).pow(2) ).sqrt(); d3:=( (x2 - x) .pow(2) + (y2 - y) .pow(2) ).sqrt(); (d1 - d2 - d3).closeTo(0,EPS);

     };
     i1,i2 := within(ux,uy), within(vx,vy);
     if(d==0) return(if(i1) T(ux,uy) else T);
     return(T( i1 and T(ux,uy), i2 and T(vx,vy) ).filter())
  }
  if(d==0)  return( T( T(ux,uy) ) );
  return( T(ux,uy), T(vx,vy) )

}</lang> <lang zkl>circle:=T( T(3,-5),3 ); p1,p2 := T(-10,11), T( 10,-9); println("Circle @ ",circle); lcpp(p1,p2,circle); p2:=T(-11,12); lcpp(p1,p2,circle,True); p1,p2 := T(3,-2), T(7,-2); lcpp(p1,p2,circle);

circle:=T( T(0,0),4 ); p1,p2 := T(0,-3), T(0,6); println("\nCircle @ ",circle); lcpp(p1,p2,circle); lcpp(p1,p2,circle,True);

circle:=T( T(4,2),5 ); p1,p2 := T(6,3), T(10,7); println("\nCircle @ ",circle); lcpp(p1,p2,circle); p1,p2 := T(7,4), T(11,8); lcpp(p1,p2,circle,True);

fcn lcpp(p1,p2,circle,segment=False){

  println("  %s %s -- %s intersects at %s"
     .fmt(segment and "Segment" or "Line   ",
          p1,p2,intersectLineCircle(p1,p2, circle,segment)));

}</lang>

Output:
Circle @ L(L(3,-5),3)
  Line    L(-10,11) -- L(10,-9) intersects at L(L(6,-5),L(3,-2))
  Segment L(-10,11) -- L(-11,12) intersects at L()
  Line    L(3,-2) -- L(7,-2) intersects at L(L(3,-2))

Circle @ L(L(0,0),4)
  Line    L(0,-3) -- L(0,6) intersects at L(L(0,4),L(0,-4))
  Segment L(0,-3) -- L(0,6) intersects at L(L(0,4))

Circle @ L(L(4,2),5)
  Line    L(6,3) -- L(10,7) intersects at L(L(8,5),L(1,-2))
  Segment L(7,4) -- L(11,8) intersects at L(L(8,5))