EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 659 (March 2022) A&A, 659 (2022) A149 Full HTML
Free Access
Issue
A&A
Volume 659, March 2022
Article Number A149
Number of page(s) 12
Section Astrophysical processes
DOI https://doi.org/10.1051/0004-6361/202142437
Published online 22 March 2022

© ESO 2022

1. Introduction

Cosmic dust is a vital ingredient of the circumstellar and interstellar media. Dust grains have an important role in several astrophysical processes, including the photoelectric heating of interstellar gas (Draine 1978) and the formation of molecular hydrogen (Gould & Salpeter 1963; Hollenbach & Salpeter 1971; Bron et al. 2014). One of the most important characteristics of cosmic dust grains is the efficient interaction with radiation: dust grains efficiently absorb ultraviolet (UV) and optical radiation and subsequently re-emit the absorbed energy as thermal radiation in the infrared. Even though the dust-to-stellar mass ratio in galaxies typically ranges between 10−5 and 10−3 (Nersesian et al. 2019), on average one-third of all starlight emitted in galaxies is absorbed by dust and transformed to thermal infrared radiation (Popescu & Tuffs 2002; Viaene et al. 2016; Bianchi et al. 2018).

Dust grains are also very efficient at scattering radiation, especially at UV and optical wavelengths. At these wavelengths, absorption and scattering have roughly the same importance (Draine 2003a; Zubko et al. 2004). This implies that if we want to infer the intrinsic properties of cosmic systems from UV and optical radiation, we need to properly account for scattering effects. Properly modelling scattering is particularly challenging as multiple scattering effects are often hard to predict or are counter-intuitive (see e.g. Witt et al. 1992a; Byun et al. 1994).

Properly modelling scattering by interstellar dust requires a good characterisation of the scattering phase function, that is, the probability density of the scattering angle after a single scattering event. Observations of reflection nebulae, dark clouds, and diffuse Galactic light have shown that scattering by interstellar grains is strongly anisotropic, especially at short wavelengths (Lillie & Witt 1976; Witt et al. 1992b; Calzetti et al. 1995; Gibson & Nordsieck 2003).

The Henyey-Greenstein (HG) phase function was originally introduced by Henyey & Greenstein (1941) to describe scattering by interstellar dust grains in the Milky Way. It has become the standard phase function for describing anisotropic scattering by dust grains in astronomy (Bruzual et al. 1988; Witt & Gordon 1996; Xilouris et al. 1999; Ferrara et al. 1999; Bailey & Kedziora-Chudczer 2012; De Geyter et al. 2014) and has been widely used in other fields as well, including atmospheric sciences (Danielson et al. 1969; Boucher 1998), oceanography (Hakim & McCormick 2003; Stramski et al. 2004; Light et al. 2003), and biomedical imaging (Tuchin 1993; Arridge 1999; Roggan et al. 1999). The HG phase function is a simple analytical formula that depends on just a single parameter, g. By varying this anisotropy parameter, the HG phase function can describe a wide variety of scattering regimes, from completely isotropic to extremely forward or backward scattering. Another interesting characteristic of the HG phase function is that it is easily integrated into the most popular radiative transfer solution methods. Unfortunately, the HG phase function has a number of disadvantages too. The main disadvantage is that it is a relatively poor approximation to the real phase function for scattering off dust grains outside the optical range, in particular at UV and near-infrared (NIR) wavelengths (Draine 2003b; Sharma & Roy 2008).

Several alternatives have been proposed for the HG phase function. One interesting alternative was proposed by Draine (2003b), a two-parameter generalisation of the HG phase function. It has a number of interesting characteristics and has been shown to provide a better approximation to the scattering function for interstellar dust at optical and NIR wavelengths compared to the classic HG phase function. At UV wavelengths, however, the Draine phase function still poorly describes the scattering properties of astrophysical dust. An alternative two-parameter extension of the HG phase function was presented by Reynolds & McCormick (1980) and is known as the Reynolds-McCormick (RM) or the Gegenbauer kernel phase function. The advantage of the RM over the classic HG phase function is that it better describes strongly anisotropic scattering. It has primarily been applied in the context of scattering in biological tissues and fluids (e.g. Yaroslavsky et al. 1997; Hammer et al. 1998; Lindbergh et al. 2009; Karlsson et al. 2012; Calabro & Cassarly 2015). The RM phase function offers more flexibility than the one-parameter HG phase function but generally fails for distributions that are both forward and backward peaked. This is particularly problematic in the UV and NIR regimes, where the scattering phase function corresponding to interstellar dust grains shows this feature. A solution to that problem can be to use a linear combination of two HG scattering functions, resulting in the two-term Henyey-Greenstein (TTHG) phase function (Irvine 1965, 1968).

Wang et al. (2019) recently presented an interesting study of the light scattering properties from a polymer layer containing sub-micron metallic and dielectric particles. They experimentally measured the scattering phase function, which exhibited strongly forward and backward peaked scattering patterns. They proposed a two-term Reynolds-McCormick (TTRM) phase function that is a linear combination of two RM phase functions, resulting in a five-parameter analytical phase function. They fitted this phase function and several other analytical phase functions to their experimental data and found that the TTRM phase function consistently gave the best fit in all cases. Very recently, Harmel et al. (2021) confirmed the superiority of the TTRM phase function proposed by Wang et al. (2019) to describe the scattering by microalgae and mineral hydrosols in water.

The main issue with this TTRM phase function is that it contains five free parameters and hence lacks the attractive simplicity and ease of use of a simpler function, such as the HG phase function. The same disadvantage applies to more numerical phase functions, which are often based on expansions of the phase function in terms of Legendre or Chebyshev polynomials (e.g. Chandrasekhar 1950; Sharma et al. 1998; Sharma & Roy 2008; Zhang et al. 2017).

In this paper we investigate alternatives for the HG phase function to describe the scattering properties of dust grains in the interstellar medium in galaxies. Our goal is to find a balance between realism and complexity: on the one hand the scattering phase function should be flexible enough to provide an accurate fit to the actual scattering properties of dust grains over a wide wavelength range, and on the other hand it should be simple enough to be easy to handle, especially in the context of radiative transfer calculations.

This paper is structured as follows. In Sect. 2 we present the dust model and the resulting empirical scattering phase function data we use for our calculations. In Sect. 3 we fit different analytical phase functions to the empirical data. In Sect. 4 we present a new three-parameter phase function and discuss its properties and benefits. Our conclusions are summed up in Sect. 5.

2. Dust model

We adopted the BARE-GR-S dust model presented by Zubko et al. (2004), one of the most popular and commonly adopted models for interstellar dust. It consists of a mixture of spherical graphite and silicate grains and polycyclic aromatic hydrocarbons (PAHs). The size distribution and relative abundances of each component are derived by simultaneously fitting the far-ultraviolet (FUV) to NIR extinction, the diffuse infrared emission in the Milky Way, and elemental abundance constraints on the dust.

