Experimental , theoretical and computational investigation of the inelastic neutron scattering spectrum of a homonuclear diatomic molecule in a nearly spherical trap : H 2 @ C 60

In this paper we report a methodology for calculating the inelastic neutron scattering spectrum of homonuclear diatomic molecules confined within nano-cavities of spherical symmetry. The method is based on the expansion of the confining potential into multipoles of the coupled rotational and translational angular variables. The Hamiltonian and the INS transition probabilities are evaluated analytically. The method affords a fast and computationally inexpensive way to simulate the inelastic neutron scattering spectrum of molecular hydrogen confined in fullerene cages. The potential energy surface is effectively parametrized in terms of few physical parameters comprising an harmonic term, anharmonic corrections and translation–rotation couplings. The parameters are refined by matching the simulations against the experiments and the excitation modes are identified for transfer energies up to 215 meV.


Introduction
The topic of confined particles lies at the foundations of quantum mechanics and textbook examples of a particle confined in a well or moving in a harmonic trap provides the starting point for a deeper understanding of many quantum systems.With the advent of ''molecular surgery'' techniques in organic synthesis, pioneered in the laboratories of Komatsu, 1 Murata 2 and more recently Whitby, 3 molecular complexes have become available in which a variety of quantum rotors may be fully enclosed and confined inside a fullerene cage.These are the small molecule endofullerenes, of which H 2 @C 60 is the primary example, representing a model system of its kind. 4ith a hydrogen molecule entrapped inside C 60 , the endofullerene H 2 @C 60 affords an actual realization of a quantum rotor confined in a three dimensional cage.The dynamics of H 2 encapsulated within a fullerene is fundamentally quantum mechanical for several reasons.Firstly, owing to its low mass, free molecular hydrogen is characterized by a large rotational constant with sparse levels whose separation is in excess of 170 K in temperature units.Secondly once hydrogen is confined inside the C 60 nano-cavity its translational degrees of freedom becomes quantised.Significantly, since the internuclear bond length in hydrogen (0.74 Å) is in the same scale as the fullerene radius, the quantised rotational and translational degrees of freedom become coupled.Thirdly, since the nuclei comprising H 2 are indistinguishable fermions, the antisymmetry principle applies in determining the allowable eigenstates and there exists ortho-and para-nuclear spin isomers of H 2 .At low temperature these remain largely decoupled from each other and the two co-exist as independent and identifiable species in the absence of time-dependent inhomogeneous magnetic fields.
The study of the quantum dynamics of H 2 @C 60 is relevant for understanding the interaction between carbon based materials and hydrogen.Molecular endofullerenes are ideal compounds in which the guest molecules are relatively isolated from each other within the homogenous environment provided by the cages.6][7][8] This reveals in fine detail the quantum dynamics of the entrapped translator-rotator and enables the potential energy surface (PES) experienced by the rotor to be probed meticulously along with the various interactions that characterise the Hamiltonian.With its unique ability to induce transitions between ortho-and para-spinisomers, inelastic neutron scattering (INS) is a prime technique for studying confined hydrogen. 9ioneering work by the group of Bac ˇic ´using a rigorous computational approach to rigorously solve the Hamiltonian [10][11][12][13][14] has been instrumental in developing an understanding of

View Article Online
6][17] Their methodology has been extended to determine the INS spectrum of hydrogen confined in nanocavities 18,19 and in fullerenes. 20,21There, the potential energy surface is modelled using modified Lennard-Jones functions to represent the atom-atom potentials, then the Hamiltonian is solved exactly from which the spectrum may be determined.
The scope of this paper is different.Starting from the experimental spectrum, the observed low energy transitions are described in terms of excitations of a rotator angularly coupled to a harmonic oscillator.The quantum dynamics of a homonuclear diatomic molecules in a weakly perturbed isotropic harmonic trap is investigated by expanding the confining potential in symmetry adapted multipoles of the coupled angular variables.The Hamiltonian and the INS transition probabilities are represented analytically in the basis of the coupled rotor and harmonic oscillator.The INS spectrum is simulated by refining few Hamiltonian parameters, including harmonic and anharmonic terms along with translation-rotation coupling constants.In this methodology, the PES is not built from specific microscopic interactions but it is determined from a direct fit of the Hamiltonian to the experimental spectrum with the advantage that the fitted parameters represent physical attributes, such as the harmonic potential, the anharmonicity and translationrotation coupling constants.
In the context of the spectroscopy of H 2 @C 60 , it is wellestablished that the quantum states of the coupled harmonic translator-rotator can be described in terms of the quantum numbers (n,l,j,l). 22The entrapped rotor possesses both rotational and translational energy.The rotational states are characterised by j, while the translational states, defined by displacements of the centre-of-mass of the particle within the 3D-cage, are characterised by n and l.The former, n, is the principle quantum number while the latter, l, determines the translational angular momentum.Rotational and translational angular momenta couple according to the usual vector addition rules and we identify the quantum number l to characterise the coupled translation-rotation states.A more complete account of the rotational and translational states will be discussed in Section 4. Finally, applying the antisymmetry principle, the Pauli Exclusion Principle (PEP) demands that the nuclear spin isomer para-H 2 is a nuclear spin singlet with total spin I N = 0 and characterised by even j, while ortho-H 2 is a nuclear spin triplet with I N = 1 and odd values of j.In this way the spatial and spin degrees of freedom are entangled.This has consequences for experiments since ortho 2 para transitions are spin-restricted and forbidden to photon spectroscopy, whereas neutron scattering interactions couple space and spin enabling ortho 2 para transitions to be directly observed in the INS spectrum with high intensity.
The paper is organised as follows.First the experimental INS spectrum is presented in Sections 2 and 3.The analytical approach used to represent the Hamiltonian and to determine the INS transition probabilities is discussed in Sections 4 and 5, respectively.The method is then applied to simulate the INS spectrum of H 2 @C 60 in Section 6.The parameters characterising the cage potential and the spectral assignments are discussed in Section 7.

