# 2D-kromme

import numpy as np

J = 1.7321
Rs = 1

h1 = 1 - 3 * Rs**2 / J**2
h2 = np.sqrt (h1)
h3 = np.sqrt (h2)

n1 = (1 + 2 * h2) / (3 * Rs)

rbegin = 0.999999999 / n1
reinde = 0.999999999
aantalstappen = 250000000
stap = (reinde - rbegin) / aantalstappen

phi0 = np.pi / 2
r = rbegin - stap

while r > reinde:

 r = r + stap
 phi = 2 / h3 * (phi0 - np.arcsin (h3 / np.sqrt (Rs / r - 1 / 3 + h2 / 3)))
 x = r * np.cos (phi)
 y = r * np.sin (phi)

 print ('{:10.5e}'.format (x), '{:10.5e}'.format (y))