Thursday, February 2, 2017

How to build a projection (Schlegel diagram) of a tesseract, and show a four-dimensional point inside it

This is an explanation about how to build a projection (Schlegel diagram) of a tesseract, and to show a point inside it (this is a mirror of my question at MSE). Regarding the methodology applied to visualize the 4D point, basically, if we want to show a point inside the tesseract, we need to project the tesseract first, and then project the desired point as well, following the same projection rules.

1. The definition of the tesseract is as follows (please see my post at MSE for the credits):

The tesseract is a four dimensional cube. It has 16 edge points $v=(a,b,c,d)$, with $a,b,c,d$ either equal to $+1$ or $-1$. Two points are connected, if their distance is $2$. Given a projection $P(x,y,z,w)=(x,y,z)$ from four dimensional space to three dimensional space, we can visualize the cube as an object in familiar space. The effect of a linear transformation like a rotation
$$
R(t)=\pmatrix{1&0&0&0\\0&1&0&0&\\0&0&\cos(t)&\sin(t)\\0&0&-\sin(t)&\cos(t)}
$$
in $4D$ space can be visualized in $3D$ by viewing the points $v(t) = P R(t) v$ in $\mathbb R^3$.

2. The definition of the projection (please see my post at MSE for the credits):

$$
P(x, y, z, w) = \frac{h}{h - w}(x, y, z).
$$

3. And finally, the definition of the distance between two four-dimensional points is calculated as follows (this is used to show the edges of the tesseract properly, making lines between the correct projected vertices):

$$d=\sqrt{(x_0-x'_0)^2+(x_1-x'_1)^2+(x_2-x'_2)^2+(x_3-x'_3)^2}$$

4. I have prepared a Python code snippet that creates the frames (jpg) of an animation of a tesseract including an internal point. In this case, the length of the edge is $1000$ so the distance between the vertices is not $2$, but $1000$. For the projection I have used a light source located at three times the length of the edge, this is $h=3000$. Finally, I have applied a rotation as defined above. The star is marking the location of the point $(\frac{3}{4} \cdot \frac{1000}{2}, \frac{3}{4} \cdot \frac{1000}{2}, \frac{3}{4} \cdot \frac{1000}{2}, \frac{3}{4} \cdot \frac{1000}{2})$ while we rotate the tesseract. Be aware that the position of the camera in the animation is lateral. The typical location of the camera is from above, which is the usual view of a square inside a square. But for visualization purposes (we want to see clearly the movement of the projection of the point due to the rotation of the tesseract) the camera in this case was located in a lateral position. Please use and modify it freely (for instance instead of one point it is possible to show a set of points and verify if there is symmetry, etc.):






from math import pi, sin , cos, sqrt
import matplotlib.pyplot as plt
import matplotlib as mpl
from mpl_toolkits.mplot3d import Axes3D

edge_length=1000
edge_half_length= int(edge_length/2)

lotuples=[]
list_of_loxt_lists=[]
list_of_loyt_lists=[]
list_of_lozt_lists=[]

