Implement the new particle spray model by Chen+24 (#670)

Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com>
Co-authored-by: Jo Bovy <bovy@astro.utoronto.ca>

Author: 3 people

HISTORY.txt

@@ -3,6 +3,11 @@ v1.10.1 (Expected around 2024-11-01)
 
 - Propagate general plotting keywords in Potential.plot/plotPotentials.
 
+- Added the generalized particle-spray model as galpy.df.basestreamspraydf with two
+  subclasses: chen24spraydf and fadal15spraydf. Enabled integrating orbits of stream
+  particles with the progenitor's potential. Deprecating the old particle-spray
+  model galpy.df.streamspraydf.
+
 v1.10.0 (2024-07-07)
 ====================
 

doc/source/images/chen24spraydf-b14-xz-wstreamdf.png

@@ -132,7 +132,8 @@ functionality is introduced in separate papers:
 * ``galpy.actionAngle.actionAngleIsochroneApprox``: please cite `Bovy (2014) <http://adsabs.harvard.edu/abs/2014ApJ...795...95B>`__.
 * ``galpy.df.streamdf``: please cite `Bovy (2014) <http://adsabs.harvard.edu/abs/2014ApJ...795...95B>`__.
 * ``galpy.df.streamgapdf``: please cite `Sanders, Bovy, & Erkal (2016) <http://adsabs.harvard.edu/abs/2016MNRAS.457.3817S>`__.
-* ``galpy.df.streamspraydf``: please cite `Fardal et al. (2015) <https://ui.adsabs.harvard.edu/abs/2015MNRAS.452..301F/abstract>`__ for the method and `Qian et al. (2022) <https://ui.adsabs.harvard.edu/abs/2022MNRAS.511.2339Q/abstract>`__ for the ``galpy`` implementation
+* ``galpy.df.chen24spraydf``: please cite `Chen et al. (2024) <https://ui.adsabs.harvard.edu/abs/2024arXiv240801496C/abstract>`__ for the method and the ``galpy`` implementation
+* ``galpy.df.fardal15spraydf``: please cite `Fardal et al. (2015) <https://ui.adsabs.harvard.edu/abs/2015MNRAS.452..301F/abstract>`__ for the method and `Qian et al. (2022) <https://ui.adsabs.harvard.edu/abs/2022MNRAS.511.2339Q/abstract>`__ for the ``galpy`` implementation
 * ``galpy.potential.ttensor`` and ``galpy.potential.rtide``: please cite `Webb et al. (2019a) <https://ui.adsabs.harvard.edu/abs/2019MNRAS.488.5748W/abstract>`__.
 
 * ``galpy.potential.to_amuse``: please cite `Webb et al. (2019b) <http://arxiv.org/abs/1910.01646>`_.

doc/source/images/chen24spraydf-fardal15spraydf-b14-xz.png

@@ -334,7 +334,9 @@ Helper routines to compute kicks
 The distribution function of a tidal stream using a particle-spray technique
 ----------------------------------------------------------------------------
 
-Model from `Fardal et al. (2015)
+Model from `Chen et al. (2024)
+<https://ui.adsabs.harvard.edu/abs/2024arXiv240801496C/abstract>`__ and
+`Fardal et al. (2015)
 <https://ui.adsabs.harvard.edu/abs/2015MNRAS.452..301F/abstract>`__ with full
 details of the ``galpy`` implementation given in `Qian et al. (2022)
 <https://ui.adsabs.harvard.edu/abs/2022MNRAS.511.2339Q/abstract>`__;
@@ -346,5 +348,13 @@ General instance routines
 .. toctree::
    :maxdepth: 1
 
-   __init__ <streamspraydf.rst>
    sample <streamspraydfsample.rst>
+
+Specific particle-spray models
++++++++++++++++++++++++++++++++
+
+.. toctree::
+   :maxdepth: 1
+
+   Chen et al. (2024) <streamspraydfchen24.rst>
+   Fardal et al. (2015) <streamspraydffardal15.rst>

doc/source/images/streamspraydf-b14-xz-wstreamdf.png

@@ -0,0 +1,5 @@
+The particle-spray stream DF using Chen et al. (2024) formalism
+===============================================================
+
+.. autoclass:: galpy.df.chen24spraydf
+   :members: __init__

doc/source/index.rst

@@ -0,0 +1,7 @@
+.. _api_streamspraydf:
+
+The particle-spray stream DF using Fardal et al. (2015) formalism
+==================================================================
+
+.. autoclass:: galpy.df.fardal15spraydf
+   :members: __init__

doc/source/reference/df.rst

@@ -1,4 +1,4 @@
 galpy.df.streamspraydf.sample
-==========================================
+=============================
 
 .. automethod:: galpy.df.streamspraydf.sample

doc/source/reference/streamspraydf.rst

@@ -509,16 +509,24 @@ for the difference in :math:`v_Y` as a function of unperturbed :math:`X`:
 
 .. _streamspray-tutorial:
 
-Particle-spray modeling of streams with ``streamspraydf``
----------------------------------------------------------
-
-``galpy`` also contains an implementation of the particle-spray method
-for generating tidal streams, roughly following the parametrization of
-`Fardal et al. (2015)
+Particle-spray modeling of streams
+----------------------------------
+
+``galpy`` also contains implementations of two particle-spray methods
+for generating tidal streams. ``chen24spraydf`` follows the method by
+`Chen et al. (2024)
+<https://ui.adsabs.harvard.edu/abs/2024arXiv240801496C/abstract>`__.
+``fardal15spraydf`` roughly follows the
+parametrization of `Fardal et al. (2015)
 <https://ui.adsabs.harvard.edu/abs/2015MNRAS.452..301F/abstract>`__. Full
 details on the ``galpy`` implementation are given in `Qian et al. (2022)
 <https://ui.adsabs.harvard.edu/abs/2022MNRAS.511.2339Q/abstract>`__. Here,
-we give a simple example of the method.
+we give a simple example of the the two methods.
+
+.. note::
+   ``fardal15spraydf`` was previously known as ``streamspraydf`` before
+   version ``v1.10.1``. While the old name is still supported  until ``v.1.11.0``
+   for backward compatibility, we recommend using the new name ``fardal15spraydf``.
 
 Like in the ``streamdf`` example above, we use the same orbit, potential, and
 cluster mass as in
@@ -531,24 +539,31 @@ as a simple ``LogarithmicHaloPotential``):
 >>> o= Orbit([1.56148083,0.35081535,-1.15481504,0.88719443,-0.47713334,0.12019596])
 >>> lp= LogarithmicHaloPotential(normalize=1.,q=0.9)
 
