NLTE Monte Carlo Radiation Transfer in Circumstellar Disks

Presented at the meeting of the Working Group on Active B Stars during the 25th IAU General Assembly in Sydney, Australia on 2003 July 16

J.E. Bjorkman and A.C. Carciofi

Ritter Observatory, MS 113, Department of Physics & Astronomy, University of Toledo, Toledo, OH 43606-3990, USA
email: jon@physics.utoledo.edu, acarcio@physics.utoledo.edu

Received: 2005 January 9; Accepted: 2005 January 24.

A popular model for the circumstellar disks of Be stars is a rotationally supported (i.e., Keplerian) viscous decretion disk (see reviews by Bjorkman 2000; Porter & Rivinius 2003). This model is essentially the same as that employed for protostellar disks, the primary difference being that Be disks are outflowing, while pre-main-sequence disks are inflowing. The essential physics that determines the geometrical structure of Keplerian disks is reasonably well understood (at least in the case of pre-main-sequence stars). The primary result is that the disks are hydrostatically supported in the vertical direction, while the radial structure is governed by the viscosity. Since the disk is pressure-supported (in the vertical direction), the geometrical structure of the disk is determined by the temperature of the disk. Consequently, to critically test Keplerian disk models of Be stars against observations, we must determine the temperature structure of the disk. To do so, we have developed a 3-D NLTE Monte Carlo radiation transfer code that self-consistently solves the radiative equilibrium temperature, vertical hydrostatic equilibrium, and steady state density of a gaseous hydrogen decretion disk.

In brief, the Monte Carlo simulation performs a full spectral synthesis by emitting stellar photons with random frequencies (sampled using a Kurucz model atmosphere for the B star). Each photon is tracked as it travels through the envelope (where it may be scattered, or absorbed and reemitted, many times) until it escapes. As the photons escape, they are binned according to their emergent direction and frequency, which gives the emergent spectrum.

During the simulation, whenever a photon scatters, it changes direction, Doppler shifts, and becomes partially polarized. Similarly, whenever a photon is absorbed, it is not destroyed; it is reemitted immediately (on the spot) with a new frequency and direction determined by the local emissivity, , of the gas. Note that we include both continuum processes and spectral lines in the opacity and emissivity of the gas. Since photons are never destroyed (absorption is always followed by reemission of an equal energy packet), our procedure automatically enforces radiative equilibrium and conserves flux exactly (see Bjorkman & Wood 2001). Since the interactions of the photons with the gas provide a direct sampling of all the radiative rates (as well as the heating and cooling of the gas), we can solve the rate equations at the end of the simulation to update the level populations and electron temperature of the gas. Statistical equilibrium then is solved by iteration. At the same time, we also solve the hydrostatic equilibrium equations for the vertical structure of the disk, as well as the radial fluid equations for steady state outflow. We then repeat the simulation until the temperature and level populations converge (this typically requires about 10 iterations).

Figure 1 shows the radiative equilibrium temperature and hydrogen level populations in the equatorial plane of the disk for our best fit model to the Be star  Tau. Initially, the temperature drops like a flat blackbody reprocessing disk (red line; see Adams, Lada, & Shu 1988) until it reaches a minimum when the disk becomes optically thin (vertically). Beyond this location, the temperature rises back up and becomes roughly isothermal (blue line) at large radii.

Figure 1. Mid-Plane Temperature and Level Populations. Shown are the NLTE radiative equilibrium temperature (left) and level populations (right) for the first five levels of hydrogen in the mid-plane of a Keplerian disk.

The temperature controls the disk scale height, so it determines the geometrical thickness of the disk. A flat reprocessing disk has a temperature T  r^-3/4, giving a scale height H  r^9/8, while an isothermal disk has a scale height H  r^3/2. Figure 2 shows the disk temperature and corresponding density. Since the temperature in the inner disk falls rapidly, there is little flaring (increase of opening angle) in the inner disk, but as the temperature rises back to the isothermal value, the disk begins to flare quite dramatically at large radii.

Figure 2. Disk Temperature and Density. Shown for three different radial scales are the NLTE radiative equilibrium temperature (top) and hydrostatic equilibrium density (bottom) for a Keplerian hydrogen disk.

As confirmed by optical interferometry (Quirrenbach et al. 1997), the disk is responsible for producing the Balmer emission lines, IR excess, and intrinsic polarization of Be stars. Since the radial dependence of the disk density, temperature, and opening angle all affect the slope of the IR excess (Wright & Barlow 1975; Cassinelli & Hartmann 1977; Waters 1986), as well as the detailed shape of the intrinsic polarization (Cassinelli, Nordsieck, & Murison 1987; Wood, Bjorkman, & Bjorkman 1997), reproducing the IR excess and detailed spectropolarimetry is a non-trivial test of the Keplerian disk model.

Figure 3 shows the comparison of the predicted SED and intrinsic polarization to observations for  Tau (Bjorkman et al. 1991; Wood et al. 1997). Not only does the predicted flux match the observed Balmer jump, which is partially filled-in by bound-free disk emission, it also matches the slight excess shortward of the Paschen jump, as well as the IR excess observed by IRAS. Similarly, the predicted polarization also agrees with the observed polarization Balmer and Paschen jumps, as well as the slope in the Paschen continuum. Finally, we note that our model fit to  Tau has only two free parameters: the inclination angle (i=70°) and the density scale of the disk (n₀=3×10¹³cm^-3), which corresponds to a disk mass loss rate ∼10^-10 Myr^-1.

Figure 3. Spectral Energy Distribution and Polarization. The top panels show the SED predicted by our NLTE Keplerian disk model viewed at an inclination of 70°, while the bottom panels show the corresponding polarization. Note the excellent agreement with the observations of  Tau.

Our results indicate that a Keplerian decretion disk model does reproduce the detailed continuum observations (both flux and polarization). From this we conclude that the temperature, density and geometry of the disk are consistent with a Keplerian decretion disk in (vertical) hydostatic equilibrium. Our next goal is to make detailed comparisons of the hydrogen emission line profiles. In particular, we hope to use the line profile shapes to test whether the disk rotation speed is in fact Keplerian as required by these models.

Acknowledgements. This work was supported by NSF grants AST-9819928, AST-0307686, and NASA grant NAG5-8794 to the University of Toledo.

References:

Adams, F.C., Lada, C., & Shu, F.H. 1988, ApJ, 326, 865

Bjorkman, J.E. 2000, "The Formation and Structure of Circumstellar Disks", in The Be Phenomenon in Early-Type Stars, IAU Colloquium 175, ASP Conference Series, Vol. 214, eds. M.A. Smith, H.F. Henrichs, & J. Fabregat (San Francisco: ASP), 435

Bjorkman, J.E., & Wood K. 2001, ApJ, 554, 615

Bjorkman, K.S. et al. 1991, ApJ, 383, L67

Cassinelli, J.P., & Hartmann, L. 1977, ApJ, 212, 488

Cassinelli, J.P., Norsdieck, K.H., & Murison, M.A. 1987, ApJ, 317, 290

Porter, J.M., & Rivinius, T. 2003, PASP, 115, 1153

Quirrenbach, A., et al. 1997, ApJ, 479, 477

Waters, L.B.F.M. 1986, A&A, 162, 121

Wood, K., Bjorkman, K.S., & Bjorkman, J.E. 1997, ApJ, 477, 926

Wright, A.E., & Barlow, M.J. 1975, MNRAS, 170, 41

back to Meeting Summary