Problem:
Tangent to a Circle Find the tangent lines between a point $p$ and a circle $c$. Input The input is given in the following format. $px \; py$ $cx \; cy \; r$ $px$ and $py$ represents the coordinate of the point $p$. $cx$, $cy$ and $r$ represents the center coordinate and radius of the circle $c$ respectively. All input values are given in integers. Output Print coordinates of the tangent points on the circle $c$ based on the following rules.  Print the coordinate with smaller $x$ first. In case of a tie, print the coordinate with smaller $y$ first.  The output values should be in a decimal fraction with an error less than 0.00001. Constraints  $-1,000 \leq px, py, cx, cy \leq 1,000$  $1 \leq r \leq 1,000$  Distance between $p$ and the center of $c$ is greater than the radius of $c$.  Sample Input 1 0 0
2 2 2 Sample Output 1 0.0000000000 2.0000000000
2.0000000000 0.0000000000 Sample Input 2 -3 0
2 2 2 Sample Output 2 0.6206896552 3.4482758621
2.0000000000 0.0000000000

Solution (Python):
px, py = map(int, input().split())
cx, cy, r = map(int, input().split())
dx, dy = px - cx, py - cy
s2 = dx ** 2 + dy ** 2
df = (s2 - r ** 2) ** 0.5
ix1, ix2 = (dx * r + dy * df) * r / s2, (dx * r - dy * df) * r / s2
iy1, iy2 = (dy * r - dx * df) * r / s2, (dy * r + dx * df) * r / s2
if (ix1, iy1) > (ix2, iy2):
    ix1, ix2 = ix2, ix1
    iy1, iy2 = iy2, iy1
print(ix1 + cx, iy1 + cy)
print(ix2 + cx, iy2 + cy)

Solution (Python):
from math import sqrt
x1, y1 = map(int, raw_input().split())
x0, y0, r0 = map(int, raw_input().split())

def solve(x1, y1, r1, x2, y2, r2):
    rr0 = (x2 - x1)**2 + (y2 - y1)**2
    xd = x2 - x1
    yd = y2 - y1
    rr1 = r1**2; rr2 = r2**2
    cv = (rr0 + rr1 - rr2)
    sv = sqrt(4*rr0*rr1 - cv**2)
    return (
        (x1 + (cv*xd - sv*yd)/(2.*rr0), y1 + (cv*yd + sv*xd)/(2.*rr0)),
        (x1 + (cv*xd + sv*yd)/(2.*rr0), y1 + (cv*yd - sv*xd)/(2.*rr0)),
    )

r1 = sqrt((x1 - x0)**2 + (y1 - y0)**2 - r0**2)

p0, p1 = sorted(solve(x1, y1, r1, x0, y0, r0))
print "%.08f %.08f" % p0
print "%.08f %.08f" % p1

Solution (Python):
from math import pi, cos, sin, atan2
EPS = 10**(-9)

def eq(value1, value2):
    return abs(value1-value2) <= EPS

class Point():
    def __init__(self, x, y):
        self.x = x
        self.y = y
        self.arg = atan2(y, x) 
    
    def __str__(self):
        return "{0:.8f} {1:.8f}".format(self.x, self.y)
    
    def __add__(self, other):
        return Point(self.x + other.x, self.y + other.y)
    
    def __sub__(self, other):
        return Point(self.x - other.x, self.y - other.y)
    
    def __mul__(self, scal):
        return Point(self.x*scal, self.y*scal)
    
    def __truediv__(self, scal):
        return Point(self.x/scal, self.y/scal)
    
    def __eq__(self, other):
        return eq(self.x, other.x) and eq(self.y, other.y)

    
    def __abs__(self):
        return (self.x**2+self.y**2)**0.5
    

def Rotation(vec: Point, rad):
    return Point(vec.x*cos(rad)-vec.y*sin(rad), vec.x*sin(rad)+vec.y*cos(rad))


class Circle():
    def __init__(self, p, r):
        self.p = p
        self.r = r


class Line():
    
    def __init__(self, a, b):
        self.a = a
        self.b = b
        self.arg = (a-b).arg % pi
    
    def __str__(self):
        return "[({0}, {1}) - ({2}, {3})]".format(self.a.x, self.a.y, self.b.x, self.b.y)

    
    def par(self, point):
        return Line(point, point+(self.a-self.b))

    
    def tan(self, point):
        return Line(point, point + Rotation(self.a-self.b, pi/2))


class Segment(Line):
    def __init__(self, a, b):
        super().__init__(a, b)



def cross(vec1: Point, vec2: Point):
    return vec1.x*vec2.y - vec1.y*vec2.x