-Then, we setup ``streamspraydf`` models for the leading and trailing arm of
-the stream:
+Then, we setup ``chen24spraydf`` and ``fardal15spraydf`` models for the leading
+and trailing arm of the stream. The ``chen24spraydf`` model requires integrating
+the orbits of stream stars with the progenitor's potential. Here, we use a
+Plummer potential for the prognenitor:
 
 >>> from astropy import units
->>> from galpy.df import streamspraydf
->>> spdf= streamspraydf(2*10.**4.*units.Msun,progenitor=o,pot=lp,tdisrupt=4.5*units.Gyr)
->>> spdft= streamspraydf(2*10.**4.*units.Msun,progenitor=o,pot=lp,leading=False,tdisrupt=4.5*units.Gyr)
+>>> from galpy.df import chen24spraydf, fardal15spraydf
+>>> progpot = PlummerPotential(2*10.**4.*units.Msun, 4.*units.pc)
+>>> spdf_c24= chen24spraydf(2*10.**4.*units.Msun,progenitor=o,pot=lp,tdisrupt=4.5*units.Gyr, progpot=progpot)
+>>> spdft_c24= chen24spraydf(2*10.**4.*units.Msun,progenitor=o,pot=lp,leading=False,tdisrupt=4.5*units.Gyr, progpot=progpot)
+>>> spdf_f15= fardal15spraydf(2*10.**4.*units.Msun,progenitor=o,pot=lp,tdisrupt=4.5*units.Gyr)
+>>> spdft_f15= fardal15spraydf(2*10.**4.*units.Msun,progenitor=o,pot=lp,leading=False,tdisrupt=4.5*units.Gyr)
 
-To sample a set of 300 stars in both arms, we do
+First, we sample a set of 300 stars in both arms:
 
->>> orbs,dt= spdf.sample(n=300,returndt=True,integrate=True)
->>> orbts,dt= spdft.sample(n=300,returndt=True,integrate=True)
+>>> orbs_c24,dt_c24= spdf_c24.sample(n=300,returndt=True,integrate=True)
+>>> orbts_c24,dt_c24= spdft_c24.sample(n=300,returndt=True,integrate=True)
+>>> orbs_f15,dt= spdf_f15.sample(n=300,returndt=True,integrate=True)
+>>> orbts_f15,dt= spdft_f15.sample(n=300,returndt=True,integrate=True)
 
 which returns a ``galpy.orbit.Orbit`` instance with all 300 stars. We can plot
