\subsection{Hydrodynamical Test Problems}

\subsubsection{Noh's shock tube}
In this shock tube problem a piston of constant speed $v = -1$ pushes a cold gas with 
$\gamma = {5\over3}$ into the origin. Courant \& Friedrichs~[1] have discussed this problem 
qualitatively and showed that an instantaneous shock of infinite Mach number is generated. 
The shock speed turns out to be ${1\over2}(\gamma-1)|v|={1\over3}$. Noh~[9] defined this 
problem as a challenging benchmark for hydrodynamics codes that can be run in cartesian, 
cylindrical, and spherical geometry. In addition, there are analytical formulae against 
which the numerical results can be tested. We follow Noh's examinations at time $t=0.6$ 
before the shock collides with the piston. 

We represent the computational domain with 100 grid points and employ the density, pressure 
and energy in logarithmic scaling as our grid parameters together with $\alpha = 1.5$ and 
$\tau = 10^{-8}$. We obtain solutions to the Noh problem in all three geometries. 
Let us first consider the planar case which we calculate with a 
fixed and the adaptive grid. The superiority of the adaptive grid is well 
established: since the shock is the only nonlinear feature in the flow it is resolved with 
some 60 zones. The time $t=0.6$ is reached after 160 time steps while constraining the time
step to $\le 10^{-2}$. The adaptive grid is locked into the shock frame and thus the code 
evolved almost statically. In contrast to that the shock has to move from one grid point to
the next in the Eulerian computation thus allowing the time step to only vary between 
$10^{-5}$ and $10^{-4}$ (note that this still orders of magnitudes larger than the 
practical CFL number!). This computation exhibits all the problems induced by the use 
of the artificial viscosity as discussed in~[9]. The adaptive grid, however, overcomes 
these difficulties very early on ($t<10^{-2}$). % after about 90 time steps. 

\psfig{file=noh.ps,width=\textwidth}
\vspace{-2.0in}
{\tenrm Fig. 1 (Noh's shock tube): density ($+$), energy ($\Diamond$), and pressure ($\triangle$) 
profiles in planar, cylindrical,}

\hspace{10mm}{\tenrm and spherical geometry.}
\vspace{5mm}

In figure~1 we plot profiles of the density, pressure, and energy at time $t=0.6$ for
planar, cylindrical, and spherical geometry. The adaptive grid behaves very similar
in all cases producing sharp steps at the shock front. The step in energy is always
${1\over2}$ because of the isothermal nature of the shock. The step in density equals to
$3\cdot4^\mu$. The density profiles in the unshocked regions satisfy the analytical 
relations: $\rho(r,t) = (1+t/r)^\mu$ which are valid for $t<{3\over4}$ and $r>{1\over3}t$. 
In the shocked regions the density plateau at $\rho(r,t) = 4^{\mu+1}$ for 
$0\le r\le {1\over3}t$ is well maintained in planar geometry and closely established in 
cylindrical and spherical geometry. We have evaluated the relative errors and found them to 
be between one tenth of a and a few percent in all cases. Since the theoretical plateau 
rises rapidly with $\mu$ the error becomes most visible in spherical geometry. We would 
like to remind that fixed grid computations usually generate errors of order unity for 
this problem! Therefore their solutions do not resemble the correct answer. The adaptive
grid allows, for the first time, to generate good numerical results with little computational
effort.

Next we simulate two different kinds of blast waves, one in planar and the other in
spherical geometry.

\subsubsection{Woodward's interacting blast waves}
First we repeat a planar problem of two interacting blast waves that was invented by 
Woodward and is described in~[14]. The computational domain is divided into three sections 
with three different but constant pressure regimes in them. In the left part we have 
$P(0\le r\le {1\over10}) = 1000$, in the middle part $P({1\over10}< r< {9\over10}) = 0$, 
and in the right part $P({9\over10}\le r\le 1) = 100$ with an ideal gas resting at unit 
density. The adiabatic index is $\gamma = 1.4$. The density is initially constant, $\rho = 1$, 
across the domain $0\le r\le 1$. We adopt again the same grid parameters as 
before but this time run the code with 200 grid points. The artificial viscosity
is kept linear at $\ell_o = 10^{-4}$.
\vspace{2mm}

{\hbox{%
{\vbox{\hsize=3.0in \strut
The large pressure steps at both ends of the cylinder create two shocks which travel 
towards the middle. The graph to the right (figure 2) shows their density profiles at an
early stage ($+$'s marking individual grid points). One can see that the left front is 
much stronger than the right one. Further, rarefaction waves begin to appear in form 
of the density dips, which connect to the unit density near the walls.\strut}}
{\vbox{\hsize=10mm}}
{\vbox{%
\psfig{figure=wood0.ps,width=3.0in}
{\tenrm Fig. 2 (Woodward's blast waves): density profile}

\hspace{10mm}{\tenrm at time t = 10$^{-4}$.}
}}
}}
\vspace{5mm}

In the following figure, density profiles at four different times $t\in($1.57, 2.94, 3.85, 
5.00$)\times 10^{-2}$ are plotted. The top left graph shows the approach of the two shocks. 
Note in particular that the left one is much faster than the right one. The fast shock has 
a resolution of $10^5$ and the slow shock one of 
several $10^4$. The density between the shocks is still at the $\rho = 1$ plateau thus 
indicating that they move independently. The shocks are followed by contact discontinuities, 
both being resolved by $10^3$. They connect to the rarefaction fans which stretch to the walls.
The high resolution of the discontinuities leaves only a crude representation of the 
rarefaction fans with resolution of about $20$. It is by choice that the adaptive grid is 
allowed to zone only into the shocks and has to neglect the rarefaction waves. At time
$t=0.028$ the shocks collide near position 0.69 and penetrate each other. When they separate
a new third contact discontinuity forms in between them (lower right plot). At $t=0.05$
the fast shock has reflected from the right wall but has not yet run into the third 
contact discontinuity.\\

\centerline{\psfig{figure=wood1.ps,width=5.5in}}
{\tenrm Fig. 3 (Woodward's blast waves): shock profiles at four different times 
($+$'s indicate grid points).}
\vspace{5mm}

A complete overview of the actions of all nonlinear waves can be obtained by investigating
the space-time diagram of the density. The respective contours in logarithmic scale are
drawn in figure 4. The shock fronts (due to their high resolutions) appear as thick lines, 
the contact discontinuities as broad lines, and the rarefaction waves as an ensemble of 
single contour lines.\\

\centerline{\psfig{figure=wood2.ps,width=5.5in}}
{\tenrm Fig. 4 (Woodward's blast waves): space-time diagram of the logarithmic density; 
25 contour lines have}

\hspace{10mm}{\tenrm been selected equidistantly.}
\vspace{5mm}

The five density profiles (figures 2 \& 3) are horizontal cross sections at the respective
times. For times $t<0.3$ there are the fast left traveling and slow right traveling shocks
(speed given by the slope of the thick contour lines) each followed by a contact
discontinuity which connects to a rarefaction fan. The fans are reflecting soon from the
walls indicating that the region near the walls are being depleted quickly. At $t=0.028$
the two shocks collide and slow down after the interaction while creating a secondary contact
discontinuity. Shortly after the fast shock runs into the primary discontinuity from the
right and accelerates again to near its original speed while the discontinuity is reflected.
At $t=0.045$ this shock reflects at the right wall and interacts with the discontinuity at
about $t=0.05$. Both nonlinear waves are slowed down considerably. The situation is comparable
Sod's shock tube. New, however, is that the shock interacts with the secondary discontinuity 
and speeds up again. Shortly before $t=0.045$ the slower right propagating shock collides
with the other primary discontinuity and is greatly accelerated while the discontinuity is
reflected. This shock now propagates uninterruptedly to the left wall. There it reflects 
at time $t=0.063$ and moves with similar speed towards the right again until it runs again
into the discontinuity. This time, however, the other shock front meets also with them leading
to a multiple interaction of nonlinear waves.\vspace{2mm}

A comparison of {\TITAN}'s computation with the original one by Woodward and Collela (using a 
special version of the PPMLR scheme on a domain of 3096 zones) tells that the shocks are 
well represented but the contact discontinuities lack sharpness. This can be accounted for 
in part to a fairly large artificial mass and energy diffusion which allowed our code to 
complete the run in a little over 1000 time steps (about a factor 5 less than the number of 
time steps PPMLR required). For comparison reasons the space-time diagram of the grid motion 
is presented in the following figure.\\ 

\centerline{\psfig{figure=wood3.ps,width=5.5in}}
{\tenrm Fig. 5 (Woodward's blast waves): space-time diagram of the grid motion; 
the trajectory of every fourth}

\hspace{10mm}{\tenrm zones is depicted.}
\vspace{5mm}

Figure~5 demonstrates well that the adaptive grid keeps track of all discontinuities. 
Note the strong correlation between the grid motion and the density contours.
Its typical behavior at the momentary absence of the shock(s) due to reflection from a wall
or shock-shock interaction is also evident.
In conclusion of this exercise we can highlight the generic behavior of the adaptive grid: 
Once the nonlinear features of the flow are resolved the adaptive grid locks into them and 
the code can evolve quasi-statically at time steps that are considerably larger than those 
of a fixed grid. When the grid needs to redistribute, as during the reflection of the shock 
from the wall or its interaction with the other shock or contact discontinuity, the time step 
drops by several orders of magnitude to accomplish this in a very short time interval. 

\subsubsection{Sedov-Taylor self-similar blast wave}
We conclude our hydrodynamics benchmarking with the Sedov-Taylor blast wave. 
An energy deposition of $10^{50}~ergs$ at the center of a ball with radius $R_{max} = 10^9~km$
generates a spherical shock that propagates in a self-similar fashion. The gas is assumed
to have a polytropic index $\gamma = {5\over3}$ and is initially at rest with a constant
density of $10^{-8}~g/cm^3$ and temperature of $50~K$. This time, density and
pressure serve as the only grid parameters. Logarithmic scaling is choosen for all quantities
involved. The grid $\alpha$ is again kept at 1.5, the shock width is $\ell_1 = 10^{-4}$
and we take 100 zones to cover the ball. \vspace{2mm}

{\hbox{%
{\vbox{\hsize=3.0in \strut
This problem is set up so that the shock wave reaches the surface of the ball 
after $t = 2\times 10^5~sec$. The code completes the calculation in about 1000 time steps. The
resolution of the shock front is about $10^6$ during the evolution. In the figure to the
right the evolution of this blast wave is illustrated. Snapshots of the density and the 
velocity at five different times $t\in($140, 2900, 26590, 93023, 197670$)~sec$ are overlayed.
The upwind density is always $10^{-8}$ $g/cm^3$ with the material at rest. At the shock front 
the density jump remains constant $4\times 10^{-8}~g/cm^3$. The velocity jump decreases
from $117$ $km/sec$ to $1.5~km/sec$ during these five instances. Both the density and the 
velocity profiles fall off to zero behind the shock.\strut}}
{\vbox{\hsize=10mm}}
{\vbox{%
\psfig{figure=sedov1a.ps,width=3.0in}
\psfig{figure=sedov1b.ps,width=3.0in}
{\tenrm Fig. 6 (Sedov-Taylor blast wave): density and}

\hspace{10mm}{\tenrm velocity profiles at five consecutive times.}}}\\
}}
\vspace{5mm}

The position, density, pressure, temperature, and velocity of the shock as a function of 
time are given in next figure. All variables vary strictly monotonically except for the 
shock density which remains at a steady value $\rho_s = 4\times 10^{-8}~g/cm^3$ as it should. 
In particular the expected relations
$$       R_s(t) \propto t^{2/5}                                   \eqno{(31)}$$ 
$$       v_s(t) \propto t^{-3/5}                                  \eqno{(32)}$$
$$       T_s(t) \propto P_s(t) \propto t^{-6/5}                   \eqno{(33)}$$
are satisfied. These are the analytical solution as discussed in the textbook of Landau 
\& Lifshitz.\\

\centerline{\psfig{figure=sedov2.ps}}
{\tenrm Fig. 7 (Sedov-Taylor blast wave): position R$_s$~(km), density $\rho_s$~(10$^{-8}$ g/cm$^3$), 
temperature T$_s$~(K),}

\hspace{10mm}{\tenrm pressure P$_s$~(dyn/cm$^2$), and speed v$_s$~(km/sec) of the shock as a 
	      function of time; notice the}

\hspace{10mm}{\tenrm power law dependencies.}
\vspace{5mm}

The self-similarity of the computed flow is demonstrated in figure~8. It shows the profiles of
density, pressure, and velocity in terms of the position, all measured with respect to the 
shock data. The graphics were generated by overlaying the actual grid values at six different 
times $t\in($140, 799, 12090, 42398, 110070, 174900$)~sec$. One can verify that the shocked
region obeys the analytical relations: 
$$          v/v_s\propto ~R/R_s~~~~~~~~~~                         \eqno{(34)}$$
$$    \rho/\rho_s\propto (R/R_s)^{3/(\gamma-1)}                   \eqno{(35)}$$
and similarly for the pressure. Their functions jump to the preshock values and stay at a 
constant level. Compare this figure with the one in Landau \& Lifshitz!\\

\centerline{\psfig{figure=sedov3.ps}}
{\tenrm Fig. 8 (Sedov-Taylor blast wave): self-similarity in the profiles of the density 
$\rho/\rho_s$, pressure P/P$_s$, and}

\hspace{10mm}{\tenrm velocity v/v$_s$ generated through overlays from six different 
times t $\in$ (140, 532, 576, 1160, 2290,}

\hspace{10mm}{\tenrm 4540) sec.}

\subsection{Radiation Test Problems}
In this section we follow L. Ensman's suggestions for radiative tests~[4]. We consider spheres 
of about $0.5~M_\odot$ and radial extension of $R=2.84\times 10^3~R_\odot$ that are filled
with pure hydrogen at a constant density $\rho = 2.9677\times 10^{-11}~g/cm^3$. Further
we assume the hydrogen to be fully ionized so that the opacity source are scattering
electrons; hence we set the mean flux opacity proportional to the Thomson scattering 
opacity: $\rho\kF = \chi_e = 1.19\times 10^{-11}~cm^{-1}$. In addition we always equate the
energy and Planck mean opacities, $\kE = \kP$, and allow them to differ from the mean flux
opacity by a constant factor, $\kE = \fr \kF$, in order to illustrate the effects of a 
scattering versus absorption dominated gas.

\subsubsection{Radiative heating tests}
We set up this sphere to be in diffusion equilibrium at a luminosity ${\cal L} = 
826~{\cal L}_\odot$. From the available data one computes the radiative energy according to
$$
  \er(r) = {{\cal L}\over{4\pi c R}}~
           \Bigl[\sqrt{3} + 3\chi_e\Bigl({R\over r}-1\Bigr)\Bigr]            \eqno{(36)}
$$

\noindent
and from that the temperature profile via $T(r) = {^4}\!\sqrt{\er(r)/\ar}$. At time $t=0$ 
we turn on a strong light source at the center with a $10^{1.5}$ times higher luminosity, 
${\cal L}_c=26,130~{\cal L}_\odot$, which illuminates the gas sphere. The sphere now has to
evolve into a new radiative equilibrium at this higher luminosity. 
\vfill

\centerline{\psfig{figure=heat0.ps,width=5.5in}}
{\tenrm Fig. 9 (Heating towards radiative equilibrium): profiles of the radiant energy ($+$) 
	       at different times}
	       
\hspace{10mm}{\tenrm during the relaxation of a scattering dominated atmosphere with $\rho\kP$ = 10$^{-5}$ $\chi_e$.}
\vspace{5mm}

This graph illustrates how the higher central luminosity changes the profiles of the 
radiative energy in time. Initially the sphere is at equilibrium at the lower luminosity, 
${\cal L}_c = 826~{\cal L}_\odot$, which is represented by curve $\er(r)$ at the bottom 
($+$'s indicating individual grid points). This curve relaxes and the radiative energy 
increases until the new equilibrium at the higher luminosity, ${\cal L}_c = 
26,130~{\cal L}_\odot$, is reached. For both the initial and the final $\er(r)$ curve the
analytical solution from equation (36) were overlayed. The agreement between the analytical
(assuming diffusion) and numerical (assuming variable Eddington factors) result is remarkable.
This should be expected because the shell is optically thick except right underneath the
surface and hence the diffusion regime is established everywhere else.\vspace{2mm}

\centerline{\psfig{figure=heat1.ps,width=5.5in}}
{\tenrm Fig. 10 (Heating towards radiative equilibrium): profiles of the gas ($+$) and 
                  radiation (---) temperature}
		  
\hspace{12mm}{\tenrm at different times during the relaxation of a scattering dominated atmosphere with $\rho\kP$ = }

\hspace{12mm}{\tenrm 10$^{-5}$ $\chi_e$.}
\vspace{5mm}

The heating of the gas ($+$) and radiative (---) temperature by the central light source is 
given in figure~10. It shows the relaxation of their profiles as a function of the optical 
depth at various times. Both bottom curves are the initial condition in which the two 
temperatures agree because of equilibrium diffusion. When the light source is turned on a 
radiation wave moves through the sphere thereby heating up the gas. Since a weak coupling 
between the gas and the radiation field is assumed, $\fr = 10^{-5}$, the gas temperature rises 
noticebly slower than the radiative temperature. As soon as the higher luminosity value has 
reached the surface both temperatures agree once again and the new radiative equilibrium is 
established. Because most of the gas is optically thick the diffusion regime prevails and the 
temperature follows again the profile from equation (36). \vspace{2mm}

This experiment is repeated with different factors $\fr$. The sphere relaxes to the same 
equilibrium state every time. Figure~11 depicts the evolution of the surface luminosity for six
different choices $0\le \fr \le 1$. The correspondence between the symbols and the ratios $f_r$
for the gray opacity is as follows: $\fr = 0$ ($\Box$), $10^{-9}$ ($\triangle$), $10^{-7}$ 
($\Diamond$), $10^{-5}$ ($\ast$), $10^{-2}$ ($+$), 1 (---). \break
The relaxation of the temperature profiles from the previous figure is represented here by
the asterisk. Note that the time the radiation wave takes to reach the surface and hence be-
{\hbox{%
{\vbox{\hsize=3.0in \noindent \strut
gins to raise the luminosity is, with- in a factor of 2, the same in all cases, namely around 
$6\times 10^6~sec$. This number is in reasonable agreement  with the classical diffusion time 
scale $\chi_e R^2/c$ $\approx$ $1.4\times 10^7~sec$. The spheres with almost pure scattering 
opacities reach the higher luminosity plateau first. Those spheres with almost pure absorptive 
opacities attain it later.\strut}}
{\vbox{\hsize=10mm}}
{\vbox{%
\psfig{figure=heat2.ps,width=3.0in}
{\tenrm Fig. 11: Relaxation of the surface luminosity for}

\hspace{13mm}{\tenrm six different opacities}}}
}}
\vspace{5mm}

The experimental setup allows us further to investigate the different radiative 
transport schemes in {\TITAN}. There are three choices: 

\begin{description}
\item``equilibrium diffusion'':\\
The gas and radiative temperature are not distinguished, $\er = \ar T^4$, and the 
radiative flux is evaluated according to equation (7). The differentiation between the 
mean opacities becomes meaningless and only the total opacity is considered. 

\item``non-equilibrium diffusion'':\\
The radiative flux is still computed according to equation (7). But now the radiation 
energy $\er$ is computed from equation (6). Two temperatures are possible and the
absorptive character of the gas can be taken into account. Both diffusion cases demand a
constant Eddington factor $\fe = {1\over3}$. 

\item``full transport'':\\
The radiation field is described by the two equations (4) \& (6). Since they necessitate 
an estimate for $\fe$ one obtains also information about the angular distribution of the 
intensity field. 
\end{description}

The consistency of the three approaches can be checked in the context of 
the hydrogen spheres. The equilibrium diffusion computation can also be simulated with
non-equilibrium diffusion by setting $\kE = \kF$. Since the gas is optically thick (aside
right underneath the surface) the corresponding full transport calculation should behave
in the same manner. Performing the actual runs demonstrates a tight agreement between the
three transport schemes.

We calculate the relaxation with the full transport scheme under the assumption that the
absorptive opacity is considerably smaller than the total opacity: $\rho\kE = 10^{-5}~\rho\kF 
= 10^{-5}~\chi_e$.
This defines the case of a scattering dominated gas. The computation is carried out with
our adaptive grid. This time we take the radiation energy and the luminosity as the grid
parameters in logarithmic scaling and set $\alpha=3$, $\tau=10^3$. 

\subsubsection{Radiative cooling tests}
We now take the equilibrium model described above as the initial condition. At $t=0$ we turn the 
light source at the center off. This time we choose a fixed (Eulerian) grid instead of the adaptive 
grid and perform the runs with the full radiative transport scheme.
\vfill

\centerline{\psfig{figure=cool0.ps,width=5.5in}}
{\tenrm Fig. 12 (Cooling away from radiative equilibrium): profiles of the radiant energy ($+$) 
	       at different times}
	       
\hspace{12mm}{\tenrm during the relaxation of a scattering dominated atmosphere with $\rho\kP$ = 10$^{-7}$ $\chi_e$.}
\vspace{5mm}

This graph illustrates how the radiative energy decreases in time once the central luminosity 
is turned off. Initially the sphere is at equilibrium at ${\cal L}(0) = 26,130~{\cal L}_\odot$
which is represented by the curve $\er(r)$ ($+$'s again marking the grid points) at the top 
with the analytical solution overlayed. As time progresses, $\er$ decreases throughout the sphere 
as radiant energy leaks out, and the gradient of $\er$ flattens, consistent with ${\cal L}(t)$
deminishing to zero. At the surface a sharp gradient continues to drive ${\cal L}(t)$
outward. For $t\rightarrow\infty$ it goes to zero, $\er\rightarrow 0$, which represents the new 
equilibrium when ${\cal L}_c = 0$.

Just like before one can consider different strengths of the gas-radiation coupling.
{\hbox{%
{\vbox{\hsize=3.0in \noindent \strut
In figure~13 the decline of the surface luminosity for various gray opacities
defined by the ratios $0\le \fr\le 1$ are plotted. Again these ratios correspond to the
symbols in the following way: $\fr = 0$ ($\Box$), $10^{-9}$ ($\triangle$), $10^{-7}$ 
($\Diamond$), $10^{-5}$ ($\ast$), $10^{-2}$ ($+$), 1 (---). The curve of the diffusion 
regime ($\fr = 1$) is representative of all regimes with $\fr > 0$. Their analyses tells 
that the long time behavior 
\strut}}
{\vbox{\hsize=10mm}}
{\vbox{%
\psfig{figure=cool2.ps,width=3.0in}
{\tenrm Fig. 13: Relaxation of the surface luminosity for}

\hspace{13mm}{\tenrm six different opacties}}}
}}
\noindent 
the luminosity decay follows a 
power law, ${\cal L}(t) \propto t^{-s}$, with the power being close to $s={1\over3}$. In 
contrast to that the regime of complete decoupling between the gas and radiative temperature, 
$\fr = 0$, is characterized by an exponential decay, ${\cal L}(t) \propto e^{-t/s}$, with the 
time scale $s\approx 1.7\times 10^7~sec$ being close to the classical  diffusion time. This 
regime behaves differently because the cooling wave alters only the radiative temperature but 
leaves the profile of the gas temperature unchanged.\vspace{5mm}

In figure~14 the profiles of the gas ($+$) and radiative (---) temperature for 
the case $\fr = 10^{-7}$ are plotted. Initially they are in equilibrium but the decoupling 
begins immediately. One observes that the profile of the radiative temperature is merely 
shifted downward in time. Simultaneously the profile of the gas temperature is modified to 
become a constant. Since the bulk of the total opacity comes from electron scattering the 
radiation cools off faster than the gas so that the gas energy contribution to the luminosity 
rises soon to domination. Note that the radiative temperature is related to the radiative
energy via equation (36) and therefore shows the same behaviour as $\er$ in figure 12.

\centerline{\psfig{figure=cool1.ps,width=5.5in}}
{\tenrm Fig. 14 (Cooling away from radiative equilibrium): profiles of the gas ($+$) and 
                  radiation (---) temper-}

\hspace{12mm}{\tenrm ature at different times for a scattering dominated
                  atmosphere with $\rho\kP$ = 10$^{-7}$ $\chi_e$.}
 
\subsection{Radiation Hydrodynamics Test Problems}

The fundamental understanding of radiating shocks was developed by Zel'dovich \& Raizer and additional 
discussion of the phenomena involved is given by Mihalas \& Mihalas. 
This section describes three radiation hydrodynamics benchmarks. All tests are carried out
in spherical geometry. The full transport scheme for the radiative field is employed. 
We follow again L. Ensman's test setups.\vspace{2mm}

We first consider a thin shell with an inner radius $R_i = 8\times 10^6~km$ and outer
radius $R_o = 8.7\times 10^6~km$ filled with a cold gas at constant density $7.78\times
10^{-10}~g/cm^3$. Further, we assume a shallow temperature profile 
$$           T(r) = 10 + 75~{{r-R_o}\over{R_i-R_o}}~K .                \eqno{(37)}
$$
The Rosseland mean opacity is taken to be constant $\rho\kR=3.115\times 10^{-10}~cm^{-1}$, and 
all of it is absorptive, $\fr=1$. The gas is initially at rest. At time $t=0$ a piston at the inner 
radius is moved outwards at a constant high speed. A shock front is thus formed which couples strongly 
to the radiation field. If the maximum temperature immediately ahead of the shock is less than the 
downstream temperature one speaks of a ``subcritical'' shock, the case of equality one speaks of a
``supercritical'' shock.\vspace{2mm}

In the following calculations we represent the sphere by 100 grid points. Density and mass in 
logarithmic scaling are selected as the sole grid parameters along with
$\alpha=1.5$ and $\tau = 10^{-20}$. The shock width is taken to be proprotional to
the radius, $\ell_1 = 10^{-5}$, to demand a high resolution. All runs are executed fully 
implicit, \ie~with time centering parameter $\theta=1$.

\subsubsection{Subcritical shock}
The piston speed is $6~km/sec$. Figure~15 shows the radiative shock front propagating 
towards the surface. The four quantities gas $T$ and radiative temperature $T_r$, 
radiative flux $\Fr$, and density $\rho$ are plotted as functions of the radius at seven 
different times. The shock front is located where the flux is largest. The temperatures 
just ahead of the shock are always smaller than those in the shocked region thus defining 
this run as a subcritical shock. 
\vfill

\centerline{\psfig{figure=sub1.ps,width=5.5in}}
{\tenrm Fig. 15 (Subcritical shock): evolution of the gas and radiative temperature, the
                  radiative flux, and the }

\hspace{12mm}{\tenrm velocity profiles while a piston of speed 
                  6 km/sec pushes the left boundary towards right.}
\vspace{5mm}

{\hbox{%
{\vbox{\hsize=3.0in \strut
It is often useful to look at the dependencies of the physical quantities on other coordinates.
The following figures illustrate this for the example of the gas temperature. The same seven
snapshots were taken but instead of showing them as a function of the radius they are now plotted
as a
\strut}}
{\vbox{%
\psfig{figure=sub1a.ps,width=3.0in}}}
}}

{\hbox{%
{\vbox{\hsize=3.0in \noindent \strut
function of the exterior mass, of the optical depth, and of the grid point (zone 
index). The dependency $T=T(\Me)$ is the Lagrangean representation of the data which emphazises 
the outermost layers of the gas sphere. The graph $T=T(\tau)$ shows the temperature profile 
from a radiation transport perspective which stresses how transparent the features are. These 
functional dependencies help the physical interpretation of the results. On the other hand, 
the profile $T=$ \{$T_k$ ; $k=1,\dots,100$\} is interesting from a numerical point of view 
because one can learn from them how well the calculation was carried out. In the given example 
one can, for instance, see that the temperature jumps in the shock are well resolved by about 
10 grid points.\strut}}
{\vbox{\hsize=10mm}}
{\vbox{%
\psfig{figure=sub1b.ps,width=3.0in}
\psfig{figure=sub1c.ps,width=3.0in}
{\tenrm Fig. 16 (Subcritical shock): three different repres-}

\hspace{12mm}{\tenrm entations of the evolution of the tempera-}

\hspace{12mm}{\tenrm ture profile}}}
}}
\vspace{5mm}

Its characteristic features, as explained in Zel'dovich \& Raizer, are an overshoot of the 
gas temperature, a nearly symmetric profile of the radiative flux, and an extended 
non-equilibrium region upstream from the shock, which is generated by the effects of radiative 
preheating from the flux. In figure~17 the density, gas and radiative temperature, flux 
and Eddington factor (scaled so that 1.0 corresponds tp $\fe = {1\over3}$), and velocity 
are plotted in dependency on the optical depth away from the front, $(\tau-\tau_s)$, at the half 
piston crossing time ($t = 56,581~sec$ to be exact). \vspace{2mm}

\centerline{\psfig{figure=sub2.ps,width=5.5in}}
{\tenrm Fig. 17 (Subcritical shock): snap shots of the density, gas ($+$) and radiative
                  ($\Diamond$) temperature, radiative}
		  
\hspace{12mm}{\tenrm flux, and gas velocity as functions
                  of the optical depth $\tau$ relative to the optical depth of the}
		  
\hspace{12mm}{\tenrm shock $\tau_s$ at half piston crossing time; overlayed over the flux plot is
                  three times the Eddington}
		  
\hspace{12mm}{\tenrm factor ($\cdots$).}
\vspace{5mm}

A non-equilibrium diffusion analysis, \cf~Zel`dovich \& Raizer, demonstrates that the density, gas
pressure, radiative flux and energy all decrease upstream depend exponentially on the optical
depth away from the shock:
$$  \rho(\tau) \propto \Pg(\tau) \propto \Fr(\tau) \propto \er(\tau) 
	       \propto e^{-\sqrt{3}|\tau-\tau_s|}                  \eqno{(38)}
$$ 
One consequence is that the radiative temperature must decrease slower than the gas temperature. 
This functional dependency on the relative optical depth holds also downstream for the flux and 
the gas temperature. The maximum flux is determined by the postshock temperature $T_2\;$:
$F_s \approx {1\over2\sqrt{3}}\ar c T_2^4$. The analytical values differ a factor less than
2 from the actual numerical data. The behavior of our simulation agrees quite well with the
analytical estimates, and the existing deviations can be explained by the use of variable
Eddington factors.

\subsubsection{Supercritical shock}
The piston speed is $20~km/sec$. Figure~18 shows the radiating shock while it propagates 
across the shell. The temperatures in the regions immediately ahead of and trailing the shock 
front are equal. Therefore this run is classified as a supercritical shock. \vspace{2mm}

\centerline{\psfig{figure=sup1.ps,width=5.5in}}
{\tenrm Fig. 18 (Supercritical shock): evolution of the gas and radiative temperature, the
                  radiative flux, and}
		  
\hspace{12mm}{\tenrm the velocity profiles while a piston of speed 
                  20 km/sec pushes the left boundary towards right.}
\vspace{5mm}

Its characteristic features are again a strong overshoot of the gas temperature, this time a 
highly asymmetric flux profile, and an extended radiative precursor whose head is out of 
equilibrium with the gas. In contrast to the subcritical case the immediate surroundings of 
the shock satisfy the condition of radiative equilibrium. The temperature overshoot itself 
is out of equilibrium and its width is narrower than the unit mean free path of the photons. 
This implies that the radiation flow from this feature cannot be treated correctly by the
diffusion approximation. In figure~24 the density, gas and radiative temperature, flux and 
Eddington factor (scaled so that 1 corresponds to $\fe = {1\over3}$), and velocity 
are plotted in dependency on the optical depth away from the front, $\tau-\tau_s$, 
at one quarter piston crossing time ($t = 8,541~sec$). \vspace{2mm}

A non-equilibrium diffusion analysis according to Zel'dovich \& Raizer predicts again that
the radiative flux and radiative energy decline exponentially and so does the gas 
temperature in the upstream part of the precursor that is out of equilibrium. In the 
equilibrium region around the shock one finds for the gas temperature a dependency on 
the optical depth away from the front like
$$  T(\tau) \propto ~{^3}\!\sqrt{1 + {3\sqrt{3}\over4}|\tau-\tau_s|}       \eqno{(39)}
$$

There is also an estimate for the maximum flux which is not as accurate as for the 
subcritical case. Radiative equilibrium is well established within $(\tau-\tau_s)<-5$ aside 
of the shock front itself. Density, flux, and velocity undergo jumps from the up- to the 
downstream regions. Our numerical results reflect 
the analytical predictions. Because the analytical results are not as reliable as in the 
previous case a close comparison does not give good results despite the fact we actually 
have $\fe = {1\over3}$.
\vfill

\centerline{\psfig{figure=sup2.ps,width=5.5in}}
{\tenrm Fig. 19 (Supercritical shock): snap shots of the density, gas ($+$) and radiative
                  ($\Diamond$) temperature, radia-}
		  
\hspace{12mm}{\tenrm tive flux, and gas velocity as functions
                  of the optical depth $\tau$ relative to the optical depth of}

\hspace{12mm}{\tenrm the shock $\tau_s$ at quater piston crossing time; overlayed over the flux plot is
                  three times the}
		  
\hspace{12mm}{\tenrm Eddington factor ($\cdots$).}
\vspace{5mm}

The spike in the gas temperature is prominent and very narrow. Its width is given by 
$(\rho\chi_e)^{-1}$ $=$ $3,210~km$. This is much larger than the rezoning limitation of 
$\ell \approx 80~km$ imposed by the artificial shock width. In figure~25 on the right gas 
($+$) and radiative ($\Diamond$) temperature are plotted in the vicinity of the shock 
front. Their dependency in the relative optical depth is almost independent on time. The 
spike exists clearly only for the gas temperature while the radiative temperature varies 
smoothly. The spike is optically 

{\hbox{%
{\vbox{\hsize=3.0in\noindent\strut
thin and hence practically transparent to the radiation (although the actual location is mostly 
at large optical depths). Note that there are many (around 25) grid points in this region, 
each zone with an optical depth of at most $5\times 10^{-3}$. The spike is therefore 
numerically well represented and represents the physical, in particular radiative, 
properties correctly.
\strut}}
{\vbox{\hsize=10mm}}
{\vbox{%
\psfig{figure=sup3.ps,width=3.0in}
{\tenrm Fig. 20: Temperature spike of supercritical shock}}}
}}

\subsubsection{Radiative blast}
We now consider a shell stretching from an inner radius $R_i = 10^7~km$ to an outer radius 
of $R_o = 10^8~km$ containing a mass of $0.16~M_\odot$. Instead of a constant density 
profile we assume a power law atmosphere: 
$$      \rho(r) = 10^{-6}~\Bigl({R_i\over r}\Bigr)^7~g/cm^3              \eqno{(40)}
$$
The gas is purely absorptive and the mean opacity is taken to be constant and equal to the 
Thomson free electron scattering opacity (hence we imply the gas to be completely ionized hydrogen). 
The adiabatic index is $\gamma={5\over3}$. The sphere is illuminated by a central light source with a 
luminosity ${\cal L} = 10^{34}~ergs/sec$. We obtain a self-consistent temperature profile for the 
initial model by relaxating the sphere's luminosity to this value using the full transport 
scheme, \cf~the discission of radiative heating.\vspace{2mm}
 
At time $t=0$ we deposit a large amount of energy inside the sphere (corresponding to a 
central temperature of $10^9~K$). This sets off a radiative blast wave (differing from the 
Sedov-Taylor blast wave by the fact that now the effects of radiation are included). Its evolution 
is followed employing the adaptive grid. As grid parameters we choose the mass, density, and 
radiative energy in logarithmic and the velocity $u$ ($1000~km/sec$) in linear scaling, together 
with $\alpha=1.5$, $\tau=10^{-20}$, and $\theta=1$. The shock width is $\ell_1=10^{-4}$. After 
initiating the blast wave it reaches the surface in 700 time steps. \vspace{2mm}

\centerline{\psfig{figure=rblast1.ps,width=5.5in}}
{\tenrm Fig. 21 (Radiative blast): snap shots of the density, gas temperature, radiative 
                  flux, and gas velocity}
		  
\hspace{12mm}{\tenrm as functions of the exterior mass log M$_{ext}$ 
                  at six different times t $\in$ (437, 2162, 7134, 10329,}
		  
\hspace{12mm}{\tenrm 15562, 29099) sec; 
		  the corresponding optical depths of the shock front are $\tau\in$
                  (3.2$\times$10$^4$,}
		  
\hspace{12mm}{\tenrm 6.25$\times$10$^3$, 184, 61.5, 15.9, 0.17).}
\vspace{5mm}

In figure~21 the evolution of the gas sphere is depicted at six different times. 
Initially the shock moves down the density and temperature gradient. There is a very prominent 
flux peak and the shock profiles for density, temperature and velocity are sharp. The shock 
looks like an adiabatic one. In time the radiation preheats an increasing region ahead of

{\hbox{%
{\vbox{\hsize=3.0in\noindent\strut
it. When the shock has penetrated to an 
optical depth of about 200 (which occurs at time $t\approx 7000~sec$) this radiative precursor 
has reached the surface. The surface luminosity rises very fast and reaches its maximum at 
about $t=1.13\times 10^4~sec$, see figure~22. The gas and radiative temperatures are out of 
equilibrium at the head of the precursor. The heating induces a quasi-\strut}}
%{\vbox{\hsize=10mm}}
{\vbox{%
\psfig{figure=rblast2.ps,width=3.0in}
{\tenrm Fig. 22: surface luminosity as a function of time;}

\hspace{12mm}{\tenrm the $+$ signs mark individual time steps.}}}
}}
\noindent
isothermal wind in the outer shell and the gas equilibrates again. Before reaching an 
optical depth of 60 the shock has dissipated through the radiative processes. The shocked gas, 
however, maintains its momentum and begins a ``snow plow'' phase where it sweeps up wind material. 
The surface luminosity declines and an almost constant flux profile becomes established. The 
sphere begins to cool off. Behind the ``snow plow'' a very dense shell is formed which is well 
established by $t\approx 1.5\times 10^4~sec$ at optical depth $\tau = 16$. While the plow steepens 
up at its head to a strong supercritical shock (with the typical temperature spike) another,
reverse shock forms at its tail. It propagates slowly into the homologously expanding inner
regions. The dense shell expands rapidly and moves at nearly the shock speed outwards.
Through the radiation the gas keeps cooling down, and, since the almost
constant flux level decreases, the surface luminosity keeps dropping.\vspace{2mm}

\centerline{\psfig{figure=rblast3.ps,width=5.5in}}
{\tenrm Fig. 23 (Radiative blast): optical depth-time diagram of the logarithmic density in
                  g/cm$^3$.}
\vspace{5mm}

The evolution of the radiative blast simulation can be illustrated via ``optical 
depth-time diagrams'' of the density, velocity, and temperature (figures~23 -- 25). One can 
recognize the adiabatic shock during the initial phase up to several 1000 $sec$. While it 
dissipates through radiative losses the radiative precursor surfaces at about $10^4~sec$ and 
the luminosity there breaks out. The shocked regions plow into the isothermal wind while 
the gas keeps cooling off. A new supercritical shock is born and runs fast towards the surface. 
Simultaneously a reverse adiabatic shock forms. It moves into the inner regions which already 
undergo a homologous expansion and a high density shell at a constant speed is created.

\centerline{\psfig{figure=rblast4.ps,width=5.0in}}
{\tenrm Fig. 24 (Radiative blast): optical depth-time diagram of the velocity in 1000 km/sec.}
\vspace{5mm}

\centerline{\psfig{figure=rblast5.ps,width=5.0in}}
{\tenrm Fig. 25 (Radiative blast): optical depth-time diagram of the logarithmic gas temperature in K.}
\vspace{5mm}

The grid motion is depicted in figure~25 in which every fourth zone is plotted 
as a function of optical depth and time. After a very short relaxation period the grid 
points begin to zoom into the forming shock, $t\approx 2000~sec$. They actually lock into 
the radiative precursor until it surfaces at $t\approx 7000~sec$. The grid then rezones 
quickly to resolve a new shock system that forms around $t\approx 10^4~sec$. It tracks a 
reverse shock and follows a forward shock all the way to the surface.\\

\centerline{\psfig{figure=rblast6.ps,width=5.5in}}
{\tenrm Fig. 26 (Radiative blast): optical depth-time diagram of the grid point motion;
                  every fourth zone is}
		  
\hspace{12mm}{\tenrm shown.}
.
