The torsion of stellar streams and the overall shape of galactic gravity’s source

Context. Flat rotation curves, v ( r ), are naturally explained by elongated (prolate) dark matter (DM) distributions, and we have provided competitive ﬁts to the SPARC database. To further probe the geometry of the halo, or the equivalent source of gravity in other formulations, one needs observables outside the galactic plane. Stellar streams, poetically analogous to airplane contrails, but caused by tidal dispersion of massive substructures such as satellite dwarf galaxies, would lie on their own plane (consistently with angular momentum conservation) should the DM-halo gravitational ﬁeld be spherically symmetric. Tracks resembling entire orbits are seldom available because their periods are commensurable with Hubble time, with streams often presenting themselves as short segments. Aims. Therefore, we aim to establish stellar stream torsion, a local observable that measures the deviation from planarity in di ﬀ erential curve geometry, as a diagnostic providing sensitivity to aspherical DM distributions and ensuring the use of even relatively short streams. Methods. We performed small-scale simulations of tidally distorted star clusters to check that indeed a central force center produces negligible torsion, while distorted halos can generate it. Turning to observational data, we identiﬁed among the known streams those that are at the largest distance from the Galactic center, and that are likely not a ﬀ ected by the Magellanic clouds, as the most promising for the study, and by means of polynomial ﬁts we extracted their di ﬀ erential torsion. Results. We ﬁnd that the torsion of the few known streams that should be sensitive to most of the Milky Way’s DM halo is much larger than expected for a central spherical bulb alone. This is consistent with the nonsphericity of the halo. Conclusions. Future studies of stellar stream torsion with larger samples and further out of the galactic plane should be able to extract the ellipticity of the halo to see whether it is just a slight distortion of a spherical shape or whether it rather resembles a more elongated cigar.

