I've written a function (in python 3) which computes if a line segment (constraint) and a circle (body) collide, and returns the point of intersection (closest point to the centre of the circle):

def constraintCollide(self, constraint):

    p1 = constraint.point1 # Sets the first point object of the constraint to p1
    p2 = constraint.point2 # Sets the second point object of the constraint to p2

    if self.p.distance_to(p1.p) - p1.p.distance_to(p2.p) > self.r \
    or self.p.distance_to(p2.p) - p1.p.distance_to(p2.p) > self.r:
        return False 
    ''' Checks if difference in distance from the radius to a point and the lenth of the line segment is greater than the radius of the circle.
        If this is so, then the line segment should be nowhere near intercepting and the function should return false'''

    #Calculating the gradient m of the line segment:
    if p2.p.x - p1.p.x == 0: # In case of no difference in x cooridnates,
        m = float('inf') # Gradient is infinity
    elif p2.p.y - p1.p.y == 0: # In as of no difference in y coordinates,
        m = 0 # Gradient is zero
    else: # Otherwise,
        m = (p2.p.y - p1.p.y)/(p2.p.x - p1.p.x) # Normal gradient camculation

    # Calculating the point of interception (if any):
    if p1.p.distance_to(self.p) <= self.r: # Checks if relative distance from point 1 is less that the radius,
        i = p1.p # In with case, 'intercection' is at the edge of the line segment,
    elif  p2.p.distance_to(self.p) <= self.r: # And likewise for point 2
        i = p2.p # "
    else: # Otherwise, solve quadratic equation for the point of inteption 
        a = self.p.x; b = self.p.y; c = p1.p.y -m*p1.p.x; r = self.r;
        # Defines centre of circle as (a, b), y-itercept of 'line' as c, and the radius as r
        x1 = (-(-a**2*m**2+2*a*b*m-2*a*c*m-b**2+2*b*c-c**2+m**2*r**2+r**2)**0.5 +a+b*m-c*m)/(m**2+1)
        x2 = ((-a**2*m**2+2*a*b*m-2*a*c*m-b**2+2*b*c-c**2+m**2*r**2+r**2)**0.5 +a+b*m-c*m)/(m**2+1)
        y1 = (-m*(-a**2*m**2+2*a*b*m-2*a*c*m-b**2+2*b*c-c**2+m**2*r**2+r**2)**0.5+a*m+b*m**2+c)/(m**2+1)
        y2 = (m*(-a**2*m**2+2*a*b*m-2*a*c*m-b**2+2*b*c-c**2+m**2*r**2+r**2)**0.5+a*m+b*m**2+c)/(m**2+1)
        ''' As this is a quadratic, there is a possiblitiy of two real roots (points of inteception), which 
            is why the x and y are computated twice to acomodate the (+/-) sign in a quadratic'''

        x = (x1.real+x2.real)/2 # Average x cooridnate of the real solutions (As python 3 returns **.5 as complex number)
        y = (y1.real+y2.real)/2 # Average y cooridantes of the real solutions ''
        i = pygame.math.Vector2(x1.real, y1.real) # Stores 'the point of interseption' as a pygame vector object 
        # In actuality, i is not  the point of intersection, but the 'closest point to the centre' (I hope)

    if i.distance_to(self.p) > self.r: # Checks if the closest point to the circle is larger than the circle radius,
        return False, # In which case the function returns false

    # Checks if the combined distance from the closest point i to each point is longer than the line segnent,
    if (i-p1.p).length() + (i-p2.p).length() > p1.p.distance_to(p2.p): 
        return False # In which case return false

    # Line segment has passed all tests so, 
    return i # Returns i

The function surprisingly works as expected, with a success rate of I'd say 99%. The problem is that it seems very expensive for its purpose. My guess for the cause would be the solving of the quadratic with the lengthy formula.

Would there be a way to optimize this algorithm, or a different approach altogether with the ability to solve for two points of interception of a circle of known radius and a line segment of known endpoints?

Any suggestions would be greatly appreciated...