rotation_accuracy=100
filled_once=False
for ratio in range(0,rotation_accuracy):
    angle= ((2*pi)*ratio)/rotation_accuracy
    loxt=[]
    loyt=[]
    lozt=[]
  
    #t=edge_half_length (positive)
    a=-edge_half_length
    b=edge_half_length
    ret0=-edge_half_length
    ret1=edge_half_length
    finala=a
    finalb=b
    finalret0=(ret0*cos(angle))+(ret1*sin(angle))
    finalret1=(ret0*(-sin(angle)))+(ret1*cos(angle))
  
    light_projection_factor = ((edge_length*3)/((edge_length*3)-(finalret1)))
    loxt.append(light_projection_factor*finala)
    loyt.append(light_projection_factor*finalb)
    lozt.append(light_projection_factor*finalret0)
    if filled_once==False:
        lotuples.append([a,b,ret0,ret1])
      
    a=-edge_half_length
    b=-edge_half_length
    ret0=-edge_half_length
    ret1=edge_half_length
    finala=a
    finalb=b
    finalret0=(ret0*cos(angle))+(ret1*sin(angle))
    finalret1=(ret0*(-sin(angle)))+(ret1*cos(angle))  

    light_projection_factor = ((edge_length*3)/((edge_length*3)-(finalret1)))
    loxt.append(light_projection_factor*finala)
    loyt.append(light_projection_factor*finalb)
    lozt.append(light_projection_factor*finalret0)
    if filled_once==False:
        lotuples.append([a,b,ret0,ret1])

    a=edge_half_length
    b=edge_half_length
    ret0=-edge_half_length
    ret1=edge_half_length
    finala=a
    finalb=b
    finalret0=(ret0*cos(angle))+(ret1*sin(angle))
    finalret1=(ret0*(-sin(angle)))+(ret1*cos(angle))  

    light_projection_factor = ((edge_length*3)/((edge_length*3)-(finalret1)))
    loxt.append(light_projection_factor*finala)
    loyt.append(light_projection_factor*finalb)
    lozt.append(light_projection_factor*finalret0)
    if filled_once==False:
        lotuples.append([a,b,ret0,ret1])

    a=edge_half_length
    b=-edge_half_length
    ret0=-edge_half_length
    ret1=edge_half_length
    finala=a
    finalb=b
    finalret0=(ret0*cos(angle))+(ret1*sin(angle))
    finalret1=(ret0*(-sin(angle)))+(ret1*cos(angle))  

    light_projection_factor = ((edge_length*3)/((edge_length*3)-(finalret1)))
    loxt.append(light_projection_factor*finala)
    loyt.append(light_projection_factor*finalb)
    lozt.append(light_projection_factor*finalret0)
    if filled_once==False:
        lotuples.append([a,b,ret0,ret1])

    a=edge_half_length
    b=edge_half_length
    ret0=edge_half_length
    ret1=edge_half_length
    finala=a
    finalb=b
    finalret0=(ret0*cos(angle))+(ret1*sin(angle))
    finalret1=(ret0*(-sin(angle)))+(ret1*cos(angle))  

    light_projection_factor = ((edge_length*3)/((edge_length*3)-(finalret1)))
    loxt.append(light_projection_factor*finala)
    loyt.append(light_projection_factor*finalb)
    lozt.append(light_projection_factor*finalret0)
    if filled_once==False:
        lotuples.append([a,b,ret0,ret1])

    a=edge_half_length
    b=-edge_half_length
    ret0=edge_half_length
    ret1=edge_half_length
    finala=a
    finalb=b
    finalret0=(ret0*cos(angle))+(ret1*sin(angle))
    finalret1=(ret0*(-sin(angle)))+(ret1*cos(angle))  

    light_projection_factor = ((edge_length*3)/((edge_length*3)-(finalret1)))
    loxt.append(light_projection_factor*finala)
    loyt.append(light_projection_factor*finalb)
    lozt.append(light_projection_factor*finalret0)
    if filled_once==False:
        lotuples.append([a,b,ret0,ret1])

    a=-edge_half_length
    b=edge_half_length
    ret0=edge_half_length
    ret1=edge_half_length
    finala=a
    finalb=b
    finalret0=(ret0*cos(angle))+(ret1*sin(angle))
    finalret1=(ret0*(-sin(angle)))+(ret1*cos(angle))  

    light_projection_factor = ((edge_length*3)/((edge_length*3)-(finalret1)))
    loxt.append(light_projection_factor*finala)
    loyt.append(light_projection_factor*finalb)
    lozt.append(light_projection_factor*finalret0)
    if filled_once==False:
        lotuples.append([a,b,ret0,ret1])

    a=-edge_half_length
    b=-edge_half_length
    ret0=edge_half_length
    ret1=edge_half_length
    finala=a
    finalb=b
    finalret0=(ret0*cos(angle))+(ret1*sin(angle))
    finalret1=(ret0*(-sin(angle)))+(ret1*cos(angle))  

    light_projection_factor = ((edge_length*3)/((edge_length*3)-(finalret1)))
    loxt.append(light_projection_factor*finala)
    loyt.append(light_projection_factor*finalb)
    lozt.append(light_projection_factor*finalret0)
    if filled_once==False:
        lotuples.append([a,b,ret0,ret1])

    #t=-edge_half_length (negative)
    a=-edge_half_length
    b=edge_half_length
    ret0=-edge_half_length
    ret1=-edge_half_length
    finala=a
    finalb=b
    finalret0=(ret0*cos(angle))+(ret1*sin(angle))
    finalret1=(ret0*(-sin(angle)))+(ret1*cos(angle))  

    light_projection_factor = ((edge_length*3)/((edge_length*3)-(finalret1)))
    loxt.append(light_projection_factor*finala)
    loyt.append(light_projection_factor*finalb)
    lozt.append(light_projection_factor*finalret0)
    if filled_once==False:
        lotuples.append([a,b,ret0,ret1])

    a=-edge_half_length
    b=-edge_half_length
    ret0=-edge_half_length
    ret1=-edge_half_length
    finala=a
    finalb=b
    finalret0=(ret0*cos(angle))+(ret1*sin(angle))
    finalret1=(ret0*(-sin(angle)))+(ret1*cos(angle))  

    light_projection_factor = ((edge_length*3)/((edge_length*3)-(finalret1)))
    loxt.append(light_projection_factor*finala)
    loyt.append(light_projection_factor*finalb)
    lozt.append(light_projection_factor*finalret0)
    if filled_once==False:
        lotuples.append([a,b,ret0,ret1])

    a=edge_half_length
    b=edge_half_length
    ret0=-edge_half_length
    ret1=-edge_half_length
    finala=a
    finalb=b
    finalret0=(ret0*cos(angle))+(ret1*sin(angle))
    finalret1=(ret0*(-sin(angle)))+(ret1*cos(angle))  

    light_projection_factor = ((edge_length*3)/((edge_length*3)-(finalret1)))
    loxt.append(light_projection_factor*finala)
    loyt.append(light_projection_factor*finalb)
    lozt.append(light_projection_factor*finalret0)
    if filled_once==False:
        lotuples.append([a,b,ret0,ret1])

    a=edge_half_length
    b=-edge_half_length
    ret0=-edge_half_length
    ret1=-edge_half_length
    finala=a
    finalb=b
    finalret0=(ret0*cos(angle))+(ret1*sin(angle))
    finalret1=(ret0*(-sin(angle)))+(ret1*cos(angle))  

    light_projection_factor = ((edge_length*3)/((edge_length*3)-(finalret1)))
    loxt.append(light_projection_factor*finala)
    loyt.append(light_projection_factor*finalb)
    lozt.append(light_projection_factor*finalret0)
    if filled_once==False:
        lotuples.append([a,b,ret0,ret1])

    a=edge_half_length
    b=edge_half_length
    ret0=edge_half_length
    ret1=-edge_half_length
    finala=a
    finalb=b
    finalret0=(ret0*cos(angle))+(ret1*sin(angle))
    finalret1=(ret0*(-sin(angle)))+(ret1*cos(angle))  

    light_projection_factor = ((edge_length*3)/((edge_length*3)-(finalret1)))
    loxt.append(light_projection_factor*finala)
    loyt.append(light_projection_factor*finalb)
    lozt.append(light_projection_factor*finalret0)
    if filled_once==False:
        lotuples.append([a,b,ret0,ret1])

    a=edge_half_length
    b=-edge_half_length
    ret0=edge_half_length
    ret1=-edge_half_length
    finala=a
    finalb=b
    finalret0=(ret0*cos(angle))+(ret1*sin(angle))
    finalret1=(ret0*(-sin(angle)))+(ret1*cos(angle))  

    light_projection_factor = ((edge_length*3)/((edge_length*3)-(finalret1)))
    loxt.append(light_projection_factor*finala)
    loyt.append(light_projection_factor*finalb)
    lozt.append(light_projection_factor*finalret0)
    if filled_once==False:
        lotuples.append([a,b,ret0,ret1])

    a=-edge_half_length
    b=edge_half_length
    ret0=edge_half_length
    ret1=-edge_half_length
    finala=a
    finalb=b
    finalret0=(ret0*cos(angle))+(ret1*sin(angle))
    finalret1=(ret0*(-sin(angle)))+(ret1*cos(angle))  

    light_projection_factor = ((edge_length*3)/((edge_length*3)-(finalret1)))
    loxt.append(light_projection_factor*finala)
    loyt.append(light_projection_factor*finalb)
    lozt.append(light_projection_factor*finalret0)
    if filled_once==False:
        lotuples.append([a,b,ret0,ret1])

    a=-edge_half_length
    b=-edge_half_length
    ret0=edge_half_length
    ret1=-edge_half_length
    finala=a
    finalb=b
    finalret0=(ret0*cos(angle))+(ret1*sin(angle))
    finalret1=(ret0*(-sin(angle)))+(ret1*cos(angle))  

    light_projection_factor = ((edge_length*3)/((edge_length*3)-(finalret1)))
    loxt.append(light_projection_factor*finala)
    loyt.append(light_projection_factor*finalb)
    lozt.append(light_projection_factor*finalret0)
    if filled_once==False:
        lotuples.append([a,b,ret0,ret1])
        filled_once=True

    list_of_loxt_lists.append(loxt)
    list_of_loyt_lists.append(loyt)
    list_of_lozt_lists.append(lozt)