The BARE-GR-S model was used as the standard dust model in the suite of TRUST radiative transfer benchmark studies (Camps et al. 2015; Gordon et al. 2017). In the frame of this benchmarking effort, data files containing the detailed optical properties of the dust mixture are available online1. In particular, size- and population-integrated effective scattering phase functions are available with a 1° angular resolution for 1201 wavelengths from 10 Å to 10 mm. Throughout this paper we use the normalisation convention

(1)

We denote these BARE-GR-S as the empirical phase functions.

In this paper we focus on the UV–NIR wavelength regime, covering the range from about 0.1 to 5 μm. In Fig. 1 we show the empirical scattering phase function Φ(θ) as a function of the scattering angle θ at the effective wavelengths of several standard broadband filters. In Fig. 2 we show two representative cases in the more familiar polar representation. At UV wavelengths, the scattering is strongly forward peaked (θ ∼ 0°) with a very minor secondary backward peak (θ ∼ 180°). As the wavelength increases into the optical regime, the phase function becomes gradually less forward peaked (i.e. more isotropic), but the relative contribution of the backward scattering peak gradually increases. At NIR wavelengths, the phase function is largely isotropic, with a forward scattering peak that is only mildly dominant over the backward scattering peak. At the longest wavelengths considered here, the empirical phase function converges to the Rayleigh phase function,

(2)

thumbnail Fig. 1.

Empirical scattering phase functions (grey dots) of the Zubko et al. (2004) BARE-GR-S dust model at the effective wavelengths of eight common UV–NIR bands. For the sake of clarity we only plotted the data points at a resolution of 5 deg, whereas the BARE-GR-S dataset contains phase function data at 1 deg resolution. The solid lines in each upper panel correspond to analytical phase function fits to these empirical data: the HG phase function (red, Sect. 3.1), the Draine phase function (purple, Sect. 3.2), the TTRM phase function (cyan, Sect. 3.4), and the new TTU2 phase function (gold, Sect. 4). Bottom panels: relative difference between the analytical model and the empirical data.

3. Analytical phase functions

In this section, we fit different analytical phase functions used in the literature to the empirical data corresponding to the BARE-GR-S dust model. For an analytical phase function Φana(θ, p) characterised by free parameters p, we determine at each wavelength the best-fitting parameters pfit by minimising the relative error, defined by Draine (2003b) as

(3)

This metric is used to determine the best-fitting parameters, and can be used to compare the global quality of the fit of different analytical phase functions. For a more detailed investigation of how well a given parameterisation describes the empirical phase function, it is also useful to consider the relative difference between the empirical phase function and the best-fitting analytical phase function,

(4)

explicitly as a function of the scattering angle θ.

3.1. The HG phase function

Scattering off dust grains is usually described using the HG phase function,

(5)

This expression contains a single free parameter, g, which in this particular case is nothing but the anisotropy parameter

(6)

In general, the anisotropy parameter indicates to which degree the scattering phase function deviates from isotropy. The wavelength dependence of the asymmetry parameter is usually one of the constraints or sanity checks for the construction of physical dust models (e.g. Zubko et al. 2004; Siebenmorgen et al. 2014; Ysard et al. 2018; Hensley & Draine 2021). The fact that the HG phase function contains the anisotropy parameter as an explicit free parameter makes it very attractive, and explains to some degree its popularity.

The solid red lines in Figs. 1 and 2 represent HG phase functions fitted to the empirical BARE-GR-S data2. Overall, the HG phase function provides a reasonable first-order approximation to the actual phase function, but there are significant deviations. These are most clearly seen in the bottom parts of each panel of Fig. 1, where we show the relative difference between the empirical phase function and the best-fitting HG phase function. The HG phase function underestimates the backward scattering peak at all wavelengths. At UV and NIR wavelengths the underestimation amounts up to 30%, while at optical wavelengths it is reduced to ∼20%. Similarly, the HG phase function underestimates the forward scattering peak of the phase function at all wavelengths, with magnitudes varying from ∼10% in the optical up to 30% in the NIR and more than 40% in the near-ultraviolet (NUV) band. The only exception is the FUV band, where the distribution function has a particular shape in the forward direction and the actual peak is underestimated by 30%. For scattering angles around 60°, the HG phase function overestimates the empirical phase function by different degrees, ranging from around 10% in the optical up to 30% in the NIR and even up to 50% in the NUV.

thumbnail Fig. 2.

Scattering phase function for two representative wavelengths (0.23 and 1.65 μm), shown in polar format. The meanings of the lines and symbols are the same as in Fig. 1.

3.2. The Draine phase function

As indicated in the Introduction, a generalisation of the HG phase function was proposed by Draine (2003b). He proposed the following analytical phase function:

(7)

This phase function was designed to bridge the cases of HG and Rayleigh scattering: for α = 0 the Draine phase function (7) reduces to the HG phase function (5), for g = 0 and α = 1 it reduces to the Rayleigh phase function (2).

The purple lines in Figs. 1 and 2 represent the best fitting Draine phase functions to the empirical BARE-GR-S data. Compared to the HG phase function, the improvement is particularly evident at the longest wavelengths, shown on the bottom row. This is not surprising as the Draine phase function is a generalisation of both the HG phase function and the Rayleigh phase function that is applicable at long wavelengths. At blue and particularly UV wavelengths there is less improvement. In the FUV the fit to the scattering phase function is slightly improved, but in the NUV band the forward scattering peak is still underestimated by about 40% and the phase function at ∼60° is overestimated by 50%.

3.3. The TTHG phase function

The TTHG phase function consists of a linear combination of two HG phase functions, one with a forward scattering and one with a backward scattering peak,

(8)

This three-parameter phase function was first proposed by Irvine (1965, 1968) and subsequently applied by several other authors (e.g. Kattawar 1975; Witt 1977; Zhang & Li 2016; Chen-Chen et al. 2021).

The pink lines in Figs. 1 and 2 represent the best fitting TTHG phase functions to the empirical BARE-GR-S data. Compared to the single HG phase function, the TTHG phase function results in a significant improvement in the fit quality over the entire wavelength range. In particular, by combining two HG phase functions with a similar relative weight and with g1 ≈ −g2, the phase function at NIR wavelengths can properly be fitted (see bottom panel of Fig. 2). At UV wavelengths, a significant discrepancy remains, as clearly illustrated in the top panel of Fig. 2.

3.4. The TTRM phase function

Inspired by the conclusions by Wang et al. (2019) and Harmel et al. (2021), our next step was to test whether the TTRM phase function provides a suitable option to describe the scattering phase function for interstellar dust. The TTRM phase function can be written as

(9)

