dyna3/notes/inversive.md

Inversive Coordinates

(proposed by Alex Kontorovich as a practical system for doing 3D geometric calculations)

These coordinates are of form I=(c, b, x, y, z) where we think of c as the co-radius, b as the "bend" (reciprocal radius), and x, y, z as the "Euclidean" part, which we abbreviate E_I. There is an underlying basic quadratic form Q(I_1,I_2) = (c_1b_2+c_2b_1)/2 - x_1x_2 -y_1y_2-z_1z_2 which aids in calculation/verification of coordinates in this representation. We have:

Entity or Relationship	Representation	Comments/questions
Sphere s with radius r>0 centered on P = (x,y,z)	I_s = (\frac1{c}, \frac1{r}, \frac{x}{r}, \frac{y}{r}, \frac{z}{r}) satisfying Q(I_s,I_s) = -1, i.e., c = r/(\|P\|^2 - r^2).	Note that 1/c = \|P\|^2/r - r, so there is no trouble if \|P\| = r; we just get first coordinate to be 0. Using the point representation I_P from below, let's orient the sphere so that its normals point into the "positive side," where Q(I_P, I_s) > 0. The vector I_s then represents a sphere with outward normals, while -I_s represents one with inward normals.
Plane p with unit normal (x,y,z) through the (Euclidean) point (sx,sy,sz)	I_p = (-2s, 0, -x, -y, -z)	Note that Q(I_p, I_p) is still -1. We orient planes using the same convention we use for spheres. For example, (-2, 0, -1/\sqrt3, -1/\sqrt3, -1/\sqrt3) and (2, 0, 1/\sqrt3, 1/\sqrt3, 1/\sqrt3) represent planes that coincide in space, which have their normals pointing away from and toward the origin, respectively. Note that the ray from (sx, sy, sz) \in p in direction (-x, -y, -z) is the ray perpendicular to the plane through the origin; since (-x, -y, -z) is a unit vector, (sx, sy, sz) and hence p is at distance s from the origin. These coordinates are essentially the limit of a sphere's coordinates as its radius goes to infinity, or equivalently, as its bend goes to 0.
Point P with Euclidean coordinates (x,y,z)	I_P = (\|P\|^2, 1, x, y, z)	Note Q(I_P,I_P) = 0. This gives us the freedom to choose a different normalization. For example, we could scale the representation shown here by (\|P\|^2+1)^{-1}, putting it on the sphere where the light cone intersects the plane where the first two coordinates sum to 1.
∞, the "point at infinity"	I_\infty = (1,0,0,0,0)	The only solution to Q(I,I) = 0 not covered by (some normalization of) the above case.
Point P lies on sphere or plane given by I	Q(I_P, I) = 0	Actually also works if I is the coordinates of a point, in which case "lies on" simply means "coincides with".
Sphere/planes represented by I and J are tangent	If I and J have the same orientation where they touch, Q(I,J) = -1. If they have opposing orientations, Q(I,J) = 1.	For example, the xy plane with normal -e_z, represented by (0,0,0,0,1), is tangent with matching orientation to the unit sphere centered at (0,0,1) with outward normals, represented by (0,1,0,0,1). Accordingly, their Q - product is -1.
Sphere/planes represented by I and J intersect (respectively, don't intersect)	\lvert Q(I,J)\rvert \le (\text{resp. }>)\; 1	Follows from the angle formula and the tangency condition, at least conceptually. One subtlety: parallel planes have Q - product \pm 1, because they intersect at infinity (and in fact, are "tangent" there)!
P is center of sphere rep'd by I	Q(I, I_P) = -r/2, where 1/r = 2Q(I_\infty, I) is the signed bend of the sphere, and I_P is normalized in the standard way, which is to set Q(I_\infty, I_P) = 1/2	This relationship is equivalent to both of the following. (1) The point P has signed distance -r from the sphere. (2) Inversion across the sphere maps \infty to P.
Distance between points P and R is d	Q(I_P, I_R) = d^2/2	If P and R are represented by non-normalized vectors V_P and V_R, the relation becomes Q(V_P, V_R) = 2\,Q(V_P, I_\infty)\,Q(V_R, I_\infty)\,d^2. This version of the relation makes it easier to see why d goes to infinity as P or R approaches the point at infinity.
Signed distance between point rep'd by V and sphere/plane rep'd by I is d	In general, \frac{Q(I, V)}{2Q(I_\infty, V)} = Q(I_\infty, I)\,d^2 + d. When V is normalized in the usual way, this simplifies to Q(I, V) = d^2/r + d for a sphere of radius r, and to Q(I, V) = d for a plane.	We can use a Euclidean motion, represented linearly by a Lorentz transformation that fixes I_\infty, to put the point on the z axis and put the nearest point on the sphere/plane at the origin with its normal pointing in the positive z direction. Then the sphere/plane is represented by I = (0, 1/r, 0, 0, -1), and the point can be represented by any multiple of I_P = (d^2, 1, 0, 0, d), giving Q(I, I_P) = d^2/2r + d. We turn this into a general expression by writing it in terms of Lorentz-invariant quantities and making it independent of the normalization of I_P.
Distance between sphere/planes rep by I and J	Note that for any two Euclidean-concentric spheres rep by I and J with radius r and s, Q(I,J) = -\frac12\left(\frac rs + \frac sr\right) depends only on the ratio of r and s. So this can't give something that determines the Euclidean distance between the two spheres, which presumably grows as the two spheres are blown up proportionally. For another example, for any two parallel planes, Q(I,J) = \pm1.	Alex had said: Q(I,J)=\cosh(d/2)^2 maybe where d is distance in usual hyperbolic metric. Or maybe \cosh(d). That may be right depending on what's meant by the hyperbolic metric there, but it seems like it won't determine a reasonable Euclidean distance between planes, which should differ between different pairs of parallel planes.
Sphere centered on point P through point R		Probably just calculate distance etc.
Plane rep'd by I goes through center of sphere rep'd by J	This is equivalent to the plane being perpendicular to the sphere: that is, Q(I, J) = 0.
Dihedral angle between planes or spheres rep by I and J	\theta = \arccos(Q(I,J))	Aaron Fenyes points out: The angle between spheres in S^3 matches the angle between the planes they bound in R^{(1,4)}, which matches the angle between the spacelike vectors perpendicular to those planes. So we should have Q(I,J) = \cos(\theta). Note that when the spheres do not intersect, we can interpret this as the "imaginary angle" between them, via \cosh(t) = \cos(it).
Points R, P, S are collinear	Maybe just cross product of two differences is 0. Or, R,P,S,\infty lie on a circle, or equivalently, I_R,I_P,I_S,I_\infty span a plane (rather than a three-space). Or we can add two planes constrained to be perpendicular with one constrained to contain the origin, and all three points constrained to lie on both. But that's a lot of auxiliary entities and constraints...	R,P,S lying on a line isn't a conformal property, but R,P,S,\infty lying on a circle is.
Plane through noncollinear R, P, S	Should be, just solve Q(I, I_R) = 0 etc.
circle	Maybe concentric sphere and the containing plane? Note it is easy to constrain the relationship between those two: they must be perpendicular.	Defn: circle is intersection of two spheres. That does cover lines. But you lose the canonicalness
line	Maybe two perpendicular containing planes? Maybe the plane perpendicular to the line and through origin, together with the point of the line on that plane? Or maybe just as a bag of collinear points?	The first is the limiting case of the possible circle rep, but it is not canonical. However, there is a distinguished "standard" choice we could make: always choose one plane to contain the origin and the line, and the other to be the perpendicular plane containing the line. That choice or Plücker coordinates might be the best we can do. If we use the standardized perpendicular planes choice, then adding a line would be equivalent to adding two planes and the two constraints that one contains the origin and the other is perpendicular to it. That doesn't seem so bad. The second convention (perpendicular plane through the origin and a point on it) appears to be canonical, but there doesn't seem to be a circle representation that tends to it in the limit.
Inversion of entity represented by v across sphere s, rep'd by I_s	v \mapsto v + 2Q(I_s, v)\,I_s	This is just an educated guess, but its behavior is consistent with inversion in at least two ways. (1) It fixes points on s and spheres perpendicular to s. (2) It preserves dihedral angles with s.

The unification of spheres/planes is indeed attractive for a project like Dyna3. The relationship between this representation and Geometric Algebras is a bit murky; likely it somehow fits under the Geometric Algebra umbrella.

Additional more disorganized notes

Discussed coordinates with Alex Kontorovich. He was suggesting "inversive coordinates" -- for a sphere, that's 1/coradius, 1/radius, center/radius (where coradius is radius of sphere inverted in the unit sphere.) The advantage is tangent to and perpendicular to are linear in these coordinates (in the sense that if one is known, the condition of being tangent to or perpendicular to that one are linear). Planes have 1/radius = 0, and in fact, you can take the coordinates to be (2s, 0, x, y, z) where s is the distance to the origin and x,y,z are the normal direction. (Note the normal direction is only determined up to a scalar multiple. So could always scale so that the first non-zero coordinate is 1, or if you like only allow x, y to vary and let z be determined as sqrt(1-x^2^-y^2^). ) Points can be given by (r^2,1,x,y,z) where x,y,z are the coordinates and r is the distance to the origin. Quadratic form that tells you if something is a sphere/plane, or in the boundary, or up in the hyperbolic plane above. There are some details, but not quite explicit for modeling R^3, at http://sites.math.rutgers.edu/~alexk/files/LetterToDuke.pdf -- all this emphasize need to be agnostic with respect to geometric model so that we can experiment. Not really sure exactly how this relates or not to conformal geometric algebra, and whether it can be combined with geometric algebra. As formulated, there are clear-ish reps for planes/spheres and for points, but not as clear for lines. Have to see how to compute distance and/or specify a given distance. To combine inversive coordinates and geometric algebra, maybe think dually; there should be a lift from a normal vector and distance from origin to the five-vector; bivectors would rep circles/lines; trivectors would rep point pairs/points. What is the signature of this algebra, i.e. how many coordinates square to +1, -1, or 0? But it doesn't seem worth it for three dimensions, because there is a natural representation of points, as follows:

The signature of Q will be (1,4), and in fact Q(I1,I2) = 1/2(ab+ba) - E1\dot E2, where a is the "first" or "coradius" coordinate, "b" is the "second" or "radius" coordinate, and E is the Euclidean part (x,y,z). Then the inversive coordinates of a sphere with center (x,y,z) and radius r will be I = (1/\hat{r},1/r,x/r,y/r,z/r) where \hat{r} = r/(|E|^2 -r^2). These coordinates satisfy Q(I,I) = -1. For this to make sense, of course r > 0, but we get planes by letting the radius of a tangent sphere to the plane go to infinity, and we get I = (2s, 0, x0, y0, z0) where (x0,y0,z0) is the unit normal to the plane and s is the perpendicular distance from the plane to the origin. Still Q(I,I) = -1. Since r>0, we can't represent individual points this way. Instead we will use some coordinates J for which Q(J,J) = 0. In particular, if you take for the Euclidean point E = (u,v,w) the coordinates J = (|E|^2,1,u,v,w) then Q(J,J) = 0 and moreover it comes out that Q(I,J) = 0 whenever E lies on the sphere or plane described by some I with Q(I,I) = -1. The condition that two spheres I1 and I2 are tangent seems to be that Q(I1,I2) = 1. So given a fixed sphere, the condition that another sphere be tangent to it is linear in the coordinates of that other sphere. This system does seem promising for encoding points, spheres, and planes, and doing basic computations with them. I guess I would just encode a circle as the intersection of the concentric sphere and the containing plane, and a line as either a pair of points or a pair of planes (modulo some equivalence relation, since I can't see any canonical choice of either two planes or two points). Or actually as described below, there is a more canonical choice. I will have to work out formulas for the Euclidean distance between two entities, and the angle between them, and especially the intersection of two lines and the condition that three points are collinear. In this vein, it seems as though if J1 and J2 are the reps of two points, then Q(J1,J2) = d^2/2. So then the sphere centered at J1 through J2 is (J1-(2Q(J1,J2),0,0,0,0))/sqrt(2Q(J1,J2)). Ugh has a sqrt in it. Similarly for sphere centered at J3 through J2, (J3-(2Q(J3,J2),0000))/sqrt(2Q(J3,J2)). J1,J2,J3 are collinear if these spheres are tangent, i.e. if those vectors have Q-inner-product 1, which is to say Q(J1,J3) - Q(J1,J2) - Q(J3,J2) = 2sqrt(Q(J1,J2)Q(J2,J3)). But maybe that's not the simplest way of putting it. After all, we can just say that the cross-product of the two differences is 0; that has no square roots in it.

