In-Class Assignment 2

In-Class Assignment 2#

Credit: These exercises follow from work by Mike Zingale’s computational lectures.

Following Exercise 1.3 from HKT, Chapter 1, Available Online.

Work in groups of 2

Learning Objectives#

  • familiarize with symbolic python tools

  • use sympy to compute properties of stellar structure

  • make estimates of stellar density and comparisons to the Sun

Note

Install sympy via conda using

conda install sympy

Documentation is available here

from sympy import init_session
init_session(use_latex="mathjax")
%matplotlib inline

IPython console for SymPy 1.13.3 (Python 3.11.8-64-bit) (ground types: python)

These commands were executed:
>>> from sympy import *
>>> x, y, z, t = symbols('x y z t')
>>> k, m, n = symbols('k m n', integer=True)
>>> f, g, h = symbols('f g h', cls=Function)
>>> init_printing()

Documentation can be found at https://docs.sympy.org/1.13.3/

1. Stellar density profile#

rho = symbols('rho', cls=Function)
rhoc = symbols('rho_c')
qc = symbols('q_c')
Pc = symbols('P_c')
G = symbols('G')
Mstar, Rstar = symbols('M_\star R_\star')
r = symbols('r')
xi = symbols('xi')
mu = symbols('mu')
k = symbols('k')
m_u = symbols('m_u')

Imagine a star has a density profile of the form:

\[\rho(r) = \rho_c \left [ 1 - \left (\frac{r}{R_\star} \right )^3 \right ]\]

where \(\rho_c\) is the central density.

# using the variables above, define this function for density and very its form by outputting the eqn
#rho = #

Find an expression for the central density in terms of \(R_\star\) and \(M_\star\).

Recall the formulation for the mass equation

\[ M(r) = \int^{r}_{0} 4 \pi r^2 \rho dr \]

  1. First, integrate the mass equation to find \(M(r)\) using sympys integrate function evaluated from 0 to \(R_\star\) (defined above).

# 1. solution here for M

# M =

  1. Next, we want to solve for \(\rho_{c}\) by setting \(M_\star\) to the \(M\) we computed in 1 using sympys Solve function.

    Hint: in sympy you can set define an equality using Eq relation function. For example to set the LHS of an equation to the RHS, you would write Eq(LHS,RHS).

# solve the equation M = Mstar for rhoc, we will need to name this variable different than we defined above,
# such as rhoc_value to distinguish it from the symbolic definition above. 
rhoc_value = #solve(Eq(# finish here
#rhoc_value

  Cell In[5], line 3
    rhoc_value = #solve(Eq(# finish here
                 ^
SyntaxError: invalid syntax

  1. Next, confirm that the average density of this star is \(\bar{\rho}=0.4\rho_{c}\), assuming \(\bar{\rho}=M/V\).

Hint: for fractions in sympy you would use the Rational function.

# substitute our found rhoc_value expression into our original rho expression here, 
# then output the function to verify
#rhobar = symbols(r"\bar{\rho}")
#rhobar = #
#rhobar

  1. Next, use substition again to plug in the solar radius and mass in cgs and evaluate \(\rho_{c}\) using evalf to make an estimate for the central density with the equation we found in 2.

msun = 2e33    # g
rsun = 6.96e10 # cm

# compute a numerical value here by substituting in solar radius and mass, then evaluate
#rhoc_value = rhoc_value.subs(###).evalf()
#rhoc_value

  1. Confirm that \(\bar{\rho}\approx 1.4 \ (\rm{g \ cm}^{-3})\) for this density profile using the numerical value found in 4. and Subs. Similar to the Sun.

# compute a numerical value for rhobar here for this density profile and the rhoc found above
Contents