Escher Knot Side View

Note: Similar to the Knot 1 example, however this object class is Surface3DCollection which has Cartesian native coordinates. Colormaps are:

import copy
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.colors as colors
from matplotlib.animation import FuncAnimation
import s3dlib.surface as s3d
import s3dlib.cmap_utilities as cmu

#.. influenced by M.C.Escher - Knots
#   https://mcescher.com/gallery/mathematical/

# 0. Define animation control parameters ............................

totalTime, f_domain, numFrames = 2, (0.0,1.0), 45   # time in seconds
frames=np.linspace(*f_domain, numFrames, endpoint=False)
interval = int(1000.0*totalTime/numFrames)          # milliseconds

# 1. Define function to examine ....................................

wdth = 0.75
twst, twstOff = 0.75, 0.25  # MC Escher controls
elev = 50

def getMesh(Tu) :
    Cu, Cv = 2, 5
    Su, Sv  = Tu*Cu, 4*Cv 

    Uarr = np.linspace(0,1.0,Su+1)
    Varr = np.linspace(0,1.0,Sv+1)
    y,x = np.meshgrid(Varr,Uarr)
    z = np.zeros_like(y)
    verts = np.reshape(np.array( [x,y,z] ) .T, (-1,3))

    fvInd = []
    for v in range(0,Sv) :
        for u in range(0,Su) :
            a = u + v*(Su+1)
            d = a + (Su+1)
            if not (u%Cu and (v+3)%Cv) : fvInd.append( [a,a+1,d+1,d] )
    return verts, fvInd

def SquareRing_XYZ(xyz) :
    width,height = wdth,wdth
    x,y,z = xyz        #domain: x,y [0,1], z=0
    r = np.ones(len(x))
    t = 2*np.pi*x
    z = 2*y - 1
    zeros = np.zeros(len(z))
    width_ar = np.full(len(z),width)
    # fold the cylinder into 4 parts..
    alpha = -2*width*z+width
    alpha = np.where( z <= 0.5, zeros ,     alpha )
    alpha = np.where( z <= 0.0, 2*width*z , alpha )
    alpha = np.where( z <= -.5, -width_ar , alpha )
    beta = height
    beta = np.where( z <= 0.5, 2*height*z,         beta)
    beta = np.where( z <= 0.0, zeros,              beta)
    beta = np.where( z <= -.5, -2*height*z-height, beta)
    R = r + alpha
    R = R + width/2
    Z = beta - height/2
    xyz = s3d.PolarSurface.coor_convert([R,t,Z],tocart=True)
    return xyz

def twistFunction_XYZ(xyz, twists=twst, toff=twstOff) :
    r,t,z = s3d.CylindricalSurface.coor_convert(xyz)
    offset = toff*np.pi
    x0 = 1-r
    y0 = z
    r0, t0, temp = s3d.PolarSurface.coor_convert([x0,y0,np.zeros(len(z))])
    t0 = t0 - t*twists + offset
    x1, y1, temp = s3d.PolarSurface.coor_convert([r0,t0,np.zeros(len(z))],True)
    R = 1 - x1
    Z = y1
    x,y,z = s3d.CylindricalSurface.coor_convert([R,t,Z],True)
    return x,y,z

def Trefoil_XYZ(xyz) :
    r,t,z = s3d.PolarSurface.coor_convert(xyz)
    rw = 1-wdth/2
    X = rw*(np.sin(t)+2*np.sin(2*t))
    Y = rw*(np.cos(t)-2*np.cos(2*t))
    R0,T,Z = s3d.PolarSurface.coor_convert([X,Y,z])
    R = R0 + r - rw
    Z = z - np.sin(3*t)
    x,y,z = s3d.CylindricalSurface.coor_convert([R,T,Z],True)
    return x,y,z

# 2. Setup and map surface .........................................
ba = colors.rgb_to_hsv( [ 0.482, 0.333, 0.267 ] )
bb = colors.rgb_to_hsv( [ 0.855, 0.584, 0.427 ] )
cmap1 = cmu.hsv_cmap_gradient(ba,bb,name='cpr')
cmap2 = cmu.hsv_cmap_gradient("darkolivegreen",'olivedrab',name='dog_od')
bcmap = cmu.stitch_cmap(cmap2,cmap1,name='stgrad')

ring = s3d.Surface3DCollection(*getMesh(80)).map_geom_from_op(SquareRing_XYZ)
orig_ring = copy.copy(ring)
ring.map_geom_from_op(twistFunction_XYZ)
ring.map_geom_from_op( Trefoil_XYZ )

# 3. Construct figure, add surface, and plot ......................
info, infocolor, facecolor = 'S3Dlib.org',  [0.898,0.843,0.800], [ 0.933, 0.902, 0.859 ] 
fig = plt.figure(figsize=plt.figaspect(1), facecolor=facecolor)
text = fig.text(0.11, 0.05, info, color=infocolor, fontsize=45, fontweight='bold'  )
ax = plt.axes(projection='3d',aspect='equal', proj_type='ortho', facecolor=facecolor )
minmax = (-1.8,1.8)
ax.set(xlim=minmax, ylim=minmax, zlim=minmax )
ax.set_axis_off()
ax.view_init(elev=elev)

ring.map_cmap_from_normals(direction=ax,cmap=bcmap)
ring.shade(0.2, ax=ax,direction=[-1,-1,1])
ring.hilite(0.5,ax=ax,direction=[-1,-1,1])
ax.add_collection3d(ring)

fig.tight_layout(pad=0)
#plt.show()

# 4. Animation ======================================================

def update_fig(frame):
    global ring
    ring.remove()

    ofst = 0.5*frame
    ring = copy.copy(orig_ring)
    ring.map_geom_from_op( lambda rtz : twistFunction_XYZ(rtz,toff=ofst) )
    ring.map_geom_from_op( Trefoil_XYZ )
    ring.map_cmap_from_normals(direction=ax,cmap=bcmap)
    ring.shade(0.3, ax=ax,direction=[-1,-1,1])
    ring.hilite(0.5,ax=ax,direction=[-1,-1,1])

    ax.add_collection3d(ring)
    return

anim = FuncAnimation(fig, update_fig, frames, interval=interval, repeat=True)
anim.save('mceknot_side.html',writer='html')

msg = "saved {} frames, values: [{:.3f} to {:.3f}] @ {} milliseconds/frane"
print(msg.format(numFrames,np.min(frames),np.max(frames),interval))