list_of_loxi_lists=[]
list_of_loyi_lists=[]
list_of_lozi_lists=[]

list_of_loxi_lists_axis=[]
list_of_loyi_lists_axis=[]
list_of_lozi_lists_axis=[]

for ratio in range(0,rotation_accuracy):
    angle= ((2*pi)*ratio)/rotation_accuracy
    loxi=[]
    loyi=[]
    lozi=[]

    finala=int((3/4)*edge_half_length)
    finalb=int((3/4)*edge_half_length)
    finalret0=(int((3/4)*edge_half_length)*cos(angle))+(int((3/4)*edge_half_length)*sin(angle))
    finalret1=(int((3/4)*edge_half_length)*(-sin(angle)))+(int((3/4)*edge_half_length)*cos(angle))

    light_projection_factor = ((edge_length*3)/((edge_length*3)-finalret1))
    loxi.append(light_projection_factor*finala)
    loyi.append(light_projection_factor*finalb)
    lozi.append(light_projection_factor*finalret0)

    list_of_loxi_lists.append(loxi)
    list_of_loyi_lists.append(loyi)
    list_of_lozi_lists.append(lozi)
  
    # Show projection of refence axes BEGIN
    loxi=[]
    loyi=[]
    lozi=[]

    finala_axis=finala
    finalb_axis=0
    finalret0_axis=0
    finalret1_axis=0
  
    light_projection_factor = ((edge_length*3)/((edge_length*3)-finalret1_axis))
    loxi.append(light_projection_factor*finala_axis)
    loyi.append(light_projection_factor*finalb_axis)
    lozi.append(light_projection_factor*finalret0_axis)
  
    finala_axis=0
    finalb_axis=finalb
    finalret0_axis=0
    finalret1_axis=0
  
    light_projection_factor = ((edge_length*3)/((edge_length*3)-finalret1_axis))
    loxi.append(light_projection_factor*finala_axis)
    loyi.append(light_projection_factor*finalb_axis)
    lozi.append(light_projection_factor*finalret0_axis)
  
    finala_axis=0
    finalb_axis=0
    finalret0_axis=finalret0
    finalret1_axis=0
  
    light_projection_factor = ((edge_length*3)/((edge_length*3)-finalret1_axis))
    loxi.append(light_projection_factor*finala_axis)
    loyi.append(light_projection_factor*finalb_axis)
    lozi.append(light_projection_factor*finalret0_axis)
  
    finala_axis=0
    finalb_axis=0
    finalret0_axis=0
    finalret1_axis=finalret1
  
    light_projection_factor = ((edge_length*3)/((edge_length*3)-finalret1_axis))
    loxi.append(light_projection_factor*finala_axis)
    loyi.append(light_projection_factor*finalb_axis)
    lozi.append(light_projection_factor*finalret0_axis)
  
    list_of_loxi_lists_axis.append(loxi)
    list_of_loyi_lists_axis.append(loyi)
    list_of_lozi_lists_axis.append(lozi)
    # Show projection of refence axes END