-these in :math:`Z` versus :math:`X` and compare to Fig. 1 in
-Bovy (2014). First, we also integrate the orbit of the progenitor forward
-and backward in time for a brief period to show its location in the area
-of the stream:
+the ``galpy.orbit.Orbit`` instance in :math:`Z` versus :math:`X`
+and compare to Fig. 1 in Bovy (2014). First, we also integrate the orbit of the
+progenitor forward and backward in time for a brief period to show its location
+in the area of the stream:
 
 >>> ts= numpy.linspace(0.,3.,301)
 >>> o.integrate(ts,lp)
@@ -559,19 +574,21 @@ Then we plot
 
 >>> o.plot(d1='x',d2='z',color='k',xrange=[0.,2.],yrange=[-0.1,1.45])
 >>> of.plot(d1='x',d2='z',overplot=True,color='k')
->>> plot(orbs.x(),orbs.z(),'r.')
->>> plot(orbts.x(),orbts.z(),'b.')
+>>> plot(orbs_c24.x(),orbs_c24.z(),'r.', alpha=0.5)
+>>> plot(orbts_c24.x(),orbts_c24.z(),'r.', alpha=0.5)
+>>> plot(orbs_f15.x(),orbs_f15.z(),'b.', alpha=0.5)
+>>> plot(orbts_f15.x(),orbts_f15.z(),'b.', alpha=0.5)
 
 which gives
 
-.. image:: images/streamspraydf-b14-xz.png
+.. image:: images/chen24spraydf-fardal15spraydf-b14-xz.png
   :width: 600
 
 We can also compare to the track for this stream as predicted by ``streamdf``.
 For this, we first setup a similar ``streamdf`` model (they are not exactly
 the same, as ``streamdf`` uses a velocity dispersion to set the progenitor's
-mass, while ``streamspraydf`` uses the mass directly); see the ``streamdf``
-documentation for a full explanation of this code:
+mass, while ``fardal15spraydf`` and ``chen15spraydf`` uses the mass directly);
+see the ``streamdf`` documentation for a full explanation of this code:
 
 >>> from galpy.actionAngle import actionAngleIsochroneApprox
 >>> from galpy.df import streamdf
@@ -586,23 +603,23 @@ Then, we can overplot the track predicted by ``streamdf``:
 
 >>> o.plot(d1='x',d2='z',color='k',xrange=[0.,2.],yrange=[-0.1,1.45])
 >>> of.plot(d1='x',d2='z',overplot=True,color='k')
->>> plot(orbs.x(),orbs.z(),'r.',alpha=0.1)
->>> plot(orbts.x(),orbts.z(),'b.',alpha=0.1)
+>>> plot(orbs_c24.x(),orbs_c24.z(),'r.',alpha=0.1)
+>>> plot(orbts_c24.x(),orbts_c24.z(),'r.',alpha=0.1)
 >>> sdf.plotTrack(d1='x',d2='z',interp=True,color='r',overplot=True,lw=1.)
 >>> sdft.plotTrack(d1='x',d2='z',interp=True,color='b',overplot=True,lw=1.)
 
 This gives then
 
-.. image:: images/streamspraydf-b14-xz-wstreamdf.png
+.. image:: images/chen24spraydf-b14-xz-wstreamdf.png
   :width: 600
 
 We see that the track from ``streamdf`` agrees very well with the location
-of the points sampled from ``streamspraydf``.
+of the points sampled from ``chen24spraydf``.
 
-The ``streamspraydf`` ``sample`` function can also return the points at
+The ``sample`` function can also return the points at
 the time of stripping, that is, not integrated to the present time
 (when using ``integrate=False``); this can be useful for visualizing where
-stars get stripped from the progenitor. When initializing ``streamspraydf``,
+stars get stripped from the progenitor. When initializing a particle-spray DF,
 you can also specify a different potential for computing the tidal radius
 and velocity distribution of the tidal debris, which can be useful when the
 overall potential contains pieces that are irrelevant for computing the tidal

doc/source/reference/streamspraydfchen24.rst

@@ -66,4 +66,6 @@
 osipkovmerrittdf = osipkovmerrittdf.osipkovmerrittdf
 osipkovmerrittNFWdf = osipkovmerrittNFWdf.osipkovmerrittNFWdf
 constantbetadf = constantbetadf.constantbetadf
+chen24spraydf = streamspraydf.chen24spraydf
+fardal15spraydf = streamspraydf.fardal15spraydf
 streamspraydf = streamspraydf.streamspraydf