One conceivable way to canonicalize lines is to use the perpendicular plane that goes through the origin, that's uniquely defined, and anyway just amounts to I = (0,0,d) where d is the ordinary direction vector of the line; and a point J in that plane that the line goes through, which just amounts to J=(r^2,1,E) with Q(I,J) = 0, i.e. E\dot d = 0. It's also the point on the line closest to the origin. The reason that we don't usually use that point as the companion to the direction vector is that the resulting set of six coordinates is not homogeneous. But here that's not an issue, since we have our standard point coordinates and plane coordinates; and for a plane through the origin, only two of the direction coordinates are really free, and then we have the one dot-product relation, so only two of the point coordinates are really free, giving us the correct dimensionality of 4 for the set of lines. So in some sense this says that we could take naively as coordinates for a line the projection of the unit direction vector to the xy plane and the projection of the line's closest point to the origin to the xy plane. That doesn't seem to have any weird gimbal locks or discontinuities or anything. And with these coordinates, you can test if the point E=x,y,z is on the line (dx,dy,cx,cy) by extending (dx,dy) to d via dz = sqrt(1-dx^2 - dy^2), extending (cx,cy) to c by determining cz via d\dot c = 0, and then checking if d\cross(E-c) = 0. And you can see if two lines are parallel just by checking if they have the same direction vector, and if not, you can see if they are coplanar by projecting both of their closest points perpendicularly onto the line in the direction of the cross product of their directions, and if the projections match they are coplanar.

