;Created for an ISM homework.  
;This is a small program to calculate the angular
;momentum of a constant density sphere that must 
;hold itself up at every point by Keplarian motion.
;I break the sphere up into cylindrical shells that
;rotate at Keplarian speeds due to the force from
;the mass below it.  In effect it is like projecting
;all of the mass onto a disk and calculating the 
;motion that way.   
;To call the program in IDL, type '.com agm.pro',
;then type 'agm' and the program will print the total
;angular momentum in the IDL prompt. Be sure that you
;have this text file saved in the directory you are
; working in as agm.pro

pro agm

r = 4.37d17            ;radius of the sphere in cm
n=250                  ;number of shells to use 
delta=r/n              ;step size in cm
pho=5.09d-20           ;density of sphere in g/cm^3
g=6.67d-8              ;gravitational constant 
speed=findgen(n+1)     ;array for the speed of each shell
rad=findgen(n+1)       ;array for the radius of each shell

d=delta
m=pho*(!pi)*(d^2)*sqrt(r^2-d^2)
l=0

for i=1, n do begin

l=l+pho*(!dpi)*(d^2-(d-delta)^2)*2*sqrt((r^2-d^2)*g*d*m)

m=m+pho*(!dpi)*(d^2-(d-delta)^2)*2*sqrt(r^2-d^2)

d=d+delta

speed[i]=sqrt(g*m/d)
rad[i]=d

endfor

print, l
;print, rad
;print, speed

end