import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
import ipywidgets as widgets
from tqdm import tqdm
def binary_expansion(x: float, nmax=10):
assert x>=0.0 and x<=1.0
coeffs = []
for i in range(1,nmax+1):
binary_term = 2**i
coeffs.append( int(x * binary_term) )
x -= coeffs[-1]/binary_term
return coeffs
# for debugging
"""def binary_approximation(coeffs):
x=0.0
for i in range(len(coeffs)):
x += coeffs[i] / 2**(i+1)
return x"""
# define function composition with one argument
def compose (*functions):
def inner(arg):
for f in reversed(functions):
arg = f(arg)
return arg
return inner
# De Rham map: [0,1] -> Complex plane
def c(d0, d1, x, a, nmax=10):
assert x>=0.0 and x<=1.0
global binary_expansion
# find the binary expansion of x
bs = binary_expansion(x, nmax=nmax)
# compose the functions d0() and d1() accordingly
functions = [d1 if b else d0 for b in bs]
composed = compose(*functions)
return composed(x)
def render(N=1000, a_real=0.5, a_imag=0.5, nmax=10):
a = a_real + 1j*a_imag
def d0(z, a=a):
assert abs(a)<1. and abs(1.-a)<1.
return a*z
def d1(z, a=a):
assert abs(a)<1. and abs(1.-a)<1.
return a + (1-a)*z
xs = np.linspace(0.0, 1.0, N)
px = []
py = []
for x in tqdm(xs, desc='Calculating points'):
C = c(d0, d1, x, a, nmax)
px.append( C.real )
py.append( C.imag )
fig,ax = plt.subplots(figsize=(8,8))
ax.axis('equal')
ax.set_title(f"Cesàro curve: N={N}, nmax={nmax}, a={a_real:.2f}+{a_imag:.2f}i")
ax.scatter(px,py, s=1, c='r', marker='.', alpha=0.7)
plt.show()
dict_args = {
"N": widgets.IntSlider(min=100, max=1e5, step=100, value=1000),
"a_real": widgets.FloatSlider(min=-1.0, max=1.0, step=0.01, value=0.5),
"a_imag": widgets.FloatSlider(min=-1.0, max=1.0, step=0.01, value=0.5),
"nmax": widgets.IntSlider(min=1, max=30, step=1, value=10)
}
widgets.interact_manual(render, **dict_args);
def render(N=1000, a_real=0.5, a_imag=0.5, nmax=10):
a = a_real + 1j*a_imag
def d0(z, a=a):
assert abs(a)<1. and abs(1.-a)<1.
return a * z.conjugate()
def d1(z, a=a):
assert abs(a)<1. and abs(1.-a)<1.
return a + (1-a)*z.conjugate()
xs = np.linspace(0.0, 1.0, N)
px = []
py = []
for x in tqdm(xs, desc='Calculating points'):
C = c(d0, d1, x, a, nmax)
px.append( C.real )
py.append( C.imag )
fig,ax = plt.subplots(figsize=(8,8))
ax.axis('equal')
ax.set_title(f"Koch-Peano curve: N={N}, nmax={nmax}, a={a_real:.2f}+{a_imag:.2f}i")
ax.scatter(px,py, s=1, c='r', marker='.', alpha=0.7)
plt.show()
dict_args = {
"N": widgets.IntSlider(min=100, max=1e5, step=100, value=1000),
"a_real": widgets.FloatSlider(min=-1.0, max=1.0, step=0.01, value=0.6),
"a_imag": widgets.FloatSlider(min=-1.0, max=1.0, step=0.01, value=0.37),
"nmax": widgets.IntSlider(min=1, max=30, step=1, value=10)
}
widgets.interact_manual(render, **dict_args);
"""
def render(N=1000, a_real=0.5, a_imag=0.5, nmax=10):
a = a_real + 1j*a_imag
# here we do not use a.
def d0(z, a=a):
return z / (z+1)
def d1(z, a=a):
return 1. / (2.-z)
xs = np.linspace(0.0, 1.0, N)
px = []
py = []
for x in tqdm(xs, desc='Calculating points'):
C = c(d0, d1, x, a, nmax)
px.append( C.real )
py.append( C.imag )
fig,ax = plt.subplots(figsize=(8,8))
ax.axis('equal')
ax.set_title(f"Minkowski's ?(x): N={N}, nmax={nmax}, a={a_real:.2f}+{a_imag:.2f}i")
ax.scatter(px,py, s=1, c='r', marker='.', alpha=0.7)
plt.show()
dict_args = {
"N": widgets.IntSlider(min=100, max=1e5, step=100, value=1000),
"a_real": widgets.FloatSlider(min=0, max=0, step=0.01, value=0),
"a_imag": widgets.FloatSlider(min=0, max=0, step=0.01, value=0),
"nmax": widgets.IntSlider(min=1, max=30, step=1, value=10)
}
widgets.interact_manual(render, **dict_args);
"""
They apply to the 2D real plane
# General affine map: [0,1] -> 2D plane
def c_affine(d0, d1, x:float, a,b,d,e,f,g, nmax=10):
global binary_expansion
# find the binary expansion of x
bs = binary_expansion(x, nmax=nmax)
# compose the functions d0() and d1() accordingly
functions = [d1 if b else d0 for b in bs]
composed = compose(*functions)
X = np.array([1.,x,x]).reshape(3,1) #column vector
return composed(X)
# rendering function, implementing the algorithm
def render_affine(N=1000, a=0.5, b=1.0, d=0.5, e=0.5, f=0.5, g=0.5, nmax=10):
def d0(x, a=a,b=b,d=d,e=e):
M = np.matrix([
[1., 0., 0.],
[0., a , d ],
[0., b , e ]
])
return M @ x
def d1(x, a=a,b=b,f=f,g=g):
M = np.matrix([
[1., 0., 0.],
[a ,1.-a,f],
[b , -b, g]
])
return M @ x
xs = np.linspace(0.0, 1.0, N)
px = []
py = []
for x in tqdm(xs, desc='Calculating'):
C = c_affine(d0, d1, x, a,b,d,e,f,g, nmax) # shape (3,1)
px.append( C[1,0] )
py.append( C[2,0] )
fig,ax = plt.subplots(figsize=(8,8))
ax.axis('equal')
ax.set_title(f"General affine map: N={N}, nmax={nmax},\n"
+"a,b,d,e,f,g=%.2f,%.2f,%.2f,%.2f,%.2f,%.2f"%(a,b,d,e,f,g))
ax.scatter(px, py, s=1, c='r', marker='.', alpha=0.9)
plt.show()
dict_args = {
"N": widgets.IntSlider(min=100, max=1e5, step=100, value=1000),
"a": widgets.FloatSlider(min=0.5, max=0.5, step=0.01, value=0.5),
"b": widgets.FloatSlider(min= 1., max=1., step=0.01, value=1.0),
"d": widgets.FloatSlider(min=-1.0, max=1.0, step=0.01, value=0.33),
"e": widgets.FloatSlider(min=-1.0, max=1.0, step=0.01, value=-0.38),
"f": widgets.FloatSlider(min=-1.0, max=1.0, step=0.01, value=-0.18),
"g": widgets.FloatSlider(min=-1.0, max=1.0, step=0.01, value=-0.42),
"nmax": widgets.IntSlider(min=1, max=30, step=1, value=10)
}
print("a,b paramers define the mid-point of the curve. We choose a=0.5, b=1.0 to be fixed.")
widgets.interact_manual(render_affine, **dict_args);