Engine Conventions

The coordinate conventions used in the engine are different from the ones used in these notes. Marking the engine vectors and coordinates with ', we have

I' = (x', y', z', b', c'),

where

\begin{align*}
x' & = x & b' & = b/2 \\
y' & = y & c' & = c/2. \\
z' & = z
\end{align*}$$
The engine uses the quadratic form $Q' = -Q$, which is expressed in engine coordinates as
$$Q'(I'_1, I'_2) = x'_1 x'_2 + y'_1 y'_2 + z'_1 z'_2 - 2(b'_1c'_2 + c'_1 b'_2).$$
In the `engine` module, the matrix of $Q'$ is encoded in the lazy static variable `Q`.

In the engine's coordinate conventions, a sphere with radius $r > 0$ centered on $P = (P_x, P_y, P_z)$ is represented by the vector
$$I'_s = \left(\frac{P_x}{r}, \frac{P_y}{r}, \frac{P_z}{r}, \frac1{2r}, \frac{\|P\|^2 - r^2}{2r}\right),$$
which has the normalization $Q'(I'_s, I'_s) = 1$. The point $P$ is represented by the vector
$$I'_P = \left(P_x, P_y, P_z, \frac{1}{2}, \frac{\|P\|^2}{2}\right).$$
In the `engine` module, these formulas are encoded in the `sphere` and `point` functions.