Sage code to plot a conformal map

triangle-conformal-map.sage

# Plot the conformal mapping taking the unit disk to the doubly ideal
# triangle having vertices 0, 1, i in the Poincaré disk model (with 1,
# i, -i being mapped to these three vertices).  See
# https://mathoverflow.net/questions/495782/riemann-mapping-to-a-specific-curvilinear-triangle/495809#495809
# for details.

k = 2 * gamma(3/4)^2 / gamma(1/4)^2
kN = N(k)
var('z')
s(z) = k*sqrt(z) * hypergeometric([3/4,3/4],[3/2],z) / hypergeometric([1/4,1/4],[1/2],z)
ss(z) = s(1-I-2/(z+I))

# A memoized version of the above function, for efficiency's sake:
ss_tab = {}
def memoize_ss(z):
    z = N(z)
    if z in ss_tab:
        return ss_tab[z]
    w = N(1-I-2/(z+I))
    ss_tab[z] = N(kN*sqrt(w) * hypergeometric([3/4,3/4],[3/2],w) / hypergeometric([1/4,1/4],[1/2],w))
    return ss_tab[z]

def plothelp(rad, ang):
    return memoize_ss(rad*exp(2*I*pi*ang))

def plothelp_x(rad, ang):
    return real(plothelp(rad,ang))

def plothelp_y(rad, ang):
    return imag(plothelp(rad,ang))

circplots = []
for i in range(20):
    print("computing circular plot at radius %d/20 (%d values computed so far)"%(i,len(ss_tab)))
    circplots.append(parametric_plot([lambda t: plothelp_x(i/20, t), lambda t: plothelp_y(i/20, t)], (0,1)))

radplots = []
for i in range(48):
    print("computing radial plot at angle %d/48 (%d values computed so far)"%(i,len(ss_tab)))
    radplots.append(parametric_plot([lambda t: plothelp_x(t, i/48), lambda t: plothelp_y(t, i/48)], (0,19/20)))

redplot = parametric_plot([lambda t: plothelp_x(99/100, t), lambda t: plothelp_y(99/100, t)], (0,1), color="red", plot_points=960)
pinkplot = parametric_plot([lambda t: plothelp_x(999/1000, t), lambda t: plothelp_y(999/1000, t)], (0,1), color="magenta", plot_points=3840)

fullplot = sum(circplots + radplots + [redplot, pinkplot])
fullplot.save(filename="plot.png", dpi=300, aspect_ratio=1)