for ratio in range(0,rotation_accuracy):  
    fig = plt.figure()
    ax = fig.gca(projection='3d')
    ax.view_init(elev=17., azim=-152) 

    for i in range(0,len(lotuples)):
        for j in range(i+1,len(lotuples)):
            distance = int(sqrt(((lotuples[i][0]-lotuples[j][0])**2)+((lotuples[i][1]-lotuples[j][1])**2)+((lotuples[i][2]-lotuples[j][2])**2)+((lotuples[i][3]-lotuples[j][3])**2)))
            if distance<=edge_length:
                ax.plot([list_of_loxt_lists[ratio][i],list_of_loxt_lists[ratio][j]],[list_of_loyt_lists[ratio][i],list_of_loyt_lists[ratio][j]],[list_of_lozt_lists[ratio][i],list_of_lozt_lists[ratio][j]],"r")


        ax.plot([-edge_length],[edge_length],[-edge_length],"w")
        ax.plot([-edge_length],[-edge_length],[-edge_length],"w")
        ax.plot([edge_length],[edge_length],[-edge_length],"w")
        ax.plot([edge_length],[-edge_length],[-edge_length],"w")
        ax.plot([edge_length],[edge_length],[edge_length],"w")
        ax.plot([edge_length],[-edge_length],[edge_length],"w")
        ax.plot([-edge_length],[edge_length],[edge_length],"w")
        ax.plot([-edge_length],[-edge_length],[edge_length],"w")

    ax.plot([list_of_loxi_lists[ratio][0]], [list_of_loyi_lists[ratio][0]], [list_of_lozi_lists[ratio][0]], "r*")

    #projection of refence axes around the point  
    ax.plot([0,list_of_loxi_lists_axis[ratio][0]],[0,list_of_loyi_lists_axis[ratio][0]],[0,list_of_lozi_lists[ratio][0]],"b")
    ax.plot([0,list_of_loxi_lists_axis[ratio][1]],[0,list_of_loyi_lists_axis[ratio][1]],[0,list_of_lozi_lists_axis[ratio][1]],"b")
    ax.plot([0,list_of_loxi_lists_axis[ratio][2]],[0,list_of_loyi_lists_axis[ratio][2]],[0,list_of_lozi_lists_axis[ratio][2]],"b")
    ax.plot([0,list_of_loxi_lists_axis[ratio][3]],[0,list_of_loyi_lists_axis[ratio][3]],[0,list_of_lozi_lists_axis[ratio][3]],"b")
  
    ax.dist=8
    mpl.pyplot.savefig("tesseract_movie_"+str(ratio)+".png")
    print("End ratio "+str(ratio))


The animated gif was generated joining the jpg files with VirtualDub.

No comments:

Post a Comment