Hopf Fibration

Extracted from Wikipedia: In mathematics, differential topology, the Hopf fibration describes a hypersphere in four-dimensional space in terms of circles and an ordinary sphere. Stereographic projection of the Hopf fibration induces a remarkable structure, in which all of 3-dimensional space, except for the z-axis, is filled with nested tori made of linking Villarceau circles.

\documentclass[border=10pt]{standalone}
\usepackage[inline]{asymptote}
\begin{document}
\begin{asy}
settings.render = 0;
import graph3;

size(20cm);

currentprojection = orthographic(-5,-4,2);

typedef triple surfaceparam(pair);
typedef pair curveparam(real);

surfaceparam revolve(curveparam F) {
  return new triple(pair uv) {
    real t = uv.x, theta = uv.y;
    pair rz = F(t);
    real r = rz.x, z = rz.y;
    return (r*cos(theta), r*sin(theta), z);
  };
}

pair circlecenter(real r) {
  return r + 1 / (1 + r);
}

curveparam circleparam(real r) {
  pair center = circlecenter(r);
  return new pair(real t) {
    return center + r * expi(2 pi * t);
  };
}

surfaceparam torusparam(real r) { return revolve(circleparam(r)); }

surfaceparam hopfparam(real r) {
  surfaceparam mytorus = torusparam(r);
  return new triple(pair uv) {
    return mytorus((uv.x, 2 pi * (uv.y - uv.x)));
  };
}

surface surface(triple f(pair z), pair a, pair b, int nu=nmesh, int nv=nu,
        bool smooth=true, pen color(pair z)) {
  surface s = smooth ? surface(f, a, b, nu, nv, Spline) : surface(f, a, b, nu, nv);
  real delta_u = (b.x - a.x)/nu;
  real delta_v = (b.y - a.y)/nv;
  for (int i = 0; i < nu; ++i) {
    for (int j = 0; j < nv; ++j) {
      patch currentpatch = s.s[s.index[i][j]];
      pair z = (interp(a.x, b.x, i/nu), interp(a.y, b.y, j/nv));
      currentpatch.colors = new pen[];
      currentpatch.colors.push(color(z));
      currentpatch.colors.push(color(z + (delta_u, 0)));
      currentpatch.colors.push(color(z + (delta_u, delta_v)));
      currentpatch.colors.push(color(z + (0, delta_v)));
    }
  }
  return s;
}

real rescale(real t) { return atan((pi/2)*t) / (pi / 2); }

for (real r : new real[] {0.2, 0.8, 2.0}) {
  pen color(pair uv) {
    real z = 2 * rescale(r) - 1;
    if (z > 1)  z = 1;
    else if (z < -1) z = -1;
    real cylradius = sqrt(1 - z^2);
    real theta = 2pi * uv.y;
    real x = cylradius * cos(theta), y = cylradius * sin(theta);
    x = (x+1)/2; y = (y+1)/2; z = (z+1)/2;
    return x*red + y*green + z*blue;
  }
  int nu = 8, nv = 8;
  surface s = surface(hopfparam(r), (0,0), (1,1/2), nu, nv, smooth=true, color=color);
  draw(s);
  for (int j = 0; j < nv; ++j)
    draw(s.vequals(j), linewidth(0.5));
  draw(s.vequals(nv-.01), linewidth(0.5));
}

\end{asy}
\end{document}

Source: math.SE

Author: Charles Staats (License)