{ "cells": [ { "cell_type": "markdown", "id": "27aa8974-0487-45af-aa90-275ce6c7935a", "metadata": {}, "source": [ "\n", "# White dwarfs" ] }, { "cell_type": "markdown", "id": "78f4bc8d-2d6d-4fd2-848b-8f7b00a7c2c6", "metadata": {}, "source": [ " White dwarfs are the final evolutionary state of stars whose mass is\n", " not high enough to become a neutron star or a black hole. After the\n", " hydrogen–fusing ($\\mathrm{H} \\to \\mathrm{He}$) lifetime of a star\n", " ends, such a star fuses helium to carbon and oxygen,\n", " ($\\mathrm{He} \\to \\mathrm{C, \\, O}$). If a star has insufficient\n", " mass to generate the core temperatures required to further fuse\n", " carbon and oxygen, an inert mass of carbon and oxygen builds up at\n", " its center. After shedding its outer layers, the star leaves behind\n", " the core, which is the white dwarf.\n", "\n", " Since the material in a white dwarf no longer undergoes fusion\n", " reactions, the star is not supported against gravitational collapse\n", " by the heat generated by fusion. It is supported only by a much\n", " weaker electron gas pressure. Therefore, the star collapses into an object\n", " very small size and extremely high density.\n", "\n", " Surprisingly, the larger is the mass of a white dwarf, the smaller\n", " is it radius. There is a characteristic mass, called\n", " *Chandrasekhar mass* or *Chandrasekhar limit*, above which\n", " electron degeneracy pressure in the star's core is insufficient to\n", " balance the star's own gravitational self-attraction. A star with a\n", " mass greater than the limit is evolving into a neutron star or black\n", " hole. Chandrasekhar limit corresponds to the point where the graph\n", " of radius of white dwarf vs mass, $r(m)$ crosses the $m$ axes." ] }, { "cell_type": "code", "execution_count": null, "id": "1a5306a5", "metadata": {}, "outputs": [], "source": [ "\n", "using OrdinaryDiffEqTsit5\n", "using PyPlot" ] }, { "cell_type": "markdown", "id": "539fb344-c8bf-493c-8573-985597eacaad", "metadata": {}, "source": [ "\n", "The system of 2 dimensionless differential equations\n", "describing the radial distribution of the density, $\\rho(r)$,\n", "and mass, $m(r)$, inside a white dwarf star:\n", "\n", "$$\\frac{\\mathrm{d}m}{\\mathrm{d} r} = \\rho r^2, $$\n", "$$\\frac{\\mathrm{d}\\rho}{\\mathrm{d} r} = - \\frac{m \\rho}{\\gamma(\\rho) r^2},$$\n", "\n", "where \n", "\n", "$$\\gamma(\\rho) = \\frac{\\rho^{2/3}}{3 \\sqrt{1 + \\rho^{2/3}}}.$$" ] }, { "cell_type": "markdown", "id": "144e52ca-4f44-4ce3-9e01-ee941b510392", "metadata": {}, "source": [ "In the equations above the density is measured in units of $\\rho_0$,\n", "\n", "\\begin{equation}\n", " \\rho_0 = \\frac{M_p \\, m_e^3 \\, c^3}{3 \\,\n", " \\pi^2 \\, \\hbar^3 Y_e} = 9.82 \\times 10^8 \\, Y_e^{-1} \\; \\mathrm{kg}\n", " \\, \\mathrm{m}^{-3}, \n", "\\end{equation}\n", "\n", "the distances are measured in units of $R_0$,\n", "\n", "\\begin{equation}\n", " R_0 = \\left( \\frac{m_e\\,c^2\\,Y_e}{4\\,\\pi\\,\\rho_0\\,G\\,M_p} \\right)^{\\frac{1}{2}} \n", " = 7.71 \\times 10^6\\,Y_e \\; \\mathrm{m},\n", "\\end{equation}\n", "\n", "and the mass is measured in units of $M_0$\n", "\n", "\\begin{equation}\n", " M_0 = 4\\,\\pi\\,R_0^3\\,\\rho_0 = 5.66 \\times 10^{30} \\, Y_e^2 \\; \\mathrm{kg}.\n", "\\end{equation}\n", "\n", "where $M_p$ is the mass of the proton, $m_e$ is the mass of of the\n", "electron, $Y_e$ is the number of electrons per nucleon, $c$ is the\n", "speed of light, $\\hbar$ is the Planck constant." ] }, { "cell_type": "markdown", "id": "310bc898-1853-4228-bd33-58f20def7476", "metadata": {}, "source": [ "\n", "Below we consider a white dwarf star consisting of $^{12}$C, a chemical\n", "element with 6 protons, six neutrons, and six electrons, then\n", "$Y_e = \\frac{1}{2}$ and $M_0 = 0.71 \\times M_{\\odot}$ and $R_0 = 0.006\n", "\\times R_{\\odot}$, where $M_{\\odot}$ and\n", "$R_{\\odot}$ are the mass and the radius of the Sun." ] }, { "cell_type": "markdown", "id": "629a9489-23aa-43cf-b8a5-df6125ddb08a", "metadata": {}, "source": [ "\n", "The pair of equations is integrated from $r =\n", "0,\\,\\rho = \\rho_c,\\,m = 0$ to the value of $r$\n", "at which $\\rho = 0$, which defines the dimensionless radius of\n", "the star $R$, and the dimensionless mass of the star is then\n", "$M = m(R)$." ] }, { "cell_type": "markdown", "id": "3ed9fea5-1770-45aa-abac-c65e22a420da", "metadata": {}, "source": [ "To avoid numerical difficulties in calculating the right hand side of\n", "the equation for $\\frac{\\mathrm{d}\\rho}{\\mathrm{d} r}$ for small values of $r$, notice that for\n", "small $r$\n", "\\begin{equation}\n", " m(r) \\approx \\frac{1}{3} r^3 \\, \\rho_c.\n", "\\end{equation}\n", "where $\\rho_c$ is the density of the material in the center of a white dwarf.\n", "\n", "Hence, for small $r$ the equation can be written in the following\n", "form:\n", "\n", "\\begin{equation}\n", " \\frac{\\mathrm{d}\\rho}{\\mathrm{d}r} =\n", " - \\, \\frac{r \\, \\rho_c^2}{3 \\gamma(\\rho_c)}. \n", "\\end{equation}\n", "\n", "which avoids diverging factor $1/r^2$.\n" ] }, { "cell_type": "code", "execution_count": null, "id": "1a391c2e", "metadata": {}, "outputs": [], "source": [ "\n", "\"\"\"\n", " white_dwarf_eqs!(dudr, u, p, r)\n", "\n", "The right hand side of the system of 2 dimensionless differential \n", "equations describing the radial distribution of the density, rho(r),\n", "and mass, m(r), inside a white dwarf star\n", "\n", " m' = rho r^2\n", " rho' = -m rho /(gamma(\\rho) r^2)\n", "\n", "where gamma(rho) = rho^(2/3)/(3 sqrt(1 + rho^(2/3)))\n", "\"\"\"\n", "function white_dwarf_eqs!(dudr, u, p, r)\n", " m = u[1]\n", " rho = u[2]\n", " rho_c = p\n", " if rho >= 0.0\n", " w = rho^(2/3)\n", " gamma = w/(3 * sqrt(1 + w))\n", " dudr[1] = rho * r * r\n", " if (r > 1.e-6)\n", " dudr[2] = -m * rho/(gamma * r * r)\n", " else\n", " dudr[2] = -rho_c/3 * r * rho/gamma\n", " end\n", " else\n", " dudr[1] = 0.0\n", " dudr[2] = 0.0\n", " end\n", "end" ] }, { "cell_type": "code", "execution_count": null, "id": "3993e7bd-eb6e-4860-bd18-148ea508cf08", "metadata": {}, "outputs": [], "source": [ "\n", "rho_c1 = 0.08 # min density in the star center\n", "rho_c2 = 100000.0 # max density in the star center\n", "np = 100\n", "\n", "rhos = exp.(range(log(rho_c1), log(rho_c2), np));" ] }, { "cell_type": "markdown", "id": "9e257e9d-fcd9-4e01-9604-07c758889190", "metadata": {}, "source": [ "Preallocate storage:" ] }, { "cell_type": "code", "execution_count": null, "id": "e1725d1e-6982-47e2-a3d2-c0b1bb95c43a", "metadata": {}, "outputs": [], "source": [ "\n", "radi = zeros(np)\n", "mass = zeros(np);" ] }, { "cell_type": "code", "execution_count": null, "id": "6afd6149-93e8-49ef-b4df-3c9269484c43", "metadata": {}, "outputs": [], "source": [ "\n", "rspan = (0.0, 10.0)" ] }, { "cell_type": "markdown", "id": "09ba7786-9e67-403c-819c-4fdc5ff8cee7", "metadata": {}, "source": [ "Main loop:" ] }, { "cell_type": "code", "execution_count": null, "id": "08191845-0cc1-4f75-ade6-569a08ab2885", "metadata": {}, "outputs": [], "source": [ "\n", "for (i, rho_c) in enumerate(rhos)\n", "\n", " # New initial conditions\n", " u0 = [0.0, rho_c]\n", "\n", " # Set up the new problem\n", " prob = ODEProblem(white_dwarf_eqs!, u0, rspan, rho_c)\n", "\n", " # Set up the callback\n", " condition(u, t, integrator) = u[2] # stop integration when rho = 0.\n", " affect!(integrator) = terminate!(integrator)\n", " cb = ContinuousCallback(condition, affect!)\n", "\n", " # Integrate the ODEs\n", " sol = solve(prob, Tsit5(), callback=cb)\n", "\n", " # Store the relevant part: the radius and the mass of the white dwarf.\n", " radi[i] = sol.t[end]\n", " mass[i] = sol(radi[i])[1]\n", " \n", "end" ] }, { "cell_type": "markdown", "id": "b40761fb", "metadata": {}, "source": [ "Download white dwarfs observational data from \\\"The VizieR database of astronomical catalogues" ] }, { "cell_type": "code", "execution_count": null, "id": "3fbd1c16-552b-4bc6-b4be-913356e6720e", "metadata": {}, "outputs": [], "source": [ "\n", "using DataFrames\n", "using CSV" ] }, { "cell_type": "code", "execution_count": null, "id": "d80fa8d4", "metadata": {}, "outputs": [], "source": [ "\n", "url = \"http://vizier.u-strasbg.fr/viz-bin/asu-tsv?-source=J/A+A/420/507/tablea1&-out=Mass&-out=Rad\"\n", "catalog = download(url)\n", "df = CSV.read(catalog, DataFrame, comment=\"#\", skipto=4, \n", " types=[Float64, Float64], delim='\\t', ignorerepeated=true, \n", " silencewarnings=true)\n", "dropmissing!(df);" ] }, { "cell_type": "code", "execution_count": null, "id": "5d605175", "metadata": {}, "outputs": [], "source": [ "\n", "# The constants convert our dimensionless units to the units of the solar mass and solar radius\n", "plot(mass*0.71, radi*0.006, label=\"theory\") \n", "plot(df.Mass, df.Rad, marker=\"o\", linestyle=\"none\", \n", " fillstyle=\"none\", markersize=4, label=\"observations\")\n", "grid(true)\n", "xlabel(L\"$M/M_{\\odot}$\")\n", "ylabel(L\"$R/R_{\\odot}$\")\n", "legend();" ] }, { "cell_type": "code", "execution_count": null, "id": "ab1fae9a", "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Julia 1.10.5", "language": "julia", "name": "julia-1.10" }, "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.10.5" } }, "nbformat": 4, "nbformat_minor": 5 }