Skip to content

Instantly share code, notes, and snippets.

@gagern
Created August 7, 2013 02:22
Show Gist options
  • Select an option

  • Save gagern/6170668 to your computer and use it in GitHub Desktop.

Select an option

Save gagern/6170668 to your computer and use it in GitHub Desktop.
Details on how the number in http://math.stackexchange.com/a/461338/35416 was obtained
sage: def angle_bisectors(g, h, p = None):
... if p is None:
... p = g.cross_product(h)
... a = (g[0]*h[1] + g[1]*h[0])
... b = 2*(g[0]*h[0] - g[1]*h[1])
... c = -a
... d = sqrt(b^2 - 4*a*c)
... r = [b + d, b - d]
... r = [vector(SR, [i, 2*a, 0]) for i in r]
... r = [i.cross_product(p) for i in r]
... return r
sage: def dist(p, q):
... x = p[0]/p[2] - q[0]/q[2]
... y = p[1]/p[2] - q[1]/q[2]
... return sqrt(x^2 + y^2)
sage: A = vector(SR, [-1, 1, 1])
sage: B = vector(SR, [0, 0, 1])
sage: C = vector(SR, [1, 0, 0])
sage: AB = A.cross_product(B)
sage: AC = A.cross_product(C)
sage: BC = B.cross_product(C)
sage: bisA = angle_bisectors(AB, AC, A)
sage: for i in bisA:
... print(N(-i[0]/i[1]))
2.41421356237309
-0.414213562373095
sage: AD = bisA[1]
sage: D = BC.cross_product(AD)
sage: bisD = angle_bisectors(BC, AD, D)
sage: for i in bisD:
... print(N(-i[0]/i[1]))
5.02733949212585
-0.198912367379658
sage: DN, DM = bisD
sage: M = AB.cross_product(DM)
sage: NN = AC.cross_product(DN)
sage: BM = dist(B, M)
sage: MN = dist(M, NN)
sage: ratio = QQbar(MN/BM)
sage: ratio
4.165617886885153?
sage: ratio.minpoly()
x^8 - 32*x^6 + 304*x^4 - 896*x^2 + 544
sage: ratios = [sqrt(i.rhs()) for i in solve(x^4 - 32*x^3 + 304*x^2 - 896*x + 544, x)]
sage: [i.N() for i in ratios]
[0.901339015602487, 1.96821733555214, 3.15615844333697, 4.16561788688515]
sage: ratios[-1].simplify()
sqrt(1/2*sqrt(16*sqrt(2) + 32) + 4*sqrt(2) + 8)
sage: (sqrt(sqrt(4*sqrt(2)+8)+4*sqrt(2)+8) - ratio).is_zero()
True
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment