Contents

• Far-infrared to millimeter astrophysical dust emission
• 1 Introduction
• 2 Basic knowledge of the FIR/mm dust emission
• 3 Temperature and spectral variations of the optical properties of dust
• 4 A model to explain temperature and spectral variation of dust optical properties
• 5 Model of the dust FIR/mm emission
• 6 Conclusion
• References
• 7 Online Material
• Appendix A: TLS optical density of states
• Appendix B: Evaluation of a simplified equation for tunneling relaxation calculation
• Appendix C: Calculation of the integral for the tunneling relaxation spectrum
• Appendix D: Calculations of the integral for hopping relaxation

A&A 468, 171-188 (2007)
DOI: 10.1051/0004-6361:20065771

Far-infrared to millimeter astrophysical dust emission

I. A model based on physical properties of amorphous solids

C. Meny¹ - V. Gromov^1,2 - N. Boudet¹ - J.-Ph. Bernard¹ - D. Paradis¹ - C. Nayral¹

1 - Centre d'Étude Spatiale des Rayonnements, CNRS, 9 avenue du Colonel Roche, 31028 Toulouse, France
2 - Space Research Institute, RAS, 84/32 Profsoyuznaya, 117810 Moscow, Russia

Received 7 June 2006 / Accepted 11 December 2006

Abstract
Aims. We propose a new description of astronomical dust emission in the spectral region from the far-infrared to millimeter wavelengths.
Methods. Unlike previous classical models, this description explicitly incorporates the effect of the disordered internal structure of amorphous dust grains. Our model is based on results from solid state physics used to interpret laboratory data. The model takes into account the effect of absorption by disordered charge distribution, as well as the effect of absorption by localized two level systems.
Results. We review constraints on the various free parameters of the model from theory and laboratory experimental data. We show that, for realistic values of the free parameters, the shape of the emission spectrum will exhibit very broad structures whose shape will change in a non trivial way with the temperature of dust grains. The spectral shape also depends upon the parameters describing the internal structure of the grains. This opens new perspectives for identifying the nature of astronomical dust from the observed shape of the FIR/mm emission spectrum. A companion paper will provide an explicit comparison of the model with astronomical data.

Key words: dust, extinction - diffuse radiation - infrared: ISM - radio continuum: ISM - scattering - radiation mechanisms: thermal

  
1 Introduction

It is now well established that the spectral energy distribution (SED) from the inter-stellar medium (ISM) emission from our Galaxy and most external galaxies is dominated by thermal emission from dust grains that spans over almost 3 orders of magnitude in wavelengths, from the near infra-red (NIR) to the millimeter wavelength range. Several dust models have been proposed to explain the observations. Recent dust models share many common characteristics (e.g. Mathis et al. 1977; Draine & Lee 1984; Draine & Anderson 1985; Weiland et al. 1986; Désert et al. 1990; Li & Greenberg 1997; Dwek et al. 1997; Draine & Li 2001). Most of them require a wide range in dust grain sizes spanning almost 3 orders of magnitude, from sub-nanometer dimensions to fractions of a micron. They also require a minimum of three components with different chemical compositions in order to explain the ultraviolet and optical extinction, along with infrared emission.

The smallest component is needed to explain the "aromatic'' emission features at 3.3, 6.2, 7.7, 8.6 and $11.3~{\rm\mu m}$. It is now routinely assigned to large aromatic molecules or polycyclic aromatic hydrocarbons (PAH) (Leger & Puget 1984; Allamandola et al. 1985; Désert et al. 1990). These PAH are transiently heated after absorbing single UV and far-UV photons, and cool down emitting NIR photons in vibrational transitions representative of their aromatic structure.

A second component, composed of very small grains (VSG) is required to explain the continuum emission in the mid infra-red (MIR). Emission in this range still requires significant temperature fluctuations of the grains, which necessitates grain sizes in the nanometer range under conditions prevailing in the inter-stellar medium (ISM). The VSGs could be composed of carbonaceous material whose absorption could also explain the $2200~\AA$ UV bump in the extinction curve (see for instance Désert et al. 1990).

A third component, composed of big grains (BG) with sizes from 10-20 nm to about 100 nm, is necessary to account for the long wavelength emission (in particular as observed in the IRAS $100~{\rm\mu m}$band and above). The BG component dominates the total dust mass and the absorption in the visible and the NIR. Observations of the "silicate'' absorption feature near $10~{\rm\mu m}$ indicate that their mass is dominated by amorphous silicates (see Kemper et al. 2004). This sort of amorphous structure for these dust grains is expected from the interaction with cosmic rays that should alter the nature of the solid (Brucato et al. 2004; Jäger et al. 2003), but also from the formation processes of these BGs by aggregation of sub-sized particles.

Recent studies attempting to explain the submillimeter observations have revealed the need for either the existence of an additional very massive component of very cold dust (e.g. Reach et al. 1995; Finkbeiner et al. 1999; Galliano et al. 2003) or for substantial modifications of dust properties in the submillimeter with respect to predictions of standard dust models (e.g. Ristorcelli et al. 1998; Bernard et al. 1999; Stepnik et al. 2003a; Dupac et al. 2001; Dupac et al. 2002, 2003a).

Similarly, laboratory spectroscopic measurements of interstellar-grain analogs reveal that noticeable variations in their optical properties in the far infra-red (FIR) and millimeter (mm) wavelength range can occur (Agladze et al. 1994; Agladze et al. 1996; Mennella et al. 1998; Boudet et al. 2005). They points out the possible role of two phonon-difference processes, disorder-induced one-phonon processes, resonant or relaxation processes in the presence of two-level systems(TLS) in amorphous dust (see Sect. 4.3). The influence of such TLS on interstellar grain absorption and emission properties was first proposed to the astronomy community by Agladze et al. (1994). More recently, a preliminary investigation by Boudet et al. (2002) has shown that optical resonant and relaxation transitions in a distribution of TLS could qualitatively explain the anticorrelation between the temperature and the spectral index of dust as observed by the PRONAOS balloon experiment (e.g. Dupac et al. 2003a).

A precise modeling of the long wavelength dust emission is important in order to accurately subtract foreground emission in the cosmological background anisotropy measurements, especially for future missions for space CMB measurements (Lamarre et al. 2003) and for surveys of foreground compact sources (see Gromov et al. 2002). A number of comprehensive foreground analysis papers (Bouchet & Gispert 1999; Tegmark et al. 2000; Bennett et al. 2003; Barreiro et al. 2006; Naselsky et al. 2005, and other) did not pay attention to dust emissivity model accuracy. In addition, precise modeling of the FIR/mm dust emission is also important to derive reliable estimates of the dust mass and to trace the structure and density of pre-stellar cold cores in molecular clouds. Indeed the dust governs the cloud opacity and influences the size of the cloud fragments. In addition, the efficiency of star formation should be related to the dust properties, which are expected to vary between the diffuse medium and the denser cold cores. It is thus important to know how some physical or chemical properties of dust can influence their FIR/mm emission.

For models of the submillimeter electromagnetic dust emission, considering the intrinsic mechanical vibrations of the grain structure is natural since they fall in the right frequency range (typically in the frequency range $\nu\sim10^{11}{-}10^{13}$ Hz, i.e. $\lambda\rm\sim30~{\rm\mu m} {-}3$ mm). Observations in this spectral region potentially provide new tools for investigating the internal structure of dust grain material. The internal mechanics of grains is also important for chemistry in the ISM. The dust vibrations supply an energy sink for newly-formed molecules that are usually formed in unstable excited states. Understanding the grain structure and intrinsic movements is necessary for considerations of grain collisions, destruction, and agglomeration, including preplanetary bodies formation in circumstellar disks.

In Sect. 2, we first recall some basics regarding the FIR/mm dust emission and the semi-classical model of light interaction with matter. In Sect. 3, we gather evidence for spectral variations in astronomical and laboratory measurements. Then, we present a model in Sect. 4, based on the intrinsic properties of amorphous materials, taking absorption by disordered charge distribution into account as well as the effect of absorption by a distribution of TLS. In Sect. 5, we propose a new model for the submillimeter and millimeter grain absorption and discuss the range of plausible values for the free parameters based on existing laboratory data. In Sect. 5.4, we discuss the implications of the model, in particular regarding the expected variations of the FIR/mm emission spectra with the dust temperature. Section 6 is devoted to conclusions. The determination of the model parameters applicable to astronomical data will be addressed in a companion paper.

2 Basic knowledge of the FIR/mm dust emission 

Basic knowledge of the FIR/mm dust emission is based on the well known semi-classical models of the interaction of light with matter. We first show that in the FIR/mm wavelength range, the intensity of interstellar thermal emission from spherical dust grains can be related in a simple way to the optical constants which characterize the dust material. Then we recapitulate these classical models for the interaction between an electromagnetic wave and the bound or free charges in the material, respectively the Lorentz and Drude models. Both models, applied to the thermal dust emission, lead to an asymptotic behavior in the long-wavelength range described by a simple power law whose exponent, named spectral index, equals 2. Thus the FIR/mm interstellar dust emission is up to now generally expected to be fully characterized by a constant and temperature-independent spectral index around 2.

2.1 Interstellar dust emission and extinction

The intensity of thermal emission from interstellar dust at temperature $T_{\rm d}$ is

$$I_{\omega}(\omega, T)=\epsilon_{\rm e}(\omega)\cdot B_{\omega}(\omega,T_{\rm d})$$ (1)

where $I_{\omega}$ is the energy flux density per unit area, angular frequency, and solid angle (or specific intensity), $\epsilon_{\rm e}$ the dust emissivity, and $B_{\omega}$ the Planck function at angular frequency $\omega$.

According to the Kirchhoff law, the emissivity $\epsilon_{\rm e}$ is equal to the absorptivity

$$\epsilon_{\rm e}(\omega)=1-{\rm e}^{-\tau(\omega)},$$ (2)

where the optical depth $\tau$ is related to the dust-mass column density on the line of sight $M_{\rm d}$ and the dust mass opacity $\kappa$ as

$$\tau(\omega)=\kappa(\omega)\cdot M_{\rm d}.$$ (3)

For spherical grains of radius a and density $\rho$, the dust opacity (effective area per mass) is given by

$$\kappa(\omega) = \frac{3}{4\rho} \frac{Q(\omega)}{a},$$ (4)

where the absorption efficiency $Q(\omega)=\sigma(\omega)/(\pi a^2)$is the ratio of the absorption cross section $\sigma(\omega)$ to the geometrical cross section of the grain $\pi a^2$. The grain equilibrium temperature $T_{\rm d}$ is deduced from the balance between the emitted and absorbed radiation from the inter-stellar radiation field (ISRF)

$$\int^{\infty}_0 B_{\omega}(\omega,T_{\rm d}) Q(\omega){\rm ... ... \int^{\infty}_0 I^{\rm ISRF}_{\omega} Q(\omega){\rm d}\omega.$$ (5)

The Mie theory, using the Maxwell equations of the electromagnetic theory, leads to an exact solution for the absorption and scattering processes by a homogeneous spherical particle of radius whose material is characterized by its complex dielectric constant $\epsilon(\omega) = \epsilon'(\omega)+{\rm i}\epsilon''(\omega)$. In that case, the absorption coefficient $Q(\omega)$ is written as an infinite series of terms related to the complex dielectric constant. For small particles compared to the wavelength such as $\vert\sqrt{\epsilon(\omega)}\vert\cdot 2\pi a/\lambda \ll 1$ (where $\lambda=2\pi c/\omega$ is the wavelength and c the speed of light in vacuum), the absorption coefficient can be approximated by

$$Q(\omega) = \frac{8\pi a}{\lambda} \cdot {\rm Im} \left[\frac{\epsilon(\omega)-1}{\epsilon(\omega)+2}\right]\cdot$$ (6)

In the FIR/mm wavelength range, the condition $\epsilon''\ll\epsilon'$ is generally satisfied for non-conducting materials like astronomical silicate-based particles, and in this case the absorption coefficient reduces to

$$\frac{Q(\omega)}{a}=\frac{\omega}{c} \frac{12}{(\epsilon'+2)^{2}}\epsilon''.$$ (7)

As a consequence, following Eq. (7) the dust opacity rewrites as

$$\kappa(\omega) = \frac{\omega}{c} \frac{9}{(\epsilon'+2)^2} \frac{\epsilon''}{\rho}\cdot$$ (8)

   
2.2 Assumption of a $\lambda$-independent emissivity spectral index

In the absence of specific knowledge about realistic dependencies for $\epsilon(\omega)$, a simple power law approximation for FIR/mm dust emission is often assumed:

$$\kappa(\omega) = \kappa(\omega_0) \left( \frac{\omega}{\omega_0} \right)^\beta,$$ (9)

where $\beta$ is referred to as the emissivity spectral index. The intensity of thermal emission from interstellar dust, assuming an optically thin medium, at a given wavelength $\lambda=2\pi/\omega$located in the FIR/mm range is

$$I_{\omega}(\omega, T)=\epsilon_{\rm e}(\lambda_0)\cdot \le... ...mega}{\omega_0}\right)^\beta\cdot B_{\omega}(\omega,T_{\rm d})$$ (10)

where $\epsilon_{\rm e}(\lambda_0)$ is the dust emissivity at a reference wavelength $\lambda_0=2\pi c/\omega_0$.

Simple semi-classical models of absorption, such as the Lorentz model for damping oscillators and the Drude model for free charge carriers (Sect. 2.4), provide an asymptotic value $\beta=2$ when $\omega \rightarrow 0$. This value of the spectral index was in satisfactory agreement with the earliest observations of the FIR/mm SED of interstellar dust emission, but not necessarily with the most recent ones.

The emissivity spectral index $\beta$ is equal to the slope of a plot of the dust opacity versus wavelength in logarithmic scale. In the general case, the power-law assumption is not valid, and the spectral index has to be considered as a wavelength and temperature-dependent parameter. Similarly, $\epsilon_{\rm e}(\lambda_0)$ in Eq. (10) should be considered a function of T.

   
2.3 Relations between absorption properties of bulk material and dust

The absorption coefficient defined as the optical depth in the bulk material per unit length is

$$\alpha =2{\rm Im}(k),$$ (11)

where k is the complex value of the wave number $k=\omega\sqrt{\epsilon}/c$ and $\epsilon$ the complex dielectric constant. Therefore for $\epsilon''\ll\epsilon'$,

$$\alpha= \frac{\omega\epsilon'' }{c\sqrt{\epsilon'}}\cdot$$ (12)

In dust and porous material measurements, a mass absorption coefficient $\alpha/\rho$ is also often used.

In the following analysis (Sect. 2.4 and other) we link the macroscopic value $\epsilon$ to the value of the local susceptibility $\chi_0=\chi_0'+{\rm i}\chi_0''$ defined on a microscopic scale. Because small astronomical particles cannot always be considered as a continuous medium, the electrodynamic equations of continuous media must be applied with care. So we explain the equations used in the following and discuss their range of application.

The equations of motion describes the response of charges to a local electric field $\vec{E}_0$, which permits a dipole moment per unit volume to be determined as $\vec{P}=\chi_0\vec{E}_0$. It should be taken into account that the local electric field $\vec{E}_0$ differs from an external field $\vec{E}$ due to an additional field produced by other parts of the dielectric. In the case of isotropic dielectric materials,

$$\vec{E}_0=\vec{E} +\frac{4\pi}{3}\vec{P}.$$

The dielectric constant $\epsilon$ links the external field $\vec{E}$ and the electric inductance $\vec{D}=\vec{E}+4\pi\vec{P}$ through $\vec{D}=\epsilon\vec{E}$. In this case, $\epsilon$ can be derived from the susceptibility $\chi_0$, using the so-called Clausius-Mossotti equation

$$\epsilon = 1 + \frac{4\pi\chi_0}{1-\frac{4}{3}\pi\chi_0}\cdot$$ (13)

When $\epsilon''\ll\epsilon'$, we have

$$\epsilon'' \approx \frac{(\epsilon'+2)^2}{9}\;4\pi\chi_0'',$$ (14)

which leads to

  

$$\alpha(\omega)= \frac{4\pi\omega}{c\sqrt{\epsilon'}}\frac{(\epsilon'+2)^2}{9}\;\chi_0''(\omega),$$ (15)

$$\kappa(\omega)=\frac{\omega}{c} \frac{4\pi\chi_0''}{\rho}\cdot$$ (16)

In Eq. (15), the term $(\epsilon'+2)^2/9$ is the local field correction factor, and the term $\sqrt{\epsilon'}$ is the refraction index.

Thus, the relation between the absorption coefficient $\alpha$deduced from transmission measurements in the laboratory and the absorption efficiency Q or opacity $\kappa$, for small spherical grains of radius a and density $\rho$ can be obtained from Eqs. (7) and (12):

$$\frac{Q(\omega)}{\alpha(\omega)}=\frac{12a\cdot\sqrt{\epsilon'}}{(\epsilon'+2)^2}$$ (17)

and, as a consequence,

$$\frac{\kappa(\omega)}{\alpha(\omega)}=\frac{9\cdot\sqrt{\epsilon'}}{\rho\cdot(\epsilon'+2)^2}\cdot$$ (18)

The condition $\epsilon''\ll\epsilon'$ is valid for dust materials in the FIR/mm wavelength range, implying that the real part of the dielectric constant $\epsilon'$ can be considered to the first order as a constant, leading to ${\rm d}\epsilon'/{\rm d}\omega=0$. This comes from the coupling between the real part $\epsilon'$ and the imaginary part $\epsilon''$ through the Kramers-Kronig relations discussed in Sect. 3.2. This implies that the absorption coefficient $\alpha$, the absorption efficiency Q and the opacity $\kappa$ are simply proportional in that spectral range. Thus the slopes of $\alpha(\omega)$, $Q(\omega)$, and $\kappa(\omega)$ are the same. As for the emissivity, an absorptivity spectral index can then be defined by the slope of the spectral plot of the absorption $\alpha(\omega)$ on a logarithmic scale. This absorptivity spectral index deduced from laboratory transmission measurements can be directly compared to the emissivity spectral index deduced from astrophysical measurements.

   
2.4 The semi-classical models

The well-known semi-classical model of light interaction with charge carriers in matter (electrons and ions) considers the motion of charged particles driven by the time-dependent electric field of a monochromatic plane wave, with magnetic forces being neglected compared with electrical forces.

The equations of motion of interacting charge carriers are

 

$$m_j\ddot{u}_j^{~\xi} + \sum_{j'\xi'}K_{jj'}^{\xi\xi'} u_{j'}^{~\xi'} + \sum_{j'\xi'}b_{jj'}^{\xi\xi'} \dot{u}_{j'}^{~\xi'} = q_j E^{\xi}(t),$$ (19)

where m_j, q_j are masses and charges of particles, $K_{jj'}^{\xi\xi'}$ are constants of elastic interaction, $E^{\xi}$is $\xi$th component of electrical field vector $\vec{E},\ u_j^{\xi}$ are components of the 3-dimensional vector of jth ion displacement. Damping parameters b_j describes dissipation of energy. Force constants K are the coefficients of power expansion of the potential energy U:

$$U=\frac{1}{2}\sum_{j'\xi'}K_{jj'}^{\xi\xi'} u_{j}^{\xi} u_{j'}^{\xi'}+...$$ (20)

where $U(\cdots u_j^{\xi} \cdots)$ depends only on relative deviations $(u_j^{\xi}-u_{j'}^{\xi})$ with corresponding symmetry of $K_{jj'}^{\xi\xi'}$, particularly

$$K_{jj'}^{\xi\xi'}=K_{j'j}^{\xi'\xi},\ \sum_{j'\xi'}K_{jj'}^{\xi\xi'}=0.$$ (21)

Linear Eq. (19) with $\vec{E}=0$ has solutions in the form of a sum of damping oscillations

$$\vec{u}= \sum a_\kappa\vec{e}_{\kappa}\cdot \exp(-\gamma_{\kappa}t+{\rm i}\omega_{\kappa}t),$$ (22)

where $\vec{u}=(\cdots u_{\kappa j}^{\ \xi}\cdots)$ is the multidimensional deviation vector, $a_\kappa$ are constants, $\vec{e}_{\kappa}=(\cdots e_{\kappa j}^{\ \xi}\cdots)$ are eigenvectors, and $\gamma_{\kappa}-{\rm i}\omega_{\kappa}$ are eigenvalues of matrix

$$\Vert\!-\!m_j\delta_{jj'}\delta^{\xi\xi'} + K_{jj'}^{\xi\xi'} + {\rm i}\ b_{jj'}^{\xi\xi'} \Vert.$$ (23)

A solution to Eq. (19) with $\vec{E}(t)=\hat{\vec{E}}{\rm e}^{-{\rm i}\omega t}$ is $\vec{u}=\hat{\vec{u}}{\rm e}^{-{\rm i}\omega t},$

$$\hat{\vec{u}}=\sum_{\kappa\ j} \frac{q_j \hat{\vec{E}}\cdot\v... ...{\omega_{\kappa}^2-\omega^2-{\rm i}\gamma_{\kappa}\omega}\cdot$$ (24)

Dielectric properties of a substance can be described by a susceptibility $\chi_0$, which is the dipole moment per unit volume generated by a unity local strenght $\vert\hat{\vec{E}}\vert=1$ electrical field,

$$\chi_0 = \sum_{\kappa} \frac{\rho_{\kappa}\omega_{\kappa}^2} {\omega_{\kappa}^2-\omega^2-{\rm i}\gamma_{\kappa}\omega}\cdot$$ (25)

For isotropic materials, constants

$$\rho_{\kappa}=\Big\vert\sum_{j\ j'} q_{j'}q_j\vec{e}_{\kappa j'}(\hat{\vec{E}}\cdot\vec{e}_{\kappa j}^*)\Big\vert$$ (26)

do not depend upon the direction of $\hat{\vec{E}}$. In the extreme case where $\omega_{\kappa}\gg \gamma_{\kappa}$ the spectrum of the dielectric absorption can be described as a collection of non-overlapping resonant Lorentz profiles

$$\chi_0'' = \frac{\rho_{\kappa}\omega_{\kappa}^2 \gamma_{\kapp... ...(\omega_{\kappa}^2-\omega^2)^2+(\gamma_{\kappa}\omega)^2}\cdot$$ (27)

This is the case for infrared absorption bands, especially in ionic crystals ("reststralen'' bands). This corresponds to the third term of Eq. (19) being small. In the opposite case ( $\omega_{\kappa}\ll \gamma_{\kappa}$) this leads to a Debye relaxation spectrum when $\omega\ll\omega_{\kappa}$:

$$\chi_0'' = \frac{\rho_\kappa\ (\omega\tau_{\kappa})}{1+(\omega\tau_{\kappa})^2},$$ (28)

where $\tau_{\kappa}=\gamma_{\kappa}/\omega_{\kappa}^2$ is a relaxation constant. This is usually the case for polar dielectric liquids. It corresponds to the first term of Eq. (19) being negligible.

This spectrum corresponds to situations where the inertial effects (first term in Eq. (19)) are comparatively small. If this condition is satisfied, the time constant in Eq. (28) can usually be determined by measuring or calculating the timescale of energy dissipation during return to the equilibrium state, without solving Eq. (19).

A third extreme case, corresponding to the second term of Eq. (28) being negligible, is usually referred to as the Drude approximation. It applies to free carriers and is considered in Sect. 2.6.

   
2.5 The Lorentz model for crystalline insulators

The Lorentz model applies to those crystalline isolators with an ionic character and describes the resonant interaction of an electromagnetic wave with the vibrations of ions that act like bound charges. Such an interaction is caused by the electric dipoles that arise when positive, and negative ions move in opposite directions in their local electric field. This kind of dust absorption can be clearly identified in astrophysical observations as emission or absorption bands. Indeed, the optical depth $\tau$ is proportional to opacity $\kappa$ (Eq. (16)):

$$% \kappa(\omega)% = \frac{4\pi}{c\rho} \sum_{\kappa} \frac{... ...(\omega_{\kappa}^2-\omega^2)^2+(\gamma_{\kappa}\omega)^2}\cdot$$ (29)

The maximum number of optically active vibrational modes in a crystal is $n_{\rm opt}=3(s-1)$, where s is the number of atoms in an elementary cell of the lattice. They are usually classified as branches in dispersion dependence of the oscillation frequency $\omega(k)$ on wave number k. All these branches are restricted by relatively high $\omega(k)$ in the infrared region. Another 3 branches with $\omega(k)\rightarrow 0$ when $k \rightarrow$0 are acoustic modes. They are optically inactive in crystal because the total cell charge is zero. Another situation for amorphous matter is considered in Sect. 4.2.

For substances of astrophysical interest, the characteristic band frequencies fall in region $\lambda \la 50~{\rm\mu m}$. Therefore, in the FIR/mm wavelength range, the asymptotic $\omega \rightarrow 0$ of Eq. (29) can be used:

$$\kappa(\omega) = \omega^2\ \frac{4\pi}{c\rho} \sum_{\kappa} \frac{\rho_{\kappa}\gamma_{\kappa}}{\omega_{\kappa}^2}\cdot$$ (30)

Thus, the Lorentz model formally applied to dust in the FIR/mm wavelength range gives a spectral index $\beta=2$, a value often used in "standard'' dust models.

Many dust models have related the long wavelength absorption by interstellar silicates and ices ( $\lambda>100~{\rm\mu m}$) to the damping wing of the infrared-active fundamental vibrational bands, visible at shorter waves (see for example Andriesse et al. 1974, and references therein).

Laboratory infrared data are satisfactorily described (see for example Rast 1968) by a sum of resonant Lorentz profiles (or reststrahlen bands) as described by Eqs. (25) and (29), in which the various constants $\rho_{\kappa}$, $\omega_{\kappa}$, $\gamma_{\kappa}$ are determined from spectroscopic measurements. This model gives a satisfactory approximation for the infrared region, not only for ionic crystals, but also for covalent ones and for their amorphous counterparts with an appropriate choice of $\rho_{\kappa}$, $\omega_{\kappa}$, $\gamma_{\kappa}$.

The situation is more complicated for the low-frequency region, where absorption is defined by damping constant $\gamma$, in contrast to IR bands, where absorption depends weakly on $\gamma$. The damping constants b_j introduced in the semi-classical model (Eq. (19)) and therefore $\gamma_{\kappa}$ in Eq. (29) have no simple microscopic explanation. The damping is linked to heating or distribution of energy to all modes of grain vibrations, so $\gamma$ characterizes intermode interactions, which are generally frequency and temperature dependent. In contrast to single phonon absorption of resonant IR photons, a wing absorption is a 2-phonon difference process. Due to the conservation of energy, the energy of the absorbed low-frequency photon is equal to the energy difference of the two phonons in the fundamental vibrational bands (Rubens & Hertz 1912; Ewald 1922). In fact in crystals, any combination of acoustical and optical phonons that satisfies the energy and momentum conservation rules can take part in this absorption. This is the main source of infrared absorption in homoatomic crystals. In the FIR wavelength range, the two-phonon difference processes are responsible for the absorption in ionic crystals. Because the efficiency of these processes is entirely determined by the difference between the phonon occupation numbers of the 2 involved phonons, the corresponding absorption decreases drastically at decreasing temperatures (Hadni 1970).

However, the optical behavior of amorphous materials is different. The disorder of the atomic arrangement leads to a breakdown in the selection rules for the wavenumber. The absorption spectra thus reflects the whole density of vibrational states. This induces a broad absorption band in the longest wavelength range due to single phonon processes, which completely dominates the multiphonon absorption (Henning & Mutschke 1997). As a consequence, lattice absorption does not vanish at low temperatures for the amorphous solids, in contrast to the crystalline ones.

   
2.6 The Drude model for conducting materials

For free charge carriers, the previous semi-classical equation of motion (Eq. (19)) can still be used, setting all restoring force coefficients K to zero, equaling all the masses m_j to the effective mass of the free carriers $m_{\rm c}$ ( $m_j=m_{\rm c}$), and using a constant value for the damping parameters $b_j/m_j={\rm const.}=\gamma_{\rm c}$. Therefore Eq. (25) for the susceptibility simplifies as $\chi_0 = - \omega_{\rm p}^2\ (\omega^2-{\rm i}\gamma_{\rm c}\omega)^{-1}$, where $\omega_{\rm p}$ is the so-called plasma frequency $\omega_{\rm p}^2 = N_{\rm c} q_{\rm c}^2/m_{\rm c}$, and $N_{\rm c}$ the number of free carriers per unit volume.

In metals, the number of free electrons per unit volume is independent of temperature, and is within the same order of magnitude as the number of atoms ( $N_{\rm c} \sim 10^{22}$ cm^-3). This leads to a frequency $\omega_{\rm p}\sim 10^{15}$ s^-1, located in ultraviolet or in the visible, and to correspondingly high conductivity $\sigma=\omega_{\rm p}^2/\gamma_{\rm c}$.

In the long wavelength range, far from the plasma frequency ( $\omega\ll\omega_{\rm p}$) and when the modulus of $\epsilon$ is large compare to 1, the asymptotic behavior of the opacity is given by

$$\kappa(\omega) = \omega^2 \frac{9\gamma_{\rm c}^2}{c\rho\ \omega_{\rm p}^4},$$ (31)

which also gives a spectral index $\beta=2$. Metallic materials and high-conductivity semiconductors are expected to follow this model. Possible dust absorbers of this type are graphite and metallic iron. Therefore, both the Lorentz and the Drude models in the low frequency limit yield $\beta=2$ in Eq. (9). These very different models have common features. Their characteristic times ( $\omega_{\kappa}^{-1}$, $\omega_{\rm p}^{-1}$) are far from the FIR/mm region. In this case a good approximation is obtained by asymptotic expansion of the absorption coefficient $\alpha(\omega)$on even powers of frequency  $\omega^{2n}$, n>0. In the general case, the first term dominates (n=1) and then we have $\beta=2$.

Other values (e.g. n=2, $\beta=4$) are possible if physical mechanisms suppressing the n=1 term are present, as in the case of screening considered in Sect. 4.2. Intermediate values of $\beta$ provide clear evidence of specific processes whose time constants are on the order of $\omega^{-1}$.

   
3 Temperature and spectral variations of the optical properties of dust

Recent astronomical observations or laboratory experiments seem to give some clues that the basic knowledge of the FIR/mm optical properties fails to give a satisfactory understanding of all data. We first review in this chapter the astronomical evidence that the FIR/mm thermal dust emission can not simply be described by a wavelength and temperature-independent spectral index. Second, we review some laboratory work carried out on analogs of interstellar grains, which clearly shows unexpected behaviour of the optical properties with wavelength and temperature that are not compatible with the classical models.

   
3.1 Astronomical evidence

Analyses of millimeter and submillimeter emission from molecular clouds have found spectral indices between $\beta=1.5$ (Walker et al. 1990) and $\beta=2$ (Gordon 1988; Wright et al. 1992). However some values in excess of 2 can also be found in the literature (Schwartz 1982; Gordon 1990). A few studies of dense clouds have yielded spectral indices around 1 (Gordon 1988; Woody et al. 1989; Oldham et al. 1994), but observations of the disks of gas and dust around young stars can indicate values sometimes less than 1 (Chandler et al. 1995). A submillimeter continuum study of the Oph 1 cloud (Andre et al. 1993) found a typical value of $\beta$ equal to 1.5 but attributed changes in observed ratios to temperature variations. Without appropriate temperature data, the authors could not be conclusive on this issue.

More recent observations have shown that the actual spectral energy distribution of dust emission could be significantly more complicated than described above. Analysis of the FIRAS results have shown that the emission spectrum is significantly flatter than expected along the galactic plane, with a slope roughly compatible with $\beta=1$ (see Reach et al. 1995, in particular their Fig. 7). They attributed the flattening compared to the $\beta=2$canonical value to an additional component peaking in the millimeter range and favored the millimeter excess being due to the existence of very cold dust at $T_{\rm d}=5{-}7$ K in our Galaxy. However, the origin of such very cold dust remained unexplained, and the poor angular resolution of the FIRAS data ( ${\simeq} 7\hbox{$^\circ$ }$) did not allow further investigation.

Finkbeiner et al. (1999) later adopted a similar model. They showed that a good fit of the overall FIRAS data could be obtained using a mixture of warm dust at $T_{\rm d}\simeq16$ K and very cold dust at $T_{\rm d}\simeq 9$ K. The temperature of the warm component was set using a fit to the DIRBE data near $\lambda=200~{\rm\mu m}$, while the temperature of the very-cold component was set assuming an independent dust component immersed in the same radiation field but with different absorption properties. They implied that these two components could be graphite (warm) and silicates (very cold). However, they gave no quantitative argument to support this hypothesis.

PRONAOS was a French balloon-borne experiment designed to determine both the FIR dust emissivity spectral index and temperatures, by measuring the dust emission in four broad spectral bands centered at 200, 260, 360, and $580~{\rm\mu m}$ (see Lamarre et al. 1994; Serra et al. 2002; Pajot et al. 2006). The analysis of the PRONAOS observations toward several regions of the sky, ranging from diffuse molecular clouds to star forming regions such as Orion or M 17, have shown an anti-correlation between the dust temperature and the emissivity spectral index (see Dupac et al. 2003a). The dust emissivity spectral index varies smoothly from high values ( $\beta\simeq2.5$) for cold dust at $T_{\rm d}\simeq 12$ K to much lower values ($\beta=1$) for dust at higher temperatures (up to $T_{\rm d}\simeq 80$ K) in star-forming regions like Orion. Dupac et al. (2003a) shows that these variations were not caused by the fit procedure used and were unlikely to be due to the presence of very cold dust. They conclude that these variations may be intrinsic properties of the dust.

Extragalactic observations have also revealed unexpected behavior of dust emission at long wavelengths. For instance, by combining JCMT data at 450 and $850~{\rm\mu m}$, IRAM 30 m data at 1.2 mm with ISO and IRAS measurements in the IR and FIR, Galliano et al. (2003, 2005) prove a strong millimeter excess towards a set of four blue-compact galaxies (NGC1569, II Zw40, He 2-10, and NGC 1140). They attribute this excess to the presence of very cold dust ( $T_{\rm d}=5{-}7$ K) with a very flat emissivity index ($\beta=1$), which would have to be concentrating into very dense clumps spread over the whole galaxy. This very cold dust would then account for about half or more of the total dust mass of those galaxies.

   
3.2 Laboratory evidence

3.2.1 Comparing laboratory and astronomical data

As shown in Sect. 2.3, laboratory data on absorption spectral index could be used for interpreting astronomical observations, taking the quantitative limitations considered here into account. Some other limitations should also be considered.

First, there is a difference between bulk laboratory samples and the cosmic dust grains. The latter probably has a high surface-to-volume ratio, and surface effects are as important for spectroscopy just as for the physics and the chemistry of dust. Second, samples used for laboratory measurements are generally synthesized in very small amounts. Along with the low values of the absorption coefficients, it requires high-sensitivity dedicated instruments, such as ³He-cooled bolometers or/and high-power sources of radiation.

In the laboratory, indirect IR/FIR/mm spectroscopy methods are also used, such as Kramers-Kronig spectroscopy. The dielectric constant counterparts $(\epsilon'-1)$ and $\epsilon''$ are Hilbert transforms of each other and are related through the Kramers-Kronig relations (see Landau & Lifshits 1982) by

 

$$\epsilon'(\omega) = \frac{1}{\pi}\int^{\infty}_{-\infty}\frac{\epsilon''(x)}{x-\omega}\ {\rm d}x +1,$$ (32)

$$\epsilon''(\omega) = -\frac{1}{\pi}\int^{\infty}_{-\infty}\frac{\epsilon'(x)}{x-\omega}\ {\rm d}x,$$

where both integrals are the main Cauchy values. A simple form of Eq. (32) uses $\omega=0$:

$$\epsilon'(0)=\frac{2}{\pi}\int^{\infty}_0 \epsilon''(x)\ {\rm d}\ln x+1.$$ (33)

These relations, together with reflection measurements $R(\omega)=\vert(\sqrt{\epsilon}-1)/(\sqrt{\epsilon}+1)\vert^2$, allows us to resolve the system of equation with respect to $\epsilon'(\omega)$ and $\epsilon''(\omega)$ and therefore to calculate an absorption spectrum. Possible method inaccuracies due to inverse problem solution should be considered carefully. In particular, it is sometimes assumed as a prior condition that the solution takes the form of a power spectrum. In practice, the most reliable approach should be to use the same model for interpretating both the astrophysical and laboratory data.

   
3.2.2 Spectral index variations in disordered solids

It has been known for more than 25 years that disordered solids (glasses and mixed crystals) have significantly higher submillimeter absorption than perfect crystals of similar chemical nature.

Mon et al. (1975) performed first absorption measurements on various amorphous dielectric materials in the millimeter wavelength range (1 mm-5 mm) and at low temperatures (in the range 0.5 K-10 K). A strong temperature-dependence of the absorption was observed, leading to a millimeter absorption excess at low temperatures. They already attributed this effect to the presence of a distribution of TLS in the studied materials.

Other experimental results with bulk materials at frequencies $\omega/2\pi c$ between 0.1 and 150 cm^-1 and T between 10 K and 300 K were discussed by Strom & Taylor (1977). The spectra follow a general empirical relation $\alpha = K\omega^2$between 15 and 150 cm^-1, where the temperature-independent constant K depends on the material considered and its internal structure. At lower frequencies $\omega/2\pi<15$ cm^-1( $\lambda > 700~{\rm\mu m}$), the slopes of the spectra are variable ($\beta$ between 0 and 3) and temperature-dependent, as shown in Strom & Taylor (1977).

Bösch (1978) performed absorption measurements for temperatures between 300 K and 1.6 K in the FIR/mm range (500 $~{\rm\mu m}$-10 mm) on amorphous silicates mainly composed of SiO₂. The studied silicate is not a typical grain analog. However, measurements performed over such a wide spectral and temperature range are particularly appropriate for studying the physics governing absorption processes. Bösch found a strong temperature and frequency dependency of the absorption coefficient in the millimeter range ( $\lambda>500~{\rm\mu m}$). The absorption coefficient was characterized by an increase in the spectral index with temperature from $\beta=1.6$ at 300 K to $\beta=3$ at 10 K and by a strong temperature dependence. To describe this temperature and frequency behavior, Bösch also referred to the existence of TLS in the material.

Therefore, laboratory spectroscopy shows deviation in the absorption spectral index from the canonical $\beta=2$ in the FIR/mm. The temperature-independent absorption process that provides $\beta=2$ in the spectral region $\lambda<700~{\rm\mu m}$ does not seem to dominate the absorption at longer wavelengths. We will see in Sect. 4 that the model proposed here includes this phenomenon in a natural manner.

3.2.3 Spectral variations in dust analogues 

FIR/mm optical properties of solids with chemical composition and structure considered as analogues of interstellar dust have been investigated in a few laboratory studies (for example Koike et al. 1980; Agladze 1996; Henning & Mutschke 1997; Mennella et al. 1998; Boudet et al. 2005). Particular interest has been given to low-temperature investigations of amorphous or other disordered materials, which reveal systematic deviations from the phonon theory of crystal solids and was interpreted using the TLS theory (see Bösch 1978; Phillips 1987; Kühn 2003).

A few laboratory measurements have been performed on amorphous solids in the FIR/mm range at low temperature. However, some of these studies were performed over a wide range of temperatures and revealed a variation in the spectral behavior of the absorption coefficient with temperature.

Agladze et al. (1996) performed absorption measurements on typical interstellar analog grains (crystalline and amorphous silicate grains) at low temperatures (1.2-30 K) in the millimeter range (700  $~{\rm\mu m}$-2.9 mm). They found an unusual behavior by the millimeter absorption of amorphous silicates: the absorption coefficient decreases with temperature down to about 20 K and then increases again with decreasing temperature. They described the frequency dependence of the absorption coefficient using a temperature-dependent spectral index. Depending on the samples, the spectral index can decrease from $\beta=2.6$ at 10 K down to $\beta=1.8$ near 30 K. They referred to the tunneling effect in the two levels systems to explain this behavior.

Mennella et al. (1998) measured the temperature dependence of the absorption coefficients between 295 K and 24 K and for wavelengths between 20 $~{\rm\mu m}$ and 2 mm on different kinds of cosmic grain analogs (amorphous carbon grains, amorphous, and crystalline silicates). They reported a significant temperature dependence of the spectral index of the absorption coefficient, which was particularly strong for their amorphous iron-silicate sample. Their measurements showed a systematic decrease in the spectral index with increasing temperature.

Boudet et al. (2005) performed measurements on different types of amorphous silicates (typical analog grains and SiO₂samples) for temperatures between $\rm 10$ K and $\rm 300$ K and wavelength from $\rm 100~{\rm\mu m}$ to $\rm 2$ mm. They found a strong temperature and frequency dependency of the absorption coefficient. They defined spectral indexes depending on the wavelength range considered. For wavelengths between $\rm 500~{\rm\mu m}$ and $\rm 1~mm$they found a pronounced decrease in the spectral index with increasing temperature whereas the spectral index was found to be constant with temperature for wavelengths between $100~{\rm\mu m}$ and $200~{\rm\mu m}$. They put some SiO₂ samples through a strong thermal treatment to remove most SiOH groups, and observed that the temperature-dependent absorption disappears totally. They thus identified the silanol groups as the defects that, in their silicate sample, are at the origin of the behavior. Considering that the OH groups cannot simply increase the coupling between the photon and the bulk phonons, they assumed it was reasonable that the defects induce closely-spaced local energy minima, which correspond to the dynamics of a "particle'' in an asymmetric double-well potential.

These laboratory studies strongly indicate that the absorption coefficients of amorphous grain analog materials vary substantially both with temperature and frequency.

   
4 A model to explain temperature and spectral variation of dust optical properties

We present in this chapter our model for the FIR/mm dust emission. It takes into account some characteristics specific to the amorphous state in general. We first introduce the way we model the disordered structure of the amorphous state; the disorder is described at two levels, using the DCD theory and the TLS theory, both coming from solid state physics and presenting a large degree of universality. Both theories are then explained in detail, and the non trivial wavelength and temperature dependence of the absorption is discussed.

4.1 Modeling disordered structure

Two mechanisms have to be considered when dealing with amorphous materials, in order to take disorder into account.

First, one should consider the acoustic oscillation excitation based upon the interaction of solids with electromagnetic waves due to disordered charge distribution (DCD, see Schlömann 1964). This effect takes place not only in amorphous media, but also in disordered crystals like mixed and polycrystals, and in some monocrystals with an inverse spinel structure, for example. This absorption is independent of temperature. The DCD model introduces a correlation length $l_{\rm c}$ that quantifies the scale over which the atomic structure of the material realizes charge neutrality. The absorption spectrum of such a structure presents two asymptotic behaviors. Towards short wavelengths, the absorption is characterized by a spectral index $\beta=2$, and in the longest wavelength range, the spectral index tends towards $\beta=4$. The location of the spectral range in which the transition between those two regimes occurs is directly related to the correlation length. The DCD model is described in Sect. 4.2.

The second mechanism describes the disorder on an atomic scale as a distribution of tunneling states, which leads to a temperature-dependent optical absorption and enables to be explained the temperature dependence of the spectral index observed in laboratory experiment. This model was originally developed to explain the unusual thermal and optical properties of amorphous material at low temperatures, and it has been described in detail by Phillips (1972), Anderson et al. (1972; see also Phillips 1987, for a review). In particular, the TLS model was developed to explain the fact that, at low temperatures, the specific heat of amorphous solids exceeds what is expected from the Debye theory for solids. This anomalous behavior is common to all amorphous materials and therefore appears independent of their detailed chemical composition and structure. The existence of these TLS has been pointed out by means of microwave and submillimeter spectroscopy experiments (e.g. Agladze et al. 1996; Bösch 1978, and references therein). The TLS model is described in Sect. 4.3.

The two mechanisms considered here for the absorption by interstellar dust probably dominate in the FIR/mm and longer wavelength domain. Both have a large degree of universality without specific signatures characterizing the chemical nature of the absorbing substance. On the other hand, they are sensitive to the internal structure of the solids, in particular, to their degree of ordering and to mechanical properties such as density and elasticity.

   
4.2 DCD absorption

4.2.1 DCD Theory

In amorphous materials, the lack of long-range order permits single phonon/photon interactions with excitation of any modes of mechanical vibrations. Low-frequency vibrations are not linked to molecular vibrational bands, which is the case for ice or silicate bands. In infinite media they correspond to traveling acoustic waves. In finite bodies like interstellar dust grains, these modes can be described as standing waves, we use the term phonon here for the quantum of vibrational motion and do not restrict it (if not stated otherwise) to periodic lattice or infinite media.

The concept of phonon quasiparticles (Tamm 1930) arises as the result of quantizing the vibrational motion in the harmonic approximation, which takes into account only the quadratic terms (on atom displacements) of energy in Eq. (20). It corresponds to Eq. (19) with b_j=0. Dissipation effects arise in the phonon concept as the result of phonon interaction, taking into account anharmonic terms in Eq. (20). Their role was first pointed out by Debye (1914).

For harmonic oscillations, $u =\hat{u}{\rm e}^{-{\rm i}\omega t}$ with electrical field $\vec{E}(t)=\hat{\vec{E}}{\rm e}^{-{\rm i}\omega t}$, Eq. (19) takes the form

$$% -\omega^2\ \hat{u}_j^{\xi} +\sum_{j'\xi'} D_{jj'}^{\xi\xi'}\ \hat{u}_{j'}^{\xi'} = % \frac{q_j}{m_j}\ \hat{E}^{\xi},$$ (34)

where the coefficients of the dynamics matrix are

$$D_{jj'}^{\xi\xi'}=K_{jj'}^{\xi\xi'}/m_j.$$

The solution of the dispersion equation

$$-\omega^2\ e_j^{\xi}+D_{jj'}^{\xi\xi'}\ e_{j'}^{\xi'}=0$$ (35)

leads to resonant frequencies $\omega_{\kappa}$ and corresponding eigenvectors $\vec{e}_{\kappa}=(\cdots e_{\kappa j}^{\ \xi}\cdots)$. In an orthonormal basis of vectors $\vec{e}_{\kappa}$, when

$$\sum_{\xi\ j} e_{\kappa j}^{\xi} e_{\kappa'j}^{*\xi}=\delta_{\kappa\kappa'},$$ (36)

Eq. (34) takes a simple form

$$% (\omega_{\kappa}^2-\omega^2 )\ \hat{u}_{\kappa} = \vec{g}_{\kappa}\cdot\vec{\hat{E}},$$ (37)

where $\hat{u}_{\kappa}$ are amplitudes of normal oscillation,

$$g_{\kappa}^{\xi} = \sum_{j} q_j m_j^{-1} e_{\kappa j}^{*\xi},$$

and ^* denotes complex conjugates.

The susceptibility corresponding to Eq. (37) is given by the tensor

$$\chi^{\xi\xi'}_0=\sum_{\kappa j'j}\frac {q_{j'} q_j m_j^{-1} ... ...'}} {\omega_{\kappa}^2-\omega^2-{\rm i}\omega\gamma_{\kappa}},$$ (38)

or in the macroscopically isotropic case, by the scalar

$$\chi_0=(\chi^{xx}_0+\chi^{yy}_0+\chi^{zz}_0)/3;$$

therefore,

$$\chi_0''={\rm Im}\ \chi_0=\frac{1}{3}\sum_{\kappa j'j\xi} \fr... ...(\omega_{\kappa}^2-\omega^2)^2+(\omega\gamma_{\kappa})^2}\cdot$$ (39)

Damping factors $\gamma_{\kappa}$ were introduced here to get a stationary solution. The $\gamma_{\kappa}$ values depend on phonon interactions, as discussed above. In the FIR/mm region, overlapping resonance curves produces a quasi continuous absorption spectrum. The integrated spectrum in Eq. (39) can be simplified by replacing the integral of the sum of weak individual profiles by the sum of their integrated values:

$$\int_{-\infty}^{+\infty}{\rm Im} \frac{{\rm d}\omega} {\omega... ...} = \frac{\pi/2}{\sqrt{\omega_\kappa^2-\gamma_\kappa^2}}\cdot$$ (40)

This leads to:

$$\left<\ \chi_0''(\omega)\ \right> =\frac{\pi}{6\omega} \sum_... ...j^{-1} e_{\kappa j}^{*\xi} e_{\kappa j'}^{\xi}}{\Delta\omega},$$ (41)

where $\Delta\omega$ is an averaging interval satisfying

$$\omega\gg\Delta\omega\gg\gamma_\kappa.$$ (42)

A number of states is proportional to $\Delta\omega$; therefore, the result does not depend on $\Delta\omega$ and $\gamma_{\kappa}$ when condition in Eq. (42) is satisfied, which also defines the validity region of the results of Vinogradov (1960) and Schlömann (1964). In macroscopically uniform and isotropic cases, the eigenvectors of Eq. (35) in the form of plane waves are

$$e_{\kappa j}^{\xi}= \tilde{e}_{\vec{k}\eta}^{\;\xi}\exp(i\vec{k}\cdot\vec{r}_j).$$ (43)

The manifold of solutions of the dispersion equation is ordered here in accordance with vector $\vert\vec{k}\vert$ and splits into several branches labeled here by index $\eta$. Each branch has a different dispersion curve $\omega=\omega_\eta(\vert\vec{k}\vert)$. Most branches are restricted to high frequencies (infrared), and only three branches extend down to the FIR/mm region. They are so-called "acoustical'' modes. For such low frequencies, these have a linear dispersion dependence

$$\omega_\eta=\upsilon_\eta\ k,$$ (44)

where $k=\vert\vec{k}\vert$, and $\ \upsilon_\eta$ is the speed of sound. One branch corresponds to longitudinal propagation ( $\eta=l$), the two others to 2 polarizations of transverse propagation ( $\eta=t$). Numerically, $\upsilon_{\rm t}<\upsilon_l$. Displacements of neighbor ions in sound oscillations are synchronous, and $\tilde{e}_{\vec{k}\eta}^{\;\xi}$does not depend on j. The vectors of displacements $\vec{e}_{\vec{k}\eta}= [\tilde{e}_{\vec{k}\eta}^{\ x},\tilde{e}_{\vec{k}\eta}^{\ y},\tilde{e}_{\vec{k}\eta}^{\ z}]$ are oriented such that $\vec{e}_{\vec{k}l}\parallel\vec{k},\ \vec{e}_{\vec{k}t}\perp\vec{k}$. Their absolute values are defined by the normalization equation, Eq. (36). The number of normal oscillations per unit volume in the frequency interval $\Delta\omega$ is $N_\omega\Delta\omega$, where

$$N_\omega=\frac{\omega^2}{4\pi^{2}\upsilon_\eta^3}$$ (45)

and therefore

$$\chi_0''(\omega)=\frac{\omega N}{24\pi} \left(\frac{2}{\upsi... ...l^3}\left < \frac{q^2}{m} \right >_{\omega/\upsilon_l}\right),$$ (46)

where N is number of ions per unit volume and

$$\left < \frac{q^2}{m} \right >_k = \frac{1}{N^2} \sum_{j,j'=... ...{q_j q_{j'}}{m_j}\exp(i\vec{k}\cdot(\vec{r}_j-\vec{r}_{j'}) ).$$ (47)

The equation derived here is a generalization of the Schlömann's (1964) expression for one type of disordered ion in a perfect crystal $\langle q^2\rangle _k^{}m^{-1}$, which was defined only by a spatial power spectrum of charge distribution $\langle q^2\rangle_k$.

4.2.2 Spectrum of DCD absorption

The absorption spectrum corresponding to Eq. (47) can be calculated using correlation functions. It is a cross-correlation function for the general case

$$\phi_{\rm c}(\Delta\vec{r})= \Big< \frac{q_j q_{j'}}{m_j}\Big > \left( \overline{q^2_{\rm m}} \right)^{-1}.$$ (48)

Here $\langle\ \rangle$ designates an ensemble average, $\Delta\vec{r}=\vec{r}_j-\vec{r}_{j'}$. Parameters $\overline{q_1^2},\ \overline{q^2_{\rm m}}$ are described below and do not depend on  $\Delta\vec{r}$.

The correlation function can be determined from a physical model of considered material. The first, simplest model has an uncorrelated charge distribution and therefore a delta-function correlation

$$\phi(\Delta\vec{r})=\delta(\Delta\vec{r})= \delta(\Delta x)\delta(\Delta y)\delta(\Delta z).$$ (49)

In the general case, a cross-correlation function, Eq. (49) gives a flat spatial spectrum with intensity not depending on k:

$$\Big < \frac{q^2}{m} \Big >_k = \overline{q^2_{\rm m}},$$ (50)

and therefore

$$\chi_0''(\omega)=\frac{\omega}{12\pi \upsilon_{\rm t}^3}\ \overline{q^2_{\rm m}}.$$ (51)

Here, we omitted the second term in Eq. (46) since $\upsilon_l^{-3}\ll 2 \upsilon_{\rm t}^{-3}$. An autocorrelation function should be used for $m=\rm const.$

$$\phi_a(\Delta\vec{r})= \langle q(\vec{r}_j) q(\vec{r}_i) \rangle \left( \overline{q_1^2}\right)^{-1},$$ (52)

with a flat spatial spectrum $\langle q^2\rangle_k = \overline{q_1^2}$ for the $\delta$ correlation (Eq. (49)) leading to Eq. (51) with a constant $\overline{q^2_{\rm m}}=\overline{q^2_1}/m$.

In perfect crystals the charge distribution is entirely correlated for long distances (long-range order) and the absorption mechanism described by Eqs. (46) and  (51) does not take place. Nevertheless, even for crystals with a perfect lattice, a disordered charge distribution is possible if the lattice permits a non unique configuration of the charges in an elementary cell, as is the case, for example, for the cubic lattice of salt (NaCl). In such cases, only stochastically distributed charges should be included in the calculation of $\overline{q^2_1}$, which can then be defined as

$$\overline{q^2_1} =\frac{1}{N_1}\sum_{j=1}^{N_1} (\Delta q_j)^2 < \frac{1}{N_1}\sum_{j=1}^{N_1} q_j^2,$$ (53)

where summing is to be done only for ions involved in the disordered distribution and $\Delta q$ denotes deviations from the mean charge value, which is different for different nodes of the crystal lattice. The ordered part of the charge distribution is linked only to the spectral components with $\omega\ \widetilde{>}\ \upsilon_{\eta}/d$, where d is the lattice period. This frequency range corresponds to mid-infrared absorption bands and are not considered here.

In amorphous media, the corresponding coefficient $\overline{q^2_{\rm m}}$ should also exclude the regular part of the charge distribution produced by the short-range order in the medium. An equation similar to Eq. (53) cannot be written in the general case. The effect of this short range screening can be expressed by the inequality

$$\overline{q^2_{\rm m}} < \frac{1}{N}\sum_{j=1}^{N} \frac{q_j^2}{m_j}\cdot$$

The second model of Schlömann (1964) can be interpreted as screening over large scales. Large-scale electrical neutrality of crystals and glasses could be the result of high charge mobility in solutions or melts from which the corresponding materials were manufactured. The corresponding correlation function should be written in generalized case as

$$\phi(\Delta\vec{r})=\delta(\Delta\vec{r})-\frac{1}{8\pi l_{\rm c}^3}\exp(\vert\Delta\vec{r}\vert/l_{\rm c}),$$ (54)

where $l_{\rm c}$ defines the correlation length. The corresponding spectra are given by

$$\langle q^2\rangle_k = \overline{q^2_1}\ g(kl_{\rm c}),$$ (55)

where

 

g(x) = 1-(1+x²)^-2 (56)

and

$$\chi_0''(\omega)= \frac{\omega}{12\pi\upsilon_{\rm t}^3}\ \overline{q^2_{\rm m}}\ g(\omega/\omega_{\rm c}),$$ (57)

where

$$\omega_{\rm c}=\upsilon_{\rm t}/l_{\rm c}.$$

The spectral dependence of Eq. (56) is

$$g\Big(\frac{\omega}{\omega_{\rm c}}\Big) \left\{ \begin{ar... ... 1 &\textrm{for }\omega\gg\omega_{\rm c}. \end{array} \right.$$ (58)

This is also true of other often-used correlation function shapes, such as step correlation functions:

$$\phi(\vec{r})= \left\{ \begin{array}{ll} \delta(\vec{r})-3/... ...trm{for }\vert\vec{r}\vert\ge l_{\rm c} . \end{array} \right.$$ (59)

The spectral dependence in Eq. (58) is common for charge distributions, which are uncorrelated over short distances and respect neutrality on larger scales.

The DCD absorption is independent of temperature. The corresponding absorption coefficient $\alpha(\omega)\sim \omega \chi_0''$(Eq. (15)) presents two asymptotics behaviors on both sides of $\omega_{\rm c}$ that depends on the correlation length. In the high-frequency range, the spectral index takes the value $\beta=2$, and the spectral index is equal to $\beta=4$ for the lower frequency range. This frequency dependence absorption is shown in Fig. 1 for various values of $l_{\rm c}$.

Figure 1: Absorption spectrum due to a DCD. The absorption coefficient $\alpha _{\rm DCD}$ is given in arbitrary units, normalized at $\lambda=100~{\rm \mu m}$ and $l_{\rm c}=10$ nm. The solid, dot, dash, dash-dot, and dash-dot-dot lines correspond to $l_{\rm c}=0.1$, 0.3, 1, 3, and 10 nm, respectively, assuming a transverse sound velocity $\upsilon _{\rm t}=3\times 10^5$ cm/s.

   
4.3 TLS absorption

4.3.1 Phenomenological TLS theory

At low temperatures, the thermal and dielectric properties of disordered solids (e.g. glasses and mixed crystals) show definite deviations from the predictions of the phonon theory developed for perfectly ordered crystals. Most of these phenomena can be described in the formalism of the TLS, which is the simplest model of tunneling states. In comparison to DCD, it could be considered as a second approximation for the description of the electromagnetic properties of solids linked to disordering. Vacancies that are unavoidable in disordered lattices produce, in first approximation, a possibility of chaotic distribution of some charges, which was considered in the previous section. This distribution is considered as frozen in the lattice due to the large height of barrier-dividing vacancies in comparison with the available thermal energy. The TLS theory takes into account transitions of charges due to quantum tunneling possible at any temperature.

The TLS theory is the result of a macroscopic analysis of phenomena as opposed to an exact microscopic description. It is based on the hypothesis of a flat distribution of tunneling states with energy-level differences that are sufficiently low compared to the usual vibrational energy. The temperature dependence of TLS properties results from deviations from the model of harmonic oscillator, which is a good approximation of the description of crystal lattice vibrations. The origin of TLS is vacancies in the lattice, which gives the possibility of having different equilibrium spatial configurations for atoms or groups of atoms. As a result, the potential curve has at least a double-well form (DWP, see Fig. 2), in contrast to the single-well potential (SWP) of harmonic oscillators in crystal without vacancies.

Figure 2: Asymmetric double-well potential (DWP). The ordinate is the potential energy. Abscissa is a coordinate in configuration space, in the simple case of a rotation angle or a linear displacement of groups of atoms in the lattice. $\Delta$ is the value of potential asymmetry. The difference in energies of the lowest levels results from the asymmetry $\Delta$ and tunnel splitting $\Delta _0$.

In the TLS model described e.g. by Phillips (1987), the energy distribution of states is defined with the assumptions that the asymmetry $\Delta$ of DWP and the height V of the DWP potential barrier are stochastic values and that their probability distribution can be assumed constant over a certain range of V and $\Delta$ values. As for the well-known tunneling probability, the energy-splitting $\Delta _0$ has an exponential dependence with barrier width d:

$$\Delta_0=\hbar\Omega\exp[-d(2mV)^{1/2}/\hbar],$$ (60)

where $\hbar\Omega=(E_1+E_2)/2$ is the average energy of levels and m the effective mass. This exponential dependence, along with the flat probability distributions of $\Delta$ and V, leads to the 2-dimensional distribution

$$P(\Delta_0,\Delta)=\overline{P}/\Delta_0,$$ (61)

where $\overline{P}$ is a constant. The integral of P over $\Delta _0$does not diverge because the lowest $\Delta _0$ values are limited by the highest energy V through Eq. (60).

The exact form of $P(\Delta_0)$ in the vicinity of $\Delta_0^{\min}$is not important in most cases because the corresponding energy ( $\Delta_0^{\min}/k^{}_{\rm\tiny B}\ll 0.01$ K) is very far from the values under consideration. A value of constant $\overline{P}$ is not defined by the TLS theory, because it varies from one material to another and can be regarded as a free parameter to be determined from experiments or observations.

In order to better understand the TLS phenomenological approach, Kühn (2003) built a microscopic model using the translational invariant Hamiltonian

  

$$\sum_j\frac{p_j^2}{2m}+U(\{u_j\}),$$                          H =

$$\frac{1}{4}\sum_{ij}J_{ij}(u_i-u_j)^2+\frac{g}{2N}\sum_{ij}(u_i-u_j)^4,$$ U = (62)

where p_j are impulses, and u_j are displacements of the jth atom from some reference points. The first term is a quadratic potential with random coefficients J_ij. The second term is a stabilizing quartic potential (not random) with a constant g>0coefficient. The stochastically simulated coefficients J_ij lead to a disordered equilibrium arrangement of atoms. Disordered potential relief (Eq. (62)) reproduces an ensemble of states both SWP and DWP, with a broad spectrum of barrier heights V and asymmetries $\Delta$.

Some TLS effects associated with excited-state transitions of TLS could be brought to the fore in some laboratory experiments (FitzGerald et al. 1994; Sievers et al. 1998) on a number of mixed fluorite crystals and soda-lime silica glass. However, these effects can only be seen in a narrow temperature range, roughly below 15 K, before relaxation processes dominate at higher temperatures. These effects will be incorporated into our model in a forthcoming paper. In the following, we restrict our consideration to the two lowest levels with energy difference  $E=E_2-E_1=\sqrt{\Delta^{2} +\Delta_0^{2}}$.

As opposed to thermal properties, the TLS absorption spectrum depends not only on the TLS density-of-state distribution P over energies $\Delta _0$, $\Delta$, but also on the values of the matrix elements of the dipole transitions between levels

$$\langle 1 \vert \mu \vert 2 \rangle=\mu_{\rm b}\frac{\Delta_0}{E},$$ (63)

where $\mu$ is an electric dipole moment,

$$\mu_{\rm b}=\frac{1}{2}\frac{\partial \Delta}{\partial \xi},$$ (64)

where $\xi$ is the value of the electrical field.

The dipole moment $\mu_{\rm b}$ is the second TLS parameter, which characterizes its interaction with the electromagnetic field. For isotropic materials, the square of dipole moment should be averaged over all directions as $\langle \mu_{\rm b}^2 \rangle=\mu_{\rm b}^2/3$. Like the parameter $\overline{P}$ discussed above, parameter $\mu_{\rm b}$ varies from one material to another and can be regarded as a free parameter to be determined by experiments or observations. In some experimental publications, values of the dipole moment are given, which are not corrected for the local field and orientation averaging, $\mu'=\mu_{\rm b}(\epsilon'+2)/3\sqrt{3}$.

   

Table 1: Parameters of disordered solids.

Material	$\rho$	$\epsilon$	$\upsilon_{\rm t}$	$l_{\rm c}$	$\langle q^2\rangle/e^2$	$\mu_{\rm b}^{ }$	$\overline{P}\mu_{\rm b}^2$	$\gamma_{\rm e}$	Refs.
	g/cm3		cm/s	nm		D		eV
Silica Glass (SiO2)	2.2	3.8	$4.1\times 10^{5}$	1-4	1.1	0.18a	$1\times 10^{-5}$	1.5	1, 2, 3
Soda-silica glassb (Na2O$\cdot$3SiO2)	2.9	2.6	$3\times 10^{5}$	3	1		$1.4\times 10^{-3}$	0.5	4
Chalogenic glasses				1-4	0.1-1	0.5-5			2
Inverse spinel crystals			$4\times 10^{5}$		1/14				5

^a A previously reported value of dipole momentum 0.6 D was overestimated due to the absence of local field correction.
^b Barrier heights' distribution parameters: maximum at $V_{\rm m}/k^{}_{\rm\tiny B}=550$ K, width $V_0/k^{}_{\rm\tiny B}=410$ K, $V_{\min}/k^{}_{\rm\tiny B}=50$ K.
References: 1: Hubbard et al. (2003); 2: Strom & Taylor (1977); 3: Phillips (1987); 4: Bösch (1978); 5: Schlömann (1964).

The TLS absorption spectrum shape can be obtained by solving the Bloch equations that describe the interaction of TLS with electromagnetic wave and lattice oscillations. The latter interaction is described by the third material dependent parameter of the TLS theory, i.e. the elastic dipole moment given by

$$\gamma_{\rm e}=\frac{1}{2}\frac{\partial\Delta}{\partial e},$$ (65)

where e is the value of the strain field. The order of magnitude of $\mu_{\rm b}$ and $\gamma_{\rm e}$ is about 1 D and 1 eV, respectively.

The solution of general equations for the TLS absorption can be described as the result of two different contributions: an absorption that has a resonant character and an attenuation that has the typical form of a relaxation. In practice it is convenient to treat these processes separately. Bösch (1978) explained his experimental results considering three processes: resonant tunneling, relaxation due to phonon-assisted tunneling, relaxation due to phonon-assisted hopping over the potential barrier. In the next sections, we will discuss each of them in more detail.

  
4.3.2 Resonant absorption

The resonant absorption of a photon of energy $\hbar \omega$concerns only those TLS characterized by the splitting energy $E_{2}-E_{1} = \hbar \omega$. Let $N(\omega)$ be the density per unit frequency of those TLS, and $N_{1}(\omega)$ and $N_2(\omega)$the densities per unit frequency of those TLS that are in the states of energy E₁ and E₂, respectively. Of course these densities are strongly temperature dependent, but only the two lower levels are considered in the TLS model and it is assumed that N=N₁ + N₂.

The resonant tunneling absorption then takes the form

$$\alpha_{\rm res}=\frac{\hbar\omega}{c\sqrt{\epsilon'}} B_{12}\Big(\frac{\epsilon'+2}{3}\Big)^2(N_1-N_2),$$ (66)

where B₁₂ is the Einstein coefficient for absorption,

$$B_{12}=\frac{4\pi^2\langle 1 \vert\mu\vert 2 \rangle^2}{3\hbar^2}\cdot$$

The population of the levels under the thermal equilibrium hypothesis gives

$$N_1-N_2=N(\omega)\tanh(\hbar\omega/2k^{}_{\rm\tiny B}T),$$

and taking into account that $E=\hbar\omega$,

$$\alpha_{\rm res}=\frac{4\pi^2\omega}{3c\sqrt{\epsilon'}} \Bi... ...2}{3}\Big)^2\tanh(\hbar\omega/2k^{}_{\rm\tiny B}T)\ G(\omega),$$ (67)

where

$$G(\omega)=\frac{\overline{P}\mu_{\rm b}^2}{(\hbar\omega)^2}\... ...elta_0{\rm d}\Delta_0}{\sqrt{1-(\Delta_0/\hbar\omega)^2}}\cdot$$ (68)

Equation (67) was derived in the work of Hubbard et al. (2003), the term $G(\omega)$ is the so-called Optical Density Of States (ODOS), which can be determined from laboratory spectroscopy. We introduced the coefficient $p^{}_{\rm\tiny G}(\Delta_0)\la 1$into Eq. (68) to describe deviation of $P(\Delta_0)$ from Eq. (61) at large $\Delta _0$. To avoid confusion in this work, we do not change the definitions of variables where it is possible, preferring transformation of the original equations taken from different sources.

In the TLS theory, confirmed by low-temperature measurements of thermal properties, $p^{}_{\rm\tiny G}=1$ for energies lower than 10^-3 eV. If $p^{}_{\rm\tiny G}=1$ for all $\Delta _0$, Eq. (68) gives $G(\omega)=\overline{P}\mu_{\rm b}^2=\rm const$. The expression for the TLS resonant-tunneling absorption, corresponding to $p^{}_{\rm\tiny G}=1$, was used by Bösch (1978) and Schickfus et al. (1975, 1976):

$$\alpha_{\rm res} = \frac{4\pi^2 \overline{P}~\mu'^2}{c\sqrt{\epsilon'}}\omega\;\tanh(\hbar\omega/2k^{}_{\rm\tiny B}T).$$ (69)

The frequency-independent ODOS was found to be in satisfactory agreement with experiments for wavenumbers $\omega/2\pi c <12$ cm^-1. Its measured values ( $\overline{P}\mu_{\rm b}^2$) are given in Table 1. Assuming $\mu'= 1D$, the density of states is $\overline{P} = 0.9\times 10^{33}$ erg^-1 cm^-3 for soda-lime-silica glass.

More accurate spectroscopic measurements in wide spectral regions show a definite drop in $G(\omega)$ to zero within the error level; and for $\omega$ higher than some frequency $\omega_{\rm m}$. Agladze & Sievers (1998) measured a profile with a cut-off at 20 cm^-1 in high-density amorphous phase ices of H₂O and D₂O. Fitzgerald et al. (2001a) detected a sharp cut-off in fluoride-mixed crystal spectra, with wavenumber $\omega_{\rm m}/2\pi c=13$ cm^-1. Hubbard et al. (2003) observed a shallow cut-off at $\omega_{\rm m}/2\pi c=15$ cm^-1, 9 cm^-1, and 6 cm^-1 for the soda lime silica, the SiO₂, and the germanium glasses, respectively.

Hubbard et al. (2003) suggests that the distribution $P(\Delta_0,\Delta)$ (Eq. (61)) has a cut-off at energy $E_{\rm m}=\hbar\omega_{\rm m}$. In our Eq. (68) it corresponds to $p^{}_{\rm\tiny G}(\Delta_0)$ in the form of a step function

$$p^{}_{\rm\tiny G}(\Delta_0) = \left\{ \begin{array}{ll} 1... ... 0& \textrm{for $\Delta_0>E_{\rm m}$ .} \end{array} \right.$$ (70)

Under this assumption, Hubbard et al. (2003) calculates a frequency-dependent expression for the ODOS by integrating the distribution function equivalent to Eq. (68). For our consideration it is convenient to write it in the form

$$G(\omega) = g(\omega/\omega_{\rm m}) \overline{P}\mu_{\rm b}^2.$$ (71)

The ODOS of (Hubbard et al. 2003) is described by the following function:

$$g(\omega/\omega_{\rm m}) = \left\{ \begin{array}{ll} 1 & ... ... \textrm{for $\omega/\omega_{\rm m}>1$ }. \end{array} \right.$$ (72)

However, the spectra corresponding to the ODOS that they propose (Eq. (72)) exhibits a prominent break that was not observed in laboratory spectra of noncrystalline materials and has not been observed in the spectra of astronomical dust emission.

We suggest that the cut-off of the distribution function should be complemented by a continuity condition and therefore continuity of $p^{}_{\rm\tiny G}(\Delta_0)$. The simplest form of continuous function for $p^{}_{\rm\tiny G}(\Delta_0)$ is considered in Appendix A and leads to the following a shape of the ODOS:

$$g(x) = \left\{ \begin{array}{ll} 1+x^2 g_1(x) & {\rm for}... ...)-g_2(x)\sqrt{1-x^{-2}} & {\rm for} \ x>1, \end{array}\right.$$ (73)

where $x=\omega/\omega_{\rm m}$, polynoms $g_1(x),\ g_2(x)$are defined in Appendix A by Eqs. (A.3) and (A.4).

This ODOS function has no break at $\omega=\omega_{\rm m}$, and it is in better agreement with laboratory measurements of amorhous materials. It also fully agrees with TLS theory at $\omega\ll\omega_{\rm m}$. For this reason it will be used in the calculations of Sect. 5.

4.4 Relaxation processes

 

Direct consideration of the relaxation processes permits the nonresonant part of TLS absorption spectra to be evaluated through analyzing the mechanism of the TLS-level population relaxations after some deviation from equilibrium. The relaxation time constant $\tau$ is determined by the rate ($\tau^{-1}$) of transitions between levels, in which transitions due to intensive interactions with the lattice thermal oscillations (phonons) dominate.

The rate of the TLS transitions as a result of strain field generated by thermal fluctuations in a lattice was evaluated by Phillips (1972) as

$$\tau(\Delta_0,E)\ ^{-1}=\tau_1^{-1}\ \Delta_0^2\ E^{-2},$$ (74)

where

$$\tau_1^{}=a\ E^{-3}\tanh(E/2k^{}_{\rm\tiny B}T);$$ (75)

and

$$a=\pi \upsilon_{\rm t}^5 \rho \hbar^4 \gamma_{\rm e}^{-2}$$ (76)

is a material constant. Parameters included in Eq. (76) have already been described above.

Phillips (1972) also calculated an equilibrium value of the TLS dipole moment $\overline{\mu}$ and therefore the resulting susceptibility $\chi(\omega)$ in the form of a Debye spectrum (Eq. (28)):

$$\chi_0''(\omega)=\frac{\mu_{\rm b}^2}{3k^{}_{\rm\tiny B} T}\... ...}_{\rm\tiny B}T}\Big) \frac{\omega\tau}{1+(\omega\tau)^2}\cdot$$ (77)

The rate of tunneling relaxation (Eq. (74)) is proportional to the tunneling probability, so it does not take into account direct transitions in TLS whose energy is higher than the barrier height. This transition rate is proportional to the number of phonons with such energy, which should follow a Boltzman distribution and depends both on temperature and on the energy barrier height V. In glasses with $V/k^{}_{\rm\tiny B}\ga 500$ K, the time of hopping relaxation $\tau$ varies with T<100 K in accordance with the Arrhenius relation (Hunklinger & Schickfus 1981):

$$\tau=\tau_0\ \exp(V/k^{}_{\rm\tiny B}T)$$ (78)

with $\tau_0\approx 10^{-13}$ s (fused silica).

Experimental values of V and $\tau_0$ were determined using the temperature dependence of the position of the spectrum maximum ( $\omega=\tau^{-1}$, see Eq. (77)). The distribution of TLS on V leads to the broadening of the corresponding spectrum, sometimes preventing definite detection of its maximum. The measurement of the spectrum broadening permits evaluation of the width and shape of the TLS distribution on V. In practice, experimentalists often use the shape of a truncated Gaussian probability distribution to fit the energy distribution of TLS:

$$P(V)=C_V\ p^{}_{\rm\tiny V}(V),$$ (79)

where

$$p^{}_{\rm\tiny V}(V)= \left\{ \begin{array}{ll}\exp(-(V-V_... ..._{\min},\\ 0 &\textrm{ for }V<V_{\min}, \end{array} \right.$$ (80)

where C_V is the usual normalization coefficient (Eq. (D.2)) depending on $V_{\rm m},\ V_0,\ V_{\min}$, which are parameters of the distribution shape. An example of the measured numerical values of these parameters is given in Table 1.

Starting from a known distribution of TLS on energies  $E,\ \Delta_0,\ V$, and spectrum (Eq. (77)) for fixed values of  $E,\ \Delta_0,\ V$, one can synthesize a final relaxation spectrum. The integrated spectra are considered in the following sections.

4.4.1 Tunneling relaxation

The absorption spectrum due to phonon-assisted tunneling relaxation can be obtained by integration of Eq. (77) on Eand $\Delta _0$. The expression for this was provided by Fitzgerald et al. (2001b) as

$$\alpha_{\rm phon}= A\cdot\omega\cdot F(\omega,T),$$ (81)

where

$$A=\frac{\overline{P}\mu_{\rm b}^2}{3c\sqrt{\epsilon}} \left( \frac{\epsilon'+2}{3} \right)^2$$ (82)

is a material constant and

  

$$F(\omega,T)=\frac{1}{2k^{}_{\rm\tiny B}T} \int^{\infty}_{0}\int^{... ...{}_{\rm\tiny B}T}\Big) \frac{\omega\ {\rm d}\tau\ {\rm d} E}{1+\omega^2\tau^2},$$

with $\tau_1$ given by Eq. (75). A change of variables made using relation between $\Delta _0$ and $\tau$ for integral calculation.

Another, simplified equation describing tunneling relaxation absorption was used by Jackle (1972) and by Bösch (1978)

 

$$\alpha_{\rm phon}= A \int^{\infty}_{0} \frac{{\rm d}E}{2k^{}_{\rm... ...g(\frac{E}{2k^{}_{\rm\tiny B}T}\Big) \frac{\omega^2\tau_1}{1+\omega^2\tau_1^2},$$ (83)

where $\tau_1$ is given by Eq. (75), and where parameters A and a (in Eq. (75)) are to be defined from experiments. The use of Eq. (83) is equivalent to evaluating the integral in Eq. (77), substituting values by that for a symmetric DWP $\Delta=0$. Formally it is equivalent to replacing $\sqrt{1-\tau_1/\tau}$ by $\tau_1\delta(\tau-\tau_1)$. Our evaluation shows large inaccuracy in this procedure for values of absorption. For details see Appendix B.

However, this approach does give satisfactory accuracy for time-constant evaluation linked to parameter a. For soda-silica glass (Bösch 1978) $a=4.2\times 10^{-56}$ erg³ s. It corresponds to $\gamma_{\rm e}=0.5$ eV (see Table 1).

Figure 3: Absorption spectrum due to resonant tunneling in TLS. The absorption coefficient is given in relative units normalized at $\lambda=100~{\rm \mu m}$ and T=10 K. The abscissa axis is normalized to the cut-off wavelength of the ODOS $\lambda_{\rm m}=700~~{\rm \mu m}$ (see text). The solid, dot, dash and dash-dot lines correspond to T=10, 20, 50 and 100 K respectively. The various curves can be scaled using the simple relation $\lambda _{\rm m} T$=const. In the low frequency wings, curves scale as $\alpha \sim \omega ^2/T$.

Figure 4: The function $F_{\rm phon}(p)$ (solid line) and the exponent of its slope $\beta_2={\rm d}\ln F_{\rm phon}/{\rm d}\ln p$ (dotted line) plotted as a function of the parameter $p=a_{\rm phon}\omega/(2k^{}_{\rm\tiny B}T)^3$ (see text). The spectral index of phonon absorption $\beta _{\rm phon}$ is related to $\beta _2$ as $\beta =1+\beta _2$ and ranges from $\beta _{\rm phon}=1$ at low frequencies (low p values) to $\beta _{\rm phon}=2$ at high frequencies (high p values). The temperature dependence p(T) generates a shift in the absorption spectra along the frequency axis.

The full Eq. (83) for calculation of phonon-relaxation spectra can be simplified without loss of accuracy. We made it using a function of one argument $F_{\rm phon}(p)$instead of a function of two arguments $F(\omega,T)$ in Eq. (81). We introduce the intermediate parameter $p(\omega,T) = a\omega/(2k^{}_{\rm\tiny B}T)^3$. This permits us to use a one-dimensional grid of argument, efficiently reducing computing time. Interpolation in precalculated values (Table C.1 in Appendix C) makes use of this function as fast as a standard one.

Therefore an accurate and relatively simple equation for the tunneling-relaxation spectrum is

$$\alpha_{\rm phon}= A\cdot\omega\cdot F_{\rm phon}\left(a\omega(2k^{}_{\rm\tiny B}T)^{-3}\right),$$ (84)

where A (Eq. (82)) is a material constant and the function $F_{\rm phon}$ is given by

 

$$F_{\rm phon}(p)=\int^1_0\int^{\pi/2}_{x_1} \sqrt{ 1-\frac{\tan x_1}{\tan x} }\ {\rm d}x {\rm d}y,$$ (85)

$$x_1 =\arctan( p\ y\textrm{ arcth}^{-3} y ).$$ (86)

A plot of function $F_{\rm phon}(p)$ is shown in Fig. 4.

The absorption spectrum for tunneling relaxation is temperature dependent in accordance with Eq. (81). As opposed to resonant absorption, the tunneling relaxation absorption increases with temperature, as shown in Fig. 5.

  
4.4.2 Hopping relaxation

The expression for the absorption spectrum due to the hopping relaxation can be obtained by integrating Phillip's Eq. (77) over a Gaussian distribution of barrier height $p^{}_{\rm\tiny V}(V)$ given in Eq. (80):

$$\alpha_{\rm hop} = B_{\rm hop}(T)\ \omega\ \int^{\infty}_{0} \frac{\omega\tau}{1+\omega^2\tau^2} P(V)\ {\rm d}V ,$$ (87)

where $\tau=\tau_0\exp(V/k^{}_{\rm\tiny B}T)$ and the coefficient $B_{\rm hop}(T)$ does not depend on $\omega$.

Examples of spectra described by Eq. (87) are shown in Fig. 6. It is known from general principles that the asymptotic behavior of $\epsilon(\omega)$ when $\omega \rightarrow 0$ leads to $\alpha(\omega)\sim\omega^2$ (spectral index $\beta=2$). Figure 6 shows that this region begins at f<10 MHz for T= 100 K and at f<100 Hz for T< 50 K.

Figure 5: Spectra of absorption due to tunneling relaxation $\alpha _{\rm phon}$ in the FIR/mm. The absorption coefficient is given in arbitrary units, normalized at $\lambda=100~{\rm \mu m}$ and T=10 K. Line styles are the same as in Fig. 3.

Figure 6: Hopping-relaxation absorption spectra in the frequency range from 1 Hz to 1 THz. The ratio of the absorption coefficient $\alpha _{\rm hop}$ to the parameter $\omega B_{\rm hop}/2$ is plotted to show the whole spectrum on comparable scales and to show variations in the peak position with temperature. The curves are computed using a barrier height of $V_{\rm m}/k^{}_{\rm\tiny B}=550$ K (see Table  1). Line styles are the same as in Fig. 3. Thin curves correspond to unbroadened spectra with $V_0\ll V_{\rm m}$. Thick curves correspond to $V_0/k^{}_{\rm\tiny B}=410$ K.

Figure 7: Spectra of absorption due to hopping relaxation $\alpha _{\rm hop}$ in the FIR/mm. The absorption coefficient is given in arbitrary units, normalized at $\lambda=100~{\rm \mu m}$ and T=10 K. Line styles are the same as in Fig. 3. In this relatively narrow spectral region, the behavior of the curves agrees with the approximation of Eq. (89).

At high frequencies (beginning from the mm-wave region), the spectra shown in Fig. 7 are determined by the low-energy part of the barrier height distribution. It can be approximated by

$$\alpha_{\rm hop} \approx B_{\rm hop}(T)\omega P_{\min} \int^{\infty}_{V_{\min}} \frac{\omega\tau}{1+\omega^2\tau^2}\ {\rm d}V,$$ (88)

where $P_{\min}=P(V_{\min})$, and therefore

$$\alpha_{\rm hop} \approx B_{\rm hop}(T)\omega P_{\min}k^{}_{... ...m e}^{-V_{\min}/k^{}_{\rm\tiny B}T}}{\omega\tau_0}\right)\cdot$$ (89)

For spectra in Fig. 7, a value of $P_{\min}$ corresponds to $P_{\min}~k^{}_{\rm\tiny B}=\exp[-(V_{\min}-V_{\rm m})^2/V_0^2]/V_0\sqrt{\pi} = 3.11\times 10^{-4}$ K^-1.

Bösch (1978) used an approximated expression $B_{\rm hop}=h_r/T$, where h_r is a constant. The approximation used by Fitzgerald et al. (2001b) gives $B_{\rm hop}=2\pi\mu^2/3c\sqrt{\epsilon'}$ which is a material constant. In Annex D, we derive a more accurate approximation, which gives

$$B_{\rm hop}(T)=\frac{8\pi}{3c\sqrt{\epsilon'}}\frac{(\epsilon'+2)^2}{9} \overline{P}\mu_{\rm b}^2(c^{}_\Delta+\ln T),$$ (90)

where T is in Kelvin and material constant $c^{}_\Delta$ is given by

$$c^{}_\Delta=\ln\frac{k^{}_{\rm\tiny B}}{\Delta_0^{\min}}+\ln 4-1+\int_0^1\ln{\rm arcth}\ x\ {\rm d}x$$ (91)

in accordance with Eq. (90) and (D.4). Therefore $c^{}_\Delta$ is an additional parameter of tunneling states. Its value corresponds to tunnel splitting in low $\Delta$ region where tunneling state density distribution becomes significantly less than that of Eq. (61). The TLS theory has proved to apply in low-temperature experiments at T down to 2 mK (Phillips 1987), therefore $\Delta_0^{\min}/k^{}_{\rm\tiny B}<2\times 10^{-3}$ K and $c^{}_\Delta$ should verify $c^{}_\Delta\ga 5.8$. The dependence of $B_{\rm hop}(T)$ is relatively weak and a tenfold variation in T from 10 K to 100 K leads to a variation in $B_{\rm hop}(T)$ less than $30\%$.

A difference in temperature dependence of considered approximations is significant, especially for absorption calculations over a wide range of temperatures. It does not manifest itself when fitting laboratory data over a limited region of temperature, leading only to a biased value of a parameter A. The main temperature dependence of $\alpha _{\rm hop}$ is defined by the exponential dependence of $\tau$ on T (Eq. (78)).

In the work of Bösch (1978), Fitzgerald et al. (2001b) and others, a comparison between experimental and calculated temperature dependencies was performed. It showed that relaxation processes should be taken into account for temperatures $T \ga 10$ K. At these temperatures, relaxation dominates resonant tunneling. Hopping relaxation becomes significant at T about 25 K and dominates at higher temperatures. This tendency has a general character and shows a similar behavior for substances as different as fluoride mixed crystals and silica glasses.

5 Model of the dust FIR/mm emission 

The theoretical considerations given in Sect. 4 provide the formalism for calculations of dust emissivity in the far-infrared and millimeter-wave regions and can be gathered into a new model for explaining the astronomical emission in this spectral range. So far, modeling the astronomical dust-emission spectra in this spectral range has been performed on purely phenomenological grounds with little connection to solid state physics. Unlike previous models, the new approach proposed here does not assume that dust is composed of crystalline material, but instead for the first time uses theoretical results applicable to disordered materials that are likely to compose astronomical dust particles.

   
5.1 Model description

The mass absorption coefficient of dust entering Eq. (1) can be expressed as

$$\kappa_{\rm dust} = \kappa_{\rm\tiny IR}^{}+\kappa_{\rm\tiny DCD}+\kappa_{\rm\tiny TLS}^{}.$$ (92)

The first term $\kappa_{\rm\tiny IR}^{}$ includes the opacity in infrared bands and other relatively short-wave absorption processes. References to dust models describing $\kappa_{\rm\tiny IR}^{}$ are given in Sect. 1. The IR bands in the dust spectrum are linked to the optical modes of lattice oscillations, like bending vibration bands of interatomic bonds such as O-Si-O or H-O-H. As shown in the previous sections, the long-wavelength wing of infrared bands is relatively weak in the FIR/mm. Our model does not provide a new description of these processes and is complementary to IR models. While $\kappa_{\rm\tiny IR}^{}$ depends mainly on the chemical composition of dust, the following terms $\kappa_{\rm\tiny DCD}$ and $\kappa_{\rm\tiny TLS}^{}$ are linked to the disordered structure of the grain material and are largely independent of the dust chemical composition.

As described in Sect. 4, $\kappa_{\rm\tiny DCD}$ and $\kappa_{\rm\tiny TLS}^{}$ are linked to the existence of disorder in the structure of the grain material. Given the fast decrease in $\kappa_{\rm\tiny IR}^{}$ with wavelength, these effects are likely to dominate absorption in the FIR/mm. Disordered media absorption takes place in all substances excluding perfect crystals. Amorphous dust, as well as partially ordered or dirty crystalline dust, can therefore be described by the model.

The first term, $\kappa_{\rm\tiny DCD}$, corresponds to the disordered charge absorption described in Sect. 4.2. It is temperature independent. The second term, $\kappa_{\rm\tiny TLS}^{}$ is connected to the transition of charges in a lattice between vacant potential minima due to tunneling or thermal activation (Sect. 4.3). This term displays a specific spectrum of absorption over a wide range of frequencies $\omega/2\pi$, from about 1 Hz to about 1 THz. It includes three components designated in the model as resonant ( $\kappa_{\rm res}$), tunneling relaxation ( $\kappa_{\rm phon}$), and hopping relaxation ( $\kappa_{\rm hop}$):

$$\kappa_{\rm\tiny TLS}^{}(\omega,T) = \kappa_{\rm res}(\omega,T)+\kappa_{\rm phon}(\omega,T) +\kappa_{\rm hop}(\omega,T).$$ (93)

The spectral index of $\kappa_{\rm\tiny TLS}^{}$ significantly deviates from $\beta=2$ and depends on temperature. Resonant opacity $\kappa_{\rm res}$ decreases with T, relaxation opacity $\kappa_{\rm phon}$ and $\kappa_{\rm hop}$ rises with T. The various opacities of the model can be summarized as,

  

$$\kappa_{\rm\tiny DCD}(\omega) K_{\rm\tiny DCD}\omega^2 \left[1-\left(1+\omega^2/\omega_{\rm c}^2\right)^{-2}\right],$$ = (94)

$$\kappa_{\rm res} (\omega,T) K_{\rm\tiny TLS}^{}\omega\ g_{\rm res}\left(\frac{\omega}{\omega_{\rm m}}\right)\tanh \left(\frac{\hbar\omega}{2 kT} \right),$$ = (95)

$$\kappa_{\rm phon}(\omega,T) \frac{1}{4 \pi^2}K_{\rm\tiny TLS}^{}\ \omega\ F_{\rm phon} (a_{\rm phon}\omega T^{-3}),$$ = (96)

$$\kappa_{\rm hop} (\omega,T) \frac{2}{\pi} K_{\rm\tiny TLS}^{}\omega (c^{}_\Delta+\ln T) \int^{\infty}_{V_{\min}} \frac{\omega\tau(V,T) P(V){\rm d} V}{1+\omega^2\tau^2(V,T)},$$ = (97)

where $\tau(V,T)=\tau_0\ \exp(V/k^{}_{\rm\tiny B}T)$.

The amplitude factors for the DCD and TLS terms can be expressed with respect to material properties as follows:

  

$$K_{\rm\tiny DCD}= \overline{q^2_{\rm m}}\ (3\upsilon_{\rm t}^3 c\rho)^{-1},$$ (98)

$$K_{\rm\tiny TLS}^{}= \frac{4}{3} \pi^2\overline{P}\mu_{\rm b}^2\ (c\rho)^{-1},$$ (99)

$$a_{\rm phon}= \pi\gamma_{\rm e}^{-2}\upsilon_{\rm t}^5 \rho \hbar^4 (2k^{}_{\rm\tiny B})^{-3}.$$ (100)

The function $g_{\rm res}$ is the normalized ODOS spectrum given by Eq. (73). The function $F_{\rm rel}$ is defined by Eq. (C.1) and can be calculated by interpolating Eq. (C.3) using precalculated values in Table C.1 of the Appendix. The Gaussian probability distribution function of the TLS potential barrier height P(V) is defined by Eq. (80) with parameters $V_{\rm m}$, V₀, and $V_{\min}$. $K_{\rm\tiny DCD}$ and $K_{\rm\tiny TLS}^{}$ define the amplitude scale of $\kappa_{\rm\tiny DCD}$ and $\kappa_{\rm\tiny TLS}^{}$. Parameters $\omega_{\rm c}$, $\omega_{\rm m}$, $a_{\rm phon}$, $\tau_o$ define the scale of $\kappa_{\rm\tiny DCD},\ \kappa_{\rm res},\ \kappa_{\rm phon},\ \kappa_{\rm hop}$, respectively, along the $\omega$ axis. The numerical values and units of the various parameters entering the model are given in Table 2. They set the frequency and temperature dependence of the model.

   

Table 2: Parameter values.

Parameter	Value	Unit
$K_{\rm\tiny DCD}$$^\dag$	$1.2\times 10^{-24}$	${\rm cm^2~g^{-1}~s}^2$
$K_{\rm\tiny TLS}^{}$$^\dag$	$1.5\times 10^{-13}$	${\rm cm^2~g^{-1}~s}$
$\omega_{\rm c}$$^\ddag$	$3\times10^{12}$	${\rm s}^{-1}$
$\omega_{\rm m}$$^\flat$	$2.7\times 10^{12}$	${\rm s}^{-1}$
$a_{\rm phon}$$^\dag$	$5.3\times 10^{-10}$	${\rm K^3~s}$
$\tau_o$ $^{\star\star}$	10-13	${\rm s}$
$V_{\rm m}$ $^\diamond$	$7.6\times 10^{-14}$	${\rm erg}$
V0 $^\diamond$	$5.7\times 10^{-14}$	${\rm erg}$
$V_{\min}$ $^\diamond$	$6.9\times 10^{-15}$	${\rm erg}$
$c^{}_\Delta$ $^\star$	5.8

$^\dag$ Values calculated using Eqs. (98)-(100) and the following approximate mean values: $\upsilon _{\rm t}=3\times 10^5$ cm/s, $\rho = 3$ g/cm³, $\overline{P}\mu_{\rm b}^2=10^{-3}$, $\gamma_{\rm e} = 1$ eV, $\overline{q^2_{\rm m}}=e^2/m_o$, where e is the electron charge and m_o is the mass of the oxygen atom. $^\ddag$ Using $l_{\rm c}=\upsilon_{\rm t}/\omega_{\rm c}=1$ nm (see Table 1). $^\flat$ Corresponds to $2\pi c/\omega_{\rm m}=0.7$ mm (see Sect. 4.3.2). $^\diamond$ From Table 1. $^{\star\star}$ $\tau_o$ from Sect. 4.4. $^\star$ Evaluated from Sect. 4.4.2.

It is important to note that the values proposed here only give an order of magnitude for the parameters and are derived only from laboratory measurements (Sects. 3.2.3 and 4) or theoretical expectations. They only provide reference values, so are not tuned to reproduce observed astronomical spectra. Derivation in the parameter values to be used for astronomical purposes will be derived in a forthcoming paper.

5.2 Model applicability

5.2.1 Wavelength range of validity

The model presented here applies to the FIR to millimeter emission from amorphous materials. In the astronomical context, it should be used to describe the emission of large dust grains that usually dominate the observed spectra in the FIR/mm. In principle, it should also apply to the emission from amorphous very small 3-D particles whose emission usually dominates the mid-IR region of astronomical spectra.

The DCD component is valid only in the frequency range dominated by acoustic oscillations. Towards high frequencies the corresponding limit is set where the sound wavelength is on the order of the interatomic distances d. For silicates this corresponds to an electromagnetic wavelength $\lambda_d \simeq 10~{\rm\mu m}$. At shorter wavelengths, IR models should be used. Towards the long wavelengths, the validity is limited to $\omega \leq \omega_a$ with

$$\omega_a \approx \pi\upsilon_{\rm t}/a,$$ (101)

where $\upsilon_{\rm t}$ is the transverse sound velocity and a is the grain size. The DCD model should therefore be valid for all wavelengths below $\lambda_{a}=3$ mm for the range of sizes usually assumed for the BG component (e.g. 15-110 nm in Désert et al. 1990), using $\upsilon _{\rm t}=3\times 10^5$ cm/s. Application of this part of the model to small transiently-heated dust particles in the long wavelength range should therefore be handled with care.

The TLS part of the model is in principle not limited to a given spectral range.

5.2.2 Applicability in dust models

In general, modeling FIR/mm dust emission from an interstellar source can be seen as a three-step process. The first step is to select the materials that constitute the grains and to specify the optical properties (complex dielectric constants or refractive index, or opacities) of the selected materials. The second step is to choose, for each dust component, a distribution of size, shape, composite structure, orientation with respect to the direction of the incoming radiation field, in order to calculate extinction, absorption, and diffusion using Mie or some multipole approximation method. The third step is to specify the dust density and radiation field distributions and to include radiative transfer to calculate the dust emission on a given line of sight. The model presented here only concerns the first step in this dust emission modeling.

The equations proposed in Sect. 5 do not include the local field correction factor given in Eq. (15). They are valid for spherical or quasi-spherical particles. Any possible effect due to anisotropic shapes, such as a polarized emission by grains, for example, is beyond the scope of this paper and refers to the second step of modeling.

The model equations do not depend on macroscopic values of the dielectric constant and average density. All the parameters involved in the model have a local, microscopic meaning. Therefore the derived emissivities are defined by the total mass of dust, and do not depend on the porosity of the grains. The model should therefore be equally valid for fluffy and for bulk grains.

It is important to note that, unlike standard modeling, our model shows that it is not sufficient for the first modeling step to select some complex dielectric constants $\epsilon(\omega,T)$ from the databank obtained for bulk mineral samples of chemical composition corresponding to interstellar element abundancies. First, the particle size should be taken into account. Indeed, if the opacity (mass emissivity) of dust (Eq. (8)) does not depend on particle size for a given dieletric constant $\epsilon(\omega)$, the dieletric constant $\epsilon(\omega)$ itself is size dependent, as shown in Sect. 5.1. Second, both DCD and TLS processes are more sensitive to the degree and the type of disorder in a given material than to the exact chemical composition of the grain.

  
5.3 Temperature and frequency dust-opacity behavior

In the general case, the model can produce a wide variety of spectral shapes, ranging from a "classical'' spectrum $\beta=2$ when only DCD is taken into account ( $K_{\rm TLS}=0$) and the correlation length characterizing DCD is assumed infinite ( $\omega_{\rm c}=0$) to spectra with much flatter behaviors when TLS effects are included, with mean spectral index values $\beta$ over the whole FIR/mm range close to one at high dust temperatures. Some extreme cases even include the possibility of having $\beta<0$ or $\beta>2$ over some limited regions of the spectrum, due to the influence of the resonant TLS absorption and to the DCD cut-off, respectively.

However, the model generates a "typical emission profile'' when a correlation length of nanometric scale is assumed and when the TLS effects are expected and taken into account. This emission profile presents some interesting characteristic deviations compared to the commonly adopted profile (a modified black-body type emission characterized by a constant spectral index $\beta=2$). To illustrate this fact, we show in Fig. 8 the dust opacity spectrum including both the DCD and TLS components ( $\kappa_{\rm dust} = +\kappa_{\rm\tiny DCD}+\kappa_{\rm\tiny TLS}^{}$), calculated using Eqs. (94)-(97)) and adopting the parameter values given in Table 2. First, the slope of the dust opacity versus wavelength, consequently the spectral index, is no longer constant over the whole FIR/mm spectral range. As the temperature increases, it begins to deviate in the longest wavelength range, starting from a local spectral index value higher than $\beta=2$ around T=10 K, towards lower values down to around $\beta=1$. Second, as the temperature increases more, this spectral behavior propagates towards the short wavelengths of the FIR/mm range. It is remarkable that only with model parameter values referenced in the solid state physics literature, the temperature-dependent spectral behavior of the modeled dust opacity is in qualitative agreement with laboratory data in that temperature range (Sects. 3.2.2 and 3.2.3).

Figure 8: Total dust opacity spectra $\kappa_{\rm\tiny DCD}+\kappa_{\rm\tiny TLS}^{}$. The spectra are given in arbitrary units, normalized at $\lambda=100~{\rm \mu m}$ and T=10 K. The dash-dot-dot line shows the temperature-independent DCD dust opacity alone computed for $l_{\rm c}=1$ nm. Line styles for the temperature dependent total effect are the same as in Fig. 3.

The temperature and frequency behaviors have been described in Sect. 4 for each individual components of the model. The DCD opacity is independent of temperature and shows a quadratic frequency behavior ( $\kappa_{\rm\tiny DCD}\sim\omega^2$) for high frequencies and drops dramatically ( $\omega<\omega_{\rm c}$) at lower frequencies. The transition between the two regimes is set by the value of the correlation frequency ( $\omega_{\rm c}$) or, equivalently, by the correlation length of the charge distribution in the grain material. We point out that, although the transition between these two regimes is very smooth, this distinctive behavior opens the possibility of measuring the correlation length from astronomical emission spectra. Note that the DCD effect is the only one with such a rapid fall-off in the absorption with wavelength, which seems to be measured toward some astronomical sources that exhibit high $\beta$ values, including those higher than 2.

Resonant tunneling opacity takes place around $\omega<\omega_{\rm m}$ and dominates at low temperatures since $\kappa_{\rm res}\sim\omega^2/T$ (Fig. 3). As shown in Sect. 4.2, the wavelength $\lambda_{m}=2\pi c/\omega_{\rm m}$ is in the spectral range from $500~{\rm\mu m}$ to 1 mm for silicates. The tunneling relaxation spectrum is closely linked to the resonant tunneling. Their added effect is shown in Fig. 9. The amplitude of the resonant feature decreases with rising temperature, while the amplitude of the relaxation wings increases.

Figure 9: Temperature dependence of the resonant tunneling and phonon-induced tunneling relaxation opacity (see Fig. 3 for comparison). The spectra are given in arbitrary units, normalized at $\lambda=100~{\rm \mu m}$ and T=10 K. Line styles are the same as in Fig. 3.

The hopping relaxation behaves as a distinct spectral feature in the case of narrow barrier-height distributions (Sect. 4.4.2, Fig. 6). The amplitude of the effect rises with frequency as $\kappa_{\rm hop}\sim\omega^2$, until it reaches a maximum at $\omega=\tau_0^{-1}\exp(-V_0/k^{}_{\rm\tiny B}T)$. At low temperatures the location of the maximum is shifted to MHz frequencies. It reaches the FIR/mm region at $T \simeq 100$ K, a value that depends upon V₀. The hopping relaxation in the case of wide barrier-height distribution, creates a very wide spectral distribution (Fig. 6). The behavior is then $\kappa_{\rm hop}\sim T\omega^\beta$ with a frequency-dependent slope value, usually $\beta<1$ (Eq. (89), Fig. 7).

At a given temperature and wavelength, the relative efficiency of the TLS effects compared to the underlying DCD effect is governed by the density and efficiency of the optically active TLS in the material, hence by the value of the parameter $\overline{P}\mu_{\rm b}$. This value should be related to the fluffiness of the material and to the types and densities of defects in its structure.

Generally, the resonant part of the TLS effects is expected to dominate the other relaxation processes at temperatures below 10 K, while the TLS hopping mechanism should fully dominate above 30 K.

  
5.4 Model relevance

Recent astronomical observations of the interstellar dust emission in the FIR/mm wavelength range show some unexpected features that can only be related with difficulty to previous knowledge about dust optical properties in that range (Sect. 3.1). All data have to be explained with a blackbody emission modified by a frequency power law, whose exponent is constant over that FIR/mm range, with a temperature-independent asymptotic value equal to 2. Indeed, very little is known about the optical properties of the dust's amorphous state in that range. Often standard synthetic dust permittivities or opacities are used, which perfectly fit the near and mid-infrared observed features on the basis of a sum of Lorentz-type resonant profiles. But in the absence of reliable data and knowledge, the FIR/mm behavior is often modeled only by the vanishing long-wavelength wings of those infrared resonances. Some laboratory data on various silicates reveal additional absorption features that are strongly temperature-dependent in the FIR/mm range.

It is clear that only specific processes with characteristic frequencies located in the FIR/mm region are needed to explain non-pure asymptotic features in this spectral domain. Two phenomena satisfying these requirements are known in solid state physics. First, acoustic oscillations with acoustic wavelengths less than grain size have frequencies in the range $3\times 10^{10}{-}3\times 10^{12}$ Hz (Sect. 5.1). These acoustic oscillations are "optically'' active (interact with electromagnetic field) in disordered media (Sect. 4.2). Second, tunneling states in the disordered solids have frequencies $\Delta_0/h$ falling in the frequency range $5\times 10^6{-}5\times 10^{11}$ Hz (Sect. 4.4.2).

These two phenomena are linked to disordered material composing amorphous solids and dirty crystalline compounds. This is also the case for interstellar grains, which dominate the FIR/mm emission of the highly amorphous interstellar dust (Sect. 2.5).

The transverse sound oscillations and tunneling between the potential minima in a disorder lattice includes rotation of groups of atoms, limited and influenced by interactions with their neighboring atoms. In the classification of motion into the electronic, vibrational, and rotational states used in fundamental quantum mechanics, this model refers to rotational motions, when the former dust-emissivity models used only asymptotic behavior of electronic and vibrational processes. In that sense, the model proposed here does not replace the former ones, but complements them by taking into account the physical phenomena that are efficient in the low frequency range.

Both processes are sensitive to the degree and the type of disorder, rather than to the exact chemical composition of the material. As a consequence, the FIR/mm emission profile is governed by a few key parameters that characterize the disorder. There is no contradiction between this and laboratory results of FIR/mm dust absorption. It is important to investigate the variations in the model parameters in laboratory measurements of various duly characterized samples. For this purpose, the model gives new avenues for such simulations, and disorder characterization. It is likely that the FIR/mm properties will be influenced more by the type and degree of disorder than by the chemical composition.

It is important to note that the parameter values adopted in this paper were chosen on purely theoretical grounds or from laboratory results obtained on materials that may not be exactly representative of ISM dust particles. We make no attempt to model the actual astronomical observations, but it will be performed in a companion paper.

It is important that only the parameters considered in that companion paper be used for astronomical purposes. As the type and degree of disorder in interstellar dust should be representative of the amorphous dust formation itself, we hope that this model will open the way to new insights into interstellar medium dust.

   
6 Conclusion

The continuum FIR/mm emission is generally attributed to an interstellar dust component of amorphous grains that are silicate-based and large enough ($\simeq$15-100 nm) to radiate at thermal equilibrium. So far, the conception of FIR/mm dust emissivity has been built on the semi-classical Lorentz or Drude models. They lead in the long wavelength limit to a dust emissivity obeying a frequency power law whose exponent (the spectral index) equals 2 and is usually considered as temperature independent. However a growing volume of FIR/mm observational data induces some doubts about the applicability of such models in that wavelength range.

In this context we present a new model for FIR/mm dust emission based on physical properties of the disordered matter. We consider the interaction between the electromagnetic field and the acoustic oscillations in a disordered charge distribution and a distribution of low-energy two level tunneling states. Both mechanisms apply to the amorphous materials with a large degree of universality, independent of the exact chemical nature of the dust.

The proposed model predicts a disordered induced FIR/mm dust emission that dominates over the weak long-wavelength wings of the infrared resonances. This emission depends strongly on temperature and leads to a millimeter emission enhancement relative to predications from classical models. The FIR/mm emission spectra can no longer be characterized by a constant spectral index over that range. The detailed shape of the emission spectra is governed by a few key parameters that characterize disorder. The exact values of these parameters for interstellar dust should be determined from observations.

Our model is complementary to the classical models, which are reliable in the infrared, and it adds emission processes due to disorder induced features.

Acknowledgements

Part of this work was supported by Egide and RFBR-CNES grant 05-02-19650, and the French national program PCMI. The Ph.D. Thesis work of D. Paradis was supported by the "Fond Social Europeen'' European grant and the Noveltis company.

References

 
• Agladze, N. I., & Sievers, A. J. 1998, , 80, 4209  
• Agladze, N. I., Sievers, A. J., Jones, S. A., Burlitch, J. M., & Beckwith, S. V. W. 1994, Nature, 372, 243  
• Agladze, N. I., Sievers, A. J., Jones, S. A., et al.  
• Allamandola, L. J., Tielens, A., & Barker, J. R. 1985, ApJ, 290, L25  
• Anderson, P. W., Halperin, B. I., & Varma, C. M. 1972, Phil. Mag., 25, 1  
• Andre, P., Ward-Thompson, D., Barsony, M., et al. 1993, ApJ, 406, 122A  
• Andriesse, C. D. 1974, A&A, 37, 257  
• Barreiro, R. B., et al., 2006, MNRAS, in press  
• Bennett, C. L., Hill, R. S., Hinshaw, G., et al. 2003, ApJS, 148, 97  
• Bernard, J.-P., Abergel, A., Ristorcelli, I., et al. 1999, A&A, 347, 640  
• Bösch, M. A. 1978, Phys. Rev., 40, 879  
• Bosomworth, D. R. 1967, Phys. Rev., 157, 709  
• Boudet, N., et al. 2002, in Chemistry as a Diagnostic of Star Formation, ed. C. Curry, & M. Fich (Ottawa, Canada: NRC Press), 257  
• Boudet, N., Mutschke, H., Nayral, C., et al. 2005, ApJ, 633, 272  
• Bouchet, F. R., & Gispert, R. 1999, New A, 4, 443  
• Brucato, J. R., Strazzulla, G., Baratta, G., et al. 2004, A&A, 413, 395  
• Chandler, C. J., Koerner, D. W., Sargent, A. I., et al. 1995, ApJ, 455L, 93  
• Debye, P. 1914, Ann. Physique, 43, 49  
• Desert, F.-X., Boulanger, F., & Puget, J. L. 1990, A&A, 237, 215  
• Draine, B. T., & Anderson, N. 1985, ApJ, 292, 494  
• Draine, B. T., & Lee, H. M. 1984, ApJ, 285, 89  
• Draine, B. T., & Li, A. 2001, ApJ, 551, 807  
• Dupac, X., Giard, M., Bernard, J.-P., et al. 2001, ApJ, 553, 604  
• Dupac, X., Giard, M., Bernard, J.-P., et al. 2002, A&A, 392, 691  
• Dupac, X., Bernard, J.-P., Boudet, N., et al. 2003a, A&A, 404, L11  
• Dupac, X., del Burgo, C., Bernard, J.-P., et al. 2003b, MNRAS, 344, 105  
• Dwek, E., Arendt, R. G., Fixsen, D. J., et al. 1997, ApJ, 475, 565  
• Ewald, P. P. 1922, Naturwiss, 10, 1057  
• Finkbeiner, D. P., Davis, M., & Schlegel, D. J. 1999, ApJ, 524, 867  
• Fitzgerald, S. A., Campbell, J. A., & Sievers, A. J. 1994, Phys. Rev. Lett., 73, 3105  
• Fitzgerald, S. A., Sievers, A. J., & Campbell, J. A. 2001a, J. Phys. Condens. Matter, 13, 2177  
• Fitzgerald, S. A., Campbell, J. A., & Sievers, A. J. 2001b, J. Phys. Condens. Matter, 13, 2095  
• Galliano, F., Madden, S. C., Jones, A. P., et al. 2003, A&A, 407, 159  
• Galliano, F., Madden, S. C., Jones, A. P., Wilson, C. D., & Bernard, J.-P. 2005, A&A, 434, 867  
• Gordon, M. A. 1988, ApJ, 331, 509  
• Gordon, M. A. 1990, ApJ, 352, 636  
• Gromov, V. D., Kardashev, N. S., Kuzmin, L. S., et al. 2002, in Experimental Cosmology at mm-waves, ed. M. De Petris, & M. Gervasi, Amer. Inst. Phys. Conf. Proc. Ser., 616, 205  
• Hadni, A. 1970, in ed. S.S. Mitra, & S. Nudelman (New York: Plenum Press), 561  
• Henning, T., & Mutschke, H. 1997, A&A, 327, 743  
• Heuer, A. 1997, , 78, 4051  
• Hubbard, B. E., et al. 2003, Phys. Rev. B, 67, 144201  
• Hunklinger, S., & Schickfus, M. V. 1981, in Amorphous Solids: Low Temperature properties, ed. W. A. Phillips (Springer-Verlag), 81  
• Jäckle, J. 1972, Z. Physik A, 257, 212  
• Jäger, C. D., Fabian, D., Schrempel, F., et al. 2003, A&A, 401, 57  
• Kemper, F., Vriend, W. J., Tielens, A. G. G. M., et al. 2004, ApJ, 609, 826  
• Koike, C., Hasegawa, H., & Manabe, A. 1980, Ap&SS, 67, 495  
• Kühn, R. 2003, Europhys. Lett., 62, 313  
• Lamarre, J.-M. 1994, Infrared Phys. and Techn., 35, 277  
• Lamarre, J. M., Puget, J. L., Bouchet, F., et al. 2003, New Astron. Rev., 47, 1017  
• Landau, L. D., & Lifshits, E. M. 1982, Course of Theoretical Physics, Vol. VIII (Moscow: Nauka)  
• Leger, A., & Puget, J. L. 1984, A&A, 137, L5  
• Li, A., & Greenberg, J. M. 1997, A&A, 323, 566  
• Mathis, J. S., Rumpl, W., & Nordsieck, K. H. 1977, ApJ, 217, 425  
• Mennella, V., Brucato, J. R., Colangeli, L., et al. 1998, ApJ, 496, 1058  
• Mon, K. K., Chabal, Y. J., & Sievers A. J. 1975, , 35, 1352  
• Naselsky, P. D., Chiang, L.-Y., Novikov, I. D., et al. 2005, Int. J. Mod. Phys. D, 14, 1273  
• Oldham, P., et al. 1994, A&A, 284, 5590  
• Pajot, F. 2006, A&A, 447, 769  
• Phillips, W. 1972, J. Low Temp. Phys., 11, 757  
• Phillips, W. 1987, Rep. Prog. Phys., 50, 1657  
• Rast, H. E., Caspers, H. H., & Miller, S. A. 1968, Phys. Rev., 171, 1051  
• Reach, W. T., Dwek, E., Fixsen, D. J., et al. 1995, ApJ, 451, 188  
• Ristorcelli, I., Serra, G., Lamarre, J. M., et al. 1998, ApJ, 496, 267  
• Rubens, H., & Hertz, G. 1912, Berlin. Ber., 14, 256  
• Schlömann, E. 1964, Phys. Rev. A, 135, 413  
• Schwartz, P. 1982, ApJ, 252, 589  
• Serra, G., Giard, M., Bouchou, F., et al. 2002, Adv. Space Res., 30, 1297  
• Sievers, A. J., Tu, J., Agladze, N., FitzGerald, S. A., & Campbell, J. A. 1998, Physica B, 244, 159  
• Stepnik, B., et al. 2003, Clim. Conf., 187, Cross-Calibration of PRONAOS & ISO, in The Calibration Legacy of the ISO Mission, ESA SP-481  
• Stepnik, B., Abergel, A., Bernard, J.-P., et al. 2003, A&A, 398, 551  
• Strom, U., & Taylor, P. C. 1977, Phys. Rev., 16, 5512  
• von Schickfus, M., & Hunklinger, S. 1976, J. Phys. C, Solid State Phys., 9, L439  
• von Schickfus, M., Hunklinger, S., & Pich, L. 1975, , 35, 876  
• Tamm, I. E. 1930, Z. Physik, 60, 345  
• Tegmark, M., Eisenstein, D. J., Hu, W., & de Oliveira-Costa, A. 2000, ApJ, 530, 133  
• Vinogradov, V. S. 1960, Fiz. Tverd. Tela, 2, 2622 (English transl. 1961, Sov. Phys. - Solid St. 2, 2338)  
• Walker, C., Adams, F. C., Lada, C. J., et al. 1990, ApJ, 349, 515  
• Weiland, J. L., Blitz, L., Dwek, E., et al. 1986, ApJ, 306, L101  
• Woody, D., Scott, S. L., Scoville, N. Z., et al. 1989, ApJ, 337, 41  
• Wright, E., Meyer, S. S., Bennett, C. L., et al. 1992, ApJ, 396, L13

 

  
7 Online Material

   
Appendix A: TLS optical density of states

The ODOS spectrum for continuous distribution of DOS is calculated by integration of Eq. (68) over the energy plitting $\Delta _0$. We use the following polynomial form of the distribution function $p^{}_{\rm\tiny G}(x)$, with $x=\Delta_0/\hbar\omega_{\rm m}$

$$p^{}_{\rm\tiny G}(x) = \left\{ \begin{array}{ll} \frac{7}... ...<1,\\ % \hspace{15mm} 0 &{\rm for} \ x>1. \end{array} \right.$$ (A.1)

Equation (A.1) is a least square polynomial fit to Eq. (70), removing the discontinuity at x=1. A similar result is obtained using the first Fourier harmonics of Eq. (70). The calculated ODOS spectrum is given by $G(\omega)=\overline{P}\mu_{\rm b}^2 g(x)$ with $x=\omega/\omega_{\rm m}$ and

$$g(x) = \left\{ \begin{array}{ll} 1+x^2 g_1(x) & {\rm for}... ...)-g_2(x)\sqrt{1-x^{-2}} & {\rm for} \ x>1, \end{array}\right.$$ (A.2)

with

  

g₁(x) = 4(5-6x²)/15, (A.3)

g₂(x) = 8(1-x²)(2+3x²)/15. (A.4)

In the low-frequency approximation ( $\omega\ll\omega_{\rm m}$) we have $p^{}_{\rm\tiny G}\approx 1$ and the ODOS function $G(\omega)$ does not depend on $\omega$, in agreement with previous results of the TLS theory.

   
Appendix B: Evaluation of a simplified equation for tunneling relaxation calculation

Using simplified single-integral expression of Bösch (Eq. (83)) instead of the full double-integral equation used here (Eq. (81)) is based on that assumption that input from asymmetric states is negligible. This approach gives the correct result for resonant absorption (Eq. (68)) where the integrand term $(1-\Delta_0/\hbar\omega)^{-1/2}$ amplifies the weight of symmetric states $\Delta_0\approx\hbar\omega$. Another situation takes place for relaxation absorption (Eq. (81)) where the integrand term $(1-\tau_1/\tau)^{+1/2}$ emphasize the weight of asymmetric states $\tau\gg\tau_1$. Direct calculations give the discrepancy of results more than a hundred times. Without detailed calculations, a comparison can be simply made by considering the terms that, within the integral over E, differ between (Eq. (83)) and (Eq. (81)). In case of (Eq. (83)) it is $f_1(\omega,\tau) = \omega^2\tau_1/(1+\omega^2\tau_1^2).$
In the case of (Eq. (81)) it is

$$f_2(\omega,\tau)=\int^{\infty}_{\tau_1}\sqrt{1-\frac{\tau_1}{\tau}} \frac{\omega^2{\rm d}\tau}{(1+\omega^2\tau^2)}\cdot$$ (B.1)

It is not difficult to define the asymptotics of both functions:

$$\begin{eqnarray*}f_1(\omega,\tau_1) &\approx& \left\{ \begin{array}{ll} % \omega... ...tau_1^{-1} & {\rm for} \ \omega\tau_1 \gg 1. \end{array} \right. \end{eqnarray*}$$

The asymptotic behaviors of f₁ and f₂ are different, and the simplified equation is also not accurate for a wide range of parameters values.

   
Appendix C: Calculation of the integral for the tunneling relaxation spectrum

The function $F_{\rm phon}(\omega,T)$ (Eq. (83)) can be simplified into a function of a single parameter $F_{\rm phon}(p)$ with $p= a\omega/(2k^{}_{\rm\tiny B}T)^3$ through a change in variables to the dimensionless variables $x=\arctan(\omega\tau)$ and $y=\tanh(E/2k^{}_{\rm\tiny B} T)$. Equation (83) then writes as

$$F_{\rm phon}(p) = \int^1_0\int^{\pi/2}_{x_1} \sqrt{ 1-\frac{\tan x_1}{\tan x} }\ {\rm d}x {\rm d}y,$$ (C.1)

where

$$x_1 = \arctan( p\ y\textrm{ arcth}^{-3} y).$$ (C.2)

The function $F_{\rm phon}(p)$ can be further approximated using a fit of the form,

$$F_{\rm phon}(p)=F_{\rm phon}(p_i)\left(\frac{p}{p_i}\right)^{\beta_{2~ i}},\ p_i<p<p_{i+1}.$$ (C.3)

Precalculated values of function $F_{\rm phon}$ and coefficients $\beta _2$entering the interpolated form are given in Table C.1.

   

Table C.1: Values of $F_{\rm phon}(p_i)$ and $\beta _2(p_i)$ entering the interpolated form of $F_{\rm phon}(p)$ in Eq. (C.3).

pi	$F_{\rm phon}(p_i)$	$\beta _2(p_i)$
0.001	1.4696	-0.0321
0.002	1.4313	-0.0446
0.005	1.3604	-0.0678
0.01	1.2875	-0.0924
0.02	1.1948	-0.1248
0.05	1.0393	-0.1826
0.1	0.8984	-0.2397
0.2	0.7433	-0.3090
0.5	0.5334	-0.4185
1	0.3866	-0.5114
2	0.2623	-0.6084
5	0.1417	-0.7343
10	0.0826	-0.8208
20	0.0455	-0.8950
50	0.0194	-0.9659
100	0.0098	-0.9897
200	0.0049	-0.9969
500	0.0020	-1.0000
1000	0.0010	...

   
Appendix D: Calculations of the integral for hopping relaxation

Following the TLS formalism, the hopping relaxation spectrum can be computed as the integral of Eq. (77) over a distribution of TLS and barrier heights, considering that, in contrast to tunneling relaxation, the time constant $\tau$ depends on barrier height V according to Eq. (87) and does not depend directly on TLS parameters $\Delta _0$ and $\Delta$:

 

$$\chi_0''(\omega)=\frac{\mu_{\rm b}^2}{3k^{}_{\rm\tiny B}T} \int_0... ... \frac{\omega\tau}{1+(\omega\tau)^2}\ {\rm d}\Delta_0\ {\rm d}\Delta\ {\rm d}V.$$ (D.1)

For calculating this expression, the convention of Hunklinger & Schickfus (1981) is used, which allows reduction of the triple integral to a single integration over distribution $P(V){\rm d}V$, where P(V) is taken to be Gaussian. This convention suggests the independence of the distribution $P(\Delta_0,\Delta)$and P(V), which permits the integrals to be separated and transforms Eq. (D.1) into Eq. (87), in a form similar to that of Fitzgerald et al. (2001b). In this case P(V) is an independent Gauss probability distribution following the standard normalization

$$\int_0^\infty P(V)\ {\rm d}V=1,$$

which defines the normalization coefficient C_V in Eq. (79) as

$$C_V=\frac{1}{V_0\sqrt{\pi}} \Big[\frac{1}{2}{\rm Erf}\Big(\frac{V_{\rm m}-V_{\min}}{V_0}\Big)+\frac{1}{2}\Big]^{-1},$$ (D.2)

where the term in brackets is usually about unity.

The coefficient $B_{\rm hop}(T)$ in Eq. (87) is equal to

 

$$B_{\rm hop}(T) \frac{4\pi}{c\sqrt{\epsilon'}}\frac{(\epsilon'+2)^2}{9} \frac{2}{3}\overline{P}\mu_{\rm b}^2$$ =

$$\times \int_{\Delta_0^{\min}}^\infty\int_{\Delta_0^{\min}}^E \tex... ...^2}{E^2}} \frac{{\rm d}\Delta_0}{\Delta_0} \frac{{\rm d}E}{k^{}_{\rm\tiny B}T}=$$

$$\frac{4\pi}{c\sqrt{\epsilon'}}\frac{(\epsilon'+2)^2}{9}\frac{2}{3... ...P}\mu_{\rm b}^2 \left(\ln\frac{k^{}_{\rm\tiny B}T}{\Delta_0^{\min}}+C_1\right),$$ = (D.3)

where the constant C₁ is given by

$$C_1=\ln 4-1+\int_0^1\ln{\rm arcth}\ x\ {\rm d}x\ =\ -0.441.$$ (D.4)

The integral (D.3) over d$\Delta _0$ has been determined precisely. Integrating over dE was performed assuming that terms of order $\Delta_0^{\min}/k^{}_{\rm\tiny B}T$ and higher orders could be neglected, taking into account that $\Delta_0^{\min}/k^{}_{\rm\tiny B}T\ll 1$.