import numpy as np
from mpl_toolkits.mplot3d.art3d import Poly3DCollection, Line3DCollection
import mpl_toolkits.mplot3d.art3d as art3d
import matplotlib.pyplot as plt
from matplotlib.patches import Rectangle
points = np.array([[-1, -1, -1],
                  [1, -1, -1 ],
                  [1, 1, -1],
                  [-1, 1, -1],
                  [-1, -1, 1],
                  [1, -1, 1 ],
                  [1, 1, 1],
                  [-1, 1, 1]])
Z = points
#Z = 10.0*Z
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
sc = [-1,1]
X, Y = np.meshgrid(sc, sc)
ax.scatter3D(Z[:, 0], Z[:, 1], Z[:, 2],color="red")
verts = [[Z[0],Z[1],Z[2],Z[3]],
 [Z[4],Z[5],Z[6],Z[7]],
 [Z[0],Z[1],Z[5],Z[4]],
 [Z[2],Z[3],Z[7],Z[6]],
 [Z[1],Z[2],Z[6],Z[5]],
 [Z[4],Z[7],Z[3],Z[0]]]
ax.add_collection3d(Poly3DCollection(verts, facecolors = 'cyan', linewidths=1, edgecolors='red', alpha=0.2))

r= np.sqrt(3)
theta_2_0 = np.linspace(0, 2*np.pi, 100) 
theta_1_0 = np.linspace(0, 2*np.pi, 100) 

alpha1 = np.pi/2-np.arcsin(1/np.sqrt(3))

theta_1, theta_2 = np.meshgrid(theta_1_0, theta_2_0)
x2 = np.cos(theta_2)*np.sin(theta_1) * r
y2 = np.sin(theta_2)*np.sin(theta_1) * r
z2 = np.cos(theta_1) * r

ax.plot_surface(x2,y2,z2, color='b',alpha=0.2)

theta_4_0 = np.linspace(0, 2*np.pi, 100)
theta_3_0 = -1/(2*np.pi)*(theta_4_0-np.pi)**2+3/4*np.pi

x_c = np.cos(theta_4_0) * np.sin(theta_3_0) * r 
y_c = np.sin(theta_4_0) * np.sin(theta_3_0) * r 
z_c = np.cos(theta_3_0) * r 

ax.plot(x_c, y_c, z_c, color="r")

ax.set_xlabel('x axis')
ax.set_ylabel('y axis')
ax.set_zlabel('z axis')
ax.set_xlim(-2,2)
ax.set_ylim(2,-2)
ax.set_zlim(-2,2)

plt.show()