where the RM phase function ΦRM(θ, g, α) is given by

(10)

As the RM phase function itself is already a generalisation of the HG phase function, it is only natural that the five-parameter TTRM phase function can cover a much richer variety in phase function shapes. In particular, by choosing g1 > 0 and g2 < 0 it can describe phase functions with a distinctive forward and backward scattering peak.

Again, we minimised the relative difference hrel to determine the parameters of the best-fitting TTRM phase function to the empirical phase functions at wavelengths between 0.1 and 5 μm. Not surprisingly, the result is a significantly improved fit. This fit is shown as the cyan lines in Figs. 1 and 2. The TTRM phase function reproduces the shape of the empirical phase functions very well at all wavelengths, and the relative difference δrel(θ) remains under 5% for all bands and all scattering angles. Only for the UV bands, δrel(θ) slightly exceeds 5% around the forward and the backward directions.

4. The TTU2 phase function

4.1. Towards a new three-parameter phase function

In principle, our analysis could stop here: in line with Wang et al. (2019) and Harmel et al. (2021) we have reached the conclusion that the TTRM phase function provides a superior fit to the scattering phase function of interstellar dust. However, this phase function contains five free parameters, which makes it rather complex. Moreover, a closer analysis of the fitting results shows some interesting results.

Figure 3 shows the model parameters of the best-fitting TTRM phase function as a function of wavelength. This figure shows a number of peculiarities of the fits. Throughout the UV and optical range, the parameter α1 tends to prefer values around α1 = 1, whereas the parameter α2 appears to prefer negative values, and usually the lowest possible value α2 = −1. For α = −1, the RM phase function reduces to the isotropic phase function Φiso(θ) = 1/4π, irrespective of the value of g. This makes the value of g2 at UV and optical values largely irrelevant, which explains the spurious jumps of g2 between −1 and 0 in the third panel of Fig. 3. In summary, at UV and optical wavelengths, the phase function is well fitted by a linear combination of a forward scattering phase function with α1 ∼ 1 and an isotropic phase function. Around λ = 1 μm, both values of α shoot to αmax = 5, the upper bound for α allowed by our fitting routine, while the g parameters converge to g1 = −g2 = 0.098. Interestingly, when the upper bound for α in the fitting routine is changed, the same qualitative behaviour persists. For example, with αmax = 10, both α1 and α2 quickly rise to this value beyond λ ∼ 1 μm, while the g parameters converge to g1 = −g2 = 0.056. The difference in hrel between these models is negligible.

thumbnail Fig. 3.

Model parameters of the best-fitting TTRM phase functions as a function of wavelength. The dots indicate the effective wavelengths of the GALEX, SDSS, 2MASS, and WISE broadband filters.

This strange behaviour clearly shows that the TTRM model contains a number of degeneracies, and suggests that a five-parameter phase function is probably overkill to describe the scattering phase function. Guided by the fact that the best-fitting value of α1 remains rather constant over the entire UV and optical range, and that the quality of the fits in the NIR seem to be relatively insensitive to the value of α1 and α2, we propose a restricted TTRM phase function, in which we fix α1 = α2 to a single fixed value α that is a wavelength-independent model parameter. The obvious choice would be , in which case we end up with the TTHG phase function. The fact that, for our general five-parameter TTRM fit, α1 always tends to prefer higher values than suggests that this might not the optimal choice.

We therefore repeated our fitting exercise using different values of α. At every wavelength, we fitted a constrained TTRM phase function with α1 = α2 = α to the empirical data, and determined the best-fitting values of the four free parameters p = (g1, g2, α, γ). We finally determined the best value of α as the one that has the lowest relative error between the model fits and the empirical phase function, averaged over 12 wavelengths in the 0.1 − 5 μm range3. Figure 4 shows this average relative error ⟨hrel⟩ as a function of α. It turns out that the optimal value for α is indeed not , the value corresponding to the TTHG phase function. Interestingly, the value for α that minimises ⟨hrel⟩ is found to be exactly equal to one.

thumbnail Fig. 4.

Relative error ⟨hrel⟩, averaged over the standard UV, optical, and NIR broadband filters, for the TTRM phase function as a function of the fixed α = α1 = α2.

4.2. Characteristics of the TTU2 phase function

Putting all pieces together, we find the following expression for a new analytical phase function,

(11)

where

(12)

We propose the name ultraspherical-2 (U2) phase function for the simple phase function (12), since it is the generating function for the set of ultraspherical polynomials of order 2 (Appendix A). In analogy with the TTRM and TTHG phase functions, we then propose the name two-term ultraspherical-2 (TTU2) phase function for the phase function represented by Eq. (11). It is easy to verify that the U2 phase function, and therefore also the TTU2 phase function, satisfies the normalisation convention (1).

The TTU2 phase function is characterised by three free parameters with the following meanings. First, the parameter g1, with values between 0 and 1, characterises the shape and the strength of the forward scattering peaking part of the phase function. The larger the value of g1, the more asymmetric and peaked the forward scattering part of the phase function. If g1 = 0, the forward scattering part is isotropic, if g1 = 1, it is completely forward scattering. Similarly, g2, with values between −1 and 0, characterises the shape and strength of the backward scattering peak. Third, the parameter γ sets the relative weight of both components and can take any value, as long as the phase function remains positive for all scattering angles. We limit its range to 0 ≤ γ ≤ 1.

It is important to realise that the parameter g that characterises the U2 phase function is not the anisotropy parameter, contrary to the case of the HG phase function. For the U2 phase function we find instead

(13)

For the TTU2 three-parameter phase function (11), expression (13) is readily generalised to

(14)

4.3. Applicability to interstellar dust

The fits of the TTU2 phase function to the empirical BARE-GR-S phase functions are shown as the golden lines in Fig. 1. Compared to the general TTRM phase function, the fits are hardly poorer, in spite of the fact that we only fit three instead of five free parameters. Especially in the UV and the optical, the TTU2 phase function provides an almost equally good fit to the empirical data, with almost identical relative differences. At NIR wavelengths the TTU2 phase function slightly overestimates the forward and backward scattering peak and also the scattering around the minimum at 90°. The relative error remains within the 5% range, however.

As a summary, Fig. 5 compares the relative error hrel of the different phase functions discussed as a function of wavelength over the wavelength range from 0.1 to 5 μm. The HG phase function, which only contains one free parameter, clearly provides the poorest fit to the empirical phase function over the entire wavelength range. Relative errors range between ≳10% at 1 μm to almost 100% at the shortest wavelengths. At UV and optical wavelengths, the two-parameter Draine phase function performs only marginally better than the HG phase function, but at NIR wavelengths, this phase function outperforms all other options. The other three phase functions shown are all two-term phase functions. The TTRM phase function, with its five free parameters, obviously has the best performance, with a relative error of a few percent over the entire wavelength range. The new TTU2 phase function, with only three free parameters, is only slightly poorer. It is characterised by a relative error below 5% over the entire optical and NIR wavelength range, and only slightly larger at the shortest wavelengths. The TTHG phase function, also with three free parameters, consistently performs more poorly than the TTU2 phase function over the entire wavelength range considered, except in a small wavelength interval between 0.6 and 0.74 μm and at the shortest UV wavelengths (below 0.12 μm). Based on this plot and the previous discussion, we argue that the TTU2 phase function is the compromise looked for: it provides an accurate fit to the empirical phase function over a wide wavelength range, and with three free parameters it is still simple enough.

thumbnail Fig. 5.

Comparison of the relative error hrel as a function of wavelength for the different phase functions considered in this paper. The dots indicate the effective wavelengths of the GALEX, SDSS, 2MASS, and WISE broadband filters.

In order to better understand the meaning of the three parameters of the TTU2 phase function, we show them explicitly as a function of wavelength in Fig. 6, and we list the values at the central wavelengths of a number of standard broadband filters in Table 1. The bottom panel of Fig. 6 shows that the forward scattering component of the TTU2 phase function always dominates over the backward scattering component. Interestingly, the parameter γ does not decrease monotonically as a function of wavelength: it first rises to a maximum of more than 0.9 at the FUV band, subsequently falls back to slightly below 0.8 in the NUV band, and then gradually increases to 0.9 again around 0.8 μm. In the NIR the importance of the first component gradually decreases towards the asymptotic value of 0.5 expected for a Rayleigh scattering phase function. This behaviour needs to be interpreted in conjunction with the parameters g1 and g2. At UV wavelengths, the forward scattering component is strongly peaked, whereas the backward scattering component is essentially isotropic. Moving towards optical wavelengths, the value of g1 gradually drops, whereas the backward component becomes more anisotropic. In the NIR regime, both components tend towards the same (opposite) value of moderate anisotropy, g1 ≈ −g2 ≈ 0.246.

thumbnail Fig. 6.

Relative error and model parameters of the best-fitting TTU2 phase functions as a function of wavelength. The golden lines correspond to our default BARE-GR-S dust model, and the green dots correspond to the alternative WD01-MW dust model.

Table 1.

Model parameters of the best-fitting TTU2 phase functions for the BARE-GR-S dust model at the central wavelengths of a number of standard broadband filters.

We repeated our analysis for an alternative dust model: the WD01-MW dust model, that is, the average Milky Way dust mixture presented by Weingartner & Draine (2001). This model also consists of a mixture of graphite, silicates, and PAHs, but it is based on different optical properties and a different size distribution. Within our wavelength range from 0.1 to 5 μm, explicit scattering phase functions with a 1 deg resolution are available for 27 different wavelengths on Bruce Draine’s web page4. Fitting the different analytical phase functions to the empirical WD01-MW phase functions, we find a similar qualitative behaviour as for our default BARE-GR-S model. The one-parameter HG model generally results in the poorest fit (relative differences of several ten percent), the five-parameter TTRM model in the best fit (relative differences of a few percent), and the three-parameter TTU2 model in an almost equally good fit. The top panel of Fig. 6 shows that the TTU2 phase function provides an excellent fit to the phase function of the WD01-MW dust mixture in the optical and NIR range, with relative errors of the order of a few percent. In the UV part of the spectrum, the fits are slightly poorer, with relative errors around 10%.

The remaining three panels of Fig. 6 show that, in spite of the differences between the BARE-GR-S and WD01-MW models, the fitted parameters of the TTU2 phase function agree fairly well. For the WD01-MW dust model too the prime component gradually becomes less forward scattering when moving to longer wavelengths, the second component is isotropic in the UV and more strongly anisotropic at optical wavelengths, and γ peaks at λ ∼ 0.8 μm. The main difference between the two models is a spurious jump for the different parameters, especially g2, between 0.17 and 0.182 μm. This apparent discontinuity reflects a sudden transition that is primarily seen in the backward scattering peak of the WD01-MW phase function.

4.4. Integration in radiative transfer codes

As stated in the Introduction, our objective was to search for analytical scattering phase function with two competing characteristics. On the one hand it should be flexible enough to accurately describe the actual scattering phase function for interstellar dust over the UV–NIR wavelength range, while on the other hand it should be simple enough to be easy to use, especially in the context of radiative transfer calculations.

For 1D radiative transfer problems in either plane-parallel or spherical geometry, there are various options to solve the radiative transfer equation. One approach is the so-called spherical harmonics or PL method, which builds on the expansion of the radiation field as a series of Legendre polynomials. This turns the dust radiative transfer equation, in general an integro-differential equation, into an infinite set of ordinary differential equations. This set is turned into a finite set of equations by imposing a closure condition, that is, by truncating the Legendre polynomial expansion. This set of equations can be reformulated as a classic eigenvalue problem and thus solved by standard matrix diagonalisation. Detailed descriptions of this method can be found in Flannery et al. (1980), Roberge (1983) and di Bartolomeo et al. (1995). In Baes & Dejonghe (2001a) we compared four different methods to solve the dust radiative transfer methods in a plane-parallel geometry, and we found that the spherical harmonics method was by far the most efficient one. Furthermore, Camps & Baes (2018) and Krieger & Wolf (2021) demonstrated that the spherical harmonics method can be applied in optical depth regimes not easily accessible to Monte Carlo codes.

In most astrophysical contexts, the geometry of the objects under study is complex on different scales, which makes 3D radiative transfer a necessity. The Monte Carlo technique, which uses a probabilistic approach to simulate the individual life cycles of a large number of individual photon packets, is currently by far the most widely used technique for 3D dust radiative transfer. There are several codes available that can solve the dust radiative transfer problem in an arbitrary 3D setting (e.g. Gordon et al. 2001; Robitaille 2011; Camps & Baes 2015, 2020; Reissl et al. 2016; Juvela 2019). Most codes that participated in the TRUST 3D dust radiative transfer benchmark effort (Gordon et al. 2017) were based on the Monte Carlo method. For an overview of Monte Carlo radiative transfer we refer to Whitney (2011), Steinacker et al. (2013), and Noebauer & Sim (2019).

In the following two subsections we discuss whether the new TTU2 phase function is easily integrated into the spherical harmonics and Monte Carlo radiative transfer frameworks. We particularly investigate whether the transition from the HG to the TTU2 phase function involves a significant increase in complexity or overhead.

4.4.1. Spherical harmonics radiative transfer

As mentioned above, the bottom-line of the spherical harmonics method is that the radiation field is expanded in a series of Legendre polynomials. The scattering phase function enters the spherical harmonics methods by means of its expansion in Legendre polynomials,

(15)

with the coefficients σℓ given by

(16)

The set of coefficients {σℓ}ℓ ≤ M, with M the chosen order of the Legendre polynomial series truncation, fully characterises the scattering phase function in this method.

The HG phase function is extremely well suited for the spherical harmonics method since its expansion in Legendre polynomials is very simple (Roberge 1983):

(17)

If we want to replace this phase function with the TTU2 phase function, we have to calculate the coefficients σℓ by inserting the expression (11) in the integral (16). The resulting integral can be evaluated analytically for each value of ℓ. For the first two coefficients we obviously find

(18)

(19)

For higher values of ℓ, the resulting expressions are similar to expression (14), but with higher-order polynomials in each of the terms. While this is clearly less attractive than the simple expression for the HG phase function, this is not a major source of overhead or increased complexity: the calculation of the Legendre coefficients σℓ can be performed numerically at the beginning of the actual radiative transfer calculations, or even beforehand.

An alternative, exploratory, option could be to use the fact that the U2 phase function, and hence also the TTU2 phase function, has a simple expansion in terms of ultraspherical polynomials of order 2 (see Appendix A),

(20)

The conventional spherical harmonics method builds on an expansion of the radiation field in terms of Legendre polynomials, but this is not the only option. Another option is an expansion in terms of Chebyshev polynomials of the first or second kind, which might result in better approximate solutions (e.g. Conkie 1959; Milgram 1991; Yaşa et al. 2006; Anlı et al. 2006; Öztürk et al. 2007; Bülbül et al. 2008). Legendre and Chebyshev polynomials are special cases of the family of ultraspherical polynomials, and it is in principle also possible to expand the radiation field in ultraspherical polynomials of order 2. Some efforts towards a general expansion in ultraspherical polynomials have been undertaken in the context of neutron and heat transport (Yılmazer & Kocar 2008; Yılmazer & Tombakoglu 2008). It remains to be seen how complex the extension towards radiative transfer with a general anisotropic phase functions will turn out to be; this is clearly beyond the scope of this paper.

4.4.2. Monte Carlo radiative transfer

In Monte Carlo radiative transfer, a new scattering angle θ has to be randomly generated from the scattering phase function after every single scattering event. As every simulation typically requires millions or billions of scattering events, it is crucial that this operation is efficient and accurate. The standard way to do this is by means of inversion sampling (Devroye 2013): one chooses a random deviate X and solves the equation X = P(θ), where P(θ) is the cumulative distribution function of the phase function

(21)

For the HG phase function, the cumulative distribution function can be calculated analytically and the resulting expression can be inverted exactly (Witt 1977),

(22)

This makes the inclusion of the HG phase function in the Monte Carlo framework very attractive.

For the U2 phase function we find that the cumulative distribution function of scattering angles can also be expressed analytically,

(23)

This expression can also be inverted analytically, resulting in the following recipe to generate a random scattering angle for the U2 phase function,

(24)

So the random generation of scattering angles can also be done in a very straightforward way for the U2 phase function.

For the TTU2 phase function, we obviously get for the cumulative distribution function a linear combination of two terms of the form (23), but, unfortunately, this expression cannot be inverted analytically anymore. This implies that we need other ways to generate random scattering angles. One option could be to use a numerical tabulation of the scattering phase function, as adopted by Witt (1977) for his Monte Carlo radiative transfer calculations using the TTHG phase function.

A more elegant, efficient and accurate approach consists of using the principle of decomposition in random number generation (Devroye 2013; Baes & Camps 2015). As a first step, we generate a first random uniform deviate U and determine which of the two U2 components will be selected. The importance of both components is determined by their relative weight: if 0 < U < γ the scattering angle will be generated from the first component, if γ < U < 1 it is the second component. Subsequently, we generate a second uniform deviate X and determine the scattering angle according to Eq. (24), where we replace g by either g1 or g2, depending on the component selected in the first step. In fact, the procedure can be made even more efficient by reusing the first random deviate in the second step. Rather than generating a second uniform deviate for the second step, we can also set

(25)

This procedure avoids the generation of a second uniform deviate and makes the algorithm to generate random scattering angles from the three-parameter TTU2 phase function essentially as efficient as for the HG phase function.

The accuracy of this algorithm is shown in Fig. 7. The thick orange line in this figure corresponds to the exact TTU2 phase function with arbitrarily chosen parameters g1 = 0.6, g2 = −0.4 and γ = 0.3. The grey line corresponds to the phase function reconstructed from the histogram of 104 randomly generated scattering angles, based on bin widths of 2 deg. The phase function is accurately reproduced, with some Poisson noise, especially around the forward and backward scattering peaks. The black line is similar to the grey line, but now corresponds to 106 random scattering angles. The Poisson noise is now reduced significantly.

thumbnail Fig. 7.

Demonstration of the method to extract random scattering angles from the TTU2 phase function. The solid orange line corresponds to the TTU2 phase function with g1 = 0.6, g2 = −0.4 and γ = 0.3. The grey and black lines are reconstructions of the phase function based on randomly generated scattering angles.

In summary, we find that the three-parameter TTU2 function can easily and efficiently be included in the Monte Carlo framework.

5. Summary

The goal of this paper was to search for alternative analytical scattering phase functions to describe the scattering properties of dust grains in the interstellar medium of galaxies. We specifically aimed at finding a balance between realism and complexity: on the one hand the scattering phase function should be flexible enough to provide an accurate fit to the actual scattering properties of dust grains over a wide wavelength range, and on the other hand it should be simple enough to be easy to handle, especially in the context of radiative transfer calculations.

We considered various analytical phase functions available in the literature and fitted them to ‘empirical’ data corresponding to the BARE-GR-S model by presented by Zubko et al. (2004). We considered the UV–NIR wavelength range, corresponding to wavelengths between 0.1 and 5 μm. Our findings are the following:

  • The standard one-parameter HG phase function provides a good fit to the empirical data at optical wavelengths but fails at describing the forward and backward scattering peak at UV and NIR wavelengths. Relative differences between the empirical data and the best-fitting HG phase function reach 30% in the NIR and 50% in the UV.

  • The two-parameter Draine phase function introduced by Draine (2003b) was designed as a generalisation of both the HG phase function and the Rayleigh phase function. It provides an excellent fit at NIR wavelengths but fails to alleviate the discrepancies in the UV.

  • The TTRM phase function was recently advocated by Wang et al. (2019) and Harmel et al. (2021) in the context of light scattering in nanoscale materials and aquatic media, respectively. We find that it provides an excellent fit to the empirical BARE-GR-S phase function data over the entire UV-NIR range. We do, however, find degeneracies between its five free parameters, suggesting that a five-parameter phase function is overkill for describing the scattering phase function.

  • We propose a simpler phase function, the two-term ultraspherical-2 (TTU2) phase function, that depends on only three free parameters. We show that it provides almost equally good fits to the empirical data as the TTRM phase function over the entire wavelength range considered, with relative differences almost always below 5%.

The new TTU2 phase function is characterised by three free parameters with a simple physical interpretation: two parameters characterise the shape of two components of the phase function, a forward and a backward scattering one, and the third parameter sets the relative weight of both components. We have investigated how these parameters change with wavelength, both for the BARE-GR-S model and an alternative dust model. In spite of the differences between the two dust models, the results are compatible:

  • The forward scattering component of the phase function systematically becomes less anisotropic as the wavelength increases.

  • The backward scattering component is isotropic in the UV and moderately anisotropic at optical and NIR wavelengths. The anisotropy peaks in the optical regime.

  • The relative strength of the backward scattering component peak is of the order of 20% at UV and blue wavelengths, has a minimal value of about 3% in the red (λ ∼ 0.8 μm), and gradually increases towards 50% at longer wavelengths.

Finally, we have investigated whether the integration of the TTU2 phase function in radiative transfer calculations implies a significant overhead or increased complexity. We have specifically focused on the spherical harmonics method and the Monte Carlo method: the former is one of the most efficient methods for 1D radiative transfer, and the latter is the most popular option for 3D radiative transfer. While the standard HG phase function is easily and almost naturally integrated into both approaches, we find that the shift from the HG to the TTU2 phase function does not imply a major source of overhead or increased complexity.

Our overall conclusion is that the TTU2 scattering phase function that we have proposed in this paper provides a nice balance between being simple enough to be easily used and realistic enough to accurately describe scattering by dust grains. We advocate its employment in astrophysical applications, in particular in radiative transfer calculations. An investigation of the effects of the choice of the phase function in radiative transfer calculations in different optical depth regimes (see e.g. Baes & Dejonghe 2001b) will be considered in future work.

2

As discussed in the first paragraph of Sect. 3, we have determined g at each wavelength by minimising the relative error hrel(g). An alternative option could have been choosing g equal to the numerically determined anisotropy parameter of the empirical phase function. Advantages of our approach are that it can be applied consistently to all analytical phase functions, and that, by definition, it selects the parameters corresponding the lowest relative error. The resulting g values of both approaches are very similar.

3

We consider the effective wavelengths of the GALEX FUV and NUV bands, the five optical SDSS bands, the 2MASS J, H, and Ks bands, and the WISE W1 and W2 bands.

Acknowledgments

The authors thank the anonymous referee for insightful and detailed comments that improved the content and presentation of this work. This paper is dedicated to the memory of Jonathan I. Davies (1955–2021). He was an inspiring researcher, mentor, leader, colleague and friend, and will be deeply missed.

