% $Id: notes.tex,v 1.1 2017/03/12 04:07:29 brandenb Exp $
%\documentclass{article}
\documentclass[twocolumn]{article}
\setlength{\textwidth}{160mm}
\setlength{\oddsidemargin}{-0mm}
\setlength{\textheight}{220mm}
\setlength{\topmargin}{-8mm}
\usepackage{graphicx,natbib,bm,url,color}
\graphicspath{{./fig/}{./png/}}
\thispagestyle{empty}
\input macros0
\input journ
\def\red{\textcolor{red}}
\def\blue{\textcolor{blue}}
\title{Test problem with Eddington approximation}
\author{}
\date{\today,~ $ $Revision: 1.1 $ $}
\begin{document}
\maketitle

\section{Radiative transfer}

The equation of radiative transfer is
\begin{equation}
\nnn\cdot\nab I=-\kappa\rho\,(I-S).
\label{dI}
\end{equation}
where $I=I(\xx,t,\nnn,\nu)$ is the intensity,
$\rho$ is the density, $\kappa$ is the intensity,
and $S=(\sigmaSB/\pi)\,T^4$ is the source function.

\EEq{dI} can be solved by taking moments.
We define 0th, 1st, and 2nd moments as follows,
\begin{equation}
J={1\over4\pi}\int_{4\pi}I\;\dd\Omega,
\end{equation}
\begin{equation}
\FF=\int_{4\pi}I\nnn\;\dd\Omega,
\end{equation}
\begin{equation}
\PPPP={1\over4\pi}\int_{4\pi}I\nnn\nnn\;\dd\Omega.
\end{equation}
Here, $J$ is the mean intensity, $\FF$ is the radiative flux, and
$\PPPP$ is the radiation pressure tensor.
Taking the 0th and 1st moment of \Eq{dI} yields
\begin{equation}
\nab\cdot\FF=-4\pi\kappa\rho\,(J-S),
\label{RadFluxDiv}
\end{equation}
\begin{equation}
\nab\cdot\PPPP=-{\kappa\rho\over4\pi}\FF.
\label{RadPressTensor}
\end{equation}
Making the closure assumption
\begin{equation}
{\sf P}_{ij}=\onethird\delta_{ij}J.
\end{equation}
\Eq{RadPressTensor} becomes
\begin{equation}
\onethird\nab J=-{\kappa\rho\over4\pi}\FF.
\label{FF}
\end{equation}
Next, to calculate $\nab\cdot\FF$, we divide first by
$\kappa\rho$ and then take the divergence, so
\begin{equation}
\nab\cdot\left({1\over3\kappa\rho}\nab J\right)=-{1\over4\pi}\nab\cdot\FF.
\label{nabFF}
\end{equation}
Using \Eq{RadFluxDiv} we obtain a closed equation for $J$.
\begin{equation}
\nab\cdot\left({1\over3\kappa\rho}\nab J\right)=\kappa\rho\,(J-S).
\label{EddingtonApprox}
\end{equation}

\section{Basic hydro equation}

We are interested in steady solutions to the time-dependent equations
\begin{equation}
{\partial\rho\over\partial t}=-\nab\cdot\rho\uu
\end{equation}
\begin{equation}
\rho {\DD u\over\DD t}=\nab P+\rho\grav+\mbox{visc force}
\label{F}
\end{equation}
\begin{equation}
\left(\rho T{\DD S\over\DD t}=\right)\quad
\rho {\DD E\over\DD t}+P\nab\cdot\uu=-\nab\cdot\FF
\label{F}
\end{equation}
where $E=\cv T$ is the internal energy with $\cv=\const$.
Using \Eq{FF}, we have
\begin{equation}
\FF=-{4\pi\over3\kappa\rho}\nab J,
\end{equation}
or, more directly using \Eq{RadFluxDiv}
\begin{equation}
\nab\cdot\FF=-\nab\cdot{4\pi\over3\kappa\rho}\nab J.
\end{equation}