def dot(vec1: Point, vec2: Point):
    return vec1.x*vec2.x + vec1.y*vec2.y


def ccw(a, b, c):
    if cross(b-a, c-a) > EPS: return +1 
    if cross(b-a, c-a) < -EPS: return -1 
    if dot(c-a, b-a) < -EPS: return +2 
    if abs(b-a) < abs(c-a): return -2 
    return 0 



def projection(l, p):
    t = dot(l.b-l.a, p-l.a) / abs(l.a-l.b)**2
    return l.a + (l.b-l.a)*t


def reflection(l, p):
    return p + (projection(l, p) - p)*2

def isPararell(l1, l2):
    return eq(cross(l1.a-l1.b, l2.a-l2.b), 0)

def isVertical(l1, l2):
    return eq(dot(l1.a-l1.b, l2.a-l2.b), 0)


def isIntersect_lp(l, p):
    return abs(ccw(l.a, l.b, p)) != 1

def isIntersect_ll(l1, l2):
    return not isPararell(l1, l2) or isIntersect_lp(l1, l2.a)

def isIntersect_sp(s, p):
    return ccw(s.a, s.b, p) == 0

def isIntersect_ss(s1, s2):
    return ccw(s1.a, s1.b, s2.a)*ccw(s1.a, s1.b, s2.b) <= 0 and ccw(s2.a, s2.b, s1.a)*ccw(s2.a, s2.b, s1.b) <= 0

def isIntersect_ls(l, s):
    return cross(l.b - l.a, s.a - l.a) * cross(l.b - l.a, s.b - l.a) < EPS

def isIntersect_cp(c, p):
    return abs(abs(c.p - p) - c.r) < EPS

def isIntersect_cl(c, l):
    return distance_lp(l, c.p) <= c.r + EPS

def isIntersect_cs(c, s):
    pass

def intersect_cc(c1, c2):
    if c1.r < c2.r:
        c1, c2 = c2, c1
    d = abs(c1.p - c2.p)
    if eq(c1.r + c2.r, d): return 3 
    if eq(c1.r - c2.r, d): return 1 
    if c1.r + c2.r < d: return 4 
    if c1.r - c2.r < d: return 2 
    return 0 


def distance_pp(p1, p2):
    return abs(p1-p2)

def distance_lp(l, p):
    return abs(projection(l,p)-p)

def distance_ll(l1, l2):
    return 0 if isIntersect_ll(l1, l2) else distance_lp(l1, l2.a)

def distance_sp(s, p):
    r = projection(s, p)
    if isIntersect_sp(s, r): return abs(r-p)
    return min(abs(s.a-p), abs(s.b-p))

def distance_ss(s1, s2):
    if isIntersect_ss(s1, s2): return 0
    return min([distance_sp(s1, s2.a), distance_sp(s1, s2.b), distance_sp(s2, s1.a), distance_sp(s2, s1.b)])

def distance_ls(l, s):
    if isIntersect_ls(l, s): return 0
    return min(distance_lp(l, s.a), distance_lp(l, s.b))


def crosspoint_ll(l1, l2):
    A = cross(l1.b - l1.a, l2.b - l2.a)
    B = cross(l1.b - l1.a, l1.b - l2.a)
    if eq(abs(A), 0) and eq(abs(B), 0): return l2.a
    return l2.a + (l2.b - l2.a) * B / A

def crosspoint_ss(s1, s2):
    return crosspoint_ll(s1, s2)

def crosspoint_lc(l, c):
    p = projection(l, c.p)
    if eq(distance_lp(l, c.p), c.r): return [p]
    e = (l.b - l.a) / abs(l.b-l.a)
    dis = (c.r**2-abs(p-c.p)**2)**0.5
    return [p + e*dis, p - e*dis]

def crosspoint_sc(s, c):
    pass

def crosspoint_cc(c1, c2):
    d = abs(c1.p-c2.p)
    if not abs(c1.r-c2.r) <= d <= c1.r+c2.r:
        return []
    mid_p = (c2.p * (c1.r**2-c2.r**2+d**2) + c1.p * (c2.r**2-c1.r**2+d**2)) / (2*d**2)
    tanvec = Rotation(c1.p-c2.p, pi/2)
    return crosspoint_lc(Line(mid_p, mid_p+tanvec), c1)



def tangent_cp(c, p):
    return crosspoint_cc(c, Circle(p, (abs(p-c.p)**2 - c.r**2)**0.5))





import sys
input = sys.stdin.readline