Experimental details
INS experiments were conducted using the IN1-Lagrange spectrometer at the Institut Laue-Langevin (ILL) in Grenoble. 23This instrument is supplied from the ''hot source'' moderator of the ILL reactor, resulting in a high-flux beam of energetic neutrons.A choice of three different single crystal monochromators is available in a primary spectrometer, where a Bragg angle defines the incident energy of the monochromatic neutron beam arriving at the sample.The monochromator crystals of Si or Cu are mounted in a double focusing geometrical arrangement, enabling selection of the Si(311), Cu(220) or Cu(331) Bragg reflections.By combining the 3 monochromators, the incident neutron energy may be scanned in the range from 5 meV to 1 eV.The scattered neutrons enter a secondary spectrometer comprising a large area array of pyrolytic graphite analyzer crystals.The Bragg scattering geometry in the secondary spectrometer ensures only neutrons of a fixed energy (4.5 meV) are detected by a 3 He detector, with neutrons being collected from a wide solid angle subtended at the sample.By scanning the incident neutron energy, INS spectra were recorded in the energy transfer range 12.7 r DE r 215 meV.Over this range of energy transfer, the momentum transfer of detected neutrons scattered by the sample varies systematically in the range 2 r q r 12 Å.For the Si monochromator the resolution is given by the resolution of the analyser, which is 0.8 meV.For the Cu monochromators the resolution function of the instrument is a constant percentage of the incident neutron energy.Identifying singlet transitions in the INS spectrum of H 2 @C 60 , the Gaussian resolution function was measured as follows (FWHM; full-width half maximum): Cu(220) 2.3%; Cu(331) 1.5%.
The sample of H 2 @C 60 was prepared according to published methods. 1 In a final stage of purification, the presence of any intercalated solvent was substantially reduced by sublimation.The powdered sample with mass 101 mg was uniformly loaded inside a rectangular Al foil sachet, folded into a cylindrical annulus of approximate diameter 1 cm and mounted in the Lagrange spectrometer.The sample temperature was maintained at 2.5 K using a closed cycle He refrigerator.To subtract background and scattering from Al and from the C 60 cage, a ''blank'' sample of empty C 60 cages with matching mass was prepared in an identical Al foil sachet.
INS spectra of H 2 @C 60 were recorded using the three monochromators and the corresponding energy transfer ranges were selected to provide the best available resolution as follows: Si(331) 12.7 o DE r 27 meV; Cu(220) 27 o DE r 63 meV and Cu(311) 63 o DE r 215 meV.The intensities of the INS peaks recorded with the three monochromators were normalized to the same scale by comparing the integrated intensities of peaks recorded in overlapping energy transfer regions.The systematic error in energy transfer is estimated to be 0.1 meV below approximately 50 meV rising to 0.5 meV at the highest energy transfers recorded.Random errors in energy transfer were small by comparison, typically 0.02 meV for singlet peaks.Unless otherwise stated, these error bars will apply in the following discussions.

