- Email: [email protected]
Monte Carlo Simulation of Lyα Line Profiles
CHINESE ASTRONOMY AND ASTROPHYSICS
ELSEVIER
Chinese Astrophysics 33 33(2009) (2009)347–364 347–364 ChineseAstronomy Astronomy and and Astrophysics
Monte Carlo Simulation of Lyα Line Profiles† WANG Yu University of Science and Technology of China, Hefei 230026
Abstract In this paper, a code for simulating the Lyα line profiles by the Monte Carlo method is introduced. This code is developed on the basis of an idea given by Anne Verhamme et al. and by combining with another simulation method. For different models, the Lyα line profiles under different conditions are simulated. The simulation method of this code is given in detail, and the simulated results are discussed as well. Key words: radiative transfer—line: profile
1. INTRODUCTION The Lyα line plays a very important role in different fields of astrophysics, especially as a tool for observing and studying the objects with high redshifts. It can be taken as a probe for studying the rate of star formation, the distribution and motion of matter, the distribution of dust, etc. Meanwhile, it is also used for studying the ionization state of intergalactic medium[1] . The Lyα photons discussed in this paper are produced in HII regions. The HII region is produced by a forming massive star which has the surrounding gas ionized. In some HII regions, the massive O, B stars emit the violet radiation, in which, a part of violet photons are absorbed by the surrounding gas and therefore lead to gas ionization, another part are absorbed by dust, and the photons of λ < 912 ˚ A make hydrogen ionized. Each photon absorbed by the hydrogen atom around the star will give rise to one Lyα photon with the wavelength of 1216 ˚ A and one Balmer photon. The absorbtion of Lyα photons is the basic origin of the dust heating in molecular clouds around HII regions[2] . Because of its value in scientific research, the Lyα line has been simulated by many people with different methods. By the Monte Carlo method, Verhamme et al. simulated the †
Supported by National Foundation for Talent Cultivation of Basic Science Received 2008–05–26; revised version 2008–07–07 A translation of Acta Astron. Sin. Vol. 50, No. 2, pp. 117–133, 2009
0275-1062/09/$-see front matter © 2009 Elsevier B.V. All rights reserved. doi:10.1016/j.chinastron.2009.09.001
348
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
radiative transfer of the Lyα photon in 2006[1]. On this basis, combining with the method suggested by Hansen et al.[4] , the author of this paper tries to simulate the Lyα line profiles with the Monte Carlo method.
2. METHOD OF SIMULATION 2.1 A Brief Introduction to Monte Carlo Method The Monte Carlo method (MC method) has been widely used in different scientific fields. Its advantage is: we can simulate various physical quantities by using the given probability distribution functions. We assume that y is evaluated in the range [a, b], and that its probability distribution function is f (y). Then the probability that y takes a value y b in the range [y, y + dy] is f (y)dy. Because that a f (x)dx/ a f (x)dx = ξ, ξ ∈ [0, 1], each y-value corresponds to a number in the range [0,1], so, from a randomly generated number in the range [0,1], the y-value corresponding to it can be derived reversely. This is the quintessence of the MC method. In practical applications, for the probability distribution function which is difficult to be integrated, we can use other methods, for example the rejection sampling method, to obtain a result as same as that obtained by integrating. 2.2 Physical Processes to be Considered The simulations are made on photons one by one, namely, after the simulation of one photon is terminated, the simulation on the next one is made. And every photon to be simulated will undergo the following processes: being emitted, interacting with hydrogen and interacting with dust. 2.3 Radiation Sources The photon emitted from a source is characterized by the radiation direction (θ, φ) and the frequency x: x = (ν − ν0 )/ΔνD , (1) 1/2
ΔνD = (12.85T4 /c)ν0 .
(2)
In this paper, the frequencies of Lyα photons are all scaled by x. The radiation source in this code is an isotropic point source. For the monochromatic source, all the emergent photons have the same frequency; for the source of flat continuum, the frequencies of the emergent photons are produced randomly; and for the sources of Gauss and Lorentz profiles, the frequencies of the emergent photons can be produced by the MC method. Since the source is isotropic, it emits particles towards different directions with the same probability. When a photon emerges from the source, its radiation direction is determined by the following formula: θ = cos−1 (2ξ1 − 1),
φ = 2πξ2 .
(3)
This formula is derived by the MC method, in which ξ1 and ξ2 are two random numbers distributed uniformly in the range [0,1]. 2.4 Photons’ Walk in Simulation Space Photons walk in space randomly[1,3] . As soon as a photon walks from one place to interact with the medium to the next place to interact with the medium, it has completed
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
349
its walk one time. The one-time walk is defined by the step length s and the walking direction (θ, φ). The walking direction is the direction that the photon emerges from, and the step length s should be obtained by calculations. As the probability P that a photon interacts with the medium, and the optical depth τint that the photon will penetrate through, satisfy the following relation[1,3]: P (τint ) = 1 − e−τint , (4) so, with the MC method, we can obtain: τint = − ln(ξ),
(5)
ξ is a random number distributed uniformly in the range [0,1]. In addition, τint is also a function of the length s: τint (s) = τx (s) + τd (s), (6) in which, the optical depth of hydrogen atoms is[1] : −1/2
τx (s) = σH (x)nH s = 1.041 × 10−13 T4 the optical depth of dust is:
τd = τa + τs =
NH
H(x, a) √ , π
(7)
s
σd nd (s)ds,
(8)
0
and H (x,a)[1] can be expressed as: 2 2 e−y a ∞ e−x , if |x| < xc H(x, a) = ≈ √a , , if |x| > xc π −∞ (y − x)2 + a2 πx2
(9)
in which τd = τa + τs , τa,s = πd2 Qa,s , Qa ≈ Qs ≈ 1, d = 10−6 cm. Therefore, to derive at first τint , then to derive reversely s, the step length of photon’s one-time walk is obtained. 2.5 Interaction between the Photon and the Medium When a photon completes its one-step walk, and arrives at a place to interact with the medium, it will undergo the process of radiative transfer. 2.5.1 The interaction with hydrogen When the photon interacts with a hydrogen atom, we assume that it is absorbed at first, then it is reradiated outward. This code assumes the reradiation being isotropic. In this case, the radiation direction is determined by the same method as the above-mentioned, namely, the direction is determined by two random numbers. As for the frequency of the reradiated photon, the Hummer’s redistribution function RII is employed[1,4,5,6]. RII (ν, n, ν , n )dνdΩdν dΩ expresses the probability that the incident photon with the frequency ν in the direction n reradiates at the frequency ν in the direction n . Specifically, RII is expressed as[5] : → → RII (x, − n ; x , − n ) =
→ → x − x 2 x + x γ γ g(− n,− n ) γ exp[−( ) csc2 ( )]H( sec( ), σ sec( )) , 2 4π sin γ 2 2 2 2 2
(10)
350
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
1 3 in which, gA ( n, n ) = 4π in the case of isotropic radiation, gB ( n, n ) = 16π (1 + cos2 γ) in the case of dipole radiation, and γ is the included angle between n and n . Since RII (x0 , x)dx indicates the probability that the incident photon with the frequency x reradiates at the frequency x, so, it can produce a random number. When it equals to 0 x R dx, by deriving reversely the value of x, the frequency of the reradiated photon can II −∞ be obtained. Considering the computing speed, this code has adopted the formula given by Hansen[4] . It is an improvement on the formula given by Neufeld[4,6] : √ x2 x 2i 3 22.5 , (11) R(xi , x) = π 22.5 x4i + (x3 − x 3i )2
˜ = xi − 2/xi . But this formula does in which xi is the frequency of the incident photon, x not suit the situations when the frequency approaches or equals to zero. It is apparent that if the incident frequency is zero, this formula will be invalid. For this reason, in this code the calculations are made by using Eq.(8) while | x |≤ ξ, and by using Eq.(9) while | x |> ξ, ξ is a constant. 2.5.2 The interaction with dust When the photon arrives at the place to interact with the medium, there will be certain probability to interact with dust. The probability of photon-hydrogen interaction is[1] : PH (x) =
nH σH (x) . nH σH (x) + nd σd
(12)
Hence, we can produce a random number in the range [0,1]. If it is less than PH , then the photon interacts with hydrogen; if it is greater than PH , then the photon interacts with dust. When photons interact with dust, a part of them are absorbed by dust, and the other part are scattered by dust. In this code, it is assumed that among the photons interacting with dust, one half of photons are absorbed and another half are scattered. So in order to judge whether a photon is absorbed or scattered, a random number distributed uniformly in the range [0,1] is produced. If it is less than 0.5, then the photon is scattered, otherwise the photon is absorbed. And this code assumes the scattering being isotropic. 2.6 The Treatment for the Medium in Motion For the medium in motion, the calculations are made by transferring photons to the same coordinate system as that of the medium. This mainly leads to a change of photon’s frequency. Because that no matter what is calculated, the step length or the frequency after being scattered by the medium, the calculations are all related with the photon’s frequency. If the frequency is failed to transformed into a motionless coordinate system relative to the medium, the calculations will be erroneous. 2.7 The Termination of Photon’s Walk When a photon walks out of the simulation space or it is absorbed by the medium, the simulation on this photon terminates, and the program starts the simulation on the next photon.
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
351
2.8 The Generation of Random Numbers The quality of the random numbers prepared for the MC method relates closely to the quality of the simulated results. The random number generator used by this code is the 16807 generator. The 16807 generator adopts a kind of linear congruence method. It is actually a function, which makes an integer transform to be an another integer, and therefore forms a number series, the numbers in which have a good independence. The specific formula of the 16807 generator is a(z mod q) − r[z/q], if ≥ 0 S(z) = , (13) a(z mod q) − r[z/q] + m, otherwise in which, a=16807, b=0, m=2147483647, q=127773, r=2836, and the symbol “[ ]” means taking the integer part. 2.9 The Statistics on Results In this code the statistics on results is realized through a chain table, in which every chain element possesses the corresponding frequency x and photon number n. When the simulation on a photon terminates, the program will have the photon number n in the corresponding chain element added with 1. When the simulations on all photons are finished, the program will have the number n in each chain element divided by the total number of photons, and obtain the photon’s probability. And the program will store the photon frequency x and the corresponding probability into a file.
3. RESULTS OF SIMULATIONS For making simulations of Lyα line profiles, this code adopts the parameters taken by Verhamme et al. 3.1 Slab Model 3.1.1 Static, dust-free and monochromatic slab model As shown in Fig.1, the static monochromatic slab model is composed of three parallel infinite planes, the intermediate plane is the radiation source, the photons generated here will enter into the medium enclosed by the other two surfaces. The two terminal surfaces are the boundaries, on these two surfaces the optical depth relative to the radiation source is τ0 . Using this model, this code has obtained a symmetrical double-peaked structure. The shape of peaks is influenced by the optical depth on the boundary relative to the source. The less the optical depth on the boundary, the sharper the peaks, and the greater the optical depth, the more “short and fat” the peaks. This is easy to understand, the greater the optical depth, the more the processes of radiative transfer a photon should undergo. Using this model, Neufeld has derived an analytic solution of line profiles[1,6] : √ 6 x2 1 , (14) J(τ0 , x) = √ 3 24 π aτ0 cosh[ π /54(|x|3 /aτ0 )] and it is compared with the simulated result obtained by this code in Fig.2. From Fig.2 we can find that the two kinds of results are not very different.
352
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
Fig. 1 Diagram of the slab model
Fig. 2 The simulated result of the static dustfree monochromatic slab model, in comparison with the analytic solution derived by Neufeld
3.1.2 The dust-free monochromatic slab model with a velocity gradient Now we consider the situation when a velocity gradient exists in the medium inside the slab. The velocity gradient can be taken as: V (r) = |r/R|Vmax , R is the source position relative to the boundary, r is the medium position relative to the source, and Vmax is the medium velocity on the boundary. Fig.3 shows the Lyα line profiles simulated under different boundary velocities. We can find that as the velocity increases gradually, the originally symmetrical double peaks become asymmetrical, and finally become a single peak. The sign of the velocity on the boundary determines whether the spectral line is redshifted or blue-shifted. Fig.4 gives the simulated results for the velocities on the boundary in opposite directions. It demonstrates that when the velocity is positive, taking the source as the center, the medium expands against the source, the spectral line is red-shifted; and that when the boundary velocity is negative, the medium in-falls toward the source, the spectral line is blue-shifted.
Fig. 3 The effect of the boundary velocity on the Lyα line profile
Fig. 4 The effect of the direction of boundary velocity on the Lyα line profile
3.1.3 Effect of dust This code takes τa , the absorbtion optical depth of dust, to characterize the amount of
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
353
dust. Fig.5 shows the simulated results under different absorbtion optical depths of dust. We can find that the greater the absorption optical depth, the more the absorbed photons, and the lower the double peaks. For the situation of existence of dust, Neufeld has derived the fraction fe of photons[7] , namely the number ratio of the Lyα photons which finally emerge from the medium with respect to the totally generated Lyα photons[1,5,7] : √ 3 fe = 1/ cosh[ 5/12 |(aτ0 )1/3 τa |1/2 ] , (15) π ξ in which ξ =0.635. Fig.6 is a comparison between the faction fe derived by this code and the result given by Neufeld. It is found that the fe − (aτ0 )1/3 τa relationship derived by this code is similar with the Neufeld’s analytic solution.
Fig. 5 The line profiles under the different absorption optical depths of dust
Fig. 6 The fraction fe as a function of (aτ0 )1/3 τa derived by this code, in comparison with Neufeld’s analytic solution
3.1.4 The case of a nonmonochromatic radiation source This code can also simulate the situation that the photons radiated from the source are not monochromatic. Fig.7 and Fig.8 are the simulated results when the frequency x of the photon radiated from the source satisfies respectively: 2 1 e−x , p(x) = √ πΔνD
p(x) =
1 Γ . 4π 2 (xΔνD )2 + (Γ/4π)2
(16) (17)
The line profiles exhibit the symmetrical double-peaked structure. In the process to produce the photons, the rejection sampling method is used for the photon’s frequency. Fig.9 gives the simulated result for the flat continuum, namely the case when the radiation source emits photons with various frequencies. For the photon of large x, it will directly penetrate through the medium, but when x approaches to zero, the probability that the photon interacts with the medium increases. As a result, the curve becomes flat lines far apart from the place of x =0 and the double peaks appear around the place of x =0.
354
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
3.2 Halo Model 3.2.1 The halo model with one point source at the center of spherical boundary The structure of the model is shown in Fig.10, photons are produced at the center of the sphere, within the sphere is the homogeneous medium, expanding or in-falling with the velocity V (r) = |r/R|Vmax .
Fig. 7 The case when the source emits photons with the Gaussian profile of Eq.(16)
Fig. 8 The case when the source emits photons with the Lorentz profile of Eq.(17)
Fig.11 shows the line profiles when the halo is under the static, expanding and in-falling three states. As is shown, when the medium expands relative to the center, the spectral line is redshifted; when the medium in-falls relative to the center, the spectral line is blueshifted; and when the medium is static, the spectral line exhibits a symmetrical double-peaked structure.
Fig. 9 The case of flat continuum
3.2.2 The model with radiation sources distributed uniformly in the medium Fig.12 is the diagram of the model with the radiation sources distributed uniformly in the medium. From this figure we can find that the radiation sources are uniformly distributed in a sphere within the simulation space, not in the whole simulation space. It is necessary, if the sources are uniformly distributed in the whole simulation space, then the
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
355
photons produced on the boundary will have a large probability to escape from the medium
Fig. 10 The structure of the expanding halo with one central point source
Fig. 11 The profiles of the static, expanding and in-falling halos with one central point source
through the boundary. As a result, the produced spectral line will have a sharp peak at the place of x =0. For this reason, this code has the photon sources distributed uniformly within the sphere with the radius of 0.85Rmax . The produced spectral lines are shown in Fig.13. When the medium is static, the spectral line exhibits a symmetrical double-peaked structure; when the medium expands relative to the center, the spectral line is redshifted; and when the medium in-falls relative to the center, the spectral line is blueshifted.
Fig. 12 The structure of the expanding halo with uniformly-distributed radiation sources
Fig. 13 The profiles of the static, expanding and in-falling halos with uniformly-distributed radiation sources
3.2.3 The effect of NHI on the spectral line When NHI becomes large, the probability of the photon-medium interaction increases. In the medium the radiative transfer of atoms onto photons will make the photon’s frequency x increasingly deviate from its original value. So, the greater the value of NHI , the greater the leftshift of the peak. In the case of relatively small NHI , many photons will leave the
356
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
medium not until they have interacted with the medium. Fig.14 and Fig.15 are the line profiles under the different values of NHI for the expanding halo models with one central point source and uniformly distributed sources.
Fig. 14 The line profiles of the expanding halo with one central point source under different values of NHI
Fig. 15 The line profiles of the expanding halo with uniformly-distributed sources under different values of NHI
3.2.4 The effect of medium expansion velocity on the spectral line Fig.16 displays the simulated result on the expanding halo model with one central point source under different expansion velocities on the boundary. From this figure we can find that as the expansion velocity increases gradually, the symmetrical double-peaked structure is gradually destroyed, and there finally remains a single peak.
Fig. 16 Line profiles of the expanding halo with one central point source under different expansion velocities on the boundary
3.2.5 The case of flat continuum When the frequencies of the photons emitted from the source form a flat continuum, the P Cygni-type spectral line will result. Fig.17 and Fig.18 are respectively the line profiles of the expanding medium with a central point source inside and the in-falling medium with the radiation sources uniformly
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
357
distributed, under different absorption optical depths of dust. As is shown, for the expanding medium, the emission line appears at the place of redshift, the absorption line appears at the place of blueshift; and for the in-falling medium, the absorption line appears at the place of redshift, the emission line appears at the place of blueshift. We can find as well that the increase of the dust content will change violently the shape of the spectral line. As the dust content increases, the spectral line gradually exhibits the shape of damped type.
Fig. 17 The profiles of the expanding halo with a central point source of flat continuum
Fig. 18 The profiles of the in-falling halo with uniformly-distributed sources of flat continuum
3.3 Expanding-shell Model The model of expanding shell is illustrated by Fig.19. The radiation source is the point source at the center of the spherical shell, the space between the source and the shell is vacuum, the medium is totally concentrated in the shell expanding with a constant velocity. In the space enclosed by the shell photons walk straightly without any obstruction, only within the shell they will interact with the medium. The spectral line produced by this model is a little peculiar. The spectral line is the superposition of the multiple peaks at nVexp /b. Fig.20 is a good example, in which the dashed lines mark the x-values of the integer multiples of Vexp /b. This result is different from the simulated result obtained by Verhamme et al.[1] , but both results have a peak of the redshift 2Vexp . In certain conditions, the composite multiple peaks of the spectral line can be separated. Fig.21 gives the simulated results under different expansion velocities of the shell. As is shown, when the shell is at rest, the spectral line exhibits the symmetrical double-peaked structure. When the velocity is increased, the spectral line is split into several peaks, especially when the velocity is quite large, the distance between two peaks is enlarged, and the intensity distribution of the peaks is also changed. The Doppler parameter will affect the separation between two peaks. As is shown by Fig.22, when the Doppler parameter is rather large, the separation between two different peaks becomes small, therefore the superposition leads to a single peak at the redshift position; and when the Doppler parameter becomes small, the separation between two different peaks is enlarged. We can find as well that the variation of the Doppler parameter
358
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
Fig. 19 Diagram of the expanding shell model
Fig. 20 An exemplified profile of the expanding-shell model. The dashed lines correspond to x = nVexp /b.
causes only the variation of the separation between two different peaks, but does not change apparently the distribution of the peak intensities.
Fig. 21 The profiles of the expanding shell with different expanding velocities
Fig. 22 The profiles of the expanding shell with different Doppler parameters
The value of NHI in the medium will affect also the shape of the spectral line. As shown in Fig.23, when NHI is relatively small, the peak closer to x = 0 tends to have a higher intensity, but it is not apparent for peaks far apart from x = 0. This is because that most protons have escaped from the medium without interacting with it, so that the spectral line has only one or two high peaks close to x = 0. When NHI is increased, the intensity distribution of different peaks varies, the peak close to x = 0 begins to decline, and the peak far apart from x = 0 begins to rise. As NHI increases continuously, the intensity distribution of peaks becomes more uniform, the highest peak shifts toward the place of larger redshift; and when NHI is quite large, the different peaks are superposed into a single “short and fat” peak. The amount of dust will change as well the shape of the spectral line. As shown in Fig.24, dust will make the brightness of the spectral line decrease.
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
Fig. 23 The profiles of the expanding shell with different NHI ’s
359
Fig. 24 The profiles of the expanding shell with different τa ’s
The value of NHI and the absorption optical depth τa of dust have a major effect on the fe -value of the spectral line. The Doppler parameter b has a less effect on the fe -curve, for the same τa , the larger b gives the smaller fe . Fig.25 compares the fe -τa relationship simulated by this code with the function exp(−τa ). We can find that the derived fe -τa relationship resembles exp(−τa ). The less the value of NHI , the closer to exp(−τa ) the relationship. This point is easy to understand. The less the value of NHI , the larger the proportion of dust in the medium, the larger the probability of the photon-dust interaction, and therefore the closer the fe -τa relationship to exp(−τa ). In reverse, the greater the value of NHI , the longer the path length that a photon walks in the medium, the greater the probability that it is absorbed, and therefore the smaller the fe value.
Fig. 25 NHI ’s
The fe -τa relationships under different
Fig. 26 The fe -τa relationships under different Doppler parameters b
From the above analysis, we can find that the line profile of the expanding shell model is composed by the equally-spaced peaks and the intensity distribution curve of peaks, as shown in Fig.27. The expansion velocity Vexp , the absorption optical depth of dust, and the Doppler parameter b determine jointly the spacing of peaks, Vexp and NHI determine jointly
360
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
the distribution curve of peak intensities, and the effect of the Doppler parameter b on the distribution curve of peak intensities is not evident. Next, we consider the case of flat continuum. From Fig.28, we can find that for the expanding shell, under the condition of flat continuum, the spectral line exhibits the shape of P Cygni-type lines, at the place of redshift occurs the emission line, and at the place of blueshift occurs the absorption line, as the amount of dust increases, the spectral line becomes gradually the damped shape.
Fig. 27 The line profile of the shell model can be decomposed into two parts.
Fig. 28 The line profile under the condition of flat continuum
3.4 Damped-type Spectral Lines As the damped-type spectral lines are observed very often, the study on the formation of the damped-type Lyα line is very meaningful. From the above results, we can find that for the P Cygni-type spectral lines, when the amount of dust increases, the spectral line tends significantly to be of the damped-type. Under the condition that the radiation source is of flat continuum, the damped-type spectral line can be constructed by increasing the amount of dust. This can be tested by the simplest model of static slab. From Fig.29, we can find that when the dust content increases to a certain extent, the spectral line becomes the damped-type. The Voigt absorption line[1] I(x) = I0 exp(−τ0 H(x, a))
(18)
is generally used for simulating the damped-type Lyα line, so we compare them in Fig.30. We can find that a difference exists in the widths of both absorption lines, but this is probably due to the different aspects being considered. The Voigt absorption line model has considered the observing angle. It assumes that some mediums exist between the radiation source and the observer. These mediums have the photons scattered out of the line of sight of the observer. The probability that the photons can pass through the mediums is proportional to exp(−τ0 H(x, a)), for the photons far apart from the place of x =0, it is easy to pass through the mediums and not to be scattered, but for the photons close to the place
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
361
of x =0, it is easy to be scattered out of the line of sight and not to be observed by the observer[1]. However, what is considered by this code is the effect of the dust absorption, only the existence of dust can make the double peaks on the spectral line flattened, and makes the spectral line exhibit the shape of damped-type lines.
Fig. 29 The line profiles of the static dusty slab with a flat-continuum source under different absorption optical depths of dust
Fig. 30 The damped-type line simulated by this code, in comparison with the Voigt absorption line
3.5 The Effect of the Velocity Gradient in the Medium on Line Profiles Now we consider the effect of the variation of the velocity gradient in the medium on spectral lines. There are many kinds of velocity distributions, here we exemplify several kinds of them to make a discussion. At first, we take the slab model and consider the case of V (r) = (r/R)n Vmax . We can find that in the case of V (r) = |r/R|10 Vmax , the symmetrical double-peaked structure is destroyed to a extent less than the case of V (r) = |r/R|Vmax . This is easy to understand, if the power index of |r/R| is very large, then the medium velocities in the inner most regions of the medium will be very small, only the medium velocities near the boundary can be effective. When a photon is far apart from the boundary, it moves as it does in a static medium, and when it is close to the boundary, because of the relatively large medium velocity, and of the relatively thin optical thickness, the photon becomes transparent. Namely the photon will not be affected by the medium near the boundary and emits directly out of the medium, thus the spectral line quite approximates to the symmetrical double-peaked spectral line in the static medium. For this point, Fig.32 provides a better description. 1 Now, we consider the case of V (r) = (r/R) n Vmax . From Fig.33, we can find that when the index n increases, the spectral line exhibits gradually the characteristics of the expanding shell model, i.e., multiple peaks are superposed together at the places of redshifts, and this is easy to understand. As shown by Fig.34, when n increases, the velocities in most regions of the medium are approximately uniform, like a shell far away from the radiation source, therefore the exhibited spectral line will have the characteristics of the spectral lines in the expanding shell model. The situation for the halo model is similar, as shown by Fig.35 and Fig.36. From the above examples, we can conclude that the change of velocity gradient has a
362
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
Fig. 31 The line profiles of the slab model with the velocity distributions of V (r) = (r/R)n Vmax
Fig. 32 The line velocity of the medium in the case of V (r) = (r/R)n Vmax
rather large effect on the shape of line profiles.
Fig. 33 The profiles of the slab model with the 1 velocity distributions of V (r) = (r/R) n Vmax
Fig. 34 The velocity of the medium in the 1 case of V (r) = (r/R) n Vmax
3.6 The Effect of the Density Gradient in the Medium on the Line Profile Now we consider the situation when a density gradient exists in the medium. For the simulation of the medium of nonuniform density, we have two difficulties: the one is the difficulty to calculate the boundary position relative to the radiation source, the another is the difficulty to calculate the optical depth. Because that the medium density depends on the position, an integration should be used for the calculation of optical depth. As a result, the time of calculations is consequently increased. Here, we discuss the expanding halo model with the velocity gradient V (r) = (r/R)Vmax . The density gradient is taken as: nH (r) =
n0 , r ≤ r0 , nb (r/R)−a , r > r0
(19)
in which, n0 , nb are the densities of hydrogen atoms at the center and on the outer boundary, respectively, and a is a parameter. At the place r0 , n0 = nb (r/R)−a , the specific values of n0 and nb are not important, what is important is the ratio between them.
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
Fig. 35 The line profiles of the halo model with the velocity distributions of V (r) = (r/R)n Vmax
363
Fig. 36 The line profiles of the halo model with the velocity distributions of V (r) = 1
(r/R) n Vmax
Fig.37 shows the simulated results of the expanding halo model under the different ratios of n0 to nb . It demonstrates that as the ratio n0 /nb increases, the destroyed symmetrical double-peaked structure is recovered gradually. This can be realized as: when n0 /nb is very large, the density is quite small for most regions in the medium, the probability that a photon is radiatively transferred is very small, large amount of photons interact with the medium mainly in the regions close to the radiation source, hence the double-peaked symmetrical structure is recovered gradually. We have considered as well the case of the expanding halo model under different values of the parameter a. From Fig.38, we can find that the variation of the parameter a will cause also the variation of the line profile, and this is demonstrated mainly by the amplitudes of the peaks at the places of blue-shifts.
Fig. 37 The line profiles of the expanding halo model with different values of the ratio n0 /nb
Fig. 38 The line profiles of the expanding halo model with different values of the parameter a
364
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
4. SUMMARY The above results can be summarized as Table 1. As is seen from the above examples, this code has approximately realized the simulations on the Lyα lines under various conditions. Besides, in technique, this code has been able to simulate the situations when hydrogen and dust have different density distributions, as well as the spatially asymmetrical models. But it has two shortages, it has not considered the observing angle of the observer, and it has not made calculations with the density distribution of the medium, which is more approximate to the reality. Table 1 Models and the simulated line profiles Line shape symmetrical double peaks asymmetrical double peaks+single peak multiple peaks P Cygni-type line damped-type line
Model slab, halo, shell models of monochromatic, static medium slab, halo, shell models of monochromatic, nonstatic medium with velocities in a certain range monochromatic shell model with expansion velocities in a certain range flat continuum + medium with a certain velocity flat continuum + a certain amount of dust
ACKNOWLEDGEMENTS The author thanks Prof. Lin Xuan-bin of the University of Science and Technology of China for his guidance and the valuable opinion on this paper, and Prof. Yuan Feng of Shanghai Astronomical Observatory for the help during the preparation of this paper. References 1
Verhamme A., Schaerer D., Maselli A., A&A, 2006, 460, 397
2
He Yun-xiang, Chen zeng-sheng, Zhou Ke-ping, Introduction of Infrared Astronomy, Beijing: Beijing
3
Rybicki G. B., Lightman A. P., Radiative Processes in Astrophysics, New York: John Wiley & Sons,
4
Hansen M., Oh S. P., New Astronomy Reviews, 2006, 50, 58
5
Hummer D. G., MNRAS, 1962, 125, 21
6
Neufeld D. A., ApJ, 1990, 350, 216
7
Neufeld D. A., ApJ, 1991, 370, 85
University Press, 1993 1979
ELSEVIER
Chinese Astrophysics 33 33(2009) (2009)347–364 347–364 ChineseAstronomy Astronomy and and Astrophysics
Monte Carlo Simulation of Lyα Line Profiles† WANG Yu University of Science and Technology of China, Hefei 230026
Abstract In this paper, a code for simulating the Lyα line profiles by the Monte Carlo method is introduced. This code is developed on the basis of an idea given by Anne Verhamme et al. and by combining with another simulation method. For different models, the Lyα line profiles under different conditions are simulated. The simulation method of this code is given in detail, and the simulated results are discussed as well. Key words: radiative transfer—line: profile
1. INTRODUCTION The Lyα line plays a very important role in different fields of astrophysics, especially as a tool for observing and studying the objects with high redshifts. It can be taken as a probe for studying the rate of star formation, the distribution and motion of matter, the distribution of dust, etc. Meanwhile, it is also used for studying the ionization state of intergalactic medium[1] . The Lyα photons discussed in this paper are produced in HII regions. The HII region is produced by a forming massive star which has the surrounding gas ionized. In some HII regions, the massive O, B stars emit the violet radiation, in which, a part of violet photons are absorbed by the surrounding gas and therefore lead to gas ionization, another part are absorbed by dust, and the photons of λ < 912 ˚ A make hydrogen ionized. Each photon absorbed by the hydrogen atom around the star will give rise to one Lyα photon with the wavelength of 1216 ˚ A and one Balmer photon. The absorbtion of Lyα photons is the basic origin of the dust heating in molecular clouds around HII regions[2] . Because of its value in scientific research, the Lyα line has been simulated by many people with different methods. By the Monte Carlo method, Verhamme et al. simulated the †
Supported by National Foundation for Talent Cultivation of Basic Science Received 2008–05–26; revised version 2008–07–07 A translation of Acta Astron. Sin. Vol. 50, No. 2, pp. 117–133, 2009
0275-1062/09/$-see front matter © 2009 Elsevier B.V. All rights reserved. doi:10.1016/j.chinastron.2009.09.001
348
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
radiative transfer of the Lyα photon in 2006[1]. On this basis, combining with the method suggested by Hansen et al.[4] , the author of this paper tries to simulate the Lyα line profiles with the Monte Carlo method.
2. METHOD OF SIMULATION 2.1 A Brief Introduction to Monte Carlo Method The Monte Carlo method (MC method) has been widely used in different scientific fields. Its advantage is: we can simulate various physical quantities by using the given probability distribution functions. We assume that y is evaluated in the range [a, b], and that its probability distribution function is f (y). Then the probability that y takes a value y b in the range [y, y + dy] is f (y)dy. Because that a f (x)dx/ a f (x)dx = ξ, ξ ∈ [0, 1], each y-value corresponds to a number in the range [0,1], so, from a randomly generated number in the range [0,1], the y-value corresponding to it can be derived reversely. This is the quintessence of the MC method. In practical applications, for the probability distribution function which is difficult to be integrated, we can use other methods, for example the rejection sampling method, to obtain a result as same as that obtained by integrating. 2.2 Physical Processes to be Considered The simulations are made on photons one by one, namely, after the simulation of one photon is terminated, the simulation on the next one is made. And every photon to be simulated will undergo the following processes: being emitted, interacting with hydrogen and interacting with dust. 2.3 Radiation Sources The photon emitted from a source is characterized by the radiation direction (θ, φ) and the frequency x: x = (ν − ν0 )/ΔνD , (1) 1/2
ΔνD = (12.85T4 /c)ν0 .
(2)
In this paper, the frequencies of Lyα photons are all scaled by x. The radiation source in this code is an isotropic point source. For the monochromatic source, all the emergent photons have the same frequency; for the source of flat continuum, the frequencies of the emergent photons are produced randomly; and for the sources of Gauss and Lorentz profiles, the frequencies of the emergent photons can be produced by the MC method. Since the source is isotropic, it emits particles towards different directions with the same probability. When a photon emerges from the source, its radiation direction is determined by the following formula: θ = cos−1 (2ξ1 − 1),
φ = 2πξ2 .
(3)
This formula is derived by the MC method, in which ξ1 and ξ2 are two random numbers distributed uniformly in the range [0,1]. 2.4 Photons’ Walk in Simulation Space Photons walk in space randomly[1,3] . As soon as a photon walks from one place to interact with the medium to the next place to interact with the medium, it has completed
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
349
its walk one time. The one-time walk is defined by the step length s and the walking direction (θ, φ). The walking direction is the direction that the photon emerges from, and the step length s should be obtained by calculations. As the probability P that a photon interacts with the medium, and the optical depth τint that the photon will penetrate through, satisfy the following relation[1,3]: P (τint ) = 1 − e−τint , (4) so, with the MC method, we can obtain: τint = − ln(ξ),
(5)
ξ is a random number distributed uniformly in the range [0,1]. In addition, τint is also a function of the length s: τint (s) = τx (s) + τd (s), (6) in which, the optical depth of hydrogen atoms is[1] : −1/2
τx (s) = σH (x)nH s = 1.041 × 10−13 T4 the optical depth of dust is:
τd = τa + τs =
NH
H(x, a) √ , π
(7)
s
σd nd (s)ds,
(8)
0
and H (x,a)[1] can be expressed as: 2 2 e−y a ∞ e−x , if |x| < xc H(x, a) = ≈ √a , , if |x| > xc π −∞ (y − x)2 + a2 πx2
(9)
in which τd = τa + τs , τa,s = πd2 Qa,s , Qa ≈ Qs ≈ 1, d = 10−6 cm. Therefore, to derive at first τint , then to derive reversely s, the step length of photon’s one-time walk is obtained. 2.5 Interaction between the Photon and the Medium When a photon completes its one-step walk, and arrives at a place to interact with the medium, it will undergo the process of radiative transfer. 2.5.1 The interaction with hydrogen When the photon interacts with a hydrogen atom, we assume that it is absorbed at first, then it is reradiated outward. This code assumes the reradiation being isotropic. In this case, the radiation direction is determined by the same method as the above-mentioned, namely, the direction is determined by two random numbers. As for the frequency of the reradiated photon, the Hummer’s redistribution function RII is employed[1,4,5,6]. RII (ν, n, ν , n )dνdΩdν dΩ expresses the probability that the incident photon with the frequency ν in the direction n reradiates at the frequency ν in the direction n . Specifically, RII is expressed as[5] : → → RII (x, − n ; x , − n ) =
→ → x − x 2 x + x γ γ g(− n,− n ) γ exp[−( ) csc2 ( )]H( sec( ), σ sec( )) , 2 4π sin γ 2 2 2 2 2
(10)
350
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
1 3 in which, gA ( n, n ) = 4π in the case of isotropic radiation, gB ( n, n ) = 16π (1 + cos2 γ) in the case of dipole radiation, and γ is the included angle between n and n . Since RII (x0 , x)dx indicates the probability that the incident photon with the frequency x reradiates at the frequency x, so, it can produce a random number. When it equals to 0 x R dx, by deriving reversely the value of x, the frequency of the reradiated photon can II −∞ be obtained. Considering the computing speed, this code has adopted the formula given by Hansen[4] . It is an improvement on the formula given by Neufeld[4,6] : √ x2 x 2i 3 22.5 , (11) R(xi , x) = π 22.5 x4i + (x3 − x 3i )2
˜ = xi − 2/xi . But this formula does in which xi is the frequency of the incident photon, x not suit the situations when the frequency approaches or equals to zero. It is apparent that if the incident frequency is zero, this formula will be invalid. For this reason, in this code the calculations are made by using Eq.(8) while | x |≤ ξ, and by using Eq.(9) while | x |> ξ, ξ is a constant. 2.5.2 The interaction with dust When the photon arrives at the place to interact with the medium, there will be certain probability to interact with dust. The probability of photon-hydrogen interaction is[1] : PH (x) =
nH σH (x) . nH σH (x) + nd σd
(12)
Hence, we can produce a random number in the range [0,1]. If it is less than PH , then the photon interacts with hydrogen; if it is greater than PH , then the photon interacts with dust. When photons interact with dust, a part of them are absorbed by dust, and the other part are scattered by dust. In this code, it is assumed that among the photons interacting with dust, one half of photons are absorbed and another half are scattered. So in order to judge whether a photon is absorbed or scattered, a random number distributed uniformly in the range [0,1] is produced. If it is less than 0.5, then the photon is scattered, otherwise the photon is absorbed. And this code assumes the scattering being isotropic. 2.6 The Treatment for the Medium in Motion For the medium in motion, the calculations are made by transferring photons to the same coordinate system as that of the medium. This mainly leads to a change of photon’s frequency. Because that no matter what is calculated, the step length or the frequency after being scattered by the medium, the calculations are all related with the photon’s frequency. If the frequency is failed to transformed into a motionless coordinate system relative to the medium, the calculations will be erroneous. 2.7 The Termination of Photon’s Walk When a photon walks out of the simulation space or it is absorbed by the medium, the simulation on this photon terminates, and the program starts the simulation on the next photon.
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
351
2.8 The Generation of Random Numbers The quality of the random numbers prepared for the MC method relates closely to the quality of the simulated results. The random number generator used by this code is the 16807 generator. The 16807 generator adopts a kind of linear congruence method. It is actually a function, which makes an integer transform to be an another integer, and therefore forms a number series, the numbers in which have a good independence. The specific formula of the 16807 generator is a(z mod q) − r[z/q], if ≥ 0 S(z) = , (13) a(z mod q) − r[z/q] + m, otherwise in which, a=16807, b=0, m=2147483647, q=127773, r=2836, and the symbol “[ ]” means taking the integer part. 2.9 The Statistics on Results In this code the statistics on results is realized through a chain table, in which every chain element possesses the corresponding frequency x and photon number n. When the simulation on a photon terminates, the program will have the photon number n in the corresponding chain element added with 1. When the simulations on all photons are finished, the program will have the number n in each chain element divided by the total number of photons, and obtain the photon’s probability. And the program will store the photon frequency x and the corresponding probability into a file.
3. RESULTS OF SIMULATIONS For making simulations of Lyα line profiles, this code adopts the parameters taken by Verhamme et al. 3.1 Slab Model 3.1.1 Static, dust-free and monochromatic slab model As shown in Fig.1, the static monochromatic slab model is composed of three parallel infinite planes, the intermediate plane is the radiation source, the photons generated here will enter into the medium enclosed by the other two surfaces. The two terminal surfaces are the boundaries, on these two surfaces the optical depth relative to the radiation source is τ0 . Using this model, this code has obtained a symmetrical double-peaked structure. The shape of peaks is influenced by the optical depth on the boundary relative to the source. The less the optical depth on the boundary, the sharper the peaks, and the greater the optical depth, the more “short and fat” the peaks. This is easy to understand, the greater the optical depth, the more the processes of radiative transfer a photon should undergo. Using this model, Neufeld has derived an analytic solution of line profiles[1,6] : √ 6 x2 1 , (14) J(τ0 , x) = √ 3 24 π aτ0 cosh[ π /54(|x|3 /aτ0 )] and it is compared with the simulated result obtained by this code in Fig.2. From Fig.2 we can find that the two kinds of results are not very different.
352
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
Fig. 1 Diagram of the slab model
Fig. 2 The simulated result of the static dustfree monochromatic slab model, in comparison with the analytic solution derived by Neufeld
3.1.2 The dust-free monochromatic slab model with a velocity gradient Now we consider the situation when a velocity gradient exists in the medium inside the slab. The velocity gradient can be taken as: V (r) = |r/R|Vmax , R is the source position relative to the boundary, r is the medium position relative to the source, and Vmax is the medium velocity on the boundary. Fig.3 shows the Lyα line profiles simulated under different boundary velocities. We can find that as the velocity increases gradually, the originally symmetrical double peaks become asymmetrical, and finally become a single peak. The sign of the velocity on the boundary determines whether the spectral line is redshifted or blue-shifted. Fig.4 gives the simulated results for the velocities on the boundary in opposite directions. It demonstrates that when the velocity is positive, taking the source as the center, the medium expands against the source, the spectral line is red-shifted; and that when the boundary velocity is negative, the medium in-falls toward the source, the spectral line is blue-shifted.
Fig. 3 The effect of the boundary velocity on the Lyα line profile
Fig. 4 The effect of the direction of boundary velocity on the Lyα line profile
3.1.3 Effect of dust This code takes τa , the absorbtion optical depth of dust, to characterize the amount of
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
353
dust. Fig.5 shows the simulated results under different absorbtion optical depths of dust. We can find that the greater the absorption optical depth, the more the absorbed photons, and the lower the double peaks. For the situation of existence of dust, Neufeld has derived the fraction fe of photons[7] , namely the number ratio of the Lyα photons which finally emerge from the medium with respect to the totally generated Lyα photons[1,5,7] : √ 3 fe = 1/ cosh[ 5/12 |(aτ0 )1/3 τa |1/2 ] , (15) π ξ in which ξ =0.635. Fig.6 is a comparison between the faction fe derived by this code and the result given by Neufeld. It is found that the fe − (aτ0 )1/3 τa relationship derived by this code is similar with the Neufeld’s analytic solution.
Fig. 5 The line profiles under the different absorption optical depths of dust
Fig. 6 The fraction fe as a function of (aτ0 )1/3 τa derived by this code, in comparison with Neufeld’s analytic solution
3.1.4 The case of a nonmonochromatic radiation source This code can also simulate the situation that the photons radiated from the source are not monochromatic. Fig.7 and Fig.8 are the simulated results when the frequency x of the photon radiated from the source satisfies respectively: 2 1 e−x , p(x) = √ πΔνD
p(x) =
1 Γ . 4π 2 (xΔνD )2 + (Γ/4π)2
(16) (17)
The line profiles exhibit the symmetrical double-peaked structure. In the process to produce the photons, the rejection sampling method is used for the photon’s frequency. Fig.9 gives the simulated result for the flat continuum, namely the case when the radiation source emits photons with various frequencies. For the photon of large x, it will directly penetrate through the medium, but when x approaches to zero, the probability that the photon interacts with the medium increases. As a result, the curve becomes flat lines far apart from the place of x =0 and the double peaks appear around the place of x =0.
354
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
3.2 Halo Model 3.2.1 The halo model with one point source at the center of spherical boundary The structure of the model is shown in Fig.10, photons are produced at the center of the sphere, within the sphere is the homogeneous medium, expanding or in-falling with the velocity V (r) = |r/R|Vmax .
Fig. 7 The case when the source emits photons with the Gaussian profile of Eq.(16)
Fig. 8 The case when the source emits photons with the Lorentz profile of Eq.(17)
Fig.11 shows the line profiles when the halo is under the static, expanding and in-falling three states. As is shown, when the medium expands relative to the center, the spectral line is redshifted; when the medium in-falls relative to the center, the spectral line is blueshifted; and when the medium is static, the spectral line exhibits a symmetrical double-peaked structure.
Fig. 9 The case of flat continuum
3.2.2 The model with radiation sources distributed uniformly in the medium Fig.12 is the diagram of the model with the radiation sources distributed uniformly in the medium. From this figure we can find that the radiation sources are uniformly distributed in a sphere within the simulation space, not in the whole simulation space. It is necessary, if the sources are uniformly distributed in the whole simulation space, then the
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
355
photons produced on the boundary will have a large probability to escape from the medium
Fig. 10 The structure of the expanding halo with one central point source
Fig. 11 The profiles of the static, expanding and in-falling halos with one central point source
through the boundary. As a result, the produced spectral line will have a sharp peak at the place of x =0. For this reason, this code has the photon sources distributed uniformly within the sphere with the radius of 0.85Rmax . The produced spectral lines are shown in Fig.13. When the medium is static, the spectral line exhibits a symmetrical double-peaked structure; when the medium expands relative to the center, the spectral line is redshifted; and when the medium in-falls relative to the center, the spectral line is blueshifted.
Fig. 12 The structure of the expanding halo with uniformly-distributed radiation sources
Fig. 13 The profiles of the static, expanding and in-falling halos with uniformly-distributed radiation sources
3.2.3 The effect of NHI on the spectral line When NHI becomes large, the probability of the photon-medium interaction increases. In the medium the radiative transfer of atoms onto photons will make the photon’s frequency x increasingly deviate from its original value. So, the greater the value of NHI , the greater the leftshift of the peak. In the case of relatively small NHI , many photons will leave the
356
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
medium not until they have interacted with the medium. Fig.14 and Fig.15 are the line profiles under the different values of NHI for the expanding halo models with one central point source and uniformly distributed sources.
Fig. 14 The line profiles of the expanding halo with one central point source under different values of NHI
Fig. 15 The line profiles of the expanding halo with uniformly-distributed sources under different values of NHI
3.2.4 The effect of medium expansion velocity on the spectral line Fig.16 displays the simulated result on the expanding halo model with one central point source under different expansion velocities on the boundary. From this figure we can find that as the expansion velocity increases gradually, the symmetrical double-peaked structure is gradually destroyed, and there finally remains a single peak.
Fig. 16 Line profiles of the expanding halo with one central point source under different expansion velocities on the boundary
3.2.5 The case of flat continuum When the frequencies of the photons emitted from the source form a flat continuum, the P Cygni-type spectral line will result. Fig.17 and Fig.18 are respectively the line profiles of the expanding medium with a central point source inside and the in-falling medium with the radiation sources uniformly
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
357
distributed, under different absorption optical depths of dust. As is shown, for the expanding medium, the emission line appears at the place of redshift, the absorption line appears at the place of blueshift; and for the in-falling medium, the absorption line appears at the place of redshift, the emission line appears at the place of blueshift. We can find as well that the increase of the dust content will change violently the shape of the spectral line. As the dust content increases, the spectral line gradually exhibits the shape of damped type.
Fig. 17 The profiles of the expanding halo with a central point source of flat continuum
Fig. 18 The profiles of the in-falling halo with uniformly-distributed sources of flat continuum
3.3 Expanding-shell Model The model of expanding shell is illustrated by Fig.19. The radiation source is the point source at the center of the spherical shell, the space between the source and the shell is vacuum, the medium is totally concentrated in the shell expanding with a constant velocity. In the space enclosed by the shell photons walk straightly without any obstruction, only within the shell they will interact with the medium. The spectral line produced by this model is a little peculiar. The spectral line is the superposition of the multiple peaks at nVexp /b. Fig.20 is a good example, in which the dashed lines mark the x-values of the integer multiples of Vexp /b. This result is different from the simulated result obtained by Verhamme et al.[1] , but both results have a peak of the redshift 2Vexp . In certain conditions, the composite multiple peaks of the spectral line can be separated. Fig.21 gives the simulated results under different expansion velocities of the shell. As is shown, when the shell is at rest, the spectral line exhibits the symmetrical double-peaked structure. When the velocity is increased, the spectral line is split into several peaks, especially when the velocity is quite large, the distance between two peaks is enlarged, and the intensity distribution of the peaks is also changed. The Doppler parameter will affect the separation between two peaks. As is shown by Fig.22, when the Doppler parameter is rather large, the separation between two different peaks becomes small, therefore the superposition leads to a single peak at the redshift position; and when the Doppler parameter becomes small, the separation between two different peaks is enlarged. We can find as well that the variation of the Doppler parameter
358
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
Fig. 19 Diagram of the expanding shell model
Fig. 20 An exemplified profile of the expanding-shell model. The dashed lines correspond to x = nVexp /b.
causes only the variation of the separation between two different peaks, but does not change apparently the distribution of the peak intensities.
Fig. 21 The profiles of the expanding shell with different expanding velocities
Fig. 22 The profiles of the expanding shell with different Doppler parameters
The value of NHI in the medium will affect also the shape of the spectral line. As shown in Fig.23, when NHI is relatively small, the peak closer to x = 0 tends to have a higher intensity, but it is not apparent for peaks far apart from x = 0. This is because that most protons have escaped from the medium without interacting with it, so that the spectral line has only one or two high peaks close to x = 0. When NHI is increased, the intensity distribution of different peaks varies, the peak close to x = 0 begins to decline, and the peak far apart from x = 0 begins to rise. As NHI increases continuously, the intensity distribution of peaks becomes more uniform, the highest peak shifts toward the place of larger redshift; and when NHI is quite large, the different peaks are superposed into a single “short and fat” peak. The amount of dust will change as well the shape of the spectral line. As shown in Fig.24, dust will make the brightness of the spectral line decrease.
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
Fig. 23 The profiles of the expanding shell with different NHI ’s
359
Fig. 24 The profiles of the expanding shell with different τa ’s
The value of NHI and the absorption optical depth τa of dust have a major effect on the fe -value of the spectral line. The Doppler parameter b has a less effect on the fe -curve, for the same τa , the larger b gives the smaller fe . Fig.25 compares the fe -τa relationship simulated by this code with the function exp(−τa ). We can find that the derived fe -τa relationship resembles exp(−τa ). The less the value of NHI , the closer to exp(−τa ) the relationship. This point is easy to understand. The less the value of NHI , the larger the proportion of dust in the medium, the larger the probability of the photon-dust interaction, and therefore the closer the fe -τa relationship to exp(−τa ). In reverse, the greater the value of NHI , the longer the path length that a photon walks in the medium, the greater the probability that it is absorbed, and therefore the smaller the fe value.
Fig. 25 NHI ’s
The fe -τa relationships under different
Fig. 26 The fe -τa relationships under different Doppler parameters b
From the above analysis, we can find that the line profile of the expanding shell model is composed by the equally-spaced peaks and the intensity distribution curve of peaks, as shown in Fig.27. The expansion velocity Vexp , the absorption optical depth of dust, and the Doppler parameter b determine jointly the spacing of peaks, Vexp and NHI determine jointly
360
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
the distribution curve of peak intensities, and the effect of the Doppler parameter b on the distribution curve of peak intensities is not evident. Next, we consider the case of flat continuum. From Fig.28, we can find that for the expanding shell, under the condition of flat continuum, the spectral line exhibits the shape of P Cygni-type lines, at the place of redshift occurs the emission line, and at the place of blueshift occurs the absorption line, as the amount of dust increases, the spectral line becomes gradually the damped shape.
Fig. 27 The line profile of the shell model can be decomposed into two parts.
Fig. 28 The line profile under the condition of flat continuum
3.4 Damped-type Spectral Lines As the damped-type spectral lines are observed very often, the study on the formation of the damped-type Lyα line is very meaningful. From the above results, we can find that for the P Cygni-type spectral lines, when the amount of dust increases, the spectral line tends significantly to be of the damped-type. Under the condition that the radiation source is of flat continuum, the damped-type spectral line can be constructed by increasing the amount of dust. This can be tested by the simplest model of static slab. From Fig.29, we can find that when the dust content increases to a certain extent, the spectral line becomes the damped-type. The Voigt absorption line[1] I(x) = I0 exp(−τ0 H(x, a))
(18)
is generally used for simulating the damped-type Lyα line, so we compare them in Fig.30. We can find that a difference exists in the widths of both absorption lines, but this is probably due to the different aspects being considered. The Voigt absorption line model has considered the observing angle. It assumes that some mediums exist between the radiation source and the observer. These mediums have the photons scattered out of the line of sight of the observer. The probability that the photons can pass through the mediums is proportional to exp(−τ0 H(x, a)), for the photons far apart from the place of x =0, it is easy to pass through the mediums and not to be scattered, but for the photons close to the place
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
361
of x =0, it is easy to be scattered out of the line of sight and not to be observed by the observer[1]. However, what is considered by this code is the effect of the dust absorption, only the existence of dust can make the double peaks on the spectral line flattened, and makes the spectral line exhibit the shape of damped-type lines.
Fig. 29 The line profiles of the static dusty slab with a flat-continuum source under different absorption optical depths of dust
Fig. 30 The damped-type line simulated by this code, in comparison with the Voigt absorption line
3.5 The Effect of the Velocity Gradient in the Medium on Line Profiles Now we consider the effect of the variation of the velocity gradient in the medium on spectral lines. There are many kinds of velocity distributions, here we exemplify several kinds of them to make a discussion. At first, we take the slab model and consider the case of V (r) = (r/R)n Vmax . We can find that in the case of V (r) = |r/R|10 Vmax , the symmetrical double-peaked structure is destroyed to a extent less than the case of V (r) = |r/R|Vmax . This is easy to understand, if the power index of |r/R| is very large, then the medium velocities in the inner most regions of the medium will be very small, only the medium velocities near the boundary can be effective. When a photon is far apart from the boundary, it moves as it does in a static medium, and when it is close to the boundary, because of the relatively large medium velocity, and of the relatively thin optical thickness, the photon becomes transparent. Namely the photon will not be affected by the medium near the boundary and emits directly out of the medium, thus the spectral line quite approximates to the symmetrical double-peaked spectral line in the static medium. For this point, Fig.32 provides a better description. 1 Now, we consider the case of V (r) = (r/R) n Vmax . From Fig.33, we can find that when the index n increases, the spectral line exhibits gradually the characteristics of the expanding shell model, i.e., multiple peaks are superposed together at the places of redshifts, and this is easy to understand. As shown by Fig.34, when n increases, the velocities in most regions of the medium are approximately uniform, like a shell far away from the radiation source, therefore the exhibited spectral line will have the characteristics of the spectral lines in the expanding shell model. The situation for the halo model is similar, as shown by Fig.35 and Fig.36. From the above examples, we can conclude that the change of velocity gradient has a
362
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
Fig. 31 The line profiles of the slab model with the velocity distributions of V (r) = (r/R)n Vmax
Fig. 32 The line velocity of the medium in the case of V (r) = (r/R)n Vmax
rather large effect on the shape of line profiles.
Fig. 33 The profiles of the slab model with the 1 velocity distributions of V (r) = (r/R) n Vmax
Fig. 34 The velocity of the medium in the 1 case of V (r) = (r/R) n Vmax
3.6 The Effect of the Density Gradient in the Medium on the Line Profile Now we consider the situation when a density gradient exists in the medium. For the simulation of the medium of nonuniform density, we have two difficulties: the one is the difficulty to calculate the boundary position relative to the radiation source, the another is the difficulty to calculate the optical depth. Because that the medium density depends on the position, an integration should be used for the calculation of optical depth. As a result, the time of calculations is consequently increased. Here, we discuss the expanding halo model with the velocity gradient V (r) = (r/R)Vmax . The density gradient is taken as: nH (r) =
n0 , r ≤ r0 , nb (r/R)−a , r > r0
(19)
in which, n0 , nb are the densities of hydrogen atoms at the center and on the outer boundary, respectively, and a is a parameter. At the place r0 , n0 = nb (r/R)−a , the specific values of n0 and nb are not important, what is important is the ratio between them.
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
Fig. 35 The line profiles of the halo model with the velocity distributions of V (r) = (r/R)n Vmax
363
Fig. 36 The line profiles of the halo model with the velocity distributions of V (r) = 1
(r/R) n Vmax
Fig.37 shows the simulated results of the expanding halo model under the different ratios of n0 to nb . It demonstrates that as the ratio n0 /nb increases, the destroyed symmetrical double-peaked structure is recovered gradually. This can be realized as: when n0 /nb is very large, the density is quite small for most regions in the medium, the probability that a photon is radiatively transferred is very small, large amount of photons interact with the medium mainly in the regions close to the radiation source, hence the double-peaked symmetrical structure is recovered gradually. We have considered as well the case of the expanding halo model under different values of the parameter a. From Fig.38, we can find that the variation of the parameter a will cause also the variation of the line profile, and this is demonstrated mainly by the amplitudes of the peaks at the places of blue-shifts.
Fig. 37 The line profiles of the expanding halo model with different values of the ratio n0 /nb
Fig. 38 The line profiles of the expanding halo model with different values of the parameter a
364
WANG Yu / Chinese Astronomy and Astrophysics 33 (2009) 347–364
4. SUMMARY The above results can be summarized as Table 1. As is seen from the above examples, this code has approximately realized the simulations on the Lyα lines under various conditions. Besides, in technique, this code has been able to simulate the situations when hydrogen and dust have different density distributions, as well as the spatially asymmetrical models. But it has two shortages, it has not considered the observing angle of the observer, and it has not made calculations with the density distribution of the medium, which is more approximate to the reality. Table 1 Models and the simulated line profiles Line shape symmetrical double peaks asymmetrical double peaks+single peak multiple peaks P Cygni-type line damped-type line
Model slab, halo, shell models of monochromatic, static medium slab, halo, shell models of monochromatic, nonstatic medium with velocities in a certain range monochromatic shell model with expansion velocities in a certain range flat continuum + medium with a certain velocity flat continuum + a certain amount of dust
ACKNOWLEDGEMENTS The author thanks Prof. Lin Xuan-bin of the University of Science and Technology of China for his guidance and the valuable opinion on this paper, and Prof. Yuan Feng of Shanghai Astronomical Observatory for the help during the preparation of this paper. References 1
Verhamme A., Schaerer D., Maselli A., A&A, 2006, 460, 397
2
He Yun-xiang, Chen zeng-sheng, Zhou Ke-ping, Introduction of Infrared Astronomy, Beijing: Beijing
3
Rybicki G. B., Lightman A. P., Radiative Processes in Astrophysics, New York: John Wiley & Sons,
4
Hansen M., Oh S. P., New Astronomy Reviews, 2006, 50, 58
5
Hummer D. G., MNRAS, 1962, 125, 21
6
Neufeld D. A., ApJ, 1990, 350, 216
7
Neufeld D. A., ApJ, 1991, 370, 85
University Press, 1993 1979