The problem of galactic rotation is the empirical statement that rotational velocity around the galactic center seems to flatten out for a large fraction of the galaxy population where this has been measured at long enough distances Rubin et al. (1978).This creates a puzzle if the data is interpreted in terms of rotation external to a spherical source around which orbital equilibrium (Kepler's third law written for the velocity v) demands falling v for stellar objects or clouds of gas, The most accepted explanation is the existence of additional Dark Matter or a modification of the law of motion under gravity.Because typical velocities in spiral galaxies are of order 200-300 km/s, v/c ∼ 10 −3 , relativity is a correction and Newtonian mechanics should get the bulk of the rotation right.Therefore, either a modification of mechanics, such as MOND Milgrom (1983), or a modification of the gravity source, typically in the form of a spherical Dark Matter halo Frenk et al. (1985), are invoked.MOND however runs into problems at larger, cosmological scales Aguirre et al. (2001); Dodelson & Liguori (2006); and a spherical DM distribution has to be fine-tuned to have very nearly an isothermal ρ(r) ∝ 1/r 2 profile to explain the flatness of the rotation curve.
If we inhabited a two-dimensional cosmos, however, the natural gravitational law would be |F| ∝ 1 r instead of ∝ 1 r 2 and the observed rotational law would be v 2D ∝ constant which is the law that the experimental data demands.We do not; but a cylindrical matter source achieves the same dimensional reduction by providing translational symmetry along the OZ symmetry axis of the cylinder Slovick (2010); Llanes-Estrada (2021).If the linear density of the cylindrical dark matter source is λ, we can write That is, the constant velocity function v(r) is natural for a filamentary source.Moreover, if the rotation curve is only measured to a finite r, obviously the case, the source does not need to be infinitely cylindrical: it suffices that it be prolate (elongated) instead of spherical, as shown by detailed fits Bariego-Quintana et al. (2022) to the SPARC database Lelli et al. (2016) and consistently with simulations of DM haloes Allgood et al. (2006).
Observables in the galactic plane alone, such as detailed rotation curves, cannot distinguish between competing models such as spherical haloes with nearly ρ DM (r) ∝ 1 r 2 profiles or elongated haloes with arbitrary profile.To lift the degeneracy be-tween shape and profile one needs to find adequate, simple observables from data outside the galactic plane.
For a while now, stellar streams Ibata & Gibson (2007) in the Milky Way galaxy have been a promising new source of information on the DM distribution Nibauer et al. (2023), as they will eventually be for other galaxies Pearson et al. (2022), or more generically, for the effective gravitational potential Sanders (2014).Streams have been proposed Erkal et al. (2016) to detect DM clumps (for example, presumed subhaloes) which may collide with a stream and leave a gap imprinted therein In the rest of this article we develop what we think is a key observable to be measured on those streams to bear on the question of the overall shape of the presumed halo.Section 2 is dedicated to reviewing the definition of torsion in differential curve geometry, because it is the parameter which locally quantifies the aplanarity of a curve in three dimensions and to showing that, around a spherical halo, orbits as well as streams are torsionless.(Note the distinction: streams can be thought of as a bundle of nearby orbits, one for each constituent body.If a whole stream is taken as a perfect orbit, biases in the full reconstruction of the gravitational potential can creep in Sanders & Binney (2013).)Section 3 then shows how we expect tidal streams around elongated gravitational sources to show torsion if there is a component of the velocity parallel to the axis of elongation of the source.Section 4 makes a reasonable selection among the known stellar streams in our galaxy and we plot the torsion calculated along each of them, showing that there seems to be a signal here.Section 5 then discusses how further studies can improve the conclusion.

Torsion quantifies separation from orbital planarity
Before explaining why we wish to propose torsion as a useful observable to probe the DM halo, to fix notation, let us recall a few concepts of differential geometry.There do Carmo (2017) the torsion of a curve measures how sharply it is twisting out of the osculating plane, instantaneously defined by the velocity and normal acceleration.
To a curve r(t) in three-dimensional space parametrized by an arbitrary variable t we can associate an arclength s(t) = t |r(t ′ )|dt ′ and the tangent vector T = dr/ds; if at a certain point P the curvatuve is non-zero, then the normal vector at P is defined by N = dT/ds (its inverse modulus giving the radius of the circumference best approximating the curve at P); and the binormal vector (that completes the Frenet-Serret trihedron) by the vector product of both, If the curve is perfectly planar the tangent and normal vectors will always lie in the same plane, and in such case the binormal vector stays parallel to itself along the curve.Any natural definition of torsion will then yield zero.
But if the curve twists out of its osculating plane (like a uniformly advancing helix, which corresponds to constant torsion), the binormal vector will acquire a rotation.Torsion will then measure the speed of that rotation of the binormal, and it is a locally defined vector at each point P along the curve r(t), as the scalar product of the intrinsic derivative of B and the normal vector (this discounts the change of the modulus of B and rather measures its twisting), If the arc length is not at hand and the arbitrary parameter t needs to be used, then a convenient formula (with the prime denoting Since up to three derivatives of the position along the curve need to be computed, several adjacent points of a discretized curve are needed to extract the torsion: but it is still quite a local observable that does not need long trajectory stretches.
We are going to demonstrate the use of this observable τ for stellar streams, particularly around the Milky Way, to determine the shape of the gravitational potential of its DM Halo.

Movement around a Newtonian spherical source
Newtonian gravity predicts, for motion around a spherical body, with M the mass inside the sphere of radius |r|.The needed third derivative can be computed in a straight-forward manner, taking into account that |r| ′ = r • r ′ is the modulus of the projection of the velocity along the visual from the origin, in terms of components along the velocity and along the position.Because of Eq. ( 6), and therefore, observing that both terms of Eq. ( 7) are proportional to either r ′ or r, we see that (r ′ × r ′′ ) ⊥ r ′′′ .Therefore, the scalar product in the numerator of Eq. ( 5) vanishes, and thus τ = 0 for motion around a spherical body.The planarity of the orbit around a central potential is, of course, a textbook consequence Goldstein et al. (2001) of the conservation of the direction of the angular momentum vector L that in this language is parallel to the binormal vector.And additionally, the Newtonian gravity law is not strictly necessary: any central potential will yield the same result.This observation is of particular interest for the MOND explanation of the galactic rotation curves in Eq. (1) since, while the intensity of the acceleration induced by matter is different from Newtonian mechanics, the central direction of the force is respected: MOND likewise predicts no torsion.

Simulation of an N-point stellar stream around a spherical gravitational source
The discussion presented in Subsec.2.2 refers to the torsion of one test body under a central field.If the body lost dust grains forming a kind of contrail, its shape through space would be a planar curve.But stellar streams are not quite of this nature, rather the result of the tidal stretching of a globular cluster or dwarf galaxy Noreña et al. (2019).Since each object in the cluster starts off at a different height z with respect to the galactic plane, its orbit around the center of force lies on a slightly different plane, the effect being that the cluster, in addition to stretching, contorts, with the upper particles passing under the center of mass (cm) to become the lower ones with each half orbit.We here show that this effect is negligible and the torsion of a stream around a central potential can safely be neglected, as the cm of the stream follows one of the trajectories of subsection 2.2, with τ = 0. Rather than entangling the discussion in detailed theory, a couple of simple simulations will illustrate the point.
We simulate a globular cluster of N (specifically 100, 200, 500, 700, and 1000 particles) pointlike stars, randomly distributed at t = 0 over a sphere of radius R 0 (of order one or few kpc, typically; in the following simulation, 2 kpc) at a distance |r| from the galactic center (of order 30 kpc in the following example), with a certain mass m * i and common initial velocity v * ∼ 220 kpc.An additional random velocity kick ∆v 0 = Gm cluster /2r 0 in a random direction is given.We then let it evolve under the gravitational force of the central source with mass M ∈ (10 9 , 10 12 )M ⊙ , standing for a galactic bulb or a spherical DM halo, and we allow for a correction due to the inner binding forces of the cluster.This is small because the random masses are taken in the interval m * ∈ (0, 20)M ⊙ and thus their mutual interactions are orders of magnitude smaller than those with the galactic center.The constant GM of the central source can conveniently be eliminated in terms of the typical velocity of circular orbits around the galactic center, from orbital equilibrium v 2 rot /r = GM/r 2 .For the Milky Way this is typically 220 km/s.
We show the simulation in Fig. 1.The concentrated orange points in the stream visualization plots (three dimensional views as well as Cartesian projections as marked in the axes) mark the initial cluster at t = 0.The evolved cluster at later times, the clouds of cyan, blue, and purple dots, is seen to stretch under tidal tensions.The stream remains near the (slightly tilted) plane containing the initial velocity, without developing motion outside that plane, thus with negligible torsion for each individual orbit.Nevertheless, a measurable τ ∼ O(10 −4 ) is visible on the lower-right plot, due to the need to fit a curve to the dispersed, criss-crossing points.

Movement around a Newtonian cylindrical source
We now move on to quickly show how torsion is expected to look for an orbit around a perfectly cylindrical source of gravity, in a discussion paralleling that of subsection 2.2.We naturally employ cylindrical coordinates (ρ, φ, z), so that, denoting the time derivative by a ′ sign as in Eq. ( 4), where in the acceleration we recognize, from left to right, the radial, centrifugal, Coriolis, azimuthal and vertical accelerations, respectively.The force law is the same as that for a line of charge in electromagnetism, except of course with the constant replaced, so that in terms of the linear mass density λ, Comparing with the general form in Eq. ( 11) we recover z ′′ = 0 =⇒ z ′ = constant (reflecting translational invariance along the OZ axis) and ρ 2 φ ′ = constant so that the third component L z /m of angular momentum per unit mass is conserved just as in the central force problem.However, now the direction of L is not conserved, so that the binormal vector changes and one expects a torsion.To obtain it, start from Eq. ( 12) and take a further derivative to obtain (valid for the Newtonian force with cylindrical symmetry only).
Calculating the cross-product of Eq. ( 10) and Eq. ( 12), while using the righthandedness of the trihedron ( ρ, φ, ẑ) to evaluate each basis vector product, yields Next we take the scalar product with r ′′′ and evaluate Eq. ( 5) to obtain the torsion, yielding that we have cast in an easier to remember form in the second expression.Clearly, for there to be a torsion we need both azimuthal and vertical velocities (so stellar streams in the galactic plane are not sensitive, as expected).Additionally, because , so its maximum magnitude is controlled by the distance from the stellar stream segment to the galactic axis.

Simulation of an N-point stellar stream around a cylindrical source
Next we proceed to repeat the exercise of subsection 2.3 with the same starting data, but replacing the central spherical Newtonian source by a cylindrical source.The force, relegated to Eq. 40 in the appendix, now distinguishes the vertical acceleration from its two components parallel to the galactic plane.Its first term is the acceleration caused by the cylindrical gravitational source (that along the OZ axis being zero by translational symmetry).Its linear mass-density λ = M/L is obtained from the typical rotation curve around a galaxy v rot = √ 2Gλ Llanes-Estrada ( 2021).The second term of Eq. 39 is, again, the correction due to the tiny binding of the stellar streams stars among themselves, together with Eq. ( 40).We show two typical computations in figures 2, a three-dimensional view, and 3, a more detailed analysis of a run.
If the starting velocity profile was perfectly set in the XY plane perpendicular to the cylinder, the torsion would still be zero as per Eq. ( 15).We give it a slight tilt and then the orbit starts behaving as a helix.This can be appreciated by studying the originally compact cluster of stars at different times up to 1 Gyr, seeing it spread and ascend in a spiral, showing a torsion which is a factor 30 larger than for the stream around a purely spherical source.
The value of the torsion extracted from the numerical simulation is τ ≃ 1.08 • 10 −2 kpc −1 .This is in excellent agreement with the analytical estimate that we now present.To obtain it, we first note that Eq. ( 15) for the torsion, with the parametrization of the stream track given in Eq. ( 27) and an almost circular planar projection (semiaxes satisfying a ≃ b) is Here, c is the vertical "velocity" v z in terms of the arbitrary curve parameter t, and ω = v φ /ρ the rotational "angular velocity" with ρ = a the radius of the orbit in the XY plane.With the setup data for the simulation in figure 3, namely c ∼ 80 kpc/Gyr, a ∼ 30 kpc (note that the curvature scale is then a −1 ∼ 0.03 kpc −1 ) and ω ∼ 7 rad/Gyr, Eq. ( 16) yields τ| helix ≃ 0.011 kpc −1 in agreement with the numerical value 0.0108 extracted from the simulation.The slight oscillation seen in figure 3 is due to the ellipticity of the helix (with the oscillation period of about half a Gyr corresponding to half of the actual period of the trajectory because of the reflection symmetry of an ellipse).
Unfortunately for observing this in a physical setup, at a fixed observation time the stream appears quite planar, and it is only in a computer where we can simulate several time intervals at the same time that we see the helix.This is because of the translation symmetry along the cylinder, which does not allow for a vertical tidal stretching of the stream.On the upside, this makes the extracted torsion very insensitive to the number of bodies in the simulation, from 10 2 to 10 3 , since the additional Fig. 3: The top row tracks the motion of a star ensemble around around a toy-model, perfectly cylindrical source with an initial v z = 90 km/s, with the same conventions as in Fig. 1.A stream forms along a spiral trajectory.From a helicoidal fit to the cloud of points, as needs to be done on real data, with χ 2 /N gdl = 1.16 we find a value for the torsion of τ ≃ 1.1 × 10 −2 kpc −1 .This is in agreement with the analytical estimate in the main text below Eq.( 16).
bodies maintain reasonably close orbits in the z-v z vertical phase space.Providing an initially larger speed dispersion along the OZ axis increases the length of the stream but also blurs it, making its reconstruction more difficult.Thus, we move onto the next case, in which both a sphere and a cylinder (a very idealized toy model of a galaxy and a prolate DM halo) respectively cause the vertical tidal stretching and the helicoidal motion.

Sphere and cylinder with initial nonvanishing v z,0
To close this section, we will combine together both types of sources, a sphere (akin to a visible-matter galactic bulb) and a cylinder (mimicking an elongated DM halo).For the sphere we take the typical mass of a galaxy M s ∈ (10 9 , 10 12 )M ⊙ Busha et al. (2011), Binney & Tremaine (1987) and for the cylinder we use the expression for the linear mass density that we obtain from the asymptotic velocity at large r in the rotation curve v rot of the Milky Way, as exposed in the previous section 3.
The updated expression for the acceleration of the stars in the stream is obtained by combining Eqs. ( 38) and (39), that is, where v is taken from the galactic rotation velocity when it has flattened out at large r, and v 0 estimated from the visible mass.
We have added a small but appreciable (5km/s) contribution to the average initial velocity in the z direction, so that v i * = (220 + ∆v i x , ∆v i y , 5 + ∆v i z ) km/s to induce sufficient vertical motion that will generate torsion as per Eq. ( 15).Additionally, there is a dispersion ∆v i which depends on the particular star and is randomized from a distribution between 1 and 1000 times the quantity G j m * j /r 0 in terms of the mass of the small cluster and the characteristic length scale.There is no particular physical meaning to this initial dispersion except to be able to handle an ensemble of stars which are not all perfectly collinear to test the fitting procedure.
In Fig. 4 we clearly observe traits of the motion around cylinder+sphere sources as described in Llanes-Estrada (2021).Along the symmetry OZ axis, a star will describe harmonic oscillations between the two hemispheres due to the Newtonian pull of the spherical part of the distribution acting towards the center (unless it is provided with escape velocity, in which case it will approach an asymptotic trajectory, a helix around the OZ axis).The orbit on the XY plane is not closed due to the additional 1/r force from the elongated DM halo.The net effect in three dimensions can be seen as a precession of the orbital plane around the OZ axis, with the trajectory creating complicated helicoidal patterns.The N-body simulation in Fig. 5 reflects this and clearly shows the appearance of nonvanishing torsion in the stellar stream.From a piecewise helicoidal fit with Eq. ( 27) to the stream at different times t = 200 My, t = 400 My, t = 600 My and t = 800 My, we obtain different values for the torsion above the value τ > 10 −4 kpc −1 that we found on streams around a purely spherical source (basically, our irreducible background).The actual value happens to be apparently smaller than that for the purely cylindrical geometry presented in the previous section, but this is a contingent effect due to the initial conditions of the stream, and bears no significance.(The χ 2 /N gdl values for each of the visible stream segments, one at each successive time, are 0.35, 0.32, 0.23 and 1.88, respectively.)The threedimensional rendering in the left bottom plot clearly shows the torsion causing an out of osculating plane movement.

Torsion in a galaxy with a spherical halo and a galactic plane
We wish to have a reference for a minimum torsion that we would consider "normal" in order of magnitude, so that if extensive studies of stellar streams show that their torsion exceeds that level, one could reject the hypothesis of a spherical halo.For this purpose we propose here a toy model in which the halo is taken spherical, but we add a disk component.This adds a vertical (not radial) velocity outside of the galactic plane that points towards it.
The simplest (and coarsest) such model takes the galactic disk as being uniform and infinite.This is a reasonable approximation only for streams that do not elevate too much along the OZ axis; otherwise it provides an upper bound to a more realistic torsion (since such additional vertical force will always be larger than that of a finite disk, whose effect will fall off with z).In that case, observed torsions above this bound would still entail an incompatibility with a spherical halo to be studied further.
Therefore, in this minimum-torsion model we take the acceleration caused by the external force as In this equation, σ is a surface mass density for the disk, in the range 50 − 100M ⊙ /parsec 2 which is a usual estimate Kuijken & Gilmore (1989), Kuijken & Gilmore (1991) at 8 kpc from the galactic center.
Fig. 6 shows the characteristic wobbling of movement near the galactic plane caused by the planar disk, which is qualitatively consistent with Cordoni et al. (2021).
We can provide an analytical estimate of the torsion following the now familiar reasoning.Since an instantaneous velocity that is parallel to any of the three coordinate vectors of the cylindrical base { ρ, φ, ẑ} will display zero torsion, we take a trajectory combining two of them, Multiplying by the acceleration in Eq. ( 19) we obtain To construct the determinant ((r ′ × r ′′ ) • r ′′′ ) necessary for the torsion, we evaluate the third derivative outside the galactic plane (over which it is undefined), that is in the plane given by position and velocity, employing A slightly tedious but straightforward calculation then yields The numerator has mechanical dimensions of a squared momentum, and the denominator of squared momentum times length, yielding the correct 1/L dimensionality of the torsion.Moreover, the structure of the numerator shows that in the presence of a spherical source (M) alone, or a galactic plane (σ) alone, the torsion vanishes as it should.Likewise, both components of the velocity have to be nonvanishing as in Eq. ( 15) for the cylindrical source; and the torsion is null both on the galactic plane (z = 0) and on its perpendicular axis through the center of the sphere (ρ = 0).We can then numerically evaluate Eq. ( 24) to obtain the floor value of the torsion that we should expect to be able to use in the galaxy.Taking into account that the galactic plane is not infinite so that the elevation z will yield a diminishing multipolar field, it may be that galactic torsions from a spherical halo plus disk are even smaller; what we mean by this estimate is that those streams that may be found with larger values need to be further investigated as they may be teaching us something about the dark matter halo or about dark matter inhomogeneities.
We then conclude that torsions of stellar streams below 10 −3 kpc −1 in our galaxy can be explained without resort to deformed dark matter haloes or exotic phenomena.Of the few streams presently known at large radial distances, some present torsions below this level and are thus of no further interest for this application of the shape of the haloes.It is those that reach τ at the percent level that deserve further scrutiny to bear on the halo shape, among the ones known and in future searches for streams.

Stellar streams in the Milky Way and their torsion
In this section we finally turn to some of the known stellar streams in the Milky Way.We select as relevant those found at distances d > 30 kPc from the galactic center, so that the internal Fig. 6: Simulation of a small tidal stream in the presence of a typical non-spherical galactic source, here the galactic plane model of Eq. ( 19).From top to bottom: three-dimensional view, XZ and YZ sections.structure of the galaxy, such as the disk and spiral arms, produces the minimum possible alteration in the stream.These streams, see Figs. 7 and 8, have been extracted from Mateu (2023).The intention in this section is to extract the value of the torsion of the parametrized stream curves with Eq. ( 5) and check for their vanishing (or not).
We have taken two of the streams out of further consideration, namely those at Orphan-Chenab and Styx.The reason is that they may be influenced by gravity sources outside the MW.Due to the proximity of the Large Magellanic Cloud (LMC) to our galaxy, the streams in its periphery in the angular direction of that cloud could suffer alterations due to this additional source of gravity Conroy et al. (2021), Lilleengen et al. (2022).Moreover, time-dependent perturbations such as those caused by the Magellanic cloud Erkal et al. (2019);Vasiliev et al. (2021) difficult the interpretation of the stream track as that of an orbit, accruing further uncertainty (though they can be used in reverse to constrain the time-dependent potential of the perturbing body).This means that even those streams not removed from the analysis could be affected by time-dependent perturbations that we do, at present, not capture.Sometimes the effect can be very pronounced, such as a visible kink in the stream, others it can be more subtly manifested by, for example, a misalignment be-Fig.7: Trajectories in galactocentric coordinates for some of the MW outer streams, extracted from Mateu (2023).
tween the stream track and the velocity vectors Vasiliev et al. (2021).
To obtain the torsion of the curves following the stream we have employed both a smooth polynomial parametrization and a (piecewise) elliptical-helicoidal parametrization that we here report.Therefore, we first fit each of the galactocentric Cartesian coordinates tracking the individual streams in the data compilation of Mateu (2023) to order-four polynomial or elliptical/helicoidal parametric curves.The parameter that describes each curve takes values in the interval t ∈ [0, 1].The generic parametrizations are relegated to the appendix, see Eq. ( 26) and following.The parameters k or t therein are arbitrary coordinates that can be converted to arc length, that has clearer geometric significance, by means of the derivatives in Eq. ( 5) are to be taken respect to the parameter t, in general, or s, if a change is variables is effected (the outcome is the same, of course) to obtain the torsion.
To work with the streams in the database1 we use the galstream library and to perform the fits we use a least-squares method with a standard python installation.
A word about the uncertainty in this extraction is warranted.The data points (x, y, z) for the extracted stream trajectories are quoted without (∆x, ∆y, ∆z) uncertainties in the original reference Mateu (2023).Possible large systematic errors, particularly due to distance tracks which are difficult to pinpoint with accuracy, can at present not be estimated.Thus, the understood uncertainty of our parametric reconstruction stems in its entirety from our own fit interpreting the stream tracks in terms of Eq. ( 27).Until uncertainties in the track data on which we lean are compiled, the systematics will remain unclear.The reconstructed tor-Fig.8: Trajectories for some of the MW outer streams in Mollweide projection, extracted from Mateu (2023).sion will only become accurate when higher-precision estimates of the distance gradient along a stream with estimated uncertainties are provided.
After the reconstruction of the parametric curves we can obtain the expression for the torsion of each curve r = (x, y, z) using once more Eq.( 5) τ = Det(r ′ , r ′′ , r ′′′ ) ||r ′ × r ′′ || 2 , where the derivative is taken respect to the t parameter, ′ = d/dt.We analytically express the derivatives in terms of the polynomial/elliptical parametrizations and then evaluate them as function of the t parameter that we used for the fit, taking values from 0 to 1.Because the torsion is parametrization independent, the torsion can also be given as function of the arclength τ(s) calculating the derivatives respect to either of t or s.
The local torsion along many of the streams in Fig. 9 takes significant values, above the expected O(10 −4 ) "noise" that we found in the simulations when the streams evolved around spherical sources, as well as the O(10 −3 ) or less from simulations that included a galactic plane.Higher values for the torsion have been found for streams such as Cetus, Willka-Yaku, Cetus-Palca and Sagittarius In fact, the two streams with largest central value of the torsion (Willka-Yaku and Cetus) do not appear, by eyeballing the actual stream, to have such important helicity.In the case of Cetus, the large fit uncertainty makes it irrelevant, and the case of Willka-Yaku, whose torsion is still near the distribution of the other streams, deserves being revisited in the future because its shape is very different from our fit functions.
In Table (1) we quantify the torsion in the stream sample: we quote the mean and median values of the torsion among the parametrized pieces of the stream, and also the standard deviation and the interquartile range of the torsion distribution.The high values τ ∼ O(10 −2 ) found for the torsion in several cases are compatible with a non-spherical gravitational source, and also perhaps with the interaction with other gravitational sources different from the overall galactic field, so the observable appears to merit future study.
Because torsion (as curvature) has dimensions of inverse length, we would expect stellar streams to perhaps show an inverse relation with respect to their distance from the galactic center, as defined in Eq. ( 5).Irrespectively, in a galaxy such as the Milky Way the torsion of galactic streams should have a characteristic scale of (10 kPc) −1 .As per the discussion around Table 1: From left to right: values for the mean, standard deviation, median and interquartile range of the absolute value of the torsion in kpc −1 , and value of the χ 2 per number of degrees of freedom of the fit.Recall that, from section 3.4 we expect that values above 10 −3 kpc −1 may signal a nonspherical DM distribution as the effect of the aspherical galactic plane is not strong enough.Eq. ( 15), where we established τ ∈ − 1 ρ , 1 ρ , our selection of streams at 30 kpc or more means that we would consider values of the torsion of order 0.01 in units of inverse kiloparsec to be sizeable and very different from zero, and table 1 shows several such.Also, as per the discussion below Eq. ( 24), those with τ above 0.001 kpc −1 could perhaps carry interesting information about the DM distribution.

Stream
Turning to the data, the torsions that we seem to observe in the MW streams show mean values above 0.001 that could give an insight into the DM distribution, but also information about the interaction of the MW with other gravitational sources.Several streams show a significant value for the torsion: only Eridanus and Jet seem to have torsion compatible with being seeded by the galactic plane, as small as τ ∼ O(10 −3 ).
From the information we obtain from the streams in the database we cannot explore specific details on the shape of the gravitational source, beyond establishing its probable nonsphericity.5).The green and red lines and bands of the same shade represent the median value with the interquartile range and the mean value plus/minus one standard deviation, respectively.Stellar streams such as Eridanus and Jet show values of the torsion of order τ ∼ O(10 −3 ), so their torsion is to be taken compatible with the influence of a source of the intensity of the galactic plane.The fit to the Sagittarius data has been performed dividing the stream into 5 slices due to its irregular shape.

Conclusions and outlook
The problem of galactic rotation curves suggests that galaxies are surrounded by significant amounts of dark matter, and the overall shape of these sources is yet to be ascertained.Whereas spherical DM distributions around galaxies have to be fine-tuned to explain the flatness of rotation curves, a cylindrical (or, generally, prolate) DM source can naturally explain that flattening.This avoids the fine-tuning of spherical DM haloes to precisely follow the 1/r 2 fall-off for a large swath of r values.Observables inside the galactic plane cannot however distinguish between spherical (though fine tuned) and cylindrical/elongated gravitational sources, but out-of-galactic-plane information could provide new strong discriminants.
The stellar streams around the Milky Way have been extensively investigated for a while now, and are still nowadays a relevant subject of research.The track left by a stream, once extracted, can be used as the best proxy we have for a long-range trajectory of size commensurate to a galaxy, with all the attending dynamical uncertainties that this identification causes Bovy (2014), and thus a tool to infer the geometry of the impacting gravitational sources.Stream tracks can be characterized by their torsion according to Eq. ( 5); around a central potential orbits move in a plane and are expected to be torsionless (see Fig. 1).In addition, test masses around cylindrical sources are expected to follow helical orbits in which the torsion is non-zero (see Fig. 3) if given a vertical velocity.Another approach is to consider an ellipsoid-shaped halo, which is not perfectly cylindrical but rather elongated.The expected orbit of the stream components would arise from the combination of the orbits around central potentials and the helical orbits around cylinders, as is seen in Fig. 4. Torsion is then an observable directly tailored to assess the prolateness of a DM halo, as opposed to a multiparameter reconstruction which, though sensitive to the halo shape Koposov et al. (2023), is affected by a degeneracy with many other features.
The streams of the Milky Way have been a subject of research for a considerable time span, and many of the objects that constitute these streams have been catalogued.From a reconstruction of the stream tracks we infer the torsion caused by the gravitational source in Fig. 9.In this work we only consider those streams that seem to be far enough away from the galactic center to (1) avoid large effects from the baryonic component of the galaxy and (2) have a bird's eye view of the DM halo from outside a large fraction thereof.From the extraction of the torsion we see that it is non-negligible in some of the streams considered.
From our evaluation of the torsion we do not dare favor one or another interpretation of the DM halo shape in view of current data; this article should be seen as a proposal for a new observable, τ, and a first exploratory study.We do find streams with significant torsion, which is encouraging and suggests that further scrutiny could eventually inform us about the shape of the halo.We have payed no attention to the sign of τ since it only distinguishes left-from right-handed tracks.That opens another interesting study in itself, as, with the predominant galactic rotation determining a sign, the preference of a given handedness would inform an asymmetry in the vertical motion.But this would require much more statistics and take us far afield at present.
The most trustworthy streams should be those stemming from kinematically cold globular cluster (low velocity dispersion), as those present the best defined tracks.In future observational work it might be interesting to actively seek streams that show both vertical motion (along the axis perpendicular to the MW plane) and also azimuthal motion around that axis, as those with large v z and v φ will be most sensitive to the torsion.Should those streams show trajectories that are compatible with lying on a plane (zero torsion), a spherical halo would be preferred.Should they however appear helicoidal, with nonnegligible torsion, they would be pointing to an elongated DM halo.
Streams that may be detected in nearby galaxies carry the same information about their respective haloes.
However, for galactic-halo population studies with streams, which will require a 10-MPc reach, three-dimensional reconstruction of the streams appears more difficult; in that case, only the projected curvature is measurable Nibauer et al. (2023), but not the torsion, and what can be learnt about the DM halo hangs on statistical analysis to reduce the impact of degeneracies due to the projection.
Finally, other observables can bear on the overall shape of the halo, and we are investigating the shape-sensitivity of the gravitational lensing of both electromagnetic and gravitational radiation.

Track fits
A simple parametrization (x(t), y(t), z(t)) that we have employed is a piecewise fourth-order polynomial fit reading x i (t) = a i t 4 + b i t 3 + c i t 2 + d i t + e i for i=1,2,3. (26) It is idle to try to relate t to a Newtonian time since, not knowing a priori the dynamics of the system, it is unknown at which time a star was at what position along its trajectory.Only the instantaneous (present) geometry of the stream is known with certainty, and therefore an arbitrary parameter t (or the arc length after computing it) should suffice.
Another simple yet useful parametrization is that of an elliptical helix, with (a, b, c, ω, ϕ, x 0 ) the parameter set such that for t = 0 the curve cuts the XY plane with z = 0 with initial phase ϕ in the x − y rotation, c controlling the advance of the helix and a, b the elliptic projection on the XY plane.Should the rotation axis also need to be chosen different from OZ an additional orthogonal matrix O(φ, θ, ξ) parametrized by three Euler angles should be introduced as x → Ox.This we do not employ here out of simplicity, assuming the OZ axis is an axis of symmetry for the DM halo and coincides with the galactic rotation axis.We need to distinguish two types of fits.First, when fitting our simulated data, there is no obstacle in choosing the coordinate axes so that OX lies along the initial visual from the galactic center to the globular cluster which will give rise to the stream.In that case, the offset x 0 in Eq. ( 27) is unnecessary.
For the purposes of the fit, the elliptical parametrization is then, eliminating time, the length parameters A and B take values proportional to the distance of the stream to the galactic center in kpc.C is a vertical "velocity" and ω an "angular velocity" in the XY plane; but because of the arbitrariness of the parameter t, they are only relevant in the combination ω/C that indicates the angle advanced per unit height in the corkscrew motion characterizing torsion.
(If one insists in thinking of the orbit of a single star instead of the track of a stream, then the angular velocity is of the order of magnitude ω ∼ v ⊥ r ⊥ ∼ 220km/s 10 4 Pc × 3.15 10 16 s/Gyr × 1 3.1 10 13 km/Pc , resulting in ω ∼ 22 rad/Gyr.) The second type of fit involves real data.Then, the osculating plane for the starting point of each piece of the orbit need not be the aligned with the galactic Cartesian axes in account for initial conditions and precession in a noncentral field.We thus need a minimum of one additional angle Θ for each fit, to orient that starting plane, Additionally, real streams may need an offset x 0 = (x 0 , y 0 , z 0 ) bringing the total number of fit parameters to nine.
Once the parametrization has been chosen, we need to decide what is the optimal track through a cloud of them in a simulation.A strategy is to employ a squared-distance minimization strategy.The function to be minimized is the following sum running over each of the stars in the cloud, D 2 (a, b, c, t 0 , ω, where t i will be the point along the curve for which its distance to the ith point will be minimized.Taking the partial derivatives the set of equations N i=1 constrain the wanted parameters a j chosen among the fittable set.The first of these equations is the condition of orthogonality of the curve's tangent to the visual to the point.Small stretches of a stellar stream can be fit with the functional form of Eq. ( 27) for a straightforward extraction of the torsion.The torsion for this elliptical helix is then This is much simplified and becomes a t-independent constant for a circular helix with a = b, as given in the main text in Eq. ( 16).

Equations of motion for N-body simulations
Finally, we specify the (simplest available) numerical method employed to simulate an N-body stellar stream in section 2.3.The positions of the stellar objects composing it are updated in Cartesian coordinates.For this we use Euler's Method with time step ∆t = t f /N t , with the velocity updated via a once-improved Euler step, x i j+1 = x i j + ∆tv i j for i=1,2,3 (36) where f i is the function yielding each component's acceleration.The acceleration is calculated at each step from standard formulae such as a i=1,2,3 = − GM x i (x 2 + y 2 + z 2 ) 3/2 − N−1 j=1 G m j x i ((x − x j ) 2 + (y − y j ) 2 + (z − z j ) 2 ) 3/2 (38) for a spherical central source.The first line of this expression is the acceleration caused by the central spherical source, and the second is the force that attempts to bind the stellar-stream stars together (and that is too weak to avoid the tidal stretching).

Fig. 1 :
Fig. 1: Top row: simulated movement of a 200 star cluster around a spherical gravitational source with parameters akin to those of the Milky Way, in the three planar projections, extracted at times t = 0 My, t = 200 My, t = 400 My, t = 600 My and t = 800 My.Lower left: three-dimensional view with the containing plane.The initial cluster has a radius of 2 kpc and is at a distance of 30 kpc from the galactic origin of coordinates.Member stars have slightly different velocities v * = 220kpc ± ∆v due to the additional random ∆v.(All axes in kpc.)Lower right: extracted torsion as function of time, from an analytical helicoidal fit to the resulting stream with χ 2 /N gdl = 0.86, which yields a value of τ ≃ 3.5 × 10 −4 kpc −1 (a one-star orbit would return exactly zero torsion).

Fig. 2 :
Fig. 2: We show the spread of a 100-body stellar cluster with initial vertical speed v z = 50km/s after t = 0, 125, 250, 375 and 500 My.The resulting stream at different times clearly shows an upward advancing helix moving around the z-axis.

Fig. 4 :
Fig. 4: The first row represents the trajectory of a point mass with, alternatively, initial vertical velocities v z =1 and 2 km/s around a source combining both a sphere and a cylinder.The second row extends this to simultaneously show the motion of five stars.

Fig. 5 :
Fig. 5: Top row: Formation of a stellar stream from an initial cluster orbiting a combined spherical+cylindrical source (a toy model representing the galaxy plus an elongated DM halo) as in Fig. 4, with positions given at the times indicated in the legend.Lower left: three-dimensional view of the stream compared with a plane containing the origin and the initial stream's position and velocity, showing that at subsequent times, the stream leaves that plane.Lower right: piecewise helicoidal fit.A nonvanishing torsion can be extracted for each of the short segments at different times, yielding a mean value of τ = 0.22 • 10 −3 kpc −1 .

Fig. 9 :
Fig. 9: Local torsion |τ| of the MW streams considered in this work, calculated using Eq.(5).The green and red lines and bands of the same shade represent the median value with the interquartile range and the mean value plus/minus one standard deviation, respectively.Stellar streams such as Eridanus and Jet show values of the torsion of order τ ∼ O(10 −3 ), so their torsion is to be taken compatible with the influence of a source of the intensity of the galactic plane.The fit to the Sagittarius data has been performed dividing the stream into 5 slices due to its irregular shape.
N gdl