Experimental results
The INS spectrum of H 2 @C 60 is presented in Fig. 1.A massmatched ''blank'' has been subtracted, ensuring the spectrum arises only from scattering by H 2 inside the fullerene cage.Compared with previous studies using time-of-flight spectrometers, this new spectrum substantially extends the range of energy transfer recorded from 60 meV to beyond 250 meV.Furthermore, in the range 30-60 meV the resolution function is typically two times better than the previously published IN4 spectrum (l n = 1.2 Å), 24 enabling TR fine structure to be resolved for the first time in the high energy modes.The sample is at low temperature so that only the ground para-H 2 and ground ortho-H 2 states are occupied.At these low temperatures and in the absence of suitable time-dependent magnetic interactions, ortho-para conversion is extremely slow with no conversion being observed on the timescale of the experiments (several days).Therefore the two nuclear spin-isomers exist as distinct species and all INS transitions originate in the ground states, either (0,0,0,0) or (0,0,1,1), using the notation (n,l,j,l) to identify the energy levels.
At low temperature, and particularly in the lower energy transfer region, the spectrum is sparse with sharply resolved peaks revealing a remarkable quantity of detail.H 2 @C 60 is clearly an extremely well-characterized and homogenous system, revealing the quantum nature of the H 2 molecule and its wavelike nature in exquisite detail.
The prominent band centered on 22.5 meV received particular interest in our earlier papers, 24 since it represents the transition from the ground state of ortho-H 2 to the first excited translational state of ortho-H 2 , (0,0,1,1) -(1,1,1,l).This is a ubiquitous ''particle-in-a-box'' mode and in the Lagrange spectrum this is observed with similar resolution to the earlier IN4 spectrum (l n = 1.6 Å).This mode is a partially resolved triplet comprising the l = 0, 1, 2 components of the excited state.This band was significant since it represented the first observation in an INS spectrum of a splitting arising from TR coupling of rotational and translational angular momentum.Moving to higher energy transfer, a doublet is observed with energies 37.4 and 38.5 meV.Considering the energy level diagram and consistent with the assignments in our earlier papers, 8,24 we infer this band is the transition (0,0,0,0) -(1,1,1,l) involving the same final state triplet referred to in the previous paragraph, but this time originating in the ground para-H 2 state.Within the resolution of the spectrometer only two components of this band can be identified while the third one is apparently missing.We shall return to this point more widely in Section 7 in the context of the discussion of the INS transition probability for diatomic molecules in spherical potentials.
The second excited translational state of ortho-H 2 , n = 2, j = 1, is revealed in a doublet with energy transfers 45.8 and 47.0 meV and in the small peak at 49.8 meV.These peaks were observed with the Cu(220) monochromator and have the resolution line width.The former doublet has a partner observed at 61.3 meV while the latter peak is closely related to the singlet at 64.2 meV observed with the Cu(331) monochromator, which also has the resolution line width.Allowing for systematic  errors that arise from two different monochromators, these two set of peaks differ in energy transfer by 2B, so we are able to make the assignment to the same final state within the manifold n = 2, j = 1; the lower energy of the two originates in ortho-H 2 , while the higher energy peak originates in para-H 2 .
The peak in the region 51.5 meV has been observed with the Cu(220) monochromator.It is slightly broader than the resolution line width and slightly asymmetric so we can infer the feature has multiple components.As with the examples above, applying an offset equal to a multiple of 2B provides an energy that coincides with a known splitting in the spectrum.In this case subtracting 4B coincides with the translational splitting so we can infer these two peaks correspond to the ortho-para transition from (0,0,1,1) -(1,1,2,l).
As the above few examples show, with good instrument resolution and a spectrum which is sufficiently sparse, substantial progress can be made in the assignment of n, l, j quantum numbers based solely on the experimental data.This enables the energy level diagram to be determined for all low lying states up to approximately 73 meV.As both para to ortho and ortho to ortho transitions are allowed, and the same final states are represented by multiple peaks, the assignments may be confirmed by correlating transition energies.However, with increasing energy transfer, the number of transitions significantly increases while anharmonicities lead to non-integer relationships between peak energies.It becomes clear that a definitive assignment of the INS peaks, including the TR multiplet components, increasingly demands the assistance of a more sophisticated theoretical framework.Therefore we shall return to consider more detailed assignments in the spectrum in Section 7 after describing the computational simulation of the H 2 @C 60 spectrum.

