The Kustaanheimo-Stiefel (KS) transform turns a gravitational two-body problem into a harmonic oscillator, by going to four dimensions. In addition to the mathematical-physics interest, the KS transform has proved very useful in N-body simulations, where it helps to handle close encounters. Yet the formalism remains somewhat arcane, with the role of the extra dimension being especially mysterious. This paper shows how the basic transformation can be interpreted as a rotation in three dimensions. For example, if we slew a telescope from zenith to a chosen star in one rotation, we can think of the rotation axis and angle as the KS transform of the star. The non-uniqueness of the rotation axis encodes the extra dimension. This geometrical interpretation becomes evident on writing KS transforms in quaternion form, which also helps to derive concise expressions for regularized equations of motio

Saha, Prasenjit

Murang'a University of Technology

MUT INSTITUTIONAL REPOSITORY

SMA 2202 Algebraic Structures

SMA2202 ALGEBRAIC STRUCTURES

RERO DOC Digital Library

Mon. Not. R. Astron. Soc. 400, 228–231 (2009) doi:10.1111/j.1365-2966.2009.15437.xInterpreting the Kustaanheimo–Stiefel transform in gravitationaldynamicsPrasenjit SahaInstitute for Theoretical Physics, University of Zu¨rich, Winterthurerstrasse 190, 8057 Zu¨rich, SwitzerlandAccepted 2009 July 21. Received 2009 July 21; in original form 2009 March 30ABSTRACTThe Kustaanheimo–Stiefel (KS) transform turns a gravitational two-body problem into aharmonic oscillator, by going to four dimensions. In addition to the mathematical-physicsinterest, the KS transform has proved very useful in N-body simulations, where it helps tohandle close encounters. Yet the formalism remains somewhat arcane, with the role of theextra dimension being especially mysterious. This paper shows how the basic transformationcan be interpreted as a rotation in three dimensions. For example, if we slew a telescope fromzenith to a chosen star in one rotation, we can think of the rotation axis and angle as the KStransform of the star. The non-uniqueness of the rotation axis encodes the extra dimension.This geometrical interpretation becomes evident on writing KS transforms in quaternion form,which also helps to derive concise expressions for regularized equations of motion.Key words: stellar dynamics – celestial mechanics.1 IN T RO D U C T I O NThe Kustaanheimo–Stiefel (KS) transform is a remarkable relationbetween the two most important elementary problems in dynamics:under a transformation of coordinates and time, a Kepler prob-lem changes into a harmonic oscillator. Especially noteworthy isthat the collision singularity in the Kepler problem is transformedinto a regular point. The name comes from the works by Kustaan-heimo (1964) and Kustaanheimo & Stiefel (1965), while the book byStiefel & Scheifele (1971), which is largely devoted to the KS trans-form and its consequences, is perhaps the best-known source. Fora very short summary, see ‘regularization’ in Binney & Tremaine(2008). An important application of the KS transformation is in nu-merical orbit integration, where the singularity removal is used togreat advantage for simulating dense stellar systems with near col-lisions (Aarseth & Zare 1974a,b; Jernigan & Porter 1989; Mikkola& Aarseth 1990, 1993). Some recent papers also re-examine theformalism itself (Bartsch 2003; Waldvogel 2006).In two dimensions there is a much simpler version of the KS trans-form going back to Levi-Civita (1920). In the Levi-Civita transform,the coordinate plane is read as the complex plane, and the complexsquare root of the coordinate becomes the transformed coordinate.The geometrical interpretation is clear: the complex phase getshalved. The KS transform is also a kind of square root, but infour dimensions. One wonders how the geometrical interpretationgeneralizes.It turns out that slewing a telescope is a convenient geometricalanalogy. Suppose the telescope is at zenith and we want to slew it toE-mail: psaha@physik.uzh.cha particular star in one rotation. Normally we would simply movealong a great circle from the zenith to the star. However, we mightprefer a different rotation (to avoid crossing the moon, say). Forexample, we could choose the mid-point on the above great circleand rotate about it by 180◦. In any case, the chosen rotation axisand the rotation angle are effectively the KS transform of the star.This idea of rotation in three dimensions about a non-unique axisgeneralizes the idea of halving the phase in a complex square root.This paper attempts to provide some new insight into the KStransform by providing some reformulations and new derivations ofknown results, and especially to make the geometrical interpretationevident.2 QUAT E R N I O N S A N D ROTAT I O NBefore considering KS theory, it is useful to review a concise alge-braic way of specifying rotations in three dimensions, not often usedin astrophysics but standard in computer graphics: quaternions.Quaternions are a generalization of complex numbers. The √−1of complex numbers is replaced by three unit quaternions i , j , k,such thati i = j j = kk = −1, i j k = −1. (1)From (1) it follows that i j = k = − j i and so on. In other words,quaternions are like a combination of dot and cross products in vec-tor algebra. (Although historically quaternions came before, havingbeen invented by none other than W. R. Hamilton of Hamilton’sequations.)A general quaternion has the formA = A0 + Ax i + Ay j + Azk, (2)C© 2009 The Author. Journal compilation C© 2009 RASThe Kustaanheimo–Stiefel transform 229where we will call A0 the real part. A quaternion with no real partis effectively a vector in three dimensions.In analogy with complex numbers, we will use the followingnotation for quaternion conjugates and absolute values:re[A] ≡ A0,A∗ ≡ A0 − A1 i − A2 j − A3k,A2 ≡ A∗ A = A20 + A21 + A22 + A23. (3)It is easy to see that re [A∗] = re[A] and (AB)∗ = B∗ A∗, and as aresultre[AB] = re[B∗ A∗] = re[B A]. (4)Rotation in quaternion notation is beautifully concise. Say wewant to rotate a vector r by angle ω about a unit vector n. Usingquaternion algebra the rotation is simplyR∗r R, (5)whereR = cos 12ω + sin 12ω n. (6)Unlike the equivalent expression using Euler angles, the expres-sion (5) has no coordinate singularities (or ‘gimbal lock’) and as aresult is numerically more stable, which explains its popularity incomputer graphics.For an arbitrary (i.e. non-unit) quaternion R, the expression (5)amounts to a rotation combined with scalar multiplication.It is possible to represent quaternions as matrices (though notnecessary, even for numerical work). A familiar representation is interms of Pauli matrices:i = iσ1 j = −iσ2 k = iσ3 (7)ori =(0 ii 0)j =(0 −11 0)k =(i 00 −i). (8)Pauli matrices are most important as operators on quantum two-statesystems (being Hermitian, whereas quaternions are anti-Hermitian).In recent years, the most exciting two-state quantum systems havebeen Qbits in quantum computing. It turns out that expressions ofthe type (5) appear in the description of quantum-computing gates(see Mermin 2007, who also provides a derivation of essentially theabove three-dimensional rotation formula).3 THE KU STA A NHEIMO–STIEFELT R A N S F O R MLetq = x i + y j + zk (9)denote a point in space. The KS transform of q is the quaternionQ = Q0 + Qx i + Qy j + Qzk, (10)the transformation formula beingq = Q∗k Q. (11)A solution for Q isQI = x i + y j + Zk√2Z, Z ≡ z +√x2 + y2 + z2, (12)as is easily verified by multiplication, following the quaternion rules.However, QI is not unique, because changing toQ = (cos ψ − sin ψ k) QI (13)leaves equation (11) invariant. Thus ψ behaves like a gauge.Everything so far is already in the literature. The new result inthis paper is that we can readily visualize Q, including its non-uniqueness.Comparing (11) and (5), it is evident that Q is a rotator that takesthe z-axis to q. To visualize Q, let us rewrite q asq = r(sin θ cos φ i + sin θ sin φ j + cos θ k), (14)where r , θ and φ are the usual polar coordinates. Rewriting QI inthe solution (12) and simplifying, we haveQI = √r(sin12θ cos φ i + sin 12θ sin φ j + cos 12θ k). (15)In other words, the zenith distance of QI is half-way along the greatcircle from k to q. From (6) we see the rotation angle ω would beπ. Now let us apply the gauge transformation (13) with ψ = π/2 toQI. This givesQII = √r(cos12θ + sin 12θ sin φ i − sin 12θ cos φ j). (16)Now the implied rotation is by θ , about an axis perpendicular toboth k and q. In general, we can writeQ = cos ψ QI − sin ψ QII, (17)which is to say, Q could be anywhere on the great circle joiningQI and QII. The telescope-slewing analogy given above is simplya description of the preceding three formulas.An interesting special case is φ = 0, which gives q = r(cos θ k +sin θ i) and QI = √r[cos(θ/2) k + sin(θ/2) i]. Then QI is effec-tively the complex square root of q (we need to read k as the realaxis and i as the imaginary axis). In other words, the planar casecan be reduced to the Levi-Civita transform by a suitable gauge.Quaternion formulations of the KS transform have been discussedby several authors: Stiefel & Scheifele (1971) mention quaternionsbut appear to dislike them, while later authors (for example Vivarelli1994; Waldvogel 2006) are more favourable. The precise definitionadopted for the transform varies, but is equivalent to equation (11).That Q represents a rotation and shrinking/stretching of q is alsoknown. Bartsch (2003) specifically notes that the rotation axis isunique in two dimensions but not in three. However, the explicitdescription of the implied rotations, as above, appears to be new.4 TH E C A N O N I C A L M O M E N T U MSo far we have just discussed geometry, but of course the real sig-nificance of the KS transform is dynamics, which we now consider.Letp = px i + py j + pzk (18)be the canonical momentum conjugate to q. We seekP = P0 + Px i + Py j + Pzk (19)that will be canonically conjugate to Q. Let us writere[ p∗ dq] = re[ p∗ d Q∗ k Q] + re[ p∗ Q∗k d Q]. (20)Using the identity (4), we can rewrite the middle term asre[ p (d Q∗k Q)∗]. Since p = − p∗ and (d Q∗k Q)∗ = Q∗(−k)d Qthe term becomes becomes re[ p∗ Q∗k d Q]. Thus we havere[ p∗dq] = 2re[ p∗ Q∗k d Q]. (21)C© 2009 The Author. Journal compilation C© 2009 RAS, MNRAS 400, 228–231230 P. SahaNow if we defineP = −2k Q p, (22)we havere[ p∗d q] = re[P∗d Q] (23)which is to say, p · dq =P · d Q. Provided the Hamiltonian de-pends on P , Q only through p, q and not on the gauge ψ , thetransformation (P , Q) → ( p, q) is canonical.To get an explicit expression for p, we multiply (22) on the leftby Q∗k, obtainingp = Q∗k P2Q2. (24)Note that while we have to be careful about the order of multipli-cation when i , j , k are involved, real numbers like Q2 commutewith everything. Since p has no real part, re[ Q∗k P] = 0 identically.We can think of it as a formal constant of motion resulting frominvariance with respect to ψ .That P (as defined in equation 22, or equivalently) completes acanonical transformation is a standard part of KS theory, but thederivation of the canonical condition using quaternion identitiesappears to be new.5 TH E T WO - B O DY P RO B L E M A N D T H EH A R M O N I C O S C I L L ATO RLet us now write the Kepler HamiltonianH = 12p2 − 1/q (25)in terms of KS variables. Multiplying each of (11) and (22) by itsquaternion conjugate, we haveq2 = Q4, P 2 = 4p2Q2 (26)and substituting these givesH = 18P 2/Q2 − 1/Q2. (27)We now use a device known in Hamiltonian dynamics as a Poincare´time transformation. This involves introducing a fictitious time vari-able s, whose relation to t we choose to bedt = Q2 ds. (28)Since Q2 is the radial distance in the Kepler problem, (28) is infact Kepler’s equation, and s is the eccentric anomaly. In the ficti-tious time variable s, the equations of motion are given by a newHamiltonian = Q2(H − E) = 18P 2 − EQ2 − 1 (29)with E being the constant initial value of H. The time-transformed Hamiltonian is zero along a trajectory, but its partial derivativesare not zero.The Hamiltonian  is remarkable indeed. For E < 0 (boundorbits) it is a harmonic oscillator. Since Q has four components, is like a mass on an isotropic spring in four Euclidean dimensions.Thus the well-known fact that the bound Kepler problem has adynamical O(4) symmetry. For the unbound case, the symmetrygroup is different: formally the Lorentz group, but with a physicalmeaning completely different from special relativity. And – perhapsmost importantly – Hamilton’s equations for  are well behavedeven at Q = 0 (a collision). This is known as regularization and wasthe original motivation for KS theory.The effect of an external force F is simple. From (22) it followsimmediately that F will add an extra contribution of −2k Q F todP/dt , which amounts to a contribution of −2Q2k Q F to dP/ds.Provided the external force is non-singular, the equations of motionin s remain regular.6 R E G U L A R I Z I N G TH E T H R E E - B O DYPROBLEMApplication of KS regularization to N-body simulations involvesexpressing the gravitating system either as a tree-like hierarchy ofcoupled two-body systems (Jernigan & Porter 1989) or as a chain(Mikkola & Aarseth 1990, 1993). The basic idea can be describedusing the three-body problem with all masses unity. Here again,quaternions enable a concise formulation.In relative coordinates, the Hamiltonian for three unit gravitatingmasses can be written (cf. equation 12 in Aarseth & Zare 1974a) asH = 12p21 +12p22 + p1 · p2 −1q1− 1q2(30)plus an additional potential V (q1, q2). Here q1, q2 expresses theposition of the first and second body relative to the zeroth body,while p1, p2 express the momenta of the first and second bodies inthe barycentric frame. Meanwhile, V (q1, q2) expresses the mutualinteraction of the first and second bodies, plus any external potential.We can regard V (q1, q2) as an external potential, and since wealready know how to deal with external forces, we set V aside andconcentrate on H.Now we introduce KS variables q1 = Q∗1 k Q1 and so on. Defining	 ≡ re [( Q∗1 k P1)∗ Q∗2 k P2] , (31)we can write p1 · p2 as 	/(4Q21Q22). Applying a Poincare´ timetransformationdt = Q21Q22ds,  = Q21Q22(H − E) (32)gives = 18P 21 Q21 +18P 22 Q22 +14	 − Q22 − Q21 − EQ21Q22, (33)where E is the value of H. The  Hamiltonian has no denomi-nators,and is thus regular for collisions with the zeroth body. (Weassume V remains regular, that is to say, the first and second bodiesdo not collide with each other or any other bodies than the zeroth.In practice, simulations redefine the relative coordinates whenevernecessary, according to who is close to whom.)For the equations of motion, we need derivatives with respectto quaternion components. First we have ∇Q1Q21 = 2 Q1. Slightlymore subtle is ∇Q1 re[ Q∗1 A] = A if A is independent of Q1. Usingthis last identity, together with the definition (31) of 	, we derive∇P1	 = −k Q1 Q∗2 k P2, ∇Q1	 = −k P1 P∗2 k Q2. (34)Hamilton’s equations are thend Q1ds= 14Q21 P1 − k P1 P∗2 k Q2,dP1ds=(2 + 2EQ22 −14P 21)Q1 + k Q1 Q∗2 k P2 (35)and similarly for P2, Q2.7 D ISCUSSIONIn dynamical astronomy the KS transformation is profound, butmay appear mysterious. This paper attempts to make it less mys-terious, and hopefully therefore more useful, by explaining itC© 2009 The Author. Journal compilation C© 2009 RAS, MNRAS 400, 228–231The Kustaanheimo–Stiefel transform 231in three-dimensional geometric terms. There are several possibledirections in which the KS transformation may turn out to beuseful.First, one can imagine new orbit integrators specialized to nearlyKeplerian problems. Work on dense stellar systems with near colli-sions has already been mentioned (for reviews see the books Aarseth2003; Heggie & Hut 2003). In the planetary regime, which differsfrom the dense-stellar case in having few bodies but many moreorbital times, time transformations reminiscent of (28) used forKS regularization have proved useful for highly eccentric orbits(Mikkola 1997; Emel’yanenko 2002), while some integration algo-rithms (Mikkola & Tanikawa 1999; Preto & Tremaine 1999) applythe time transformation (28) implicitly. Could the KS transforma-tion itself be exploited here? Fukushima (2005) has some furtherideas.Secondly, it is conceivable that KS variables could simplify per-turbation theory. Perturbation theory in classical celestial mechanics(see for example Murray & Dermott 2000) is algebraically fright-eningly complicated, basically because the natural variables for theunperturbed and perturbed parts (being the Keplerian action-anglesand real-space coordinate) are related through an implicit equation.On the other hand, the action-angles of the KS-transformed Keplerproblem are explicitly related to space coordinates – the implicitequation is transferred to the time variable. Could some major sim-plification be achieved through KS variables? Some progress hasbeen made by Vrbik (2006).Thirdly, the KS transformation might provide new insight intoanalogous quantum problem. Bander & Itzykson (1966a,b) derivethe symmetry groups of the bound and unbound Coulomb problems.These turn out to be the same four-dimensional symmetries as inKS theory. Is the KS transformation implicit in that work?AC K N OW L E D G M E N T SI am grateful to Seppo Mikkola for introducing me to KS theory, andto Marcel Zemp, Scott Tremaine and the referee, Jo¨rg Waldvogel,for suggesting improvements in the manuscript.REFERENCESAarseth S. J., 2003, Gravitational N-Body Simulations. Cambridge Univ.Press, CambridgeAarseth S. J., Zare K., 1974a, Celest. Mech., 10, 185Aarseth S. J., Zare K., 1974b, Celest. Mech., 10, 516Bander M., Itzykson C., 1966a, Rev. Modern Phys., 38, 330Bander M., Itzykson C., 1966b, Rev. Modern Phys., 38, 346Bartsch T., 2003, J. Phys. A: Math. Gen., 36, 6963Binney J., Tremaine S., 2008, Galactic Dynamics. Princeton Univ. Press,Princeton, NJEmel’yanenko V., 2002, Celest. Mech. Dynamical Astron., 84, 331Fukushima T., 2005, AJ, 129, 2496Heggie D., Hut P., 2003, The Gravitational Million-Body Problem: A Multi-disciplinary Approach to Star Cluster Dynamics. Cambridge Univ. Press,CambridgeJernigan J. G., Porter D. H., 1989, ApJS, 71, 871Kustaanheimo P., 1964, in Stiefel E., ed., Mathematische Methoden derHimmelsmechanik und Astronautik Mathematisches Forschungsinsti-tut Oberwolfach, Berichte 1. Bibliographisches Institut Mannheim, DieSpinordarstellung der energetischen Identita¨ten der Keplerbewegung,p. 330Kustaanheimo P., Stiefel E., 1965, J. Reine Angew. Math., 218, 204Levi-Civita T., 1920, Acta Math., 42, 99Mermin N. D., 2007, Quantum Computer Science: An Introduction. Cam-bridge Univ. Press, New YorkMikkola S., 1997, Celest. Mech. Dynamical Astron., 67, 145Mikkola S., Aarseth S. J., 1990, Celest. Mech. Dynamical Astron., 47, 375Mikkola S., Aarseth S. J., 1993, Celest. Mech. Dynamical Astron., 57, 439Mikkola S., Tanikawa K., 1999, Celest. Mech. Dynamical Astron., 74, 287Murray C. D., Dermott S. F., 2000, Solar System Dynamics. CambridgeUniv. Press, CambridgePreto M., Tremaine S., 1999, AJ, 118, 2532Stiefel E. L., Scheifele G., 1971, Linear and Regular Celestial Mechanics.Springer-Verlag, BerlinVivarelli M. D., 1994, Celest. Mech. Dynamical Astron., 60, 291Vrbik J., 2006, New Astron., 11, 366Waldvogel J., 2006, Celest. Mech. Dynamical Astron., 95, 201This paper has been typeset from a TEX/LATEX file prepared by the author.C© 2009 The Author. Journal compilation C© 2009 RAS, MNRAS 400, 228–231

Interpreting the Kustaanheimo-Stiefel transform in gravitational dynamics

Abstract

Similar works

Full text

Available Versions

MUT INSTITUTIONAL REPOSITORY

MUT INSTITUTIONAL REPOSITORY

RERO DOC Digital Library

MUT INSTITUTIONAL REPOSITORY