def verify_1A():
    p1x, p1y, p2x, p2y = map(int, input().split())
    l = Line(Point(p1x, p1y), Point(p2x, p2y))
    Q = int(input())
    Query = [list(map(int, input().split())) for _ in range(Q)]
    for px, py in Query:
        p = Point(px, py)
        print(projection(l, p))

def verify_1B():
    p1x, p1y, p2x, p2y = map(int, input().split())
    l = Line(Point(p1x, p1y), Point(p2x, p2y))
    Q = int(input())
    Query = [list(map(int, input().split())) for _ in range(Q)]
    for px, py in Query:
        p = Point(px, py)
        print(reflection(l, p))

def verify_1C():
    p1x, p1y, p2x, p2y = map(int, input().split())
    p1 = Point(p1x, p1y); p2 = Point(p2x, p2y)
    Q = int(input())
    Query = [list(map(int, input().split())) for _ in range(Q)]
    for px, py in Query:
        p = Point(px, py)
        result = ccw(p1, p2, p)
        if result == 1:
            print("COUNTER_CLOCKWISE")
        elif result == -1:
            print("CLOCKWISE")
        elif result == 2:
            print("ONLINE_BACK")
        elif result == -2:
            print("ONLINE_FRONT")
        else:
            print("ON_SEGMENT")

def verify_2A():
    Q = int(input())
    Query = [list(map(int, input().split())) for _ in range(Q)]
    for p0x, p0y, p1x, p1y, p2x, p2y, p3x, p3y in Query:
        l1 = Line(Point(p0x, p0y), Point(p1x, p1y))
        l2 = Line(Point(p2x, p2y), Point(p3x, p3y))
        if isPararell(l1, l2):
            print(2)
        elif isVertical(l1, l2):
            print(1)
        else:
            print(0)

def verify_2B():
    Q = int(input())
    Query = [list(map(int, input().split())) for _ in range(Q)]
    for p0x, p0y, p1x, p1y, p2x, p2y, p3x, p3y in Query:
        s1 = Segment(Point(p0x, p0y), Point(p1x, p1y))
        s2 = Segment(Point(p2x, p2y), Point(p3x, p3y))
        if isIntersect_ss(s1, s2):
            print(1)
        else:
            print(0)


def verify_2C():
    Q = int(input())
    Query = [list(map(int, input().split())) for _ in range(Q)]
    for p0x, p0y, p1x, p1y, p2x, p2y, p3x, p3y in Query:
        s1 = Segment(Point(p0x, p0y), Point(p1x, p1y))
        s2 = Segment(Point(p2x, p2y), Point(p3x, p3y))
        print(crosspoint_ss(s1, s2))  

def verify_2D():
    Q = int(input())
    Query = [list(map(int, input().split())) for _ in range(Q)]
    for p0x, p0y, p1x, p1y, p2x, p2y, p3x, p3y in Query:
        s1 = Segment(Point(p0x, p0y), Point(p1x, p1y))
        s2 = Segment(Point(p2x, p2y), Point(p3x, p3y))
        print("{:.8f}".format(distance_ss(s1, s2)))

def verify_7A():
    c1x, c1y, c1r = map(int, input().split())
    c2x, c2y, c2r = map(int, input().split())
    print(intersect_cc(Circle(Point(c1x, c1y), c1r), Circle(Point(c2x, c2y), c2r)))

def verify_7D():
    cx, cy, cr = map(int, input().split())
    c = Circle(Point(cx, cy), cr)
    Q = int(input())
    Query = [list(map(int, input().split())) for _ in range(Q)]
    for x1, y1, x2, y2 in Query:
        Points = crosspoint_lc(Line(Point(x1, y1), Point(x2,y2)), c)
        if len(Points) == 1:
            Points.append(Points[-1])
        Points.sort(key=lambda p: p.y)
        Points.sort(key=lambda p: p.x)
        print(*Points)

def verify_7E():
    c1x, c1y, c1r = map(int, input().split())
    c1 = Circle(Point(c1x, c1y), c1r)
    c2x, c2y, c2r = map(int, input().split())
    c2 = Circle(Point(c2x, c2y), c2r)

    Points = crosspoint_cc(c1, c2)
    if len(Points) == 1:
        Points.append(Points[-1])
    Points.sort(key=lambda p: p.y)
    Points.sort(key=lambda p: p.x)
    print(*Points)  

def verify_7F():
    px, py = map(int, input().split())
    cx, cy, cr = map(int, input().split())
    Points = tangent_cp(Circle(Point(cx, cy), cr), Point(px, py))
    if len(Points) == 1:
        Points.append(Points[-1])
    Points.sort(key=lambda p: p.y)
    Points.sort(key=lambda p: p.x)
    print(*Points, sep="\n")    

verify_7F()