Transport in Astrophysics: II. Diffusion with Advection in Expanding Nebulae

The structure across an expanding shell in which drift and diffusion redistri-bute material requires careful consideration in order to correlate the surface brightness with the physical parameters. A new solution in terms of Fourier series is suggested for 1D diffusion in the presence of a drift velocity. The astrophysical parameters are chosen in agreement with the astronomical data for the planetary nebula A39. The new solution is then inserted into the existing theory for the astrophysical image which allows dealing with the intensity of radiation emitted in a spherical layer.


Introduction
In order to model the images of planetary nebulae (PN) and supernova remnants (SNR) we should parametrize the thickness of the expanding layer as well the involved physics. As an example, the number density of particles in the expanding emitting layer can be assumed constant or variable due to some diffusive process. A previous analysis of diffusion with drift has covered the stationary case [1] leaving the transient state as the subject of further research. A 1D diffusion with drift is governed by Fick's second equation, which is a partial differential Equation (PDE) that, at the moment of writing, is a subject of research.
We report two examples of the ongoing research on the above equation: an analytical and numerical way has been suggested by [2] and four different methods have been implemented by [3]. The intensity of the emitted radiation versus distance from the centre in astrophysical phenomena such PN and SNR presents a characteristic "U" shape which can be explained as the solution of the transport

The Existing Solutions of Diffusion with Drift
We review 1D diffusion with drift in the framework of the mathematical theory of diffusion and asymmetric random walk with drift as developed in [1]. In order to avoid duplicates we adopt the following astrophysical units: pc for length and yr for time. In these units, the advection velocity v is expressed in pc/yr and the diffusion coefficient in pc 2 /yr.

1D Diffusion with Drift, Stationary State
In one dimension and in the presence of a drift velocity, v, along the x-direction, a diffusion is governed by Fick's second equation for the concentration, ( ) C x , see Equation (4.5) in [6], 2 2 , where v can take positive or negative values and D is the diffusion coefficient. The number density rises from 0 at x = a to a maximum value C m at x = b and then falls again to 0 at r = c. The general solution to Equation (1) for a steady state is We now assume that v is negative: the solution is The boundary conditions give

1D Diffusion with Drift, Random Walk
Given a 1D segment of length side we can implement a random walk with International Journal of Astronomy and Astrophysics step-length λ by introducing the numerical parameter side NDIM λ = . We now report the adopted rules when the injection is in the middle of the grid: 1) The first of NPART particles is chosen.
2) The random walk of the particle starts in the middle of the grid. The probability of taking one step in the negative direction (downstream), is 3) When the particle reaches one of the two absorbing points, the motion starts again from (2) with a different diffusion pattern.
4) The number of visits is recorded on  , a 1D grid.
5) The random walk terminates when all the NPART particles have been processed.
6) For the sake of normalization, the one-dimensional visitation or number density grid  is divided by NPART.
There is a systematic change of the average particle position along the x-direction: for each time step. If the time step is where tr v is the transport velocity, the asymmetry, µ , that characterizes the random walk is

1D Diffusion with Drift, Transient State
In 1D and in the presence of a drift velocity, v, along the radial direction the diffusion is governed by Fick's second equation, see Equation (4.5) in [6], On assuming that the boundary conditions are ( ) where ( ) xv (12) where N 0 is the number of particles at 0, 0 x t = = . We now analyse three types of initial conditions. The first case to be analysed is a Gaussian density profile of the type ( ) where b is an adjustable parameter. The integral (12) where 1 I = − . In the case of a Gaussian profile of the density Figure 1 presents the solution for different values of the velocity. The second case to be analysed is a parabolic density profile of the type The integral (12) In the case of a parabolic profile of the density Figure 2 presents the solution for different values of the drift velocity.
The third case to be analysed is a trigonometric density profile of the type The integral (12) The trigonometric solution as a function of x and t is presented as a surface plot in Figure 3.

Astrophysical Applications
We present the astronomical data on PN A39.

A39
The PN A39 is extremely round and therefore can be considered an example of spherical symmetry, see for example Figure 1 in [7] or Figure 4.
The radius of the shell in A39, shell R , is where 77 Θ is the angular radius in units of 77" and 21 D the distance in units of 2.1 kpc, see [7]. International Journal of Astronomy and Astrophysics  The expansion velocity has a range km 32 37 s   ↔     according to [8] and the age of the free expansion is 23,000 yr, see [7]. In our astrophysical units we will where 10 Θ is the thickness in units of 10.1" and the height above the galactic plane is 1.42 kpc, see [7]. For the sake of comparison with the observations, the length of some geometrical models will be expressed in arcsec. An occasional reader may ask whether A39 can be approximated by a sphere or not? One way to answer this question is to introduce the observational percentage of sphericity, obs ε , over the whole range of the polar angle θ , obs num obs obs, 100 1 , where num r is the theoretical averaged radius of the nebula, obs r is the observed radius of the nebula, and the index j varies from 1 to the number of available observations. In the case of Figure 4, we have obs 97% ε = over 67 directions when the central star is assumed to be the source of the expansion.  (4) and (5).

