

In [ ]:

    
import numpy as np
import matplotlib.pyplot as plt

%matplotlib inline



In [ ]:

    
from newdust import graindist



In [ ]:

    
# Some basic magic numbers

AMAX = 0.3 # maximum grain size, in um
RHO  = 3.0 # grain material density, in g cm^-3
MD   = 1.e-5 # dust mass column, g cm^-2

Making a GrainDist object from scratch



In [ ]:

    
MRN  = graindist.sizedist.Powerlaw(amax=AMAX)
ECUT = graindist.sizedist.ExpCutoff(acut=AMAX)
SIL  = graindist.composition.CmSilicate(rho=RHO)



In [ ]:

    
gd_mrn  = graindist.GrainDist(MRN, SIL, md=MD, custom=True)
gd_ecut = graindist.GrainDist(ECUT, SIL, md=MD, custom=True)



In [ ]:

    
ax = plt.subplot(111)
gd_mrn.plot(ax, color='k', lw=2, label='Powerlaw')
gd_ecut.plot(ax, color='b', lw=2, alpha=0.8, label='Exponential cut-off')
plt.legend(loc='upper left', frameon=False)

Making a GrainDist object with helper function

This is a shortcut for making some common grain distributions

You can change AMAX and RHO from here, too

For the different grain size distributions, amax acts in a different way:

graindist.sizedist.Grain: amax sets the singular grain size
graindist.sizedist.Powerlaw: amax sets the maximum grain size in the distribution
graindist.sizedist.ExpCutoff: amax sets the "acut" value



In [ ]:

    
gd_mrn2 = graindist.GrainDist('Powerlaw', 'Silicate', amax=AMAX, rho=RHO, md=MD)
gd_ecut2 = graindist.GrainDist('ExpCutoff', 'Silicate', amax=AMAX, rho=RHO, md=MD)



In [ ]:

    
ax = plt.subplot(111)
gd_mrn.plot(ax, color='k', lw=2, label='')
gd_mrn2.plot(ax, color='r', lw=2, ls='--', label='Powerlaw')

gd_ecut.plot(ax, color='k', lw=2, label='')
gd_ecut2.plot(ax, color='b', lw=2, ls='--', label='Exponential Cutoff')
plt.legend(loc='upper left', frameon=False)

Note that the values I've used in this example are different from the defaults

Silicate has a default grain material density of 3.8 g cm^-3



In [ ]:

    
gd_mrn3 = graindist.GrainDist('Powerlaw', 'Silicate', md=MD)
gd_ecut3 = graindist.GrainDist('ExpCutoff', 'Silicate', md=MD)



In [ ]:

    
ax = plt.subplot(111)
gd_mrn.plot(ax, color='k', lw=2, label='')
gd_mrn3.plot(ax, color='r', lw=2, ls='--', label='Powerlaw')

gd_ecut.plot(ax, color='k', lw=2, label='')
gd_ecut3.plot(ax, color='b', lw=2, ls='--', label='Exponential Cutoff')
plt.legend(loc='upper left', frameon=False)

You can also get the mass density of dust grains

But the GrainDist.plot function only does number density



In [ ]:

    
plt.plot(gd_mrn.a, gd_mrn.mdens, label='Powerlaw')
plt.plot(gd_ecut.a, gd_ecut.mdens, label='ExpCutoff')
plt.loglog()
plt.xlabel('Radius (um)')
plt.ylabel('Mass density (g cm$^{-2}$ um$^{-1}$)')
plt.legend(loc='upper right', frameon=False)

You can't plot single grain sizes

It prints the number density instead



In [ ]:

    
gd_single = graindist.GrainDist('Grain', 'Silicate', amax=AMAX)



In [ ]:

    
gd_single.plot()



In [ ]: