Olivier Delestre
email: delestre@unice.fr

Stéphane Cordier

Frédéric Darboux

Francois James

A limitation of the hydrostatic reconstruction technique for Shallow Water equations

ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Derived from Navier-Stokes equations, Shallow Water equations describe the water flow properties as follows

∂ t h + ∂ x (hu) = 0, ∂ t (hu) + ∂ x hu 2 + gh 2 /2 = -gh∂ x z, (1)
where the unknowns are the water height h(t, x), the velocity u(t, x), and the topography z(x) is a given function. In what follows, we note the discharge q = hu, the vector of conservative variables U = (h hu) t and the flux F (U ) = (hu gh 2 /2 + hu 2 ) t . The steady state of a lake at rest, or a puddle, (h+ z = Cst and u = q = 0) is a particular solution to [START_REF] Audusse | A fast and stable well-balanced scheme with hydrostatic reconstruction for shallow water flows[END_REF]. Since [START_REF] Bermudez | Upwind methods for hyperbolic conservation laws with source terms[END_REF], it is well known that the topography source term needs a special treatment in order to preserve at least this equilibrium. Such schemes are said to be well-balanced (since [START_REF] Greenberg | A well-balanced scheme for the numerical processing of source terms in hyperbolic equation[END_REF]).

In the following, we present briefly the so-called hydrostatic reconstruction method, which permits, when coupled to a positive numerical flux, to obtain a family of well-balanced schemes that can preserve the water height nonnegativity and deal with dry zones. We show that this method, presented in [START_REF] Audusse | A fast and stable well-balanced scheme with hydrostatic reconstruction for shallow water flows[END_REF][START_REF] Bouchut | Nonlinear stability of finite volume methods for hyperbolic conservation laws, and well-balanced schemes for sources[END_REF] and widely used, fails for some combinations of slope, mesh size and water height. We give the criteria that ensures the accuracy of the results.

The numerical method

The hydrostatic reconstruction follows the general principle of reconstruction methods. We start from a first order finite volume scheme for the homogeneous form of system (1): choosing a numerical flux F(U L , U R ) (e.g. Rusanov, HLL, VFRoe-ncv, kinetic), a finite volume scheme takes the general form

U * i = U n i - ∆t ∆x [F(U i , U i+1 ) -F(U i-1 , U i )] , (2)
where ∆t is the time step and ∆x the space step. Now the idea is to modify this scheme by applying the flux to reconstructed variables. Reconstruction can be used to get higher order schemes, in that case higher order in time is achieved through TVD-Runge-Kutta methods. The aim of the hydrostatic reconstruction, which is described in the next section, is to be well-balanced, in the sense that it is designed to preserve at least steady states at rest (u = 0). When directly applied on the initial scheme, it leads to order one methods, while coupling it with high order accuracy reconstruction increases the order.

The hydrostatic reconstruction

We describe now the implementation of this method for high order accuracy. The first step is to perform a high order reconstruction (MUSCL, ENO, UNO, WENO). To deal properly with the topography source term ∂ x z, this reconstruction has to be performed on u, h and h + z, see [START_REF] Bouchut | Nonlinear stability of finite volume methods for hyperbolic conservation laws, and well-balanced schemes for sources[END_REF]. This gives us the reconstructed values (U -, z -) and (U + , z + ), on which we apply the hydrostatic reconstruction [START_REF] Audusse | A fast and stable well-balanced scheme with hydrostatic reconstruction for shallow water flows[END_REF][START_REF] Bouchut | Nonlinear stability of finite volume methods for hyperbolic conservation laws, and well-balanced schemes for sources[END_REF] on the water height, namely

       h i+1/2L = max(h i+1/2-+ z i+1/2--max(z i+1/2-, z i+1/2+ ), 0), U i+1/2L = (h i+1/2L , h i+1/2L u i+1/2-), h i+1/2R = max(h i+1/2+ + z i+1/2+ -max(z i+1/2-, z i+1/2+ ), 0), U i+1/2R = (h i+1/2R , h i+1/2R u i+1/2+ ). (3)
With the hydrostatic reconstruction, the finite volume scheme ( 2) is modified as follows

U * i = U n i -∆tΦ(U n ) = U n i - ∆t ∆x F n i+1/2L -F n i-1/2R -F c n i , (4)
where

F n i+1/2L = F n i+1/2 + S n i+1/2L , F n i-1/2R = F n i-1/2 + S n i-1/2R (5)
are left (respectively right) modifications of the initial numerical flux for the homogeneous problem. In this formula the flux is now applied to reconstructed variables:

F n i+1/2 = F(U n i+1/2L , U n i+1/2R
), and we have introduced

S n i+1/2L = 0 g 2 (h 2 i+1/2--h 2 i+1/2L ) , S n i-1/2R = 0 g 2 (h 2 i-1/2+ -h 2 i-1/2R ) . (6)
Finally, a centered source term has to be added to preserve consistency and well-balancing (see [START_REF] Audusse | A fast and stable well-balanced scheme with hydrostatic reconstruction for shallow water flows[END_REF][START_REF] Bouchut | Nonlinear stability of finite volume methods for hyperbolic conservation laws, and well-balanced schemes for sources[END_REF]):

F c i = 0 -g h i-1/2+ + h i+1/2- 2 (z i+1/2--z i-1/2+ ) . (7)
Formula ( 7) preserves the second order accuracy, but has to be modified for higher order approximations.

The hydrostatic reconstruction rewritten

We define ∆z i+1/2 = z i+1/2+ -z i+1/2-. Once the space step and the high order reconstruction chosen, this is a fixed sequence. Now, the reconstructed variables (3) write

h i+1/2L = max(h i+1/2--max(0, ∆z i+1/2 ), 0), h i+1/2R = max(h i+1/2+ + min(0, ∆z i+1/2 ), 0). (8)
With [START_REF] Delestre | SWASHES: a library of Shallow Water Analytic Solutions for Hydraulic and Environmental Studies[END_REF], the defect of the hydrostatic reconstruction becomes apparent. To fix the ideas, suppose the local slope is positive, hence ∆z i+1/2 > 0, as in Figure 1.

Then, for all ∆z i+1/2 such that ∆z i+1/2 ≥ h i+1/2-≥ 0, the reconstruction gives h i+1/2L = 0, while h i+1/2R remains unchanged. In that case, the reconstruction

∆z i+1/2 z i+1/2+ h i+1/2+ + z i+1/2+ h i+1/2R + z i+1/2+ z i+1/2- h i+1/2-+ z i+1/2- h i+1/2- i i + 1 i + 1/2 h i+1/2L + z i+1/2-
Figure 1: Default of the hydrostatic reconstruction. Left: threshold non activated (limit case). Right: threshold activated prevents an unphysical negative value for h i+1/2L , the counterpart of this being an underestimated difference h i+1/2L -h i+1/2R . Therefore, there is a lack in the numerical flux computed from the modified Riemann problem, which gives an underestimated velocity and consequently an overestimated height. This can be interpreted in terms of reconstructed slope as well, which is underestimated.

In the general case, we can write the following local criterion of "non validity" for the hydrostatic reconstruction.

Proposition 2.1 For a fixed discretization, if for some i 0 ≤ i ≤ i 1 one has

∆z i+1/2 ≥ h i+1/2-≥ 0, or -∆z i+1/2 ≥ h i+1/2+ ≥ 0, (9)
then the hydrostatic reconstruction will overestimate (resp. underestimate) the water height (resp. velocity).

Notice that from a theoretical viewpoint, this is not limiting. Indeed, since this class of schemes is consistent with the system of partial differential equations (see [START_REF] Bouchut | Nonlinear stability of finite volume methods for hyperbolic conservation laws, and well-balanced schemes for sources[END_REF]), the problem disappears when refining the discretization. However, it has to be taken into account for practical computations, with a fixed discretization. It is particularly apparent for order 1 schemes, but also remains present for order 2.

Numerical illustration

For numerical illustration purpose, we introduce an explicit solution which consists in a supercritical steady flow on an inclined plane with constant slope ∂ x z = α (it is referenced in the SWASHES library [START_REF] Delestre | SWASHES: a library of Shallow Water Analytic Solutions for Hydraulic and Environmental Studies[END_REF]). Steady states are solutions to

q(x, t) = q 0 = Cst, ∂ x hu 2 + g h 2 2 = -gh∂ x z,
and the height profile h(x) must be a solution to Bernoulli's law rewritten as a third order equation in h:

h 3 + h 2 αx - q 2 0 2gh 2 0 -h 0 + q 2 0 2g = 0, ( 10
)
where h 0 and q 0 = h 0 u 0 are the height and discharge at x = 0, which completely determine the profile since the flow is supercritical. A careful study of the roots of the polynomial shows that the supercritical height profile is decreasing in x: h(x) ≤ h 0 for all x ≥ 0.

The numerical strategy we choose consists in the HLL flux and a modified MUSCL reconstruction. In [START_REF] Delestre | Rain water overland flow on agricultural fields simulation/Simulation du ruissellement d'eau de pluie sur des surfaces agricoles[END_REF], this combination of flux and linear reconstruction has shown to be the best compromise between accuracy, stability and CPU time cost. We refer to [START_REF] Bouchut | Nonlinear stability of finite volume methods for hyperbolic conservation laws, and well-balanced schemes for sources[END_REF][START_REF] Delestre | Rain water overland flow on agricultural fields simulation/Simulation du ruissellement d'eau de pluie sur des surfaces agricoles[END_REF] for a presentation of the HLL flux. The MUSCL reconstruction of a scalar function s ∈ R is

s i-1/2+ = s i - ∆x 2 Ds i , s i+1/2-= s i + ∆x 2 Ds i , (11)
where the "minmod" operator D is given by

Ds i = minmod s i -s i-1 ∆x , s i+1 -s i ∆x , minmod(x, y) =      min(x, y) if x, y ≥ 0, max(x, y) if x, y ≤ 0, 0 else. ( 12
) In order to keep the discharge conservation, the reconstruction of the velocity has to be modified as

u i-1/2+ = u i - h i+1/2- h i ∆x 2 Du i u i+1/2-= u i + h i-1/2+ h i ∆x 2 Du i .
Notice that if we take Ds i = 0 in [START_REF] Popinet | Quadtree-adaptative tsunami modelling[END_REF], then z i+1/2-= z i+1/2-= z i so that the centered term ( 7) disappears, and we recover the first order scheme in space.

Second order in time is achieved through a classical Heun predictor-corrector method.

We turn now to the specific data for simulations. The domain length is L = 10 m, and we choose initial data h 0 = h(x = 0, t) = 0.02 m, q 0 = q(x = 0, t) = 0.01 m 2 /s, which is indeed supercritical. All simulations are performed with a space step ∆x = 0.1 m, the time step is fixed in order to satisfy the CFL condition. The analytical solution is computed with 7 negative slopes α = 5%, 13%, 14%, 15%, 16%, 17% and 18%, both at first and second order. Numerical results are compared to the analytical results in Figure 2, where a part of the domain is displayed (x between 1 m and 3 m).

With this set of data, a domain of admissible slopes can be estimated for the order 1 hydrostatic reconstruction. Indeed since ∆z j+1/2 = |α|∆x, inserting a characteristic height of the flow h * in (9) gives a bound for the slope. One can use h * = h 0 , the incoming height. Since the height profile is decreasing, the whole profile will be wrong if h 0 ≤ |α|∆x, that is here |α| ≥ 20%. But actually, since the height decreases quite rapidly, a more accurate estimate is obtained for x = 1.5 m with h * = 0.005 m (see Fig. 2), which leads to slopes above 5%.

We observe on Figure 2(Left) that with the first order scheme, the effect of the hydrostatic reconstruction is so important that all curves for slopes between 13% and 18% are superposed. In that case, the apparent result is the simulation of a single slope, namely 13%. For a 5% slope, we still observe a slightly overestimated water height, as anticipated since 5% is a limit case as observed above. With the second order, the water heights are still overestimated, but in a very slighter way, and the different curves are no longer identical (Figure 2(Right)).

Conclusions

The hydrostatic reconstruction may fail for certain combinations of water height, slope and mesh size, namely in regions where (9) holds. The defaults are particularly apparent for order 1 schemes, leading to wrongly estimated slopes, but still remain at order 2, with some overestimated water heights. We emphasize that the problem disappears when refining the mesh, but has to be taken into account for a given discretization. The generalization of the hydrostatic reconstruction proposed in [START_REF] Castro | Well-balanced numerical schemes based on a generalized hydrostatic reconstruction technique[END_REF] does not exhibit the limitation discussed here, even for first order schemes, but positivity is not ensured. Other schemes involving threshold values (e.g. [START_REF] Bouchut | A subsonic-well-balanced reconstruction scheme for shallow water flows[END_REF][START_REF] Liang | Numerical resolution of well-balanced shallow water equations with complex source terms[END_REF]) very likely encounter the same kind of problem. Alternatively, the scheme recently introduced in [START_REF] Berthon | Efficient well-balanced hydrostatic upwind schemes for shallow-water equations[END_REF] preserves the water height positivity and does not suffer from this problem. Notice finally that criterion (9) may be of some utility for adaptive mesh schemes, such as the ones used in Gerris [START_REF] Popinet | Quadtree-adaptative tsunami modelling[END_REF].

Acknowledgments

This study is part of the ANR project METHODE #ANR-07-BLAN-0232 granted by the French National Agency for Research. The authors wishes to thank M.-O. Bristeau, Ch. Berthon and M.-H. Le for fruitful discussions.


Luca Dore
email: luca.dore@unibo.it

Luca Bizzocchi

Claudio Degli Esposti

Filippo Tamassia

Fine and Hyperfine Structure of the N = 1 ← 0 Transition of ND(X 3 Σ -) in Vibrational Excited States

Keywords: rotational transitions, fine structure, hyperfine structure, free radical, spectroscopic parameters rotational transitions, fine structure, hyperfine structure, free radical, glow discharge, spectroscopic parameters

The deuterated radical ND was produced in a DC discharge cell cooled at liquid nitrogen temperature. The discharge proved to be vibrationally hot, therefore the transient species could be detected in its vibrational excited states up to v = 6. By scanning in the 431-531 GHz frequency region, several fine-structure components of the N = 1 ← 0 transition in vibrational excited states were observed, each of them showing a complex hyperfine structure.

A global analysis, including the measured frequencies and previous submillimeter-wave and infrared data, allowed an accurate determination of the equilibrium spectroscopic parameters of the ND radical including fine and hyperfine constants.

A very precise determination of the equilibrium bond length re was obtained. This value is not consistent with the value reported in literature from NH data. This incongruity was discussed in terms of the breakdown of the Born-Oppenheimer approximation.

In view of the recent detection of ND in a solar-mass protostar (A. Bacmann et al., Astron. Astrophys. 521, L42 (2010)), an extended spectroscopic characterization of this deuterated isotopologue of the NH species may prove useful, considering the large deuterium enhancement observed in molecular clouds.

Very recently, it has been reported the first astronomical detection of the radical ND towards a young solar-mass protostar [1]; this is one among the wealth of results which are being obtained by the high resolution heterodyne instrument HIFI on board the Herschel Space Observatory. The joint observation of NH hyperfine transitions allowed to estimate deuterium fractionation of imidogen, which results to be very high, with a ratio [ND]/[NH] between 30 and 100%.

The laboratory detection of the pure rotational spectrum of the radical ND in its ground vibronic state X 3 Σ -was achieved by Saito and Goto [2], who observed the three fine-structure components of the N = 1 ← 0 rotational transition and resolved their hyperfine structure. Five years later, Takano et al. [3] reported the observation of the N = 2 ← 1 transition, with its complex fine and hyperfine structure: thanks to the determined molecular parameters, they were able to list accurate frequency predictions for ND rotational transition up to 4 THz. The first seven vibrational levels of ND(X 3 Σ -) were investigated by Ram and Bernath [4] 1 by observing six vibration-rotation bands at 0.01 cm -1 resolution using a Fourier transform spectrometer. For the low N lines the characteristic triplet pattern due to the three spin components J = N + 1, N , and N -1 of each N level was resolved, while no hyperfine structure was resolved.

The present work reports measurements of a number of fine-structure components of the N = 1 ← 0 rotational transition of ND in a number of vibrational excited states; ∆J = 0, ±1 components for v = 1, 2 and the ∆J = +1 component for v = 3, 4, 6 2 2 . Each fine-structure line was split in several hyperfine components, so that a total of 69 measured line frequencies are reported.

All the line frequencies reported in references [2][3][4] have been considered in a global fit with the present measurements in order to obtain an extensive set of equilibrium spectroscopic parameters.

Experimental details

Measurements were carried out in the frequency range of 431-531 GHz, employing phase-locked Gunn oscillators (Radiometer Physics GmbH, J. E. Carlstrom Co) as primary radiation source working in the range of 75-115 GHz; power at higher frequencies was obtained using harmonic multiplication. Two phase-lock loops allowed the stabilization of the Gunn oscillator with respect to a frequency synthesizer, which was driven by a 5-MHz rubidium frequency standard. The frequency modulation of the radiation was obtained by a sine-wave at 16.66 kHz modulating the reference signal of the wide-band Gunn synchronizer; the signal, detected by a liquid-helium-cooled InSb hot electron bolometer (QMC Instr. Ltd. type QFI/2), was demodulated at 2-f by a lock-in amplifier.

ND radical was produced in a DC discharge with a current ranging between 20 and 60 mA by flowing N 2 (ca. 5 mTorr) and D 2 (1-2 mTorr) in Argon buffer gas for a total pressure of 15-30 mTorr (2-4 Pa). The Pyrex cell, 3.25 m long and 5 cm in diameter, was equipped with two cylindrical hollow electrodes 25 cm in length at either end, and was wound with a plastic pipe for liquid nitrogen circulation to cool down to ca. 90 K.

We could detect transitions in excited vibrational states as high as v = 6, that is the discharge was vibrationally hot. Figure 1 shows a portion of the hyperfine spectrum of the ∆J = +1 fine component for the v = 0, 1, and 2 vibrational states: from a series of such spectra, recorded with a controlled power of the radiation source, we could roughly estimate a vibrational temperature of 2000-3000 K. As for the mechanism of the excitation, we could observe that discharging ND 3 or a mixture NH 3 /D 2 led to a better production of ND in the ground vibrational state, while spectra in the first excited states resulted barely observable. Also, using N 2 as buffer gas instead of Ar, the vibrational temperature was enhanced. Therefore we can argue that energy transfer from excited molecular nitrogen is responsible of the vibrational excitation.

Results and Discussion

ND is a diatomic free radical with an X 3 Σ -ground electronic state. From a theoretical point of view, this molecule is well described by a Hund's case (b) basis set. Consequently, the effective Hamiltonian including rotation, fine structure and hyperfine structure terms can be written as follows [5]:

H ef f = H rot + H f s + H hf s (N ) + H hf s (D). ( 1
)
where

H rot = BN 2 -DN 4 + HN 6 + LN 8 , (2a)
H f s = 2 3 λ(3S 2 z -S 2 ) + γN •S + 1 2 γ D [N •S, N 2 ] + + 1 2 γ H [N •S, N 4 ] + , (2b)
H hf s (N ) = b F (N )I N •S + c(N ) 3 (3I N z S z -I N •S) + eq 0 Q(N ) 4I N (2I N -1) (3I 2 N z -I 2 N ) + C(N )I N •N , (2c)
H hf s (D) = b F (D)I D •S + c(D) 3 (3I Dz S z -I D •S). (2d)
According to Hund's case (b) formalism, N is the rotational angular momentum and J = N + S, where S is the electron spin vector. Both nuclei have a non-vanishing spin angular momentum with quantum numbers I N = I D = 1. The hyperfine terms in the Hamiltonian have been split in two terms, relative to the nitrogen and deuterium nuclei, for clarity. While the magnetic hyperfine terms due to the interaction with S could be determined for both nuclei, the electric quadrupole and nuclear spin-rotation interaction terms were included only for nitrogen. No centrifugal distortion corrections of the hyperfine parameters were considered in the Hamiltonian.

The global fit was performed including our data, the submillimeter-wave data of Ref. [2], the far infrared data of Ref. [3] and the infrared data of Ref. [4]. Because of the great number of transitions involving vibrational excited states, it was necessary to expand the pure vibrational energy to the fourth anharmonicity correction:

G(v) = ω e (v + 1 2 ) -ω e x e (v + 1 2 ) 2 + ω e y e (v + 1 2 ) 3 + ω e z e (v + 1 2 ) 4 + ω e w e (v + 1 2 ) 5 . (3)
Also, the rotational, fine-structure and hyperfine structure terms were expanded in order to include the appropriate vibrational contributions. For each parameter, the equilibrium value and some vibrational corrections, depending on the case, were determined in the fit. They are reported in the following equations.

B = B e -α e (v + 1 2 ) + γ e (v + 1 2 ) 2 + δ e (v + 1 2 ) 3 + e (v + 1 2 ) 4 + ζ e (v + 1 2 ) 5 , (4) D = D e -β D e (v + 1 2 ) + η D e (v + 1 2 ) 2 + θ D e (v + 1 2 ) 3 + κ D e (v + 1 2 ) 4 , (5)
H = H e -β H e (v + 1 2 ) + η H e (v + 1 2 ) 2 , (6)
L = L e -β L e (v + 1 2 ), (7)
λ = λ e + λ v (v + 1 2 ) + λ v 2 (v + 1 2 ) 2 + λ v 3 (v + 1 2 ) 3 , (8)
γ = γ e + γ v (v + 1 2 ) + γ v 2 (v + 1 2 ) 2 , (9)
γ D = γ D,e + γ D,v (v + 1 2 ) + γ D,v 2 (v + 1 2 ) 2 , ( 10
)
γ H = γ H,e , (11)
b F = b F,e + b F,v (v + 1 2 ) + b F,v 2 (v + 1 2 ) 2 , (12)
c = c e + c v (v + 1 2 ) + c v 2 (v + 1 2 ) 2 , (13)
C = C e + C v (v + 1 2 ), ( 14
)
eq 0 Q = eq 0 Q e + eq 0 Q v (v + 1 2 ). ( 15
)
Table 1 reports the values of the parameters derived from the global analysis, where a total of 677 transition frequencies 1 were fitted to the Hamiltonian of equations ( 1) and ( 2) in a weighted-least-squares procedure, as implemented in Pickett's SPCAT/SPFIT program suite [6]. The line frequencies measured in the present work are listed in Table 2 along with the residuals from the fit. The rms error of residuals resulted to be 51 kHz for the pure rotational lines, and 0.0025 cm -1 for the IR lines; the fit standard deviation is 0.915 2 .

The fit is overall very satisfactory, and the most important parameters have been determined with great accuracy. The derived parameters are consistent with the previous ones, however a direct comparison is not possible because here are reported equilibrium values obtained with more terms in the vibrational expansions. The vibrational dependence of the hyperfine terms has been obtained for the first time.

The magnetic hyperfine structure parameters derived for a free radical molecule provide direct information on the electronic distribution and on the composition of the electronic wavefunction [7]. In the spectroscopic approximation, namely where only the unpaired electrons are considered, the Fermi contact term b F is directly linked to the spin density at the nucleus i Ψ(0 i ) 2 . Also, the anisotropic dipolar interaction term c is related to the angular distribution of the electronic wavefunc- 1 Frequencies and uncertainties of the ground-state rotational transitions included in the fit are listed in Table 1 of reference [3]; the IR lines of the vibration-rotation bands from 1-0 up to 6-5 reported in Table 1 of reference [4] were fitted with an uncertainty of 0.003 cm -1 . 2 rmsres = ( residual 2 )/(N observations), σ fit = ( (residual/uncertainty) 2 )/(degrees of freedom).

tion i (3 cos 2 θ i -1)/r 3 i . Namely: b F /M Hz = 10 -6 h 2 3 g e g N µ B µ N µ 0 i Ψ(0 i ) 2 , (16)
c/M Hz = 10 -6 h

3 8π g e g N µ B µ N µ 0 i (3 cos 2 θ i -1) r 3 i , (17)
where g e is the g-value for the free electron, g N is the nuclear g-factor, µ B is the Bohr magneton, µ N is the nuclear magneton and µ 0 the permeability of vacuum.

The sums are over the unpaired electrons. The values for the electron densities derived from the fitted equilibrium parameters b F e and c e are summarized in Table 3. In the same table the value of i 1/r 3 i relative to the nitrogen nucleus, calculated assuming π|3 cos 2 θ -1|π =-2/5, is also reported. This assumption is correct for atomic p π orbitals and a good approximation for molecular π orbitals [8]. The same approximation is not expected to be valid for the deuterium atom. A comparison with the same values for NH, although not at equilibrium, is shown. Kristiansen and Veseth [9] calculated the equilibrium values for the hyperfine constants of NH. Since the equilibrium parameters are independent on the isotopic substitution, these values can be compared to our experimental constants, provided the correction for the different nuclear magnetic moments and spins is applied. In fact, the value of g N = µ I /(µ N I) depends on the nucleus and is equal to 0.857438230 for the deuterium nucleus and 5.585694772 for hydrogen. Therefore, we can calculate the equilibrium values of b F e (H) and c e (H) from our experimental data for deuterium by applying the multiplicative factor (5.585694772/0.857438230). The comparison with Kristiansen's values, reported in Table 4, shows an excellent agreement between computed and experimental parameters.

An equilibrium bond distance of 1.03666599(70) Å was derived from the fitted equilibrium value B e . The associated error was calculated with the standard linear error propagation equations and considering the errors on B e , the atomic masses, the Planck constant and the amu/kg conversion factor as reported in Ref. [START_REF] Mills | Quantities, Units and Symbols in Physical Chemistry[END_REF]. The best determination of r e available in literature from ND data is r e =1.036651(9) [4]. Our result is rather consistent with the previous one but an improvement of more than one order of magnitude in the precision of this parameter has been achieved in this work.

It is interesting to point out that the equilibrium distance derived from NH data is different from the one obtained from ND data well beyond the parameter's standard deviation. Such a discrepancy is reproduced by several experimental determinations of r e for NH: 1.0371860 (19) [START_REF] Ram | [END_REF], 1.03756(6) [12], 1.03722(2) [13] Å. We believe that this apparent incongruity arises from the breakdown of the Born-Oppenheimer (BO) approximation. According to Watson [14], the experimentally determined equilibrium bond length in a diatomic molecule AB can be related to the isotopically invariant Born-Oppenheimer bond length through the following equation: where m e is the electron mass and d A and d B are dimensionless parameters.

r e = r BO e 1 + m e d A M A + d B M B , (18)
Strictly speaking, M A and M B are the masses of the nuclei, but a negligible error is introduced if atomic masses are used instead [14]. The Born-Oppenheimer equilibrium bond length r BO e , along with d A and d B , can be determined from the experimental data, provided that there are at least three different evaluations of r e . In this case, where the equilibrium distance has been derived only for two isotopologues, we approximate Eq. ( 18) to the following pair of equations:

r (NH) e = r BO e 1 + m e d A M H (19) r (ND) e = r BO e 1 + m e d A M D , ( 20
)
where the contribution due to the nitrogen nucleus has been neglected. Using our result for r

(ND) e
, the value of r (NH) e from. Ref. [START_REF] Ram | [END_REF] and solving for r BO e and d A , we obtain r BO e =1.036145 Å and d A =1.8454. Alternatively, we can use a more sophisticated approach, which has the advantage to provide an estimate of the error associated to r BO e . Considering that in the Dunham's expansion Y 01 B e [15], the Born-Oppenheimer breakdown can be taken into account by using the following relationship [14,16,17]:

Y α 01 = U 01 µ -1 α 1 + m e ∆ H 01 M α H + ∆ N 01 M α N . (21)
Here, α labels the isotopologue, ∆ H 01 and ∆ N 01 are the electronic BO breakdown parameters for H and N atoms, respectively. M α H and M α N , represent the atomic masses of hydrogen and nitrogen isotopes in the α isotopologue, and µ α the reduced mass of the α isotopologue. As before, we can derive U 01 and ∆ H 01 from Y , whereas ∆ N 01 is experimentally indeterminate. However, since the ∆ 01 coefficients are usually comparable in magnitude for the two atoms [17], the value of ∆ N 01 can be assumed equal to that of ∆ H 01 . A more realistic assumption is ∆ N 01 = ∆ H 01 ± 1, as proposed by Cooke and Gerry [18]. In such a case the results are: U 01 = 470801.7(268) MHz•amu and ∆ H 01 = -3.68 (13). The errors on U 01 and ∆ H has been taken into account too. The Born-Oppenheimer bond lengths r BO e for NH can be determined using the equation [19]:

r BO e = 4πU 01 m u . (22)
In Eq. ( 22) m u =1.6605401(10)•10 -27 kg [START_REF] Mills | Quantities, Units and Symbols in Physical Chemistry[END_REF]. Eventually, we obtained r BO e = 1.036071(29) Å, very close to the result obtained with the previous derivation.

Whichever procedure one wishes to follow, the lack of experimental data prevents from a simultaneous determination of the BO bond length and the two correction parameters, either ∆'s or d's, and some kind of assumption needs to be made. However, the values obtained for ∆ H 01 (-3.68( 13)) and d A (1.8454), although quite different, suggest a rather large deviation from the Born-Oppenheimer approximation. This is consistent with the sizable difference r NH e -r ND e =0.00052 Å. Further insight into the nature of the BO approximation breakdown, concerning the non-adiabatic contribution due to the interaction with excited electronic states, can be achieved if one knows the value of the g J factors for several isotopologues [14]. Since precise g J factors for the ground electronic state are known only for 14 NH in excited vibrational states [5] and are very small, we decided not to consider this contribution explicitly in our derivation. In order to obtain a calculation of r BO e without approximations, it would be useful to have further experimental determinations of r e , from 15 NH and possibly from 15 ND. This experiment is being currently set up in our laboratory.

For a free radical molecule in a X 3 Σ -electronic state, similar arguments hold also for the spin-spin and electron spin-rotation parameters λ e and γ e [20]. However, the d and ∆ coefficients are specific of the spectroscopic constant considered. Since in the case of NH no equilibrium values are readily available for these parameters, the determination of λ BO e and γ BO e can not be performed. It would be possible anyway to derive λ e and γ e for NH from the parameters reported in Ref. [START_REF] Ram | [END_REF].

Conclusions

The pure rotational spectrum of the free radical ND in its ground electronic state X 3 Σ -has been detected in the frequency range 431-531 GHz. Several fine and hyperfine components of the N = 1 ← 0 transition have been observed in the ground and in excited vibrational states with v up to 6. A global fit including our data, infrared and far infrared literature data has been performed, leading to the determination of an extensive set of vibrational, rotational, fine structure and hyperfine structure parameters. The hyperfine constants b F (Fermi contact) and c (dipolar) have been discussed in terms of the expectation values of the appropriate operators over the electronic wavefunction. The precise determination of the equilibrium rotational constant B e has allowed the derivation of the equilibrium bond length r e . The comparison with r e obtained from NH data has shown a sizable discrepancy, which has been attributed to the breakdown of the Born-Oppenheimer approximation. The isotopically invariant Born-Oppenheimer bond length r BO e has been tentatively derived and a value of r BO e = 1.036071(29) Å has been obtained. Further investigations of the rare isotopologues 15 NH and 15 ND will be carried out in the near future. b Ref. [21].

c For calculation see text.

d Ref. [8].

  Fine and Hyperfine Structure of the N = 1 ← 0 Transition of ND(X 3 Σ ∼ 1. Introduction

= 2 Figure 1 .

 21 Figure 1. Portion of the hyperfine structure of the N = 1 ← 0, J = 2 ← 1 fine structure component for different vibrational states.

01

  have been derived by calculating the difference in their values, as obtained in the two limit cases ∆ N 01 = ∆ H 01 -1 and ∆ N 01 = ∆ H 01 + 1. The contribution from the errors associated to Y

i (3 cos 2 θ i - 1 )/r 3 ia

 13 This work. These values are derived from the equilibrium hyperfine constants b F,e (N), b F,e (D), ce(N), ce(D).

Table 1 .

 1 Spectroscopic parameters derived from a global fit for the vibrational states from v = 0 to v = 6 of ND(X 3 Σ -).

	Parameter	Value a	Parameter b		Value a
	ωe	/ cm -1	2 398.9948(25)	λe	/MHz	27 541.76(16)
	ωexe	/ cm -1	41.9800(22)	λv	/MHz	9.80(50)
	ωeye	/ cm -1	0.07032(79)	λ v 2	/MHz	-5.70(39)
	ωeze	/ cm -1	-0.00440(12)	λ v 3	/MHz	-0.468(86)
	ωewe	/ cm -1 ×10 3	-0.5260(72)	γe	/MHz	-900.370(25)
	Be	/ MHz	267 067.288(35)	γv	/MHz	33.796(27)
	αe	/ MHz	7 613.20(10)	γ v 2	/kHz	-33.2(75)
	γe	/ MHz	20.326(96)	γ D,e	/kHz	114.4(37)
	δe	/ MHz	-2.253(38)	γ D,v	/kHz	48.1(34)
	e	/ kHz	-62.1(64)	γ D,v 2	/kHz	-2.94(30)
	ζe	/ kHz	-9.0(40)	γ H,e	/kHz	-0.475(19)
	De	/ MHz	14.7151(16)	b F,e (N)	/MHz	19.096(25)
	β D e	/ MHz	0.16714(78)	b F,v (N)	/MHz	-0.328(28)
	η D e	/ kHz	8.73(26)	b F,v 2 (N)	/kHz	-38.6(78)
	θ D e	/ kHz	-0.126(54)	ce(N)	/MHz	-68.067(43)
	κ D e	/ kHz	0.1087(31)	cv(N)	/MHz	0.294(25)
	He	/ Hz	576.(11)	b F,e (D)	/MHz	-9.844(16)
	β H e	/ Hz	-14.9(43)	b F,v (D)	/MHz	-0.4622(96)
	η H e	/ Hz	-5.94(48)	ce(D)	/MHz	14.06(11)
	Le	/ mHz	-202.(22)	cv(D)	/MHz	-0.192(96)
	β L e	/ mHz	-15.3(22)	c v 2 (D)	/kHz	-64.(21)
				eq 0 Qe(N)	/MHz	-3.402(64)
				eq 0 Qv(N) /MHz	0.290(32)
				Ce(N)	/kHz	93.5(83)
				Cv(N)	/kHz	-13.0(42)

a The quantity in parentheses is one standard deviation in units of the last quoted digit.

b The symbol with v n as subscript refers to the (v +

1

2 ) n dependence of the corresponding parameter. See text.

Table 2 .

 2 : Observed frequencies of the N = 1 ← 0 transition in vibrational excited states of ND(X 3 Σ -).

	1 2 3	J	F N	F	J	F N	F a Obs.-calc. Molecular Physics Frequency b
	4
	5
	6
	7
	8
	9
	10
	11
	12
	13
	14
	15
	16
	17
	18
	19
	20
	21
	22
	23
	24
	25
	26
	27
	28
	29
	30
	31
	32
	33
	34
	35
	36
	37
	38
	39
	40
	41
	42
	43
	44
	45
	46
	47
	48
	49
	50
	51
	52
	53
	54
	55
	56
	57
	58
	59

c Rel. Int. d MHz kHz URL: http://mc.manuscriptcentral.com/tandf/tmph

Table 2 -

 2 continued from previous page

	J	F N	F	J	F N	F a	Frequency b	Obs.-calc. c	Rel. Int. d
							MHz	kHz

Table 3 .

 3 Expectation values of distribution functions over the electronic wave function (in units of 10 30 m -3 ) of N, D and H for ND and NH.

Page 12 of 13 URL

	June 28, 2011 18:17 Molecular Physics Dore˙ND˙rev
	Molecular Physics
	REFERENCES	13
	F o r
	P
	e
	e r
	R
	e
	v i e
	w
	O n l
	y

: http://mc.manuscriptcentral.com/tandf/tmph

Table 4 .

 4 Experimental and calculated equilibrium hyperfine constants (in MHz) for ND and NH.

		Exp. a	Calc. b
	b F,e (N)	19.096	18.94
	ce(N)	-68.067	-67.88
	b F,e (D)	-9.844
	ce(D)	14.06
	b F,e (H)	-64.13	-66.35
	ce(H)	91.59	91.96

a This work. See text. b Ref. [9]. Page 13 of 13 URL: http://mc.manuscriptcentral.com/tandf/tmph Molecular Physics

Dore et al.

See references in the article for a full list of the spectroscopic observations concerning ND.

Failure of observing lines in v = 5 state is due to poor radiation power in the proper frequency region.

Acknowledgements

Financial support from MIUR (PRIN 2007 funds, project "Trasferimenti di energia, carica e molecole in sistemi complessi") and from the University of Bologna (RFO funds) is gratefully acknowledged. LB acknowledges support from the Science and Technology Foundation (FCT, Portugal) through the grant SFRH/BPD/62966/2009.


Luca Bizzocchi
email: bizzocchi@oal.ul.pt

Filippo Tamassia

Claudio Degli Esposti

Luciano Fusina

Elisabetta Canè

Luca Dore

Claudio Degli Esposti

Prof Di Lonardo Gianfranco

High-resolution infrared spectroscopy of diacetylene below 1000 cm -1, †

Keywords: Vibration-rotation bands, Diacetylene, Hot bands, Anharmonicity constants, Resonances Vibration-rotation bands, Diacetylene, Hot bands, Anharmonicity constants, Resonances

come

High-resolution infrared spectroscopy of diacetylene

below 1000 cm-1

Introduction

Diacetylene (C 4 H 2 ) is the simplest member of the family of poly-acetylene compounds and is a molecule of relevant astrophysical interest. It plays an important role in the stratospheres of the giant planets and their moons (see e.g., Ref. [1,2]) acting as a UV shield analogously to O 3 in the Earth's atmosphere [3]. Like ozone, C 4 H 2 is photochemically reactive and it is thought to be a major source of the larger hydrocarbons, which then form the organic aerosols (tholines) characteristic of the reducing atmospheres of these solar system bodies [4]. Diacetylene has also been detected in the dusty envelopes of the C-rich proto-planetary nebula CRL 618 [5], and even outside the galaxy in a similar object embedded in the Large Magellanic Cloud [6]. Being centrosymmetric, C 4 H 2 has no permanent dipole moment, thus searches for this species in astrophysical environments rely on its mid infrared spectrum, in particular the strong ν 8 perpendicular band at 16 µm, even if detection of the weaker bending fundamental ν 9 [1] and of the combination ν 6 + ν 8 [5] was also reported. The vibrational spectrum of diacetylene was the subject of extensive infrared and Raman investigations (see Ref. [7] and references therein). Several high-resolution studies of the infrared active bands have been performed. The ν 4 and ν 5 band systems were investigated by Guelachvili et al. [8], whereas the unusually strong ν 6 + ν 8 combination band at ca. 1250 cm -1 has been investigated by Matsumura et al. [9] by Stark modulation infrared diode laser spectroscopy and, later, by McNaughton et al. [10] with a Fourier transform infrared (FTIR) spectrometer. Below 1000 cm -1 the infrared spectrum of diacetylene is characterised by the strong ν 9 and ν 8 bands, which were the subject of the high resolution study by Arié and Johns [11]. They reported band constants for the ν 9 and ν 8 fundamentals, and for several hot bands originating from the low lying v 9 bending mode and its overtones. However, they were unable to obtain a satisfactory analysis for the ν 8 + ν 9 -ν 9 hot band due to an anharmonic interaction coupling the v 8 = v 9 = 1 bending combination state and the v 3 = 1 stretching state.

In diacetylene, the dipole moment second derivatives with respect to the C≡C-H bending curvilinear coordinates are large, and this produces enhanced transition moments for the overtone, combination, and difference bands involving the v 6 and v 8 quanta [12,13]. For instance, the intensity of the ν 8 -ν 6 band is substantially enhanced and, since its origin is located below 100 GHz, it provides a mean of detecting diacetylene through millimetre-wave (mm-wave) spectroscopy [14,15]. Very recently, we have extended to the submillimetre-wave (submm-wave) region the study of this band and of its ν 8 +ν 9 -(ν 6 +ν 9 ) hot band, performing a complete analysis of the l-type resonance effects and of the above mentioned anharmonic interaction affecting the v 8 = v 9 = 1 state [16]. However, no direct information on the v 3 = 1 state could be obtained from mm-wave spectra, thus the spectroscopic constants of this level had to be estimated from recent theoretical calculations [17,18].

In this paper we present the results obtained from the analysis of the infrared spectrum of diacetylene in the 500-1000 cm -1 range, recorded using an FTIR spectrometer. Besides the reinvestigation of the ν 8 band system, we could identify the previously unreported ν 7 + ν 9 combination band at ca. 700 cm -1 , together with the strong ν 9 -associated hot bands ν 8 + ν 9 -ν 9 and ν 7 + 2ν 9 -ν 9 . In addition, the very weak and previously unreported ν 3 -ν 9 difference band -which provides direct information on the v 3 = 1 state -was also identified and analysed, leading to a comprehensive analysis of the lowest diacetylene resonance system.

The model adopted for the analysis includes the rotational and vibrational l-type resonances active between the various l-sublevels of the multiple bending excited states, and also considers the two anharmonic interactions which couple v 3 = 1 with the v 8 = v 9 = 1 combination and the v 7 = 2 overtone states through the normal coordinates cubic potential constants φ 389 and φ 377 . The previously recorded mmwave transition frequencies [16] were also taken into account in order to obtain the best set of spectroscopic parameters for the ground state and for the v 6 = 1 and v 8 = 1 bending states.

The results of the spectral analyses also provide a comprehensive set of anharmonicity constants for the three bending-bending combination states of diacetylene located below 1000 cm -1 , namely v 8 = v 9 = 1, v 6 = v 9 = 1, and v 7 = v 9 = 1. The accurate treatment of the anharmonic resonances has enabled the determination of de-perturbed molecular parameters, whose values can be compared with the results of the recent high-level ab initio calculations [17,18].

Molecular Physics diac˙v8˙rev

Infrared spectrum of diacetylene 3

Experimental details

Diacetylene was prepared by reacting commercially available 1,4-dichloro-2butyne with acqueous NaOH using the procedure described by Armitage et al. [19]. High resolution IR spectra were recorded in the 500-1000 cm -1 range using a BOMEM DA3.002 Fourier transform spectrometer at the University of Bologna. The interferometer was equipped with a KBr beam splitter, a Globar source, and a HgCdTe detector. An optical path length of 0.18 m was adopted with sample pressures of 270 Pa and 27 Pa. The attained unapodised resolution was 0.004 cm -1 . A large number of scans were co-added in order to improve the signal-to-noise ratio: in total, 780 scans were recorded for the low pressure spectrum whereas 880 scans were co-added for the high pressure spectrum. An overview of the recorded IR spectrum of diacetylene between 590 and 730 cm -1 is shown in Figure 1 together with the assignments of the most prominent spectral features. Ro-vibrational transitions of H 2 O [20] and CO 2 were used for calibration. The most accurate determination of the CO 2 transitions used as reference (ν 2 fundamental band) was performed by J.W.C. Johns, M. Noël and T.L. Tan. The line list was received as a private communication but, to our knowledge, the paper has not been published yet. The average precision of our measurements was estimated to be around 0.004 cm -1 .

Observed spectra and analysis

In the present work, we have investigated the ro-vibrational bands ν 8 , ν 7 + ν 9 , ν 8 + ν 9 -ν 9 , ν 7 + 2ν 9 -ν 9 , and ν 3 -ν 9 . The IR data have been analysed together with the previous mm-wave measurements of the ν 8 -ν 6 and ν 8 + ν 9 -(ν 6 + ν 9 ) bands [16], so that our analyses deal with a manifold of states involving excitation of the v 3 stretch (Σ + g ) and all the four bending quanta (v 6 , v 7 , v 8 , v 9 ). A list of the normal mode fundamentals is reported in Ref. [8]. The reader is also referred to Ref. [18] for a comprehensive description of the internal symmetry coordinates of diacetylene.

Diacetylene has two equivalent hydrogen nuclei and two pairs of equivalent carbon nuclei which are interchanged by the C 2 symmetry operation. As a result, the ro-vibrational lines show a 3:1 intensity alternation due to the spin statistical weights. Table 1 summarises the spin multiplicity rules which apply for the vibrational states involved in the investigated bands.

The spectra were analysed by expressing each ro-vibrational energy term using the formalism originally developed by Yamada and co-workers [21,22], and employed previously to fit the excited state rotational spectra of several carbon chains (see for example Ref. [23]). Briefly, the ro-vibrational Hamiltonian is first represented using the symmetric top basis functions |v 3 , v l6 6 , v l7 7 , v l8 8 , v l9 9 ; J, k , with k = l 6 + l 7 + l 8 + l 9 . Employing the simplified notation |l 6 , l 7 , l 8 , l 9 ; k the elements of the Hamiltonian matrix which are diagonal in the quantum numbers v t and l t are: where G v is the "pure" vibrational energy of the |v 6 , v 7 , v 8 , v 9 state and includes the appropriate anharmonic corrections. The origin of a ro-vibrational band, ν 0 , is then given as

l 6 , l 7 , l 8 , l 9 ; k| Ĥ|l 6 , l 7 , l 8 , l 9 ; k = G v + t=6,9 x L(tt) l 2 t + t ′ >t=6,9 x L(tt ′ ) l t l t ′ +    B v + t=6,9 d JL(tt) l 2 t + t ′ >t=6,9 d JL(tt ′ ) l t l t ′    {J(J + 1) -k 2 } -D v {J(J + 1) -k 2 } 2 + H v {J(J + 1) -k 2 } 3 , (1)
G ′ v -G ′′ v .
The off-diagonal rotational l-type doubling terms (∆l t ± 2, |∆k| = 2) have the general form:

l t ± 2,l t ′ ; k ± 2| Ĥ|l t , l t ′ ; k = 1 4 {q t + q tJ J(J + 1)} × (v t ∓ l t )(v t ± l t + 2) × [J(J + 1) -k(k ± 1)][J(J + 1) -(k ± 1)(k ± 2)] ; (2)
whereas the off-diagonal vibrational l-type doubling terms (∆l t ±2, ∆l ′ t ∓2, ∆k = 0) are given by the expression:

l t ± 2, l t ′ ∓ 2; k| Ĥ|l t , l t ′ ; k = 1 4 {r tt ′ + r tt ′ J J(J + 1)} × (v t ∓ l t )(v t ± l t + 2)(v t ′ ∓ l t ′ + 2)(v t ′ ± l t ′ ) . (3)
The resulting energy matrix is factorised into symmetric and antisymmetric blocks adopting the following Wang-type linear combinations of wavefunctions [21].

|v 3 ,v l6 6 , v l7 7 , v l8 8 , v l9 9 ; k ± = 1 √ 2 {|v 3 , v l6 6 , v l7 7 , v l8 8 , v l9 9 ; k ± |v 3 , v -l6 6 , v -l7 7 , v -l8 8 , v -l9 9 ; -k } (4a) |v 3 ,v 0 6 , v 0 7 , v 0 8 , v 0 9 ; 0 + = |v 3 , v 0 6 , v 0 7 , v 0 8 , v 0 9 ; 0 (4b)
Any sublevel belonging to a given vibrational state can thus be labelled through its k value and by the '+' or '-' subscripts which designate the sign of the symmetrised linear combination or, alternatively, with the e, f parity labels [24] using the correspondence + ↔ e, -↔ f for even k and vice versa for odd k.

As mentioned in the Introduction, the three-state resonance system v 3 = 1 ∼ v 8 = v 9 = 1 ∼ v 7 = 2 has been considered in the present work, and thus two Fermitype (∆k = 0) anharmonic interaction terms were included in the Hamiltonian. The corresponding off-diagonal matrix elements in the symmetrised basis are expressed by [25]:

e 1, 0 0 , 0 0 , 0 0 , 0 0 ; J, 0| Ĥ30 |0, 0 0 , 0 0 , 1 1 , 1 -1 ; J, 0 e = W 389 + W 389J J(J + 1) = 1 2 φ 389 + W 389J J(J + 1) . ( 5
) e 1, 0 0 , 0 0 , 0 0 , 0 0 ; J, 0| Ĥ30 |0, 0 0 , 2 0 , 0 0 , 0 0 ; J, 0 e = W 377 + W 377J J(J + 1) = √ 2 4 φ 377 + W 377J J(J + 1) . ( 6
)
Over 1200 ro-vibrational transitions have been measured, and the complete list of analysed data is available as electronic supplementary material. Some details dealing with the analyses performed for the various bands are given in the following subsections. The v 8 = 1 (Π u ) state has two |k| = 1 ro-vibrational sublevels, thus ν 8 is a perpendicular band whose transitions connect Π e u ← Σ +, e g levels in the P ,R branches and Π f u ← Σ +, e g levels in the Q branch. Making use of the spectroscopic constants of Ref. [11], the assignment of the IR lines was straightforward: we measured a total of 261 transitions, up to P (105), R(109), and Q(81). These data were analysed together with our previous mm-wave ro-vibrational frequencies of the ν 8 -ν 6 difference band, allowing to simultaneously determine of the spectroscopic parameters of the ground and of the v 6 = 1 and v 8 = 1 states.

The analysis was performed by expressing the ro-vibrational energy levels through Eqs. ( 1) and (2). A proper weighting factor (w = 1/σ 2 ) was assigned to each datum, in order to take into account the different measurement precisions. The mm-wave lines were included in the data sets retaining the same weighting scheme already used in the previous work [16], whereas σ = 5 × 10 -4 cm -1 was adopted for the present IR wavenumbers. The results of the least-squares fit are reported in Table 2 The inclusion in the data set of the very precise mm-wave data was beneficial to reduce significantly the correlations between the lower and upper state parameters of the ν 8 band. Improved values of the ground state rotational and quartic centrifugal distortion constants of diacetylene were obtained. Compared with the previous best literature values of Ref. [11], i.e. B 0 = 4389.3019(38) MHz and D 0 = 0.47004(54) kHz, our precision is five and ten times better, respectively. This band was first studied by Arié and Johns [11]. However, they could not attain a satisfactory fit of the experimental wavenumbers and, even using a limited data set (J ≤ 50), the root mean square of the residuals remained substantially larger than the experimental precision (rms IR ∼ 0.007 cm -1 ). According to Matsumura and Tanaka [15] the anomaly was attributed to an anharmonic resonance between the 0 e (Σ + g ) sublevel of the v 8 = v 9 = 1 state and the nearby v 3 = 1 (Σ + g ) C-C stretching state. This interaction was actually considered in our recent mmwave investigation [16] where, with the help of the recent theoretical results [17,18], we achieved a satisfactory analysis of the ν 8 + ν 9 -(ν 6 + ν 9 ) band.

In the present work we have also identified the weak spectrum of the ν 3 -ν 9 difference band thus gaining direct information on the v 3 = 1 state. Using the values of the ν 3 and ν 9 fundamentals reported in Ref. [8], ν 3 -ν 9 is calculated at ca. 652 cm -1 . Fortunately, in that region the intensity of the ν 8 band system is rapidly decreasing. Its transitions vanish almost completely at wavenumbers higher than 660 cm -1 , where the ν 3 -ν 9 R-branch lines should have their maximum intensity. A regular sequence of lines showing the characteristic intensity alternation was readily found. The assignment of the R-branch allowed to locate the Q branch, whose lines are widely spaced toward lower wavenumbers, due to the large negative value of the B v (v 3 = 1) -B v (v 9 = 1) difference. A portion of the Q branch of the ν 3 -ν 9 band spectrum is shown in Figure 2. Despite the spectral congestion we could assign 60 transitions up to R(60), Q(57), and P (32). They were analysed simultaneously with 567 lines assigned to the ν 8 + ν 9 -ν 9 hot band (spanning J interval 5-80) and with 158 mm-wave transitions of the ν 8 + ν 9 -(ν 6 + ν 9 ) difference band [16].

The weighting scheme was similar to that adopted for the simultaneous fit of ν 8 and ν 8 -ν 6 described in the previous subsection: experimental uncertainties (σ) for mm-wave lines were as in the original paper [16], while σ = 1 × 10 -3 cm -1 was assigned to each infrared line belonging to the ν 8 + ν 9 -ν 9 and ν 3 -ν 9 bands.

The ro-vibrational level energies were expressed using Eqs. ( 1)-( 3), where the l-type resonance term of Eq. ( 3) applies for the v 8 = v 9 = 1 and v 6 = v 9 = 1 combination states only. Initially, the analysis of the resonance followed the same scheme adopted for the ν 8 + ν 9 -(ν 6 + ν 9 ) band [16], i.e. we included in the calculations the matrix elements of Eq. ( 5) coupling the v 3 = 1 state and the (v 8 = v 9 = 1) 0e sublevel. The interaction constant of W 389 was fixed at the value derived from the ab initio φ 389 normal coordinate cubic force constant (13.49 cm -1 , Ref. [26]), while the G v (v 3 = 1) vibrational energy was adjusted in the least-squares fit in order to reproduce the observed ν 3 -ν 9 line positions.

The fit was very satisfactory, but the value of the B v rotational constant of the v 3 = 1 state was anomalously large. The corresponding α 3 was ca. 3.1 MHz, much smaller than the ab initio computed value of 9.4 MHz [18]. This suggested that the de-perturbation of the v 3 = 1 state energy was incomplete and that it was still affected by a mixing with a vibrational state characterised by a larger rotational constant. The only state which can be coupled to v 3 = 1 through a first-order (cubic) anharmonic resonance is v 7 = 2, k = 0 e (Σ + g ), which is at ca. 75 cm -1 higher energy and interacts through the normal coordinate force constant φ 377 , whose theoretical value is 113.87 cm -1 [26]. The same scheme of interactions has previously been found for other similar carbon chains, such as HCCCN [27], HCCNC [28] and HCCCP [29], where the lowest-energy stretching state is coupled, through cubic interactions, with nearby bending overtone and combination states.

The analysis was thus performed considering the v 3 = 1, v 8 = v 9 = 1, and v 7 = 2 states simultaneously in order to take into account the resonance effects as completely as possible. The energy level diagram for the three interacting states is illustrated in Figure 3. We adjusted the upper state parameters, keeping the constants of the lower v 9 = 1 state fixed at the values derived from the analysis of the ν 9 band [11]. The rotational parameters B v , D v , and q 7 of the unobserved v 7 = 2 state were evaluated from the experimental results for the ν 7 fundamental [8]; whereas a suitable assumption for the anharmonicity constant x L(77) was obtained by removing from the theoretical value [18] the relevant near-singular term (see Discussion).

The main off-diagonal resonance parameters W 389 and W 377 were constrained to the values derived from the corresponding cubic force constants [17,26], whereas the unperturbed energy difference between the v 3 = 1 and the v 7 = 2 states was assumed equal to the quantity 2ω 7 -ω 3 [18]. Both the centrifugal corrections W 389J and W 377J (see Eqs. ( 5) and ( 6)) were found to be important in reproducing the rotational mixing caused by the v 3 ∼ v 8 + v 9 and v 3 ∼ 2v 7 resonances. While W 389J could be determined, W 377J had to be constrained to a suitable value. As already pointed out, the energy level structure of the lowest resonance system of diacetylene is analogue to that of the isoelectronic molecule HCCCN, so it is likely that their molecular parameters show some similarities. From Ref. [30] we derived a scaled W 377J value of -1.41 MHz (see Discussion for the details), which was constrained in the final analysis. The spectroscopic constants determined by the least-squares fit are collected in Table 3. The overall fit quality (rms IR = 1.4 × 10 -3 cm -1 ) is satisfactory considering that the data set includes the very weak ν 3 -ν 9 band transitions, whose positions are clearly affected by larger measurements errors. The derived unperturbed value of the B v (v 3 = 1) rotational constant yields α 3 = 9.325 MHz which compares very well with the ab initio computed value of 9.395 MHz [18]. A similar agreement (within 1%) is observed for the α 9 vibration-rotation interaction constants, which in turn derives from the analysis of unperturbed vibrational levels [11]. The spectral survey carried out in the present investigation allowed us to identify the ν 7 +ν 9 , Σ + u ← Σ + g parallel band, which appears as a weak but nevertheless clear feature in the high pressure spectra (see Figure 1). Making use of the values of B 0 and ∆B = -(α 7 +α 9 ) (from this work and Refs. [8,11]), we assigned 137 transitions with J ′ up to 75 for both P -and R-branches.

A rather high root mean square deviation, rms IR ∼ 0.008 cm -1 , and an anomaloulsy large value for the upper state D v constant was obtained from the fit of these data. This unsatisfactory result was due to the neglect of the l-type resonance between the k = 0 e and k = 2 e sublevels of the v 7 = v 9 = 1 state. A correct treatment of this resonance requires the knowledge of the x L(79) and r 79 anharmonicity constants, as the observed ∆k = 0 transitions of the ν 7 + ν 9 band cannot provide information on the upper state manifold. These parameters were obtained from the analysis of the associated ν 7 + 2ν 9 -ν 9 hot band (see next subsection), thus allowing to characterise the ro-vibrational levels of the v 7 = v 9 = 1 state as described in the following.

The energy terms were expressed using Eqs. ( 1)-( 3), and the fit was carried out with the ground state constants fixed at the values derived in the present work (see subsection 3.1). As for the l-type doubling constants q 7 , q 9 and q 9J the literature values of the v 7 = 1 [8] and v 9 = 1 [11] states were adopted. The x L(79) and r 79 parameters were from the fit of the ν 7 + 2ν 9 -ν 9 hot band. A new fit, carried out adjusting the band origin, B v , and D v yielded a significantly smaller root mean square deviation (rms IR ∼ 0.001 cm -1 ). We then refined also x L(79) in the final analysis to reproduce accurately the magnitude of the l-type resonance effects. The results of the final least-squares fit are gathered in Table 4. Besides the mere improvement of the fit quality (rms IR = 4 × 10 -4 cm -1 ), the effectiveness of the treatment is supported by the value of D v for the v 7 = v 9 = 1 state, 0.5035 kHz, in good agreeement with that extrapolated from the singly excited bending states. The ν 7 + 2ν 9 -ν 9 hot band appears as a feature superimposed to the ν 7 + ν 9 band. The intensity of its lines is ca. one third of that of the corresponding cold band transitions, owing to the low energy of the v 9 = 1 state (∼ 220 cm -1 [11]). The (v 7 = 1, v 9 = 2) state is composed of three pairs of ro-vibrational sublevels, which can be labelled as k = 3 e,f (Φ g ), k = 1 e,f (Π g ), and k = -1 e,f (Π g ). The k = 3 e,f and k = -1 e,f sublevels have |l 9 | = 2 and are separated from the k = 1 e,f sublevels by the anharmonic contribution x L(99) l 2 9 . On the other hand, the k = 3 e,f and k = -1 e,f sublevels are split symmetrically by the term x L(79) l 7 l 9 , and finally the degeneracy between each pair of e,f levels is lifted by vibrational and rotational l-type resonance effects. Allowed ∆k = 0 transitions are observed between v 9 = 1, k = 1 e,f and (v 7 = 1, v 9 = 2), k = 1 e,f sub-states. The spectrum consists of two sequences of e and f lines rather than close doublets, owing to the large separation of e and f components. We identified a total of 218 lines up to P (57) and R(67).

The analysis was carried out by expressing the energies of the ro-vibrational 1)-( 3), keeping the spectroscopic constants of the v 9 = 1 state fixed at the literature values [11]. The excited-state constants q 7 and x L(99) were also assumed from the results of previous analyses [8,11], whereas both the anharmonicity constants x L(79) and r 79 were adjusted. In fact, our experimental data are sensitive to both these parameters, since x L(79) contributes to the energy difference between the various k sub-states, and r 79 is responsible for the vibrational l-type resonance between the k = 1 and k = -1 sub-states. The results of the leastsquares fit are reported in Table 5.

Discussion

In the present study we have determined several spectroscopic constants for a number of states involved in vibrational and ro-vibrational resonances, due both to the l-sublevel structure of the bending excited state and to accidental neardegeneracies. In our spectral analyses, we succeeded in retrieving molecular parameters with unambiguous physical meaning. This is particularly relevant if one aims at a meaningful comparison with the results of theoretical calculations. Due to correlations between the model Hamiltonian parameters, not all the quantities involved in the two anharmonic interactions analysed, i.e. v 3 ∼ v 8 +v 9 and v 3 ∼ 2v 7 , could be determined from the least-squares fits, so that a number of assumptions had to be made. Whenever possible, the values of the parameters have been derived from the experimental results for vibrational states at a lower level of excitation (e.g. B v and D v for the unobserved v 7 = 2 state) or from the most recent ab initio anharmonic force field calculations [17,18]. The only exception was the centrifugal correction parameters W 377J , whose value was estimated empirically from the results of a global ro-vibrational analysis performed for the isoelectronic carbon chain HCCC 15 N [30]. The resonance between v 4 = 1 and v 6 = 2, present in HCCC 15 N, is analogous to that between v 3 = 1 and v 7 = 2 which has been considered for diacetylene: Ref. [30] gives W 466 = 43.03 cm -1 and W 466J = -1.84 MHz. According to Aliev and Watson [START_REF] Aliev | Molecular Spectroscopy: Modern Research[END_REF], it is possible to estimate the order of magnitude of an Hamiltonian coefficient C mn using

C mn ≈ κ m+2n-2 ω vib , (7)
where an approximate value of the Born-Oppenheimer expansion parameter κ ≈ 0.1 can be adopted for most vibrational problems. The ratio between W stt and W sttJ is of the order κ -4 , thus it is appropriate to derive an estimate for an unknown W 377J of diacetylene using the scaling equality:

W 377 (HC 4 H) W 466 (HC 3 N) 4 = W 377J (HC 4 H) W 466J (HC 3 N) . (8)
Inserting in this relation the above mentioned HCCC 15 N values and W 377 = 40.255 cm -1 of diacetylene derived from the force field, one obtains W 377J = -1.41 MHz used in the present analysis. It is to be noted that we could have omitted W 377J without degrading significantly the quality of the fit, but then we would have been unable to model the correct rotational mixing produced by the strongest anharmonic interaction v 3 = 1 ∼ v 7 = 2, as already pointed out in a number of previous studies (e.g. Ref.

[32]). On the other hand, the effectiveness of such a treatment can be assessed by the inspection of the converged values of the rotational and quartic centrifugal distortion constants for the v 8 = v 9 = 1 and

Molecular Physics diac˙v8˙rev

Infrared spectrum of diacetylene v 6 = v 9 = 1 states, which compare nicely with those obtained from the ground and singly excited state parameters. The derived α 3 exhibits a remarkable agreement with the ab initio value [18] (within 1%), and also D v of the v 3 = 1 state converges to a value slightly smaller than the ground state one, as expected for a stretching excited state of a well-behaved semi-rigid linear molecule.

We have also determined the x L(tt ′ ) anharmonicity constants and the r tt ′ vibrational l-type constants for the three bending-bending combination states of diacetylene below 1000 cm -1 , i.e. v 8 = v 9 = 1, v 6 = v 9 = 1, and v 7 = v 9 = 1. These parameters can be compared to the theoretically computed ones in Ref. [18] if one removes the near-singular terms entering in the second-order perturbation expressions in cases where 2ω t ≈ ω s or ω t + ω t ′ ≈ ω s [33]. The relevant terms can be factored as follows:

x L(tt) = x * L(tt) - 1 32 φ 2 stt (2ω t -ω s ) -1 , (9)
x L(tt

′ ) = x * L(tt ′ ) - 1 16 φ 2 stt ′ (ω t + ω t ′ -ω s ) -1 , (10)
r tt ′ = r * tt ′ + 1 8 φ 2 stt ′ (ω t + ω t ′ -ω s ) -1 , (11)
where the asterisk at the superscript labels "de-perturbed" anharmonicity constants. Eq. ( 9) was used to compute the de-perturbed x * L(77) parameter from the theoretical x L(77) prediction [18]. x * L(77) was then fixed in the analysis of the polyad involving the v 7 = 2 state (see subsection 3.2 and Table 3). The comparison between the experimentally derived and ab initio computed values of these constants is presented in Table 6. As can be seen, the agreement is satisfactory for x L(tt ′ ) , but rather poor for the vibrational l-type parameters. This different behaviour was somewhat expected since the expression of r tt ′ contains a semi-diagonal quartic potential constant of the type φ ttt ′ t ′ , whereas the x L(tt ′ ) are defined in terms of the quadratic and cubic potential constants only. It is likely that high-order potential constants are affected by larger uncertainties, (usually not estimated in theoretical works), which then lead to less accurate predictions of the derived spectroscopic parameters.

This hypotesis seems to be supported by the energy difference between the k = 0 f and k = 2 f sublevels of the v 8 = v 9 = 1 state, which are both unaffected by the anharmonic resonance. A spectral analysis limited to (v 8 = v 9 = 1) 0f,2f data gives E 2f -E 0f = 0.62 cm -1 . This difference equals, to a very good approximation, the sum 2x L(89) + r 89 -4B v , in which the two near-singular terms included in x L(89) and r 89 cancel each other out. Using the ab initio calculated constants of [18] it results E 2f -E 0f = -0.09 cm -1 , which differs from the experimental value by 0.71 cm -1 (i.e. 21.3 GHz), an amount close to the discrepancy between experimental and theoretical r 89 values.

Conclusions

The present paper reports on the results of the accurate analysis of the rovibrational spectrum of diacetylene in the 500-1000 cm -1 region. New infrared measurements allowed for the identification of the previously unreported ν 7 + ν 9 , ν 7 + 2ν 9 -ν 9 , and ν 3 -ν 9 bands. Also, we reinvestigated the ν 8 fundamental and its associated ν 8 + ν 9 -ν 9 hot band, analysing these infrared data with the mm-vawe transitions of the ν 8 -ν 6 and ν 8 + ν 9 -(ν 6 + ν 9 ) difference bands [16]. The ro-vibrational treatment of this enlarged data set was performed using the linear molecule Hamiltonian developed by Yamada et al. [21], resulting in a single set of spectroscopic constants for each of the vibrational state involved in the investigated bands. The anharmonic resonances which couples the v 3 = 1 (C-C stretch, Σ + g ) state with the Σ + g (k = 0 e ) sublevel of the v 8 = v 9 = 1 and v 7 = 2 states have been explicitly taken into account in the least-squares fit in order to obtain de-perturbed spectroscopic parameters having clear physical meaning.

This work also achieves the complete characterisation of the lowest bend-bend combination states v 8 = v 9 = 1, v 6 = v 9 = 1, and v 7 = v 9 = 1, through the determination of their anharmonicity and vibrational l-type doubling constants a The quantity in parentheses is one standard deviation in units of the last quoted digit.

b Fixed at the value of Ref. [11].

c Constrained.

d Fixed at the value of the v 7 = 1 state, Ref. [8].

Infrared spectrum of diacetylene 5 3. 1 .

 51 The ν 8 fundamental band

3. 2 .

 2 The ν 8 + ν 9 -ν 9 hot band and ν 3 -ν 9 difference bandThe v 8 = v 9 = 1 state consists of three sub-states: Σ + g and Σ - g having k = 0 (species e and f , respectively) and the doubly degenerate Π g with |k| = 2, split by ro-vibrational l-type resonances into e, f doublets. The ν 8 + ν 9 -ν 9 hot band has a perpendicular structure with ∆k = ±1, and selection rules e ← e, f ← f for Pand R-type transitions, and e ↔ f for Q branch lines.

3. 3 .

 3 The ν 7 + ν 9 combination band

3. 4 .

 4 The ν 7 + 2ν 9 -ν 9 hot band

  Taylor & Francis and I.T. Consultant

Figure 1 .Figure 2 .Figure 3 .

 123 Figure 1. Overview of the infrared spectrum of diacetylene in the ν 8 band region. The assignment of some prominent features is indicated. Spectrum recorded at 270 Pa, 880 scans coadded for a final unapodised resolution of 0.004 cm -1 .

Table 1 .

 1 x * L(89) , x L(69) , x L(79) , r * 89 , r 79 , and r 69 . Bizzocchi, C. Degli Esposti and C. Puzzarini, Mol. Phys. 104, 2627 (2006). [33]D. Papoušek and M.R. Aliev, Molecular vibrational-rotational spectra(Elsevier, 1982). Spin statistical weights of the ro-vibrational levels of diacetylene.

	F o r
	P
	e
	e r
	R
	e
	v i e
	w
	O n l
	y

URL: http://mc.manuscriptcentral.com/tandf/tmph

Table 2 .

 2 Spectroscopic parameters determined by the analysis of the ν8 fundamental and of the ν8 -ν6 difference bands of diacetylene a .

			Molecular Physics
	Parameter		ground state	v8 = 1	v6 = 1
	Gv	/ cm -1	0.0 b	628.040776(36)	625.643507(36)
	Bv	/ MHz	4389.29687(60)	4391.201952(52)	4391.327830(53)
	Dv	/ kHz	0.470150(61)	0.472111(19)	0.472085(20)
	Hv	/ mHz		0.0326(24) c	0.0326(24) c
	qt	/ MHz		2.409552(64)	2.485126(64)
	qtJ	/ Hz		-1.309(11)	-1.187(11)
	rms IR	/ cm -1		2.5 × 10 -4
	rms mmw	/ kHz		22.5

a The quantity in parentheses is one standard deviation in units of the last quoted digit. b Constrained. c Fitted as a mean value for both states.

Page 13 of 18 URL: http://mc.manuscriptcentral.com/tandf/tmph

Table 3 .

 3 Spectroscopic parameters determined by the analysis of the ν8 + ν9 -ν9 hot band and of the ν3 -ν9, ν8 + ν9 -(ν6 + ν9) difference bands of diacetylene a .

					Molecular Physics
	Parameter		v9 = 1	v8 = v9 = 1	v6 = v9 = 1	v3 = 1	v7 = 2
	Gv	/ cm -1	219.974 b	848.365918(53)	845.655513(54)	888.36279(30)	962.96 d
	x * L(77)	/GHz					-50.33 c,d
	x L(99)	/ GHz	17.702 b	17.702 b	17.702 b
	x L(t9)	/ GHz		12.19037(72) c	15.87245(71)
	Bv	/ MHz	4401.8473 b	4403.70993(22)	4403.86433(20)	4379.9715(59)	4405.205 d
	d JL(99)	/ kHz	-17.1 b	-17.1 b	-17.1 b
	d JL(t9)	/ kHz		0.0 d	-4.96(13)
	Dv	/ kHz	0.4883 b	0.490137(80)	0.490189(74)	0.46159(41)	0.47 d
	q7	/ MHz					3.131 d
	qt	/ MHz		2.50214(33)	2.50544(21)
	qtJ	/ Hz		-1.309 d	-1.187 d
	q9	/ MHz	6.155 b	6.21994(19)	6.155 b
	q9J	/ Hz	-12.86 b	-12.86 b	-12.86 b
	rt9	/ GHz		11.8226(16) c	0.2970(16)
	rt9J	/ kHz		-27.98(27)	-5.61(11)
	W389	/ cm -1			6.746 d
	W389J	/ MHz			-1.1569(11)
	W377	/ cm -1			40.255 d
	W377J	/ MHz			-1.41 d
	rms IR	/ cm -1			1.4 × 10 -3
	rms mmw	/ kHz			33.4

a The quantity in parentheses is one standard deviation in units of the last quoted digit. b Fixed at the value of Ref.

[11]

c Deperturbed constant, see text. d Constrained, see text. Page 14 of 18 URL: http://mc.manuscriptcentral.com/tandf/tmph

Table 4 .

 4 Spectroscopic parameters determined by the analysis of the ν7 + ν9 combination band of diacetylene a . The quantity in parentheses is one standard deviation in units of the last quoted digit. Fixed at the value of Ref.[11]. Fixed at the value of the v 7 = 1 state, Ref.[8].

	Parameter		ground state	v 7 = v 9 = 1
	Gv	/ cm -1	0.0 b	701.8939(20)
	x L(77)	/ GHz		0.0 b
	x L(99)	/ GHz		17.702 c
	x L(79)	/ GHz		13.457(57)
	Bv	/ MHz	4389.2969 b	4409.9697(85)
	d J L(99)	/ kHz		-17.1 c
	Dv	/ kHz	0.4702 b	0.50351(49)
	q 7	/ MHz		3.131 d
	q 9	/ MHz		6.155 c
	q 9J	/ Hz		-12.86 c
	r 79	/ GHz		-13.193 b
	rms IR	/ cm -1		4.0 × 10 -4
	a

b Constrained, see text. c d

Table 5 .

 5 Spectroscopic parameters determined by the analysis of the ν7 + 2ν9 -ν9 hot band of diacetylene a .

	Parameter		v 9 = 1	v 7 = 1, v 9 = 2
	Gv	/ cm -1	219.974 b	920.9939(11)
	x L(77)	/ GHz	0.0 c	0.0 c
	x L(99)	/ GHz	17.702 b	17.702 b
	x L(79)	/ GHz		15.71(11)
	Bv	/ MHz	4401.8473 b	4422.5999(34)
	d J L(99)	/ kHz	-17.1 b	-17.1 b
	Dv	/ kHz	0.4883 b	0.54224(94)
	q 7	/ MHz		3.131 d
	q 9	/ MHz	6.155 b	6.155 b
	q 9J	/ Hz	-12.86 b	-12.86 b
	r 79	/ GHz		-13.193(63)
	rms IR	/ cm -1		5.5 × 10 -4

Table 6 .

 6 Comparison between experimental and computed anharmonicity constants for diacetylene a . CO 2 , ν 2 ν 3 -ν 9 ν 7 + ν 9 ν 7 + 2ν 9 -ν 9

	Parameter	experimental	ab initio
	x * L(89)	/ GHz	12.19	14.4
	x L(79)	/ GHz	13.46	15.0
	x L(69)	/ GHz	15.87	14.7
	r * 89	/ GHz	11.82	-13.8
	r 69	/ GHz	0.297	3.9
	r 79	/ GHz	-13.19	-32.8

a Asterisks at the superscript label deperturbed constants (see text). Page 16 of 18 URL: http://mc.manuscriptcentral.com/tandf/tmph

URL: http://mc.manuscriptcentral.com/tandf/tmph Molecular Physics

Taylor & Francis and I.T. Consultant

Acknowledgements

The authors wish to thank M. Di Lauro for the preparation of the sample. Financial support from MIUR (PRIN 2007 funds, project "Trasferimenti di energia, carica e molecole in sistemi complessi") and from the University of Bologna (RFO funds) is gratefully acknowledged. LB acknowledges support from the Science and Technology Foundation (FCT, Portugal) through the grant SFRH/BPD/62966/2009.