diff --git a/001.jl b/001.jl
new file mode 100644
index 0000000..ab9da01
--- /dev/null
+++ b/001.jl
@@ -0,0 +1,3 @@
+f(c) = c % 3 == 0 || c % 5 == 0
+
+println(sum(filter(f, 0:999)))
diff --git a/727.py b/727.py
index 76f9edb..88d8fdd 100644
--- a/727.py
+++ b/727.py
@@ -1,7 +1,7 @@
-import sympy
 import joblib
 import tqdm
 import multiprocessing
+from math import *
 
 def gcd(a, b):
     if b == 0: return a
@@ -21,56 +21,59 @@ class Point:
     def distance_to(self, other):
         dx = other.x - self.x
         dy = other.y - self.y
-        return sympy.sqrt(dx * dx + dy * dy)
+        return sqrt(dx * dx + dy * dy)
 
 def circumcenter(A: Point, B: Point, C: Point) -> Point:
-    D = sympy.S(2) * (A.x * (B.y - C.y) + B.x * (C.y - A.y) + C.x * (A.y - B.y))
+    D = 2.0 * (A.x * (B.y - C.y) + B.x * (C.y - A.y) + C.x * (A.y - B.y))
     Ux = (A.magnitude() * (B.y - C.y) + B.magnitude() * (C.y - A.y) + C.magnitude() * (A.y - B.y)) / D
     Uy = (A.magnitude() * (C.x - B.x) + B.magnitude() * (A.x - C.x) + C.magnitude() * (B.x - A.x)) / D
     return Point(Ux, Uy)
 
 def lerp(src: Point, dst: Point, dist: float) -> Point:
     total = src.distance_to(dst)
-    ratio = sympy.Rational(dist, total)
+    ratio = dist / total
     dx = src.x + (dst.x - src.x) * ratio
     dy = src.y + (dst.y - src.y) * ratio
     return Point(dx, dy)
 
 def heron(a, b, c):
-    a, b, c = sympy.S(a), sympy.S(b), sympy.S(c)
-    p = (a + b + c) / sympy.S(2)
-    return sympy.sqrt(p * (p - a) * (p - b) * (p - c))
+    p = (a + b + c) / 2.0
+    return sqrt(p * (p - a) * (p - b) * (p - c))
 
 def inv_law_cosines(a, b, c):
     # given a, b, c return C
-    return sympy.acos((a * a + b * b - c * c) / (sympy.S(2) * a * b))
+    return acos((a * a + b * b - c * c) / (2.0 * a * b))
 
 def compute_d(A: Point, B: Point, C: Point, r1: float, r2: float, r3: float) -> float:
-    r1, r2, r3 = sympy.S(r1), sympy.S(r2), sympy.S(r3)
-
     ABmid = lerp(A, B, r1)
     BCmid = lerp(B, C, r2)
     CAmid = lerp(C, A, r3)
     D = circumcenter(ABmid, BCmid, CAmid)
 
-    result = r1 + r2 + r3 + sympy.S(2) * sympy.sqrt(r1 * r2 + r2 * r3 + r1 * r3)
+    k1, k2, k3 = 1.0 / r1, 1.0 / r2, 1.0 / r3
+    e1 = k1 + k2 + k3 + 2.0 * sqrt(k1 * k2 + k2 * k3 + k1 * k3)
+    e2 = k1 + k2 + k3 - 2.0 * sqrt(k1 * k2 + k2 * k3 + k1 * k3)
+    result = 1.0 / max(e1, e2)
     # re = sympy.symbols("re")
+
+    # re = result
     # t1 = heron(r1 + r3, r1 + re, r3 + re)
     # t2 = heron(r1 + r2, r1 + re, r2 + re)
     # t3 = heron(r2 + r3, r2 + re, r3 + re)
     # tfull = heron(r1 + r2, r1 + r3, r2 + r3)
+    # print(t1 + t2 + t3, tfull, "=", tfull - t1 - t2 - t3)
     # eq = sympy.Eq(t1 + t2 + t3, tfull)
     # result = sympy.solve(eq, re)
     # result = result[0].evalf()
     # print(f"radius of E is {result}")
 
     EAB = inv_law_cosines(r1 + result, r1 + r2, r2 + result)
-    Ex = (r1 + result) * sympy.cos(EAB)
-    Ey = (r1 + result) * sympy.sin(EAB)
+    Ex = (r1 + result) * cos(EAB)
+    Ey = (r1 + result) * sin(EAB)
     E = Point(Ex, Ey)
 
-    # print(D, E)
-    return D.distance_to(E).evalf()
+    # print(r1, r2, r3, D, E)
+    return D.distance_to(E)
 
 total = 0
 count = 0
@@ -88,22 +91,23 @@ def calc(x):
     # B is at (c, 0), b = r1 + r3
     # c = r1 + r2
     # C is at ()
-    a = sympy.S(r2 + r3)
-    b = sympy.S(r1 + r3)
-    c = sympy.S(r1 + r2)
+    a = r2 + r3
+    b = r1 + r3
+    c = r1 + r2
 
-    A = Point(sympy.S(0), sympy.S(0))
-    B = Point(c, sympy.S(0))
+    A = Point(0, 0)
+    B = Point(c, 0)
     
-    CAB = sympy.acos(sympy.Rational(b * b + c * c - a * a, 2.0 * b * c))
-    Cx = b * sympy.cos(CAB)
-    Cy = b * sympy.sin(CAB)
+    CAB = inv_law_cosines(b, c, a)
+    Cx = b * cos(CAB)
+    Cy = b * sin(CAB)
     C = Point(Cx, Cy)
 
     d = compute_d(A, B, C, r1, r2, r3)
     return d
-    print(f"d({r1}, {r2}, {r3}) = {d}")
+    # return (f"d({r1}, {r2}, {r3}) = {d}")
 
 total = len(list(create_loop()))
 L = joblib.Parallel(n_jobs=multiprocessing.cpu_count())(joblib.delayed(calc)(tup) for tup in tqdm.tqdm(create_loop(), total=total))
-print(sum(L) / len(L))
+# print(L)
+print(round(sum(L) / len(L), 8))