References

  1. Anlı, F., Yaşa, F., Güngör, S., & Öztürk, H. 2006, J. Quant. Spectr. Rad. Transf., 101, 129 [CrossRef] [Google Scholar]
  2. Arridge, S. R. 1999, Inverse Probl., 15, R41 [NASA ADS] [CrossRef] [Google Scholar]
  3. Baes, M., & Camps, P. 2015, Astron. Comput., 12, 33 [NASA ADS] [CrossRef] [Google Scholar]
  4. Baes, M., & Dejonghe, H. 2001a, MNRAS, 326, 722 [NASA ADS] [CrossRef] [Google Scholar]
  5. Baes, M., & Dejonghe, H. 2001b, MNRAS, 326, 733 [NASA ADS] [CrossRef] [Google Scholar]
  6. Bailey, J., & Kedziora-Chudczer, L. 2012, MNRAS, 419, 1913 [NASA ADS] [CrossRef] [Google Scholar]
  7. Bianchi, S., De Vis, P., Viaene, S., et al. 2018, A&A, 620, A112 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  8. Boucher, O. 1998, J. Atmos. Sci., 55, 128 [NASA ADS] [CrossRef] [Google Scholar]
  9. Bron, E., Le Bourlot, J., & Le Petit, F. 2014, A&A, 569, A100 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  10. Bruzual, G. A., Magris, G., & Calvet, N. 1988, ApJ, 333, 673 [NASA ADS] [CrossRef] [Google Scholar]
  11. Bülbül, A., Ulutas, M., & Anlı, F. 2008, Kerntechnik, 73, 61 [CrossRef] [Google Scholar]
  12. Byun, Y. I., Freeman, K. C., & Kylafis, N. D. 1994, ApJ, 432, 114 [NASA ADS] [CrossRef] [Google Scholar]
  13. Calabro, K. W., & Cassarly, W. 2015, in Biomedical Applications of Light Scattering IX, eds. A. Wax, & V. Backman (SPIE), Int. Soc. Opt. Photon., 9333, 22 [NASA ADS] [Google Scholar]
  14. Calzetti, D., Bohlin, R. C., Gordon, K. D., Witt, A. N., & Bianchi, L. 1995, ApJ, 446, L97 [NASA ADS] [CrossRef] [Google Scholar]
  15. Camps, P., & Baes, M. 2015, Astron. Comput., 9, 20 [Google Scholar]
  16. Camps, P., & Baes, M. 2018, ApJ, 861, 80 [NASA ADS] [CrossRef] [Google Scholar]
  17. Camps, P., & Baes, M. 2020, Astron. Comput., 31, 100381 [NASA ADS] [CrossRef] [Google Scholar]
  18. Camps, P., Misselt, K., Bianchi, S., et al. 2015, A&A, 580, A87 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  19. Chandrasekhar, S. 1950, Radiative Transfer (London: Oxford University Press) [Google Scholar]
  20. Chen-Chen, H., Pérez-Hoyos, S., & Sánchez-Lavega, A. 2021, Icarus, 354, 114021 [NASA ADS] [CrossRef] [Google Scholar]
  21. Conkie, W. R. 1959, Nucl. Sci. Eng., 6, 260 [CrossRef] [Google Scholar]
  22. Danielson, R. E., Moore, D. R., & van de Hulst, H. C. 1969, J. Atmos. Sci., 26, 1078 [NASA ADS] [CrossRef] [Google Scholar]
  23. De Geyter, G., Baes, M., Camps, P., et al. 2014, MNRAS, 441, 869 [CrossRef] [Google Scholar]
  24. De Micheli, E., & Viano, G. A. 2013, J. Comput. Phys., 239, 112 [NASA ADS] [CrossRef] [Google Scholar]
  25. Devroye, L. 2013, Non-Uniform Random Variate Generation (New York: Springer-Verlag) [Google Scholar]
  26. di Bartolomeo, A., Barbaro, G., & Perinotto, M. 1995, MNRAS, 277, 1279 [NASA ADS] [CrossRef] [Google Scholar]
  27. Doman, B. G. S. 2015, The Classical Orthogonal Polynomials (Singapore: World Scientific) [CrossRef] [Google Scholar]
  28. Draine, B. T. 1978, ApJS, 36, 595 [Google Scholar]
  29. Draine, B. T. 2003a, ARA&A, 41, 241 [NASA ADS] [CrossRef] [Google Scholar]
  30. Draine, B. T. 2003b, ApJ, 598, 1017 [NASA ADS] [CrossRef] [Google Scholar]
  31. Ferrara, A., Bianchi, S., Cimatti, A., & Giovanardi, C. 1999, ApJS, 123, 437 [CrossRef] [Google Scholar]
  32. Flannery, B. P., Roberge, W., & Rybicki, G. B. 1980, ApJ, 236, 598 [NASA ADS] [CrossRef] [Google Scholar]
  33. Gibson, S. J., & Nordsieck, K. H. 2003, ApJ, 589, 362 [CrossRef] [Google Scholar]
  34. Gordon, K. D., Misselt, K. A., Witt, A. N., & Clayton, G. C. 2001, ApJ, 551, 269 [NASA ADS] [CrossRef] [Google Scholar]
  35. Gordon, K. D., Baes, M., Bianchi, S., et al. 2017, A&A, 603, A114 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  36. Gould, R. J., & Salpeter, E. E. 1963, ApJ, 138, 393 [NASA ADS] [CrossRef] [Google Scholar]
  37. Hakim, A. H., & McCormick, N. J. 2003, Appl. Opt., 42, 931 [NASA ADS] [CrossRef] [Google Scholar]
  38. Hammer, M., Schweitzer, D., Michel, B., Thamm, E., & Kolb, A. 1998, Appl. Opt., 37, 7410 [NASA ADS] [CrossRef] [Google Scholar]
  39. Harmel, T., Agagliate, J., Hieronymi, M., & Gernez, P. 2021, Opt. Lett., 46, 1860 [NASA ADS] [CrossRef] [Google Scholar]
  40. Hensley, B. S., & Draine, B. T. 2021, ApJ, 906, 73 [NASA ADS] [CrossRef] [Google Scholar]
  41. Henyey, L. G., & Greenstein, J. L. 1941, ApJ, 93, 70 [Google Scholar]
  42. Hollenbach, D., & Salpeter, E. E. 1971, ApJ, 163, 155 [NASA ADS] [CrossRef] [Google Scholar]
  43. Hussaini, M. Y., & Zang, T. A. 1987, Annu. Rev. Fluid Mech., 19, 339 [NASA ADS] [CrossRef] [Google Scholar]
  44. Irvine, W. M. 1965, ApJ, 142, 1563 [NASA ADS] [CrossRef] [Google Scholar]
  45. Irvine, W. M. 1968, ApJ, 152, 823 [NASA ADS] [CrossRef] [Google Scholar]
  46. Juvela, M. 2019, A&A, 622, A79 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  47. Karlsson, H., Fredriksson, I., Larsson, M., & Strömberg, T. 2012, Opt. Expr., 20, 12233 [NASA ADS] [CrossRef] [Google Scholar]
  48. Kattawar, G. W. 1975, J. Quant. Spectr. Rad. Transf., 15, 839 [NASA ADS] [CrossRef] [Google Scholar]
  49. Kim, D. S., Kim, T., & Rim, S. H. 2012, Adv. Differ. Equ., 2012, 219 [CrossRef] [Google Scholar]
  50. Krieger, A., & Wolf, S. 2021, A&A, 645, A143 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  51. Light, B., Maykut, G. A., & Grenfell, T. C. 2003, J. Geophys. Res.: Oceans, 108 [CrossRef] [Google Scholar]
  52. Lillie, C. F., & Witt, A. N. 1976, ApJ, 208, 64 [NASA ADS] [CrossRef] [Google Scholar]
  53. Lindbergh, T., Fredriksson, I., Larsson, M., & Strömberg, T. 2009, Opt. Expr., 17, 1610 [NASA ADS] [CrossRef] [Google Scholar]
  54. Milgram, M. 1991, Ann. Nucl. Energy, 18, 87 [CrossRef] [Google Scholar]
  55. Nersesian, A., Xilouris, E. M., Bianchi, S., et al. 2019, A&A, 624, A80 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  56. Noebauer, U. M., & Sim, S. A. 2019, Liv. Rev. Comput. Astrophys., 5, 1 [CrossRef] [Google Scholar]
  57. Öztürk, H., Anlı, F., & Güngör, S. 2007, J. Quant. Spectr. Rad. Transf., 105, 211 [CrossRef] [Google Scholar]
  58. Popescu, C. C., & Tuffs, R. J. 2002, MNRAS, 335, L41 [NASA ADS] [CrossRef] [Google Scholar]
  59. Reissl, S., Wolf, S., & Brauer, R. 2016, A&A, 593, A87 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  60. Reynolds, L. M., & McCormick, N. J. 1980, J. Opt. Soc. Am., 70, 1206 [NASA ADS] [CrossRef] [Google Scholar]
  61. Roberge, W. G. 1983, ApJ, 275, 292 [CrossRef] [Google Scholar]
  62. Robitaille, T. P. 2011, A&A, 536, A79 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  63. Roggan, A., Friebel, M., Doerschel, K., Hahn, A., & Mueller, G. J. 1999, J. Biomed. Opt., 4, 36 [NASA ADS] [CrossRef] [Google Scholar]
  64. Sharma, S. K., & Roy, A. K. 2008, ApJS, 177, 546 [NASA ADS] [CrossRef] [Google Scholar]
  65. Sharma, S. K., Roy, A. K., & Somerford, D. J. 1998, J. Quant. Spectr. Rad. Transf., 60, 1001 [NASA ADS] [CrossRef] [Google Scholar]
  66. Siebenmorgen, R., Voshchinnikov, N. V., & Bagnulo, S. 2014, A&A, 561, A82 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  67. Steinacker, J., Baes, M., & Gordon, K. D. 2013, ARA&A, 51, 63 [CrossRef] [Google Scholar]
  68. Stramski, D., Boss, E., Bogucki, D., & Voss, K. J. 2004, Progr. Oceanogr., 61, 27 [NASA ADS] [CrossRef] [Google Scholar]
  69. Szegö, G. 1950, Proc. Am. Math. Soc., 1, 731 [CrossRef] [Google Scholar]
  70. Tuchin, V. V. 1993, in Static and Dynamic Light Scattering in Medicine and Biology, eds. R. J. Nossal, R. Pecora, & A. V. Priezzhev (SPIE), Int. Soc. Opt. Photon., 1884, 234 [CrossRef] [Google Scholar]
  71. Viaene, S., Baes, M., Bendo, G., et al. 2016, A&A, 586, A13 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  72. Wang, J., Xu, C., Nilsson, A. M., Fernandes, D. L. A., & Niklasson, G. A. 2019, Nanoscale, 11, 7404 [CrossRef] [Google Scholar]
  73. Weingartner, J. C., & Draine, B. T. 2001, ApJ, 548, 296 [Google Scholar]
  74. Whitney, B. A. 2011, Bull. Astron. Soc. India, 39, 101 [NASA ADS] [Google Scholar]
  75. Witt, A. N. 1977, ApJS, 35, 1 [NASA ADS] [CrossRef] [Google Scholar]
  76. Witt, A. N., & Gordon, K. D. 1996, ApJ, 463, 681 [NASA ADS] [CrossRef] [Google Scholar]
  77. Witt, A. N., Thronson, H. A., Jr., & Capuano, J. M., Jr. 1992a, ApJ, 393, 611 [NASA ADS] [CrossRef] [Google Scholar]
  78. Witt, A. N., Petersohn, J. K., Bohlin, R. C., et al. 1992b, ApJ, 395, L5 [NASA ADS] [CrossRef] [Google Scholar]
  79. Xilouris, E. M., Byun, Y. I., Kylafis, N. D., Paleologou, E. V., & Papamastorakis, J. 1999, A&A, 344, 868 [NASA ADS] [Google Scholar]
  80. Yaroslavsky, A. N., Yaroslavsky, I. V., & Goldbach, T. 1997, in Optical Diagnostics of Biological Fluids and Advanced Techniques in Analytical Cytology, eds. R. C. Leif, A. V. Priezzhev, T. Asakura, & R. C. Leif (SPIE), Int. Soc. Opt. Photon., 2982, 324 [NASA ADS] [CrossRef] [Google Scholar]
  81. Yaşa, F., Anlı, F., & Güngör, S. 2006, J. Quant. Spectr. Rad. Transf., 97, 51 [CrossRef] [Google Scholar]
  82. Yılmazer, A., & Kocar, C. 2008, Int. J. Therm. Sci., 47, 112 [CrossRef] [Google Scholar]
  83. Yılmazer, A., & Tombakoglu, M. 2008, Kerntechnik, 73, 260 [CrossRef] [Google Scholar]
  84. Ysard, N., Jones, A. P., Demyk, K., Boutéraon, T., & Koehler, M. 2018, A&A, 617, A124 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  85. Zhang, F., & Li, J. 2016, J. Quant. Spectr. Rad. Transf., 184, 40 [NASA ADS] [CrossRef] [Google Scholar]
  86. Zhang, F., Liu, K., Yang, Q., & Zhao, J.-Q. 2017, Adv. Meteorol., 2017, 1835169 [Google Scholar]
  87. Zubko, V., Dwek, E., & Arendt, R. G. 2004, ApJS, 152, 211 [NASA ADS] [CrossRef] [Google Scholar]