\begin{figure*}[t!]\begin{center}
\includegraphics[width=\textwidth]{psolution}
\end{center}\caption[]{
Solution to \Eqs{ToverT0}{z}.
}\label{psolution}\end{figure*}

\section{Equilibrium solution}

Using \Eq{nabFF}, the equilibrium solution, $\nab\cdot\FF$, satisfies
$J=S=(\sigmaSB/\pi)\,T^4$, and therefore
\begin{equation}
0=\nab\cdot\FF=-\nab\cdot{4\sigmaSB\over3\kappa\rho}\nab\,T^4
=-\nab\cdot K\nab\,T
\end{equation}
where
\begin{equation}
K={16\sigmaSB T^3\over3\kappa\rho},
\end{equation}
so $K\dd T/\dd z=\const$ for the equilibrium solution.
In the following, we use
\EQ
\kappa=\kappa_0 (\rho/\rho_0)^a (T/T_0)^b,
\label{Kramers}
\EN
where $a$ and $b$ are adjustable parameters, $\rho_0$ and $T_0$ are
reference values for density and temperature, respectively, and $\kappa_0$
gives the overall magnitude of the opacity.

An analytic solution that also satisfies hydrostatic equilibrium
is given by [see Equation~(68) of \cite{Bra16}]
\EQ
{T\over T_0}
=\left[(n+1)\nabrad^{(0)}\left({P\over P_0}\right)^{1+a}
+\left({T_{\rm top}\over T_0}\right)^m\right]^{1/m}\!\!\!\!\!\!\!,\;
\label{ToverT0}
\EN
where $m=4+a-b$,
$\nabrad^{(0)}=F_{\rm rad}P_0/(K_0T_0\rho_0 g)$,
$K_0=16\sigmaSB T_0^3/(3\kappa_0\rho_0)$ are constants,
and $T_{\rm top}=\Teff/2^{1/4}$ is an integration constant
that is specified such that $T\to T_{\rm top}$ as $P\to0$.
Note also that $4+a-b=(n+1)(1+a)$, where
\EQ
n=(3-b)/(1+a)
\label{n_Def}
\EN
is the polytropic index, so the ratio of $4+a-b$ to $1+a$ is just $n+1$,
which enters in front of the $\nabrad^{(0)}$ term in \Eq{ToverT0}.
Since $K\propto T^{3-b}/\rho^{1+a}\propto T^{4+a-b}/P^{1+a}$,
we have $K\to\const=K_0$ for $T\gg T_{\rm top}$.
Finally, the $z$ coordinate is obtained by numerical integration,
\EQ
z=-\int {\dd P\over\rho g}
=-\int {P\over\rho g}\,\dd\ln P.
\label{z}
\EN

\section{Example of BB14}

Consider, as a specific example, Run~B6 of \cite{BB14} with $a=1$, $b=0$,
$\kapz=10^6\Mm^{-1}\cm^3\g^{-1}$ and $\Teff=9300\K$.
where $\rho_0=4\;10^{-4}\,\g\,\cm^{-3}$ and $T_0=38,968\K$ are used.
The result is plotted in \Fig{psolution}.

%r e f
\begin{thebibliography}{}

\bibitem[Barekat \& Brandenburg(2014)]{BB14}
Barekat, A., \& Brandenburg, A.\yanaN{2014}{571}{A68}
{Near-polytropic stellar simulations with a radiative surface}

\bibitem[Brandenburg(2016)]{Bra16}
Brandenburg A.\yapjN{2016}{832}{6}
{Stellar mixing length theory with entropy rain}

\end{thebibliography}

\vfill\bigskip\noindent\tiny\begin{verbatim}
$Header: /var/cvs/brandenb/tex/notes/eddington/notes.tex,v 1.1 2017/03/12 04:07:29 brandenb Exp $
\end{verbatim}

\end{document}