Quantum dynamics
In this section we discuss the quantum dynamics of a H 2 molecule confined in a C 60 cage.The C 60 cage is assumed to be rigid, fixed and non-rotating.The laboratory reference frame is taken with origin in the centre of the cage.The geometric configuration of the cage is completely identified by the orientation of a frame fixed with the cage C with respect to the laboratory frame L. Mathematically the orientation of the cage frame with respect to the laboratory frame is described by the set of three Euler angles O LC .For example, in the case of C 60 , the cage fixed frame can be chosen so that the Z axis is directed along a C 5 axis with the X axis along one of the 5 orthogonal C 2 axes and the Y axis determined by the right-handed condition.The configuration of the diatomic molecule is determined by the centre of mass vector R and the internuclear vector r.
The molecule-cage interaction is described by the potential V(R,r;O LC ).The function V can be thought as the PES of an H 2 molecule confined inside the fullerene cage.It is worth noting that any symmetry operation O ˆof the icosahedral point group I h acting on all the position coordinates leave the potential unaffected: and the last equality follows from the fact that the Euler angles are unchanged under O ˆ.As a consequence the potential is invariant under operation of the symmetry group acting only on the variables R and r, see also ref. 25 and 26.In the following the dependence of V on O LC will be understood.To discuss the effect of symmetry on the potential at fixed cage orientation, it is convenient to expand out the angular dependence in multipoles where the bipolar spherical harmonics are a complete set of functions for the angular variables defined by 27 and where C is a Clebsch-Gordan coefficient and {y R ,f R } and {y r ,f r } are polar angles of the centre of mass and internuclear vector, respectively.L is the rank and M L is the projection of the bipolar harmonic.For sake of simplicity, the dependence on the ordered set of angular variables will be understood henceforth.
The symmetry of the confinement is reflected in the number and types of multipoles present in the potential expansion, eqn (2).In icosahedral symmetry it is possible to prove that the rank of the multipoles that form invariant functions is restricted to L = 0, 6, 10, 12, 16. 28 Since the centre of the potential is a point of inversion symmetry, only multipoles with L + J = even are allowed.For multipoles with rank L = 0, only terms with L = 0, M L = 0 and L = J are allowed.The next lowest non-zero rank terms for an icosahedral potential is: 29 ffiffiffiffiffi 11 p In treating the quantum dynamics of H 2 @C 60 , the effect of multipoles with L 4 0 can be neglected in first instance and eventually reconsidered perturbatively. 25Retaining only terms with L = 0 is equivalent to treating C 60 as spherical.The internal symmetry of the homonuclear molecule determines the invariance of the potential under the exchange of r in Àr which in turn implies that only even-J bipolar harmonics are present in eqn (2).The vibrationally averaged effective five-dimensional Hamiltonian for a diatomic molecules moving in a spherical potential can be written as where P ˆis the linear momentum operator associated with the centre of mass vector, J ˆis the angular momentum operator associated with end to end rotations, B(n) is the effective rotational constant in the vibrational state n and n V LL;k 00 are coefficients of the expansion of the averaged molecule-cage interaction potential in bipolar harmonics of the angular variables.For a homonuclear molecule, only even k terms enter in the expansion of the potential in eqn (5) because of the symmetry of the potential under inversion.In this work the attention is restricted to H 2 @C 60 in the ground vibrational state n = 0 and the label n will be dropped henceforth.
It is a general statement of quantum mechanics that the total angular momentum of a system with spherical symmetry commutes with the Hamiltonian, so providing ''good'' quantum numbers for labelling the energy levels.For a diatomic molecule in a spherical trap, the total angular momentum K = L ˆ+ J ˆis given by the sum of the orbital angular momentum of the centre of mass L ˆand of the end to end rotational angular momentum J ˆ.The bipolar spherical harmonics F lj lm l are eigenfunctions of the angular momentum operators K2 and Lz with eigenvalues l(l + 1) and m l , respectively.They provide a complete basis set for the coupled angular motion.
A complete basis set for the radial motion is provided by the set of radial eigenfunctions for the Schro ¨dinger equation of a particle in a central potential.Early ab initio determinations of the potential energy surface of H 2 @C 60 suggest to use the eigenfunctions of a three-dimensional isotropic harmonic trap.These are discussed in the appendix, Section A.
where the coefficients c do not depend on m l .In the following we will use the notation E l to denote the energy level comprising the set of 2l + 1 degenerate states of eqn (8).The symmetry under inversion implies that parity is a good quantum number.Since the parity under inversion for the bipolar harmonics states is (À1) l+j , it follows that either even-or odd-sum l, j pairs appear in the expansion (8).

INS transition probability
In a typical neutron scattering ''spectrum'' the intensity of a transition between an initial energy level i and a final energy level f is determined by the double differential cross section: 30 where P i is the population of the energy level i, @s @O f i represents the transition probability for the scattering of neutrons and the delta function d enforces the energy conservation in the scattering process.As shown with greater detail in Section A.4, the INS transition probability for homonuclear molecules in unpolarised experiments can be written as a product of a spin s fi and a space part S fi : Specifically in the case of molecular hydrogen, I = 1/2 and I N = 0 for para-hydrogen and I N = 1 for ortho-hydrogen, so that the spin dependent part of the INS cross section becomes where o and p stands for ortho and para, respectively and where b c and b i are the bound coherent and incoherent scattering lengths of bound 1 H, respectively.The space dependent part S fi for transitions between an initial energy level i = E l and a final energy level f ¼ E l 0 0 is obtained by tilting the reduced matrix element of the INS transition probability from the spherical basis to the appropriate eigen system basis: Âhn 0 ; l 0 ; j 0 ; l 0 h g n; l; j; li 2 where the INS reduced matrix element is given by: hn 0 ; l 0 ; j 0 ; l 0 h g n; l; j; li ¼ X l 0 þl a¼jl 0 Àlj This journal is © the Owner Societies 2016 Interestingly, the derivation of the transition probability via the Wigner-Eckart theorem, as discussed in Section A.4, allows us to recover the selection rule for INS transitions in diatomic molecules confined in potentials of spherical symmetry, originally found by Xu et al. 20,31 and later extended and explained by Poirier in the context of group theory. 26Indeed in the evaluation of the reduced matrix element of eqn (13) the terms with are exactly 0 since both the Clebsch-Gordan term C g0 a0b0 and the reduced matrix elements of bipolar harmonics involve 3-j symbols with 0 projections, see eqn ( 29) and (30).These 3-j symbols are zero if the sum of the arguments are odd, see eqn (20).When the initial level is characterized by l = 0 total angular momentum, l = j and g = l 0 because of the triangular condition for the rows and column of the 9-j symbols implied by the reduced matrix element of the bipolar harmonic, eqn (29).By subtracting the last two identities out of the first one in 14, it follows that the reduced matrix element hn 0 ,l 0 ,j 0 ,l 0 8h AE g 8n,l,j,li is zero when the initial level has l = 0 and the final level angular numbers satisfy l 0 À l 0 À j 0 = odd.The selection rule for the transition probability follows from the observation that l 0 + j 0 is either even or odd for all the coefficients in eqn ( 12) as a consequence of the inversion symmetry for diatomic molecules in a spherical trap, see Section 4. The equivalent selection rule for final states with l 0 = 0 follows exchanging the role of primed and nonprimed labels in the discussion above.

Calculations
The methodology described in Sections 4 and 5 can be applied to determine the INS spectrum of H 2 @C 60 according to the following prescription.Firstly a finite coupled basis set is chosen by fixing the maximum number of translational and rotational excitations to n max and j max , respectively.The para (even j) and the ortho (odd j ) coupled states |n,l, j,l,m l i are gathered separately in blocks with the same angular momentum l.For each spin isomer and each l-block, the finite matrix basis representation of the Hamiltonian ( 5) is built in the basis |n, j,l,l,m l = 0i using eqn ( 7), ( 22) and ( 23).The rotational constant in the ground vibrational state is written as B = B 0 À D 0 j( j + 1) where the coefficient D 0 accounts for centrifugal corrections to the rotational energy. 32he potential expansion is truncated to include the harmonic term V 00;2 00 , two anharmonic terms, V 00;4 00 and V 00;6 00 and one translational-rotational coupling term V 22;2 00 .The blocks with different angular momentum l are diagonalized separately.The eigenvalues are used to determine the thermal populations, the levels and the transition energies that fall within the observation range.Assuming no spin-isomer interconversion, the ortho to para ratio can be fixed to 3 and the thermal equilibrium populations can be determined for the two spin isomer independently.The eigenvector coefficients, eqn (8), are used to calculate the total intensity of the INS transitions as a function of q, by evaluating the double differential cross section, eqn (9), via eqn (11) and (12).For numerical simulation the values b c = À3.74Â 10 À15 m and b i = 25.27Â 10 À15 m have been used as scattering lengths for the coherent and incoherent scattering length of bound 1 H, respectively. 33o achieve a more direct visual comparison between experiments and simulations, the transition probability, eqn (10), is integrated over the q-range appropriate to the instrumental apparatus at energy DE fi = E f À E i .In an inverse geometry neutron spectrometer the final energy of the scattered neutron is fixed to E n , while the incident neutron energy is scanned through.A detector at angle y with respect the direction of the incident neutron beam records scattered neutrons at transferred momentum so that integrating over q across the polar angular opening of the instrument The integration can be effectively performed averaging over an appropriate angular distribution describing the geometry of the apparatus.This procedure gives rise to a stick peak for each possible transition.To emulate finite width transitions the peaks are finally convoluted with the instrumental resolution, assuming negligible intrinsic linewidth in H 2 @C 60 .
In the simulation of the IN1-Lagrange spectrum of H 2 @C 60 the number of scattered neutrons was determined by averaging over a uniform distribution of 20 polar angles in the nominal angular ranges 341 r y r 701 at fixed E n = 4.5 meV and using the instrumental resolution discussed in Section 2.

Discussion
The INS spectrum of H 2 @C 60 was simulated using a set of purpose-built packages in Mathematica. 34In the actual calculations, the basis was defined by n max = 11 and j max = 7, giving rise to 10 192 para states and 13 104 ortho states, of which only 800 distinct reduced para states and 972 reduced ortho states were necessary for the evaluation of the Hamiltonian and of the INS transition probability.The initial rotational parameters were fixed to the values of free H 2 : B 0 = 7.357 meV and D 0 = 0.0056 meV. 35,36 2MV 00;2 00 1 4 , V ˜00;4 00 = r 4 V 00;4 00 , V ˜00;6 00 = r 6 V 00;6 00 and V ˜22;2 00 = r 2 V 22;2 00 were varied until satisfactory agreement with the experimental spectrum was found.The optimized values are reported in Table 1.
Fig. 3 shows the comparison between the experimental and calculated spectrum of H 2 @C 60 at 2.5 K.
The calculated spectra confirm the peak assignment up to 73 meV given in Section 3, that were mainly based on the approximation of uncoupled rotation and translation, see eqn (7).However even at low energy the presence of couplings has a major impact on the appearance of the spectrum.By providing intensities and transition energies as well, the calculations establish a firm framework for the analysis of the spectrum.
The structure underlying multi-component peaks is clearly resolved in the stick spectrum.The experimental peak at 22 meV is resolved into three main components corresponding to the intra-ortho translational transitions (0,0,1,1) -(1,1,1,l) with l = 1, 2, 0 in order of increasing energy.As shown in ref. 37  using group theoretical methods and successively in ref. 25 by using perturbation theory, the order of the l levels in the (1,1,1,l) multiplet is dictated by the positive sign of the translational-rotational coupling term V ˜22;2 00 .As observed previously the l = 0 component is much weaker in intensity than the other two and, importantly, these observations are in agreement with the simulation of the H 2 @C 60 spectrum by Xu et al. 15,20 In a similar manner, the experimental peaks at 37.4 and 38.5 meV are assigned to the l = 2 and l = 0 components of the paraortho translational transition (0,0,0,0) -(1,1,1,l), respectively.There is no evidence of the l = 1 component of the triplet.This highlights a striking prediction from the earlier papers, 20,21 namely the existence of a selection rule for INS transitions in diatomic molecules in a spherical confinement. 26,31ransitions towards levels with two quanta of translational excitations have been identified in the experimental discussion.The peaks at 32.3 and 34.3 meV are assigned to the transitions (0,0,1,1) -(2,2,0,2) and (0,0,1,1) -(2,0,0,0), respectively.Transitions towards states with n = 2 and j = 1 reveal a larger degree of mixing between different translational-rotational states.The calculations identifies the transitions towards (2,2,1,2) and (2,2,1,3) with the peaks at 45.8 and 47.0 meV, respectively.Concealed under the latter peak, there is a transition to an energy level comprising a superposition of the (2,0,1,1) and (2,2,1,1) states with a ratio close to 1 : 2. Similarly the peak at 49.8 meV corresponds to a transition to an energy level comprising a superposition of the (2,0,1,1) and (2,2,1,1) states in a ratio close to 2 : 1.The partner transitions to these superposition states originating from the para ground state are found at 61.3 and 64.2 meV, respectively, while the para transition towards (2,2,1,2) is forbidden and the para transition towards (2,2,1,3) has negligible intensity.
The peak centered at 51.5 meV is resolved in three transitions starting from the ortho ground state towards (1,2,1,l) state with l = 2, 3, 1 in order of increasing energy.
The stick spectrum demonstrates the existence of para-para transitions that are easily overlooked in the experimental spectrum.Owing to the small ratio between the coherent and incoherent cross section for 1 H, para-para peaks generate a very low scattered neutron count and are difficult to observe.For example the fundamental translational transition in the para manifold (0,0,0,0) -(1,1,0,1) is shadowed by the ortho transitions (0,0,1,1) -(1,1,1,l).
At higher energy transfer the absolute resolution becomes worse, the number of accessible transitions increases and the peaks become relatively less sparse and less intense.This is particularly evident in the region above 96 meV.The calculations identify 162 transitions originating from the ground para state and 191 transitions originating ground ortho state in the energy range between 96 and 215 meV.In contrast to the low energy case, where most of the levels are represented by the nominal (n,l,j,l) state with at least 95% probability except where noted, states with different n, l, j become more and more mixed at high energies, so hampering their use as labels for the energy levels.A semi-empirical approach that tries to fit more prominent features to a minimal number of Gaussian lines may still work where few high intensity transitions are clustered. 8However in other regions where transitions are more spread in energy and there is no dominant component, it is practically impossible to produce any meaningful assignment of this type.Simulations remains the only method to make sense of the complexity of the spectrum at high energy transfer.Indeed the simulated spectrum provide a realistic representation of the experimental observations over all the accessed energy range.The simulated transitions follow closely but not exactly the experimental ones.Higher order spherical terms may improve the match at the expense of increasing the parameter space.However improved accuracy appears to be not worthwhile with the current experimental resolution also considering that shifts and splittings of the m l components in the energy levels of order of fractions of meV are introduced by terms of icosahedral symmetry 38 as well as intermolecular effects. 39ble 1 Summary of the optimized values of the Hamiltonian parameters used for calculating the INS spectrum of H 2 @C 60 in Fig. 3 Hamiltonian parameter Optimized value B 0 /meV 7.375 D 0 /meV 0.012 r/Å 0.311 V ˜00;4 00 /meV À3.90 V ˜00;6 00 /meV 1.52 V ˜22;2 00 /meV 8.9

Conclusions
In this paper the low-temperature INS spectrum of H 2 @C 60 covering excitations up to 215 meV was presented.The assignment of the experimental peaks was aided by developing a model which makes use of the nearly spherical symmetry of the confinement.The guest-host interaction was expanded in terms of coupled angular spherical bipolar harmonic functions.The effective potential was then described by few physical parameters comprising a harmonic term and terms describing anharmonic corrections and translational-rotational couplings.This approach affords a good flexibility in representing the overall non-bonding interaction between the hydrogen molecule and the fullerene cage when compared with other interaction models based on interactions between atomic units.The effective potential then represents a benchmark for any ab initio determination of the PES in H 2 @C 60 .The model does not involve any numerical integration: matrix elements and transition probabilities are shown to integrate into closed algebraic expressions in the coupled rotational translational oscillator basis.The analytical expressions may be useful in the analysis of the q-dependence of the INS transitions probability in similar systems.More importantly, the numerical evaluation of these expressions comes at a relatively inexpensive computational cost.Starting from the uncoupled harmonic system, the effect of each term on the spectrum can be varied independently until the spectrum is reproduced correctly.The agreement between experiments and calculations is overall very good.The INS transitions are well reproduced indicating that the icosahedral non spherical terms give rise to minor splittings and shifts within the instrumental resolution.
The theoretical and computational model developed in this paper may be useful for investigating the quantum dynamics and the INS spectrum of confined molecular hydrogen in systems such as clathrates hydrates, 40,41 zeolites 42,43 and MOFs. 9,44,45We anticipate its extension to treat heteronuclear diatomic molecules as well as water molecules confined in C 60 .In systems with high symmetry or when the effect of non spherical multipoles are small, such as H 2 @C 60 , realistic INS spectra can be generated promptly.For systems with lower symmetry the method can be used to provide a starting approximate description useful to identify the main features of the INS spectrum and it can be eventually augmented to include non-spherical terms perturbatively.

Appendix A: mathematical appendix
In this section we provide an extended presentation of the mathematical aspects that lie at the basis of the algebraic and numeric treatment developed for the analysis of the INS spectroscopy of H 2 @C 60 .The solution of the Schro ¨dinger equation for the spherical oscillator provides the radial wave-functions for the coupled translational-rotational states.The matrix elements for local operators appearing in the Hamiltonian and in the INS probability are separated into a radial and an angular part.The radial matrix elements are integrated into closed algebraic forms.The angular matrix elements are evaluated using the Wigner-Eckart theorem.Finally the INS transition probability for diatomic molecules in spherically symmetric confinement is derived.
A.1 The isotropic three-dimensional harmonic oscillator Following the treatment presented in ref. 46 and the relations between the hypergeometric confluent function 1 F 1 and Laguerre polynomials, 47 the normalized eigenfunctions of the Schro ¨dinger equation for the three-dimensional spherical oscillator with mass M and angular frequency o are written as the product of a radial and an angular term: is the scale length associated to the harmonic oscillator, L (a) n is a generalized Laguerre polynomial 47 and is a normalization constant.The principal quantum number n is a non-negative integer that determines the energy, of the eigenfunction in eqn (17).The orbital angular momentum quantum number l is a non-negative integer.For a given n the angular momentum quantum number is limited to the range l = n, n À 2,. .., 3, 1 for odd n and l = n, n À 2,. .., 2, 0 for even n.
The degeneracy of the level with energy E n is g n = n(n + 1)/2.

A.2 Clebsch-Gordan coefficients and 3-j symbols
Clebsch-Gordan coefficients appear naturally in quantum mechanics when adding the angular momentum of two subsystems as in formula (3).3-j symbols are a more symmetrized version of the Clebsch-Gordan coefficients.The two symbols are related by The 3-j symbols are invariant under even transpositions of the columns and take the sign (À1) a+b+g under odd transposition.3-j symbols are 0 if m g a m a + m b and if the triangular condition a À b r g r a + b is not satisfied.It is worth to note the following parity property of the 3-j symbol a b g Similarly 6-j symbols and 9-j symbols appear when combining the angular momentum for three subsystems and four subsystems.9-j symbols appear in the evaluation of the reduced matrix element of two coupled subsystem as well, see eqn (29).They are defined in terms of 6-j symbols as: where max(|a À j, |b À f|, |d À h|) r x r min(a + j, b + f, d + h).9-j symbols are zero if the triangular conditions is not satisfied by any of the rows or of the columns.3-j symbols and 6-j symbols can be calculated from generating formulae 27 and are available as built-in objects in computational softwares such as Mathematica. 34

A.3 Evaluation of matrix elements
The evaluation of matrix elements for a function which is a direct product of a radial and a coupled angular part, such as the potential term in eqn (2) and the INS cross section in eqn (10), are evaluated as direct product of a radial and an angular part: A.3.1 Radial matrix elements.The evaluation of radial matrix elements for non-negative power of local operator f (R) = R k over the radial translational wave-functions, |n,l;ri of eqn (17), is found by integration where G is the gamma function 47 and are the coefficients for the Laguerre expansion of L (a) n (x) in power of x m .
Similarly the radial matrix elements of the spherical Bessel function 47 For sake of completeness, we note that only matrix elements of spherical Bessel functions j b with |l 0 À l| r b r l 0 + l and l 0 + l + b = even needs to be evaluated.In such case the confluent hypergeometric function 1 F 1 reduces to: 47 where (a) n = a(a + 1)(a + 2). ..(a + n À 1), (a) 0 = a is the Pochammer symbol. 47.3.2 Angular matrix elements.The evaluation of the matrix elements of bipolar harmonics over the coupled spherical basis is easily achieved via the use of the Wigner-Eckart theorem 48 where the round brackets represent a 3j symbol and the double bars represents the reduced matrix elements for the bipolar spherical harmonics.The reduced matrix elements in the coupled angular momentum scheme 48 is given in terms of reduced matrix elements for spherical harmonics 48 The angular matrix elements in eqn (28) are diagonal in l and m l for bipolar harmonics with L = m L = 0:

A.4 INS transition probability
In neutron scattering experiments involving transitions between two energy levels in a homonuclear diatomic molecule, labelled here as i and f respectively, the main quantity of interest is the differential cross section in Born approximation 30 where the sum is over the nuclear species x in the system at position x x = R + (À1) x (r/2), k and k 0 are the wave-vectors labelling the incident and scattered neutron wave-function, q = k À k 0 is the transferred momentum, b ˆis the scattering length operator and For nuclei with non zero spin I, the scattering length operator is defined by 33 where s ˆand I ˆare the neutron and nuclear spin operators, respectively and b c and b i are the bound coherent and incoherent scattering length.When evaluating the neutron cross-section for homonuclear diatomic molecules, only the symmetric combination intervenes in transitions where the total molecular nuclear spin does not change, DI N = 0, while the antisymmetric combination intervenes in transitions where the total molecular nuclear spin changes of one unit, |DI N | = 1.In experiments with non polarised neutrons and non polarised molecular system, the cross section can be written as product of a spin dependent part s fi and space dependent part S fi : (35)   when I f N ' I i N and  (36)   when I f N AE 1 ' I i N .The spatial dependent part of the INS cross section involves the evaluation of the matrix elements of the local operator g AE (q) where q = k À k 0 is the transferred momentum.The plus and minus sign are alternatively required for direct (ortho-ortho and para-para) and cross spin transitions (ortho-para and para-ortho), respectively.For our scopes, this operator can be conveniently written as g AE ðqÞ ¼ ð4pÞ Ã y q ; f q h AE gmg (37)   where * stands for complex conjugate and i a j a ðqrÞY ama y R ; f R ð Þ Y ama Ã y q ; f q (39)   and the expansion of the product of two spherical harmonics into a single spherical harmonic (Clebsch-Gordan series) 27 Y ama y q ; f q Y bm b y q ; f q ¼ X aþb g¼jaÀbj gmg amabm b Y gmg y q ; f q (40)   together with the definition of bipolar harmonics, eqn (3).The shorthand notation j b,AE (x) = [j b (x) AE j b (Àx)]/2 has been introduced in eqn (38).
In systems with spherical symmetry, the space dependent part of the INS cross section for transitions between an initial energy level i = E l and a final energy level f ¼ E l 0 0 is given by the transition probability of the operator g AE (q) summed over the final states m l 0 and averaged over the 2l + 1 initial states m l : X l 0 m l 0 ¼Àl 0 E; l; m l jg AE ðqÞjE 0 ; l 0 ; m l 0 h i j j X l 0 m l 0 ¼Àl 0 c E 0 ;l 0 n 0 ;l 0 ;j 0 ðrÞ h i Ã c E;l n;l;j ðrÞ Â X 1 g¼0 X g mg¼Àg Y gmg Ã y q ; f q Â n 0 ; l 0 ; j 0 ; l 0 ; m l 0 ; r h AE X g mg¼Àg Y gmg Ã y q ; f q Y gmg y q ; f q The condition 20 together with the identity j b,AE (x) = [1 AE (À1) b ] j b (x)/2 valid for integer b values allows to drop the AE sign from eqn (43) and (44).

Fig. 1
Fig.1The INS spectrum of H 2 @C 60 at 2.5 K as recorded on the IN1-Lagrange spectrometer at ILL, Grenoble (France).The vertical axis represents neutron counts in arbitrary units.

Fig. 2
Fig. 2 Expansion of the INS spectrum of H 2 @C 60 at 2.5 K of Fig. 1 in the range 10-75 meV.Each peak is assigned to an INS transitions for which the labels indicate the value of n, l, j in the final states.Transitions originating from the para-H 2 ground state are in blue and transitions originating from the ortho-H 2 ground state are in red, respectively.The vertical axis represents neutron counts in arbitrary units.
1.The coupled harmonic translationrotation wave-function is written in ket notation as n; l; j; l; m l ; r j i ¼ c n;l ðR; rÞ |fflfflfflfflffl ffl{zfflfflfflfflffl ffl} jn;l;ri F lm l lj |ffl{zffl} jl;j;l;m l i (6) using the notation of Section A.1 and the definition of bipolar harmonic of eqn (3).Here r ¼ ffiffiffiffiffiffiffiffi ffi h Mo r is the scale length associated to the harmonic oscillator with mass M and frequency o.The ket in eqn (6) is an eigenstate of the spherical three-dimensional oscillator represented by the Hamiltonian in eqn (5) when terms with L = 0 and up to k = 2 are considered.In this case o ¼ ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi V 00;2 00 .ð2pMÞ r .The zero-order energy of the undecoupled translator and rotator is E 0 n,j = h o(n + 3/2) + Bj( j + 1) (7) Anharmonic terms with L = 0, k 4 2 mix states with different n while translation-rotation couplings are characterized by L = J 4 0 and mix states with different n, l, j.In presence of anharmonicities and translation-rotation couplings the eigenfunctions of the Hamiltonian given in eqn (5) are characterized by definite angular momentum quantum number l and m l and can be written as E; l; m l j i¼ X n;l;j c E;l n;l;j ðrÞ n; l; j; l; m l ; r j i

Fig. 3
Fig.3Comparison between the experimental and calculated spectrum of H 2 @C 60 .For clarity the plot is split into two equally wide panels spanning the transfer energy ranges from 10 meV to 215 meV.The black line with error bars represent the spectrum of H 2 @C 60 recorded at 2.5 K on IN1-Lagrange at ILL.The bars represent the calculated stick spectrum of H 2 @C 60 at 2.5 K, where the colors blue and red are used to distinguish transitions originating in the para and ortho ground states, respectively.The orange line overlayed to the experimental spectrum represent the calculated spectrum obtained by convolving the stick spectrum with the instrumental resolution of IN1-Lagrange.For clarity the peaks in the stick spectrum are scaled by a factor 0.5 with respect to the numerical values and vertically shifted.
j b (qR) for integer b Z 0, which enter in the space part of the INS transition matrix, can be integrated in a closed form n 0 ; l 0 ; r j b ðqRÞ n; l; r

a0b0 Â j a ðqRÞj b;AE qr 2 27 e iqÁr ¼ 4p X 1 a¼0
expansion of a plane wave into spherical harmonics and spherical Bessel functions j a