Galactic and gravitational dynamics in Python

I was wondering whether there's a way to incorporate mass evolution or mass accretion over time of a halo into the orbital calculations of gala. I know it only takes a single halo mass numerical value and spits out an orbit, but is it possible to instead include a function for mass rather than a numerical value?

This makes use of the representation differential classes in astropy/astropy#5871 to clean up a lot of the code. A natural byproduct of this is that CartesianPhaseSpacePosition and CartesianOrbit are no longer needed, since there is now a unified interface to any representations and their respective differentials.

TODO:

• [x] support <3D positions so nonlinear integrations work
• [x] remove the velocity_coord_transforms.py and use the Differential classes instead
• [x] figure out how to handle velocity frame transforms with the Differential classes
• [x] update documentation and docstrings
• [x] clean up all documentation that mentions Orbit or PhaseSpacePosition, check repr's (especially the orbits-in-detail.rst file
• [x] make sure all mention of Cartesian* is gone
• [x] make sure all code and doc tests run and don't use the old Cartesian* classes

API-breaking changes:

• Velocity frame transforms now return Differential classes
• Velocity coord transforms are gone
• CartesianPhaseSpacePosition and CartesianOrbit are deprecated

Describe your changes

This adds a .to_sympy() classmethod to the potential classes. I've also then used this method with sympy to compute all of the Hessians, and implemented these using C code generated by sympy.

Checklist

• [x] Did you add tests?
• [x] Did you add documentation for your changes?
• [x] Did you add a changelog entry? (see CHANGES.rst)
• [x] Are the CI tests passing?
• [x] Is the milestone set?

Amazingly, this closes #159, closes #56, closes #5, and closes #85 !!

running plot_contours() on an agama GalaPotential object requires a time object that is supposed to be optional, and forcing the variable to be a single value does not resolve the issue.wasn't encountering this issue until I updated agama and gala to their most recent versions. error stack below:

grid = np.linspace(-15,15,64)
fig,ax = plt.subplots(1, 1, figsize=(5,5))
fig = galapot.plot_contours(grid=(grid,grid,0), cmap='Greys', ax=ax,time=1)

/srv/conda/envs/notebook/lib/python3.8/site-packages/gala/potential/potential/core.py in plot_contours(self, grid, filled, ax, labels, subplots_kw, **kwargs)
    530                 r[ii] = slc
    531 
--> 532             Z = self.energy(r*self.units['length']).value
    533 
    534             # make default colormap not suck

/srv/conda/envs/notebook/lib/python3.8/site-packages/gala/potential/potential/core.py in energy(self, q, t)
    228         ret_unit = self.units['energy'] / self.units['mass']
    229 
--> 230         return self._energy(q, t=t).T.reshape(orig_shape[1:]) * ret_unit
    231 
    232     def gradient(self, q, t=0.):

/srv/conda/envs/notebook/lib/python3.8/site-packages/agama/py/pygama.py in <lambda>(q, t)
    898             except TypeError: PotentialBase.__init__(self, dict(), units=units)
    899             _agama.Potential.__init__(self, *args, **kwargs)
--> 900             self._energy  = lambda q,t=0.: self.potential(q, t=t)
    901             self._density = lambda q,t=0.: _agama.Potential.density(self, q, t=t)
    902             self._gradient= lambda q,t=0.: -self.force(q, t=t)

RuntimeError: Argument 'time', if provided, must be a single number or an array of the same length as points

Hi, thank you for the great package. I realized the below code in the LeapfrogIntegrator may change the sign of the velocity when _dt is negative, which may change the result of the force function because it could depend on the velocity. In my case I add a dynamical friction term in my force function and LeapfrogIntegrator will give me incorrect results.

https://github.com/adrn/gala/blob/782a8b1a19c8546d553b7c2122505e6ee82a93db/gala/integrate/pyintegrators/leapfrog.py#L146-L150

Describe your changes

Added a check of the current axis limits when plotting Orbits with autolim=True to prevent Gala from making the axis limits too small to see everything already plotted.

Checklist for contributor:

• [x] Did you add tests?
• [x] Did you add documentation for your changes?
• [x] Did you reference any relevant issues?
• [x] Did you add a changelog entry? (see CHANGES.rst)

Checklist for maintainers:

• [x] Are the CI tests passing?
• [x] Is the milestone set?

Dear Community I´m having some when installing issues. When running the recommended way to install gala, !python -m pip install gala, I get various errors and warnings.

I suspect thsis is the reason why I get errors when running imports such as

Thanks in advance

Describe your changes

Checklist

• [ ] Did you add tests?
• [ ] Did you add documentation for your changes?
• [ ] Did you reference any relevant issues?
• [ ] Did you add a changelog entry? (see CHANGES.rst)
• [ ] Are the CI tests passing?
• [ ] Is the milestone set?

so that transformation works on skyoffset frames.

Signed-off-by: Nathaniel Starkman [email protected]

Describe your changes

Checklist

• [x] Did you add tests? There currently are no tests for get_transform_matrix
• [x] Did you add documentation for your changes? yes
• [x] Did you reference any relevant issues? yes
• [ ] Did you add a changelog entry? (see CHANGES.rst)
• [x] Are the CI tests passing? yes
• [x] Is the milestone set?

python 3.7 gala, numpy versions: 1.0, 1.16.4

I get core dump / seg fault from doing the follwing:

import astropy.units as u
from gala.potential import BovyMWPotential2014

BovyMWPotential2014().hessian([[0, 8, 0]]*u.kpc)

Tracebacks: https://gist.github.com/smoh/0e803684b14d87bba4a97e5bb3e33bb0 1.out: beginning of errors from Red Hat 2.out: full traceback from Mac OSX 10.14.4

In both I did fresh install with conda

conda create -n gala-test
conda install -c conda-forge astro-gala

Any ideas?

Right now, .pericenter() and .apocenter() are slow because they do interpolation to figure out a precise value. There should be a fast=True option that skips the interpolation.

We also need to support these methods for multiple orbits in the same object.

Describe your changes

Fixed the bug described in this issue where the OrphanKoposov19 coordinate transformation called the wrong function.

Checklist

• [ ] Did you add tests?
• [ ] Did you add documentation for your changes?
• [x] Did you reference any relevant issues?
• [x] Did you add a changelog entry? (see CHANGES.rst)
• [x] Are the CI tests passing?
• [x] Is the milestone set?

The OprhanKoposov19 oph19_to_icrs() fct returns OrphanNewberg10 matrix instead of Koposov19 matrix. In this line, galactic_to_orp() should be replaced by icrs_to_orp19().

Bug report from @abonaca: Not able to use the CylSpline potential with MockStream functionality because:

ValueError: Input potential must be C-enabled: one or more components in the input external potential are Python-only.

Right now it is very slow because it must construct a spline object with the input grids each time the energy/gradient/density functions are called. It might be possible to store this object on the Wrapper class and pass a pointer in to C. To do this, we need to add functionality to the potential classes (actually the struct types) to support having an array of pointers that get passed around.