Appendix A: Ultraspherical polynomials

The ultraspherical or Gegenbauer polynomials are orthogonal polynomials on the interval [ − 1, 1] with respect to the weight function w(x) = (1 − x2)λ − 1/2. Their explicit form can be obtained in different ways (Szegö 1950; Kim et al. 2012; Doman 2015). One of the most convenient ways is by means of the recurrence relation

(A.1)

(A.2)

(A.3)

Particularly important special cases of ultraspherical polynomials are obtained for particular values of the order λ: gives the Legendre polynomials Pn(x), λ = 1 gives the Chebyshev polynomials of the second kind Un(x), and λ = 0 gives, in a suitable limiting form, the Chebyshev polynomials of the first kind Tn(x).

The generating function for the ultraspherical polynomials is a simple rational function (Szegö 1950),

(A.4)

Comparing this expression to the U2 phase function (12) we immediately find

(A.5)

The lowest-order ultraspherical polynomials of order 2 are

(A.6)

(A.7)

(A.8)

(A.9)

(A.10)

Due to their importance in spectral methods for the numerical solution of differential equations and computational fluid dynamics (Hussaini & Zang 1987), efficient methods for the calculation of the coefficients of ultraspherical polynomial expansions have been developed (De Micheli & Viano 2013).

Appendix B: TTU2 phase function parameters

The data files BARE-GS-S-TTU2.txt and WD01-MW-TTU2.txt contain the model parameters of the best-fitting TTU2 phase function to the empirical phase function for the BARE-GR-S and WD01-MW dust models, respectively. The former file contains fitting parameters for 1201 different wavelengths between 1 nm and 10 mm, the latter file for 36 different wavelengths between 60 nm and 10.2 μm.

All Tables

Table 1.

Model parameters of the best-fitting TTU2 phase functions for the BARE-GR-S dust model at the central wavelengths of a number of standard broadband filters.

In the text

All Figures

thumbnail Fig. 1.

Empirical scattering phase functions (grey dots) of the Zubko et al. (2004) BARE-GR-S dust model at the effective wavelengths of eight common UV–NIR bands. For the sake of clarity we only plotted the data points at a resolution of 5 deg, whereas the BARE-GR-S dataset contains phase function data at 1 deg resolution. The solid lines in each upper panel correspond to analytical phase function fits to these empirical data: the HG phase function (red, Sect. 3.1), the Draine phase function (purple, Sect. 3.2), the TTRM phase function (cyan, Sect. 3.4), and the new TTU2 phase function (gold, Sect. 4). Bottom panels: relative difference between the analytical model and the empirical data.
In the text
thumbnail Fig. 2.

Scattering phase function for two representative wavelengths (0.23 and 1.65 μm), shown in polar format. The meanings of the lines and symbols are the same as in Fig. 1.
In the text
thumbnail Fig. 3.

Model parameters of the best-fitting TTRM phase functions as a function of wavelength. The dots indicate the effective wavelengths of the GALEX, SDSS, 2MASS, and WISE broadband filters.
In the text
thumbnail Fig. 4.

Relative error ⟨hrel⟩, averaged over the standard UV, optical, and NIR broadband filters, for the TTRM phase function as a function of the fixed α = α1 = α2.
In the text
thumbnail Fig. 5.

Comparison of the relative error hrel as a function of wavelength for the different phase functions considered in this paper. The dots indicate the effective wavelengths of the GALEX, SDSS, 2MASS, and WISE broadband filters.
In the text
thumbnail Fig. 6.

Relative error and model parameters of the best-fitting TTU2 phase functions as a function of wavelength. The golden lines correspond to our default BARE-GR-S dust model, and the green dots correspond to the alternative WD01-MW dust model.
In the text
thumbnail Fig. 7.

Demonstration of the method to extract random scattering angles from the TTU2 phase function. The solid orange line corresponds to the TTU2 phase function with g1 = 0.6, g2 = −0.4 and γ = 0.3. The grey and black lines are reconstructions of the phase function based on randomly generated scattering angles.
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.