IMF

a = 1.99812; (* Normalization is chosen to give ∫ξ(m)dm = 1 *)

(* IMF parameters below come from Kroupa 2007 *)

a1 = 0.3;
a2 = 1.3;
a3 = 2.35;

m0 = 0.01;

m1 = 0.08;

m2 = 0.5;

m3 = 150;

(* Function below demonstrates the use of Piecewise *)

IMF[m_]=Piecewise[{{a*(m^(-a1)),m0<=m<=m1},{a*((m1^(-a1))/(m1^(-a2)))*(m^(-a2)),m1<=m<=m2},{a*((m1^(-a1))/(m1^(-a2)))*((m2^(-a2))/(m2^(-a3)))*(m^(-a3)),m2<=m<=m3}}]

aSalpeter = 0.002694;

Salpeter[m_] = aSalpeter*(m^(-2.35));

(* You can use the structure below as a plot template *)

LogLogPlot[{IMF[m],Salpeter[m]},{m,m0,m3},PlotStyle->{{Thick,Black},{Thick,Black,Dotted}},Frame->True,FrameLabel -> {"M/(Subscript[M, ⊙])","ξ(m)"},LabelStyle->Medium]

(* Plotted above are Salpeter (dotted) and Kroupa (solid) IMF both normalized to 1. *)

HW 3 2009

chi1 = 24.587*evcgs;

chi2 = 54.416*evcgs;

rho0 = 10^-6;

Eq1[T_] = (4/(NAcgs*rho0))*4* (((2*Pi*mecgs*kcgs*T)/(hcgs^2)))^(3/2)*Exp[((-chi1)/(kcgs*T))];

Eq2[T_]= (4/(NAcgs*rho0))*(((2*Pi*mecgs*kcgs*T)/(hcgs^2)))^(3/2)*Exp[((-chi2)/(kcgs*T))];

N1=25;

Ttab = Table[ 10^4 + ((10^5)/30)*(i-1),{i,N1}]

Dataz1=Table[{Ttab[[i]],Evaluate[z1/.FindInstance[1>z1>0&&1>z2>0&&(((z1+2*z2))/z1)*z2==Eq2[Ttab[[i]]]&&(((z1+2*z2))/(1-z1- z2))*z1==Eq1[Ttab[[i]]],{z1,z2},Reals]][[1]]},{i,N1}];

Dataz2=Table[{Ttab[[i]],Evaluate[z2/.FindInstance[1>z1>0&&1>z2>0&&(((z1+2*z2))/z1)*z2==Eq2[Ttab[[i]]]&&(((z1+2*z2))/(1-z1- z2))*z1==Eq1[Ttab[[i]]],{z1,z2},Reals]][[1]]},{i,N1}];

Dataz0 = Table[{Ttab[[i]],1-Dataz1[[i,2]]-Dataz2[[i,2]]},{i,N1}];

Dataze = Table[{Ttab[[i]],Dataz1[[i,2]]+2*Dataz2[[i,2]]},{i,N1}];

Needs["PlotLegends`"]

ListPlot[{Dataz0,Dataz1,Dataz2},PlotRange->{0,1},Joined->True,InterpolationOrder->2,PlotStyle->{{Thick,Dotted,Black},{Thick,Dashed,Black},{Thick,Black}},Frame->True,FrameLabel->{"T, K","Ionization Fraction"},LabelStyle->{Medium}, PlotLegend->{"z0", "z1","z2"},LegendPosition->{0.8,-0.25},PlotLabel->"ρ = (10^(-6))"]

ListPlot[Dataze,PlotRange->{0,2},Joined->True,InterpolationOrder->2,PlotStyle->{Thick,Black},Frame->True,FrameLabel->{"T, K","(Subscript[z, e]) = (Subscript[n, e])/n"},LabelStyle->{Medium},PlotLabel->"ρ = (10^(-6))"]