Can't combine Irrgang13 model with dynamical friction

[edan0314]

As the topic of discussion sugests, I just can't integrate the combined potential ChandrasekharDynamicalFrictionForce with Irrgang13

The code:

from galpy.potential import plotRotcurve, MWPotential as mwp, MWPotential2014 as mwp14, ChandrasekharDynamicalFrictionForce
import matplotlib.pyplot as plt
from galpy.potential.mwpotentials import Irrgang13I as asi
import numpy
from galpy.util.conversion import get_physical
from galpy.util import conversion
from astropy import units
from galpy.util.coords import XYZ_to_galcencyl, vxvyvz_to_galcencyl, cyl_to_rect, rect_to_cyl,rect_to_cyl_vec
from galpy.orbit import Orbit

#tempo
ts= numpy.linspace(0,-13.5,10001)*units.Gyr
# phase-space
o1=Orbit.from_name('Sgr')
o2=Orbit.from_name('Sgr')
o3=Orbit.from_name('Sgr')

cdf1=ChandrasekharDynamicalFrictionForce(GMs=4.*10.**8.*units.Msun, rhm=0.55*units.kpc, dens=mwp)
cdf2=ChandrasekharDynamicalFrictionForce(GMs=4.*10.**8.*units.Msun, rhm=0.55*units.kpc, dens=mwp14)
cdf3=ChandrasekharDynamicalFrictionForce(GMs=4.*10.**8.*units.Msun, rhm=0.55*units.kpc, dens=asi,ro=8.4,vo=242)

#mw's potentials with cdf
mw1_sgr=mwp+cdf1
mw2_sgr=mwp14+cdf2
asi_sgr=asi+cdf3

#integrations and plots
o1.integrate(ts,mw1_sgr,method='rk4_c')
o2.integrate(ts,mw2_sgr,method='rk4_c')
o3.integrate(ts,asi_sgr,method='rk4_c')
o1.plot(label="Milky-Way")
o2.plot(label="Milky-Way 2014", overplot=True)
o3.plot(label="ASI")
plt.legend(bbox_to_anchor=(1.05,1), loc='upper left', borderaxespad=0)
plt.xlim([0,100])
plt.ylim([-100,100])

and the output is:

File C:\ProgramData\Anaconda3\lib\site-packages\spyder_kernels\py3compat.py:356 in compat_exec
    exec(code, globals, locals)

  File c:\users\user\desktop\lsst\python\untitled0.py:39
    o3.integrate(ts,asi_sgr)

  File C:\ProgramData\Anaconda3\lib\site-packages\galpy\orbit\Orbits.py:1082 in integrate
    _check_consistent_units(self,pot)

  File C:\ProgramData\Anaconda3\lib\site-packages\galpy\orbit\Orbits.py:6067 in _check_consistent_units
    assert physical_compatible(orb,pot), 'Physical conversion for the Orbit object is not consistent with that of the Potential given to it'

AssertionError: Physical conversion for the Orbit object is not consistent with that of the Potential given to it

And have 2 questions:
1 - Why integrating ChandrasekharDynamicalFrictionForce with Irrgang13 doesn't work? Even if I give the input ro=8.4, vo=242 for the physical conversion, and I checked the orbits are in physical units.
2 - Does the same could occour with potentials with intrinsec physical units like McMillan17 or the other models of Irrgang?

[jobovy]

The issue here is that you have to be careful to use the unit-system of Irrgang13I for all objects with which you use Irrgang13I. The easiest way to do this is with get_physical from galpy.util.conversion:

from galpy.util.conversion import get_physical
o3=Orbit.from_name('Sgr',**get_physical(asi))
cdf3=ChandrasekharDynamicalFrictionForce(GMs=4.*10.**8.*units.Msun, rhm=0.55*units.kpc, dens=asi,**get_physical(asi))

then integrating the orbit should work. McMillan17 would require a similar procedure. This isn't necessary for MWPotential or MWPotential2014, because they don't have a fixed unit-system attached, that is, they can use any ro and vo. But for Irrgang13I, McMillan17, and some other MW potentials, the unit system is fixed.

[edan0314]

That worked fine, thank you!