Astrophysical Diffusion
In the case of a transient diffusion, we present, in Figure 6, the parabolic solution that will be used in the theory of formation of the image.   (4) and (5) is presented when

Transfer Equation
The transfer equation in the presence of emission only, see for example [9] or We now continue analysing the case of an optically thin layer in which ν τ is very small (or k ν very small) and the density ζ is substituted with ( ) C s , our number density of particles. We now assume that the emissivity is proportional to the number density: where ν is the frequency and ( ) f α γ is a slowly varying function of f γ which is of the order of unity and is given by 2   7 3 3  1  3  7  2  ,  1 12 12 for 1 2 f γ ≥ , see Formula (1.175) in [11]. Synchrotron emission is widely used to explain the radiation observed in SNR, see [12]- [17]. This non-thermal radiation continuum emission was also detected in a PN associated with a very long-period OH/IR variable star (V1018 Sco) [18], in the young PN IRAS 15103-5754 [19] and in IRAS 18062 + 2410 [20]. A discussion of thermal/nonthermal emission from PN's can be found in [21] [22].
where Δs is the spatial interval between the various values and the sum is performed over the interval of existence of the index k. The theoretical intensity is then obtained by integrating the intensity at a given frequency over the solid angle of the source.

Image Theory at Constant Density
The simplest model for the image is characterized by a radiative process with constant number density of emitting particles in a thin layer around the advancing sphere. We therefore assume that the number density C is constant in the spherical thin layer, and in particular rises from 0 at r = a to a maximum value C m , remains constant up to r = b and then falls again to 0. The line of sight, when the observer is situated at the infinity of the x-axis, is the locus parallel to the x-axis which crosses the position y in the Cartesian x-y plane and terminates at the external circle of radius b, see Figure 7. The length of this locus is ( ) The number density C m is constant between two spheres of radii a and b and therefore the intensity of radiation is ( )   .
The ratio between the theoretical intensity at the maximum ( y a = ) and at the minimum ( 0 y = ) is given by  Figure 9 shows a spherical shell source of radius b between a spherical absorber of radius a and a spherical absorber of radius c. The number density rises from 0 at r = a to a maximum value C m at r = b and then falls again to 0 at r = c. In order to continue, the 3D spherical diffusion is approximated by the 1D steady state diffusion with drift and the numbers density to be used are Formulas (4) and (5)

Image Theory in a Geometrical Model
We assume that the number density rises from 0 at r = a to a maximum value C m at r = b and then falls again to 0 at r = c following a parabolic behavior with vertex at r = b. The parabola that satisfies the above requirements is and we continue inserting In the case of the geometrical model the theoretical intensity of A39 is presented as a cut, see Figure 11 and as a surface, see Figure 12. Figure 9 can model the formation of the image in the case of a transient state developed in Section 3 once a is replaced by zero, b is replaced by L/2 and c is replaced by L. The solution considered is Equation (9) for a parabolic density profile, see Equation (15). In this case the operation of integration which leads to the image is done numerically and is presented as a cut, see Figure 13 and as an image, see Figure 14.

Image Theory in a Transient State
The effect of the insertion of a threshold intensity, tr I , given by the observational techniques, will now be analysed. The threshold intensity can be parametrized by max I , the maximum value of intensity characterizing the ring: a typical image with a hole is visible in Figure 15 when max 2 International Journal of Astronomy and Astrophysics

PDE & Boundary Conditions
A new solution for the diffusion with drift has been derived as a Fourier series, see Equation (9). The above solution has been particularized for three types of density profiles at t = 0: Gaussian, parabolic and trigonometric.

Astrophysical Numbers
Once the time has been fixed in years and the velocity as 3.521 × 10 −5 pc/yr we consequently choose the diffusion coefficient around 10 −7 pc 2 /yr.

Image Theory
The image theory of a spherical emitter such as PN A39 has been modeled by an advancing layer with a constant density, a stationary advection with drift, an existing geometrical profile for the emitting layer and a time dependent diffusive model in presence of drift.
Plans for the Future Formula (21) fixes the degree of sphericity of A39 at 97% but some asymmetries in the spatial dimensions and in the emissivity are present. An explanation of these anomalies requires introducing the concept of velocity as a function of the position angle, which was begun to be used in [1].