In [1]:

    

import os,sys
import numpy as np
import pylab as pl
import scipy.constants as cons
%pylab inline
from matplotlib import rcParams
import matplotlib.animation as animation
rcParams["savefig.dpi"] = 150
rcParams["savefig.dpi"] = 150
rcParams["axes.facecolor"] = '#000000'
rcParams["axes.linewidth"]=5
rcParams["axes.edgecolor"]='#dddddd'
rcParams["figure.facecolor"]   = '#000000'
rcParams["figure.edgecolor"]   = '#000000'
rcParams["axes.labelcolor"] = '#ffffff'
rcParams["xtick.color"] = '#ffffff'
rcParams["ytick.color"] = '#ffffff'
kpc2ly=3261.63344
ly2m=1.0/1.05702341e-16
H=2.1e-18
kpc2m=3e19
Msun=2e30 #kg


    

    




Welcome to pylab, a matplotlib-based Python environment [backend: module://IPython.kernel.zmq.pylab.backend_inline].
For more information, type 'help(pylab)'.

In [2]:

    

def newton(M,R):
    return np.sqrt(cons.G*M/R/ly2m)


    

In [3]:

    

def iso(Rc,Vinf,R):
    tmp=1.0-Rc/R/ly2m*np.arctan(R*ly2m/Rc)
    return Vinf*np.sqrt(tmp)


    

In [4]:

    

def nfw(Vv,cv,Rv,R):
    gcinv=np.log(1.0+cv)-1.0/(1.0/cv+1.0)
    s=R/Rv/kpc2ly
    cvs=cv*s
    vsq=Vv*Vv/s/gcinv*(np.log(1.0+cvs)-1.0/(1.0/cvs+1))
    return np.sqrt(vsq)


    

In [15]:

    

#finite radii
gal='M33'
Mv=2e11*Msun
tmp=4.0/3.0*np.pi*2e-24
Rv=np.power(Mv/tmp,0.33)*5
M={'M33':1e40, 'Vinf':105400, 'Rc':1.39*kpc2m, 'Vv':H*10*Rv, 'cv':4.0,'Rv':Rv/kpc2m}
R=np.arange(1,40)*1000.
model={}



model['nodm']={'v':newton(M[gal],R)}

model['ISO']={'v':iso(M['Rc'],M['Vinf'],R)}

model['NFW']={'v':nfw(M['Vv'],M['cv'],M['Rv'],R)}


    

In [16]:

    

pl.plot(R,model['nodm']['v'])
pl.xlabel("distance to galaxy center (ly)")
pl.ylabel("rotational velocity (km/s)")


    

    
Out[16]:





<matplotlib.text.Text at 0x10be1ca50>






    

In [17]:

    

pl.plot(R,model['nodm']['v'])
pl.plot(R,model['ISO']['v'])
pl.xlabel("distance to galaxy center (ly)")
pl.ylabel("rotational velocity (km/s)")


    

    
Out[17]:





<matplotlib.text.Text at 0x10be93950>






    

In [17]:

    




    

In [24]:

    

pl.plot(R,model['nodm']['v'],label='NO DM')
pl.plot(R,model['ISO']['v'],label='ISO')
pl.plot(R,model['NFW']['v'],label='NFW')
l=pl.legend()

for text in l.get_texts():
    text.set_color("white")
pl.xlim(4000,39000)
pl.xlabel("distance to galaxy center (ly)")
pl.ylabel("rotational velocity (km/s)")


    

    
Out[24]:





<matplotlib.text.Text at 0x10be00590>






    

In [241]:

    




    

In [241]:

    




    

In [ ]: