Absorption and Emission Spectra of a Flexible Dye in Solution: a Computational Time-Dependent Approach
The spectroscopic properties of the organic chromophore 4-naphthoyloxy-1-methoxy-2,2,6,6-tetramethylpiperidine (NfO-TEMPO-Me) in toluene solution are explored through an integrated computational strategy combining a classical dynamic sampling with a quantum mechanical description within the framework of the time-dependent density functional theory (TDDFT) approach. The atomistic simulations are based on an accurately parametrized force field, specifically designed to represent the conformational behavior of the molecule in its ground and bright excited states, whereas TDDFT calculations are performed through a selected combination of hybrid functionals and basis sets to obtain optical spectra closely matching the experimental findings. Solvent effects, crucial to obtain good accuracy, are taken into account through explicit molecules and polarizable continuum descriptions. Although, in the case of toluene, specific solvation is not fundamental, the detailed conformational sampling in solution has confirmed the importance of a dynamic description of the molecular geometry for a reliable description of the photophysical properties of the dye. The agreement between theoretical and experimental data is established and a robust protocol for the prediction of the optical behaviour of flexible fluorophores in solution is set.
1 Introduction<!>2.1 General approach<!><!>2.2 QM calculations<!>2.3 FF parameterization<!>2.4 MM optimizations and MD simulations<!>2.5 Statistical calculation of the UV-VIS spectra<!>3.1 NfO-TEMPO FF parameterization<!>3.2 Simulations<!>3.3 Absorption spectra<!>3.4 Emission spectra<!>4 Conclusions
<p>Organic dyes1 have recently attracted a growing attention, thanks to their wide range of applications, that covers different fields in both materials2–11 and life sciences.12–15 In particular, the high selectivity and specificity of some fluorescent molecules in detecting particular molecular targets and/or in modulating the responses to different external stimuli are often exploited to investigate the structural and dynamic properties of a wide variety of complex systems like, e.g., polymeric materials, DNA structures, protein conformations or lipid aggregates.13 In this context, the strong influence of the environment on the fluorescence properties turns out to be a fundamental key in several applications aimed to the design of new technological devices and biomimetic systems. For example, the encapsulation of a dye in a nanoparticle can enhance its fluorescence by altering the dye's flexibility and protecting it from external quenching agents,5 whereas in mechanochromic materials the optical properties of the dyes are transferred to the whole system, through dispersion of the doping agent in matrices or by grafting it to polymer chains. Indeed, the interaction between the two species (i.e. the polymer and the dye) can be tuned by mechanical stress, in order to properly control both the absorption and emission responses.8,9,11</p><p>Spectroscopic techniques are becoming the methods of choice for deeper investigations of the photophysical properties of composite materials.16–18 Unfortunately, often experimental measures alone are not sufficient to unravel the subtle interplay of several effects in determining the overall result. In such circumstances, their synergistic combination with computational approaches can be effectively used both to validate experimental assignments and to provide interpretations of the results in terms of chemical structure, conformational re-adaptation, and molecular motion. Although computational methods can be nowadays applied directly to the study of very complex systems,18 their profitable use also by non-specialists requires the definition, implementation, and validation of robust and general computational protocols. As a matter of fact, quantum mechanical (QM) calculations are unavoidable for spectroscopic characterizations and, in the specific field of electronic spectroscopy (both absorption and emission), Time Dependent Density Functional Theory (TDDFT) has undoubtedly become the most widespread method.19 Although some problems have yet to be fully solved for multiple excitations, charge transfer and Rydberg states, for conventional valence excitations TDDFT based on hybrid (possibly range-separated) density functionals shows a nearly optimal balance between accuracy (related also to the continuous improvement of functionals) and computational inexpensiveness.19–22 In this framework, the energy difference (or vertical energy, VE) between the ground (GS) and the electronically excited state (EES), which leads to the the VE absorption (λVEabs, computed at the GS equilibrium geometry) or emission (λVEemi calculated at the EES equilibrium geometry) wavelengths, can be calculated rather straightforwardly.19 Unfortunately, these quantities alone are not always sufficient to reproduce the experimental band shapes and they cannot be directly compared to the measured maximal absorption (λmaxabs) or emission (λmaxemi) wavelengths.19,23 In order to reproduce the overall line shapes and their complex features and to reach a more reliable comparison with experimental observations, it is mandatory to take nuclear motions (both of the chromophore and of its environment) into the proper account.23–27 In the case of flexible molecules in condensed phases, some sampling of large parts of the potential energy surfaces (PES) becomes mandatory5,15,28–30 and the most suitable approach is based, in our opinion, on molecular dynamics (MD) simulations followed by proper averaging (and smoothing) of TDDFT results for a sufficient number of representative snapshots. Implicit solvation methods, such as the Polarizable Continuum Model (PCM31), offer the undeniable advantage of correctly describing the average polarization of the environment and its effect on the chromophore without any explicit sampling of solvent coordinates, but they are, of course, unable to capture specific solvent effects especially in the cybotactic region.26,30 This has lead to the proposal of an integrated discrete/continuum model able to describe both short-range and bulk solvent effects in a fully dynamical framework. The so called GLOB model enforces non periodic boundary conditions (NPBC) around a rather large cluster of explicit solvent molecules by means of both PCM reaction field and specific short-range potentials minimizing boundary effects.32,33 Several applications in recent years confirm the reliability and robustness of the GLOB and other related approaches.34–42 Although the classical MD step can be performed also employing conventional periodic conditions in order to use standard codes (here Gromacs), NPBCs are more suitable for QM computations based on localized basis sets. Then the extraction of representative frames can be carried out following different schemes. Within the simple "mechanical embedding" (ME) scheme, for instance, only the dye coordinates are extracted from each snapshot, and the environmental effects are retained only in the solute geometrical variations induced by the presence of the surrounding medium. Conversely, in the "electrostatic embedding" (EE) scheme, the direct electrostatic interaction between the dye and the surrounding environment is accounted for also in the TDDFT calculations. Indeed the EE snapshots are obtained by extracting both solute and a relevant region of the surrounding medium from MD trajectories and the TDDFT calculation is thereafter performed by adopting a full QM description for the dye, a simple (possibly polarizable43–45) point charge description for the retained surrounding environment, and a PCM description of the rest of the solvent. The spectra obtained with "dynamic" approaches are then suitable for a direct comparison with the experimental line shapes: they retain information not only on the peak wavelengths, but also on the broadening and shape of the absorption and emission profiles. Moreover, the simulation of resolved spectra of dyes that are optically active within the visible region makes it possible to predict their actual color, as perceived by the human eye.30,46</p><p>In order to reproduce reliably the spectroscopic properties of composite systems, MD simulations should be based on very accurate and specific force fields (FFs). Although a lot of work has been done to describe GS structures and reliable FFs are available in the literature, much less effort has been devoted to the development of "ad hoc" FFs for the EES conformations.35,47–50 Recently, the Joyce51 parameterization protocol, devised to obtain accurate and specific FFs from QM computed data, has been extended to the automated treatment of both GS and EES.50</p><p>In the present work, the above mentioned "dynamic" approach is applied to compute the electronic spectra of a flexible dye, namely the 4-naphthoyloxy-1-methoxy-2,2,6,6-tetramethylpiperidine (NfO-TEMPO-Me) molecule (see Figure 2), in toluene solution, where solvent molecules are modeled explicitly. This compound is interesting because of its potential application in the design of "smart" polymers with mechanochromic attributes as proven by the experimental characterization obtained through spectroscopic techniques and reported in a recent paper.52–54 More specifically, the 2,2,6,6-tetramethyl-piperidine-1-oxyl (TEMPO) unit was used, in combination with functionalizing groups such as hydroxyl, benzyl, naphthyl (HO-, BzO- and NfO- TEMPO), for grafting the fluorescent probe to polyethylene derivatives via post-reactor modification.52,53 The NfO-TEMPO derivative, substituted with a methyl group at the grafting site, has been studied from both experimental and theoretical points of view in toluene solution.53,54 The theoretical analysis54 depicted a complex conformational landscape, containing six stable conformers, characterized by different values of the three flexible dihedral angles describing the relative orientation of the naphthyl moiety and the TEMPO ring system (see Figure 2). Although this flexibility would suggest that several conformers contribute to the overall spectral shape, previous results54 indicate that a conventional approach, employing vertical excitation calculations at the TDDFT/PCM level, performs a good job in reproducing experimental results at least concerning the positions of peak maxima in both absorption and emission spectra. Here the previous computational study is extended towards a statistical approach, with the two-fold aim of achieving a more realistic sampling of the NfO-TEMPO-Me dye and solvent conformations, and possibly gaining a deeper insight into the structural and dynamic features that tune the spectroscopic characteristics of the system. From a more general point of view, the present study represents a further validation of a multi-level computational approach. The good result obtained on this system encourages to prove the reliability of the estabilished method on more complex materials such as functionalized polymers52,53 and other types of composite systems.</p><!><p>The multi-level approach adopted in the present work is sketched in Figure 11 and can be sumarized as follows.</p><!><p>QM calculations are performed on the target molecule, to sample the conformational space of both GS and EES, identify their minima and compute both the energy and its first and second derivatives at stationary points. Further geometry optimizations are performed at selected values of soft variables in order to better define specific FF terms (vide infra).</p><p>QM data are used by Joyce in two separate FF parameterizations of GS and EES, respectively.</p><p>GS and EES parameterized FFs are employed in MD simulations of the target dye surrounded by explicit solvent molecules (toluene) employing periodic boundary conditions. The MD trajectories are sampled during the dynamics and representative solute geometries are extracted, with (EE) and without (ME) explicit solvent molecules.</p><p>All the resulting mechanically and electrostatically embedded structures are employed in TDDFT QM calculations to compute the transition energies to the first EESs adding to EE structures only a PCM descritpion of bulk solvent effects.</p><p>The transition energies obtained from the GS and EES snapshots are finally combined to yield the absorption and emission spectra, respectively.</p><!><p>The conformational search of NfO-TEMPO-Me in toluene solution was performed in a previous study.54 Starting from the most stable conformers identified in that investigation, DFT and TDDFT geometry optimizations in the gas phase, for the GS and EES respectively, were carried out using the PBE0 hybrid exchange-correlation functional and the N07D basis set.55,56 For these optimized structures, gradients, Hessian matrices and harmonic vibrational frequencies were calculated at the same level of theory. Moreover, QM energy scans along the δ1, δ3 and δ4 flexible dihedrals (see Figure 2) were performed by varying each dihedral angle from 0° to 360°, in steps of 30°, minimizing the energy without any constraint, but the scanned coordinate. All electronic transitions were computed by means of the TDDFT method on the ME and EE snapshots extracted from the MD simulations, considering the seven lowest transition for each frame. As far as the EE snapshots are concerned, toluene molecules extracted from the MD frames were replaced by fixed atomic point charges, whereas bulk solvent effects were accounted for through the PCM approach.31 The direct electrostatic effect of the solvent on the electronic property is entirely neglected in the ME scheme, by switching off both the implicit and the explicit solvation schemes. The performances of different functionals, namely PBE0 and CAM-B3LYP, in combination with N07D57 and cc-pvDz basis sets, were explored.</p><p>The Gaussian09 suite of programs58 has been used to carry out all QM calculations.</p><!><p>The intramolecular FF employed for the description of NfO-TEMPO-Me, Eintra, takes the standard form:
(1)Eintra=Estretch+Ebend+ERtors+EFtors+Enb
where the first three terms are represented by harmonic expressions
(2)Estretch=12∑μNskμs(rμ−rμ0)2;Ebend=12∑μNbkμb(θμ−θμ0)2;ERtors=12∑μNRtkμt(ϕμ−ϕμ0)2
In equation (2)
kμs, kμb, kμt and rμ0, θμ0, ϕμ0 are the force constants and equilibrium values for stretching, bending and stiff torsional internal coordinates,50 respectively. The potential energies governing rotations around dihedrals δ1-δ4 are described by sums of cosines, namely
(3)EFtors=∑μNFdihedrals∑j=1Nμkjμd[1+cos(nμjδμ−γμj)]
where kjÎ¼d is the force constant, Î´Î¼ the dihedral, nÎ¼j and Î³Î¼j the multiplicity and a phase factor for the jth cosine. NÎ¼ is the number of cosine functions employed for dihedral Î¼. The last term of equation (1) is computed as
(4)Enb=∑i=1Nsites∑j=1Nsites[Eijnbintra]
where the interaction between the ith and jth atoms is described through a 12-6 Lennard Jones (LJ) plus Coulomb charge-charge potential (Enb intra).</p><p>Both the solute GS and EES FFs were parameterized by the Joyce package,50,51 through a least square minimization of the functional Iintra
(5)Iintra=∑g=0NgeomWg[Ug−Egintra]2+∑K≤L3N−62WKL″(3N−6)(3N−5)[HKL−(∂2Eintra∂QK∂QL)]g=02
where Ngeom is the number of the sampled conformations, QK is the Kth normal coordinate and Ug is the QM computed energy in the gth geometry, with respect to the absolute minimum (g = 0). The QM Hessian matrix HKL and the FF Hessian are evaluated at g = 0. All the geometry weights Wg are set to the same unitary value; the diagonal and off diagonal weights of the Hessian matrix elements WKL″ are set to 5000 and 2500 respectively. More details on the minimization procedure can be found in previous works.50,51</p><p>As far as the interaction of the NfO-TEMPO molecule with the solvent is concerned, the standard 12-6 LJ plus charge-charge potential was adopted. For the solute, LJ intermolecular parameters were taken from the OPLS FF59,60 and the atomic point charges were obtained from the equilibrium GS and EES QM geometries through the CM5 procedure.61 Conversely OPLS literature parameters59 were employed for the solvent in all simulations, i.e. with the solute in its ground or excited state.</p><!><p>Full or partial energy optimizations at the MM level were performed using the conjugate gradient algorithm available within the Gromacs 4.5.4 package.62 Structures were minimized with a convergence treshold on the root-mean-square forces of 1 kJ·mol−1 · nm−1. Harmonic dihedral restraints of 1000 kJ·mol−1 were imposed on the scanned internal dihedrals when needed.</p><p>All the MD simulations were carried out with Gromacs4.5.4.62 The GS and EES FFs of the dye were employed for both simulations in vacuo and in solution. The system simulated in solution consisted of one NfO-TEMPO-Me molecule inserted in a cubic box and surrounded by 991 toluene molecules. Periodic boundary conditions were applied in all directions and the simulation time step was set to 0.5 fs. Long-range electrostatic interactions were treated with the particle mesh Ewald (PME) method, whereas a 14 Å cuto was applied to the van der Waals interactions. The production run was preceded by an equilibration phase during which the starting configuration was firstly energy minimized, to remove unfavourable steric interactions. The resulting configuration was subjected to a further equilibration protocol coupling the system to a thermal bath at temperature T = 300 K (with a coupling constant of 0.1 ps) and a pressure bath with pressure P = 1 atm (with a coupling constant of 1.0 ps) for a total period of 1 ns. Starting from the final structure, a 3 ns-production dynamics was carried out in the NPT ensemble. During both the equilibration and production phases Berendsen's weak coupling schemes63 were employed. Configurations were stored every 0.5 ps for further analysis. Simulations in the gas phase were performed in the NVT ensemble at 300 K and the total sampling time was 1.2 ns. Selected snapshots were extracted from the MD trajectories and used for the simulations of the spectra. The samples consisted of 80 and 200 configurations for the simulation in the gas phase and in solution, respectively. In the latter case, two sets of snapshots were extracted, for the ME and EE schemes, containing respectively only the solute or the solute plus all toluene molecules within a radius Rcut from its center of mass.</p><!><p>The transition energies computed on each set of snapshots were convoluted with gaussian functions in the energy domain with half width at half maximum (HWHM Δν) of 0.05 eV. Thus, the spectrum computed for the c-th frame is
(6)ϵc(ν)∝∑i∈statescfc,iΔνexp[(ν−νc,i0σν)2]
where hνc,i0 and fc,i are the energy and the oscillator strength of the i-th excited state, respectively, whereas σν=[22ln(2)]−1⋅Δν. For a vis-à-vis comparison with the experimental data, the energy spectrum is then plotted in the wavelength domain (ϵc(λ)). The signals originating from each snapshot are then averaged according to the equation:
(7)ϵ¯(λ)=∑c∈snapϵc(λ)Nsnaps
to achieve the final "statistical" UV-VIS spectrum for a given functional-basis choice.</p><!><p>As indicated in Section 2, NfO-TEMPO intermolecular parameters were transferred from OPLS59,60 FF (LJ) or obtained from QM data (CM5 point charges). The atom types chosen to model the NfO-TEMPO molecule, for both GS and EES, are reported in Figure 2. All intermolecular parameters can be found in the supporting information.</p><p>The intramolecular FF parameterization was performed by minimizing the functional (5), achieving standard deviations of 0.172 kJ/mol and 0.165 kJ/mol for GS and EES, respectively. Inspection of Table 1 (where a selected set of the optimal parameters is reported) shows that the parameterized FF is able to capture even subtle differences between similar internal coordinates, that could not be distinguished by general purpose FFs. In the GS parameterization, for instance, a non negligible difference is found between the C1-C2 and the C2-C3 aromatic bonds in the naphthalene moiety, which are described by slightly different equilibrium distances (1.42 and 1.38 Å) and remarkably different force constants (2954 and 3897 kJ/mol Å−2, respectively). Structure rearrangements taking place upon electronic transition are also accounted for, as indicated by the changes shown by FF parameters between the GS and EES FF. For instance the CC-CO equilibrium bond distance decreases from 1.49 Å to 1.45 Å and the relative force constant increases by ~ 500 kJ/mol Å−2, paralleling the partial π character that the bond assumes in the EES. It is worth noticing that in the EES, as a result of the increased delocalization of the naphthalene π cloud toward the carboxyl region, the naphthalene structure reduces its conjugate character and consequently its rigidity, as can be seen from the decrease by ~ 50 kJ/mol of the harmonic improper dihedral terms describing the out-of-plane deviations in the region.</p><p>Broadly speaking, the force constants relative to the internal coordinates of the carboxylnaphtalene group show rather substantial differences between GS and EES: changes as large as 22%, 47% and 31% are found for the stretching, bending and torsion harmonic terms, respectively. Instead, negligible modifications were found in the parameters relative to the TEMPO-Me region, with deviations on the stretching and bending harmonic terms lower than 2% and 8%, respectively. A tentative explanation lies in the fact that NfO-TEMPO-Me HOMO and LUMO orbitals are essentially localized on the carboxyl-naphtalene portion of the molecule,54 that can thus be considered the most relevant part in the description of the optical behavior of this compound.</p><p>The capability of the new FF to reproduce correctly the vibrational behavior of the molecule is quantified by comparing the normal mode frequencies computed at the QM and MM levels (see Figure 3). A good agreement is achieved, resulting in a root mean square deviation of 38 and 37 cm−1, for GS and EES, respectively. It may be worth noticing that the achieved accuracy is in line with the results obtained in previous parameterizations,50 hence a rather reliable representation of the NfO-TEMPO vibrational behavior is expected during MD simulations.</p><p>As far as the the two flexible dihedral angles, δ1 and δ3, are concerned, they are expected to affect the conformation of the NfO moiety, and therefore, for each dihedral, the reference QM torsional energy scan was performed separately for GS and EES. Conversely, for computational convenience, the QM δ4 scan was performed only for GS, and the resulting FF parameters were simply transferred from the GS to the EES FF. For all dihedrals, both FFs accurately match their QM reference energy profile. In the top panel of Figure 4 the QM energies calculated for the optimized geometries at a given torsional angle are compared to their MM counterparts, showing a very good agreement for both GS and EES. As could be expected considering the localization of the HOMO and LUMO orbitals,54 major differences between GS and EES profiles were found, at both QM and MM level, only for δ1, whereas the δ3 torsional curve does not undergo noticeable changes upon electronic transition. The different steepness around the δ1 minimum is also in agreement with the partial π character of the CC-CO bond revealed by the force constants analysis.</p><p>On the grounds of earlier studies,34,35 it is expected that the dihedral description strongly affects the optical properties of the molecule. The absorption and emission wavelength λVEabs∕emi have been calculated, at the VE TD-DFT level, for all the QM and MM geometries, obtained in the torsional energy scans. The results are shown in the bottom panel of Figure 4. As a matter of fact, both absorption and emission wavelengths for QM and MM optimized geometries show a strong dependence on the δ1 angle, which is responsible for the relative orientation of the naphthyl and carboxyl moieties within the NfO portion of the molecule. The values computed for the MM geometries are in good agreement with those identified by the QM calculations, with slightly larger differences (+ 6.0 nm for absorption, and + 11.1 nm for emission) in the case of δ1 = 180°. Nonetheless, the latter conformation is expected to be scarcely populated at room temperature, considering its relative stability (≃ 12 kJ/mol) and high interconversion barrier (≃ 60 kJ/mol) with respect to the absolute energy minimum (at δ1 = 0 deg). As long as δ3 is concerned, its effect on λVEabs∕emi is less pronounced, with some minor variations found in emission for δ3 ~ 30°. This scarce sensitivity is not surprising, considering that the HOMO and LUMO orbitals are localized in molecular regions not involved in the δ3 torsion. As for δ1, the FF mimics with rather good accuracy this behavior, with slightly larger deviations for absorption, found in high energy conformations (δ3 < −60°).</p><!><p>The GS and EES parameterized FFs were employed in two sets of MD simulations, performed on the isolated dye at 300 K and on the solvated NfO-TEMPO at 1 atm and 300 K. All four systems (i.e. GS and EES in vacuo and solvated GS and EES) were equilibrated for 1 ns, and subsequently simulated for 1.2 and 3 ns (for the in vacuo and solvated phase, respectively), saving trajectories every 0.5 ps.</p><p>The dye internal structure, along the four different simulations has been first evaluated by looking at the distributions of the two main flexible dihedrals, reported in Figure 5. As regards the GS simulation in gas phase (black dotted line) δ1 remains constrained to the main potential well, and its distribution is symmetrically broadened around 0°, the QM reference value being within this well, at about −12°. The solvated system shows similar distributions within δ1 ≃ 0° region, but the interaction with the solvent allows for exploring also the local minimum at about 180°. Conversely, due to the increased steepness around the δ1 minimum region found in the torsional profile, the EES population distribution is much narrower for both isolated and solvated systems, and the 180° region is never explored. This indicates that any possibility for δ1 to be trapped in the above-mentioned region (either by vertical excitation from the GS or internal convertion from higher EES) is neglected at this level. The δ3 distribution shows, instead, a double peak for both GS and EES simulations, with the first one corresponding to the QM value of about 85° and the second one located at δ3 ≃ 150° i.e. close to the second local minimum found in the torsional energy scan reported in Figure 4.</p><p>Furthermore, as the solvated systems are concerned, the local solvent structure around the dye has been monitored by computing the pair correlation function between the NfO-TEMPO and toluene centers of mass, for both GS and EES. Inspection of Figure 6 shows that some differences arise in the first solvation shell for GS end EES, the number of first neighbors being larger in the latter case. On the other hand, the differences become negligible for the second solvation shell and at distances larger than 13 Å both pair correlation functions converge to an homogeneous density. Clearly, a continuum model is unable to account for these local features, and for this reason a portion of the surrounding solvent was included in the snapshots extracted for the EE calculations. In particular, all solvent molecules within a radius Rcut of 15 Å were included explicitly. Indeed, beyond this distance, the pair correlation function has reached its asymptotic value, and the PCM can be confidently applied, as shown pictorially in the right panel of Figure 6. The included solvent molecules amount on average to 30, with a standard deviation of 2, for both GS and EES simulations.</p><!><p>Snapshots spaced by constant time intervals are then extracted from the MD simulations and used to calculate the UV/VIS spectra with an integrated TD-DFT/charges/PCM computational setup. In Figure 7, top panel, a typical unconvoluted absorption spectrum from a random frame is shown in black bars. A Gaussian broadening was applied according to equation (6) with two different values for Δν: 0.25 eV is the value used in the previous time-independent study of NfO-TEMPO-Me,54 whereas 0.05 eV is the value chosen for the present work. In the bottom panel, the absorption spectrum is obtained by averaging, over an increasing number of snapshots, the signals broadened with Δν = 0.05 eV for each frame. It should be noted that the overall spectrum results now from the sum of hundreds of Gaussian functions, centered at different positions, that converges to a smooth, averaged line when the sampling is reliable. The inclusion of one frame every 15 ps of MD is chosen for all the computations in the following, in order to achieve a smooth statistical description without artifacts arising from the sampling itself. The dynamic approach proposed in the present study red-shifts the absorption wavelength, λmaxabs, by about 2 nm, with respect to the static approach employed in a previuos work.54 Unfortunately, the static value (315 nm) did already overestimate the experimental λmaxabs by 16 nm, so that dynamical effects seem to worsen agreement with experiment. Since the choice of the functional / basis set could play a role in determining this trend, other functional/basis set pairs have been explored, eventually obtaining remarkably accurate results by combining the CAM-B3LYP functional with the cc-pvDz basis set (sse Figure 8). In this case, the static approach (black placemark in the graph) underestimates the experimental peak wavelength by 10 nm, but this error is corrected by the MD treatment with a 3 nm red-shift. The thermal broadening of the absorption line, obtained by including different snapshots, allows a more direct comparison of the spectra, resulting in a good agreement between the computed and the experimental54
λmaxabs (292 nm and 299 nm, respectively). Concerning solvent effects, it should be noted that solvatochromic shifts are quite negligible, leading to very similar spectra in the gas phase and in solution.</p><!><p>Emission spectra were computed at the CAM-B3LYP/cc-pvDz level with the same approach adopted for absorption and the resulting spectral lines are reported in Figure 9. Inspection of the figure shows that the emission profile of NfO-TEMPO-Me in toluene solution reproduces fairly well the experimental data, especially the line-shape and its decay, although the peak is blue-shifted by about 15 nm. At variance with absorption, the effect of the solvent on the spectral shape is quite significant. In this case the spectrum in the gas phase is less broadened and centered at a shorter wavelength (by about 17 nm) so that a non negligible solvatochromic shift is observed. It is worth noticing that the increased sensitivity to the solvent, found for EES with respect to GS, is consistent with the augmented flexibility of the napthalene skeleton (i.e. where HOMO and LUMO orbitals are essentially localized) evidenced by the FF and with the larger number of solvent molecule found in the first solvation shell by pair correlation function analysis. A deeper insight on this solvent effect can nonetheless be gained by comparing the spectra obtained with ME and EE methods. As it can be seen in Figure 9, the main contribution to the broadening is accounted for alreadty at ME level, whereas the inclusion of the direct electrostatic interaction in the TDDFT calculations, achieved with EE, only tunes the signal, evidencing two shoulders at ~ 355 and 370 nm.</p><p>Results obtained previously with a vertical excitation approach54 and with the present time dependent protocol are summarized in Table 2, whereas the final computed absorption and emission spectra are compared with their experimental counterparts in Figure 10. The dynamic approach has the clear advantage of providing a detailed band shape that can be directly compared with the experimental signals, moreover it can reproduce true maximal absorption and emission wavelength λmaxabs∕emi, instead of vertical energies λVEabs∕emi, which are closely related to the experimental values. Indeed, in the present case, a slightly better agreement with the experiments is also achieved. More important, the Stokes shift, which is defined as the difference between λmaxemi and λmaxabs, can now be computed more rigourosly, yielding a good agreement (8 nm) with the experimental one.</p><!><p>The time dependent protocol used in this work for the computation of absorption and emission spectra of the NfO-TEMPO dye in toluene solution was able to reproduce peak wavelengths, spectral line shapes and Stokes shifts in good agreement with the experimental data. With respect to the static approach, which was previously54 adopted by our group to study the NfO-TEMPO system, the description of the spectroscopic behaviour of the solvated dye is significantly improved, and a more reliable comparison with the measured spectra is now possible, due to the availability of the band shapes and the consequent correct definiton of λmaxabs∕emi.</p><p>Besides the results achieved, the choice of adopting this time dependent protocol can be a posterori validated in consideration of two important features. First, an accurate representation of the molecular flexibility is found to be important. This is because there is a strong dependence of the spectral response on the molecular conformation. As a consequence, the use of a specific and reliable FF is essential when classical MD techniques are employed for the statistical sampling. Furthermore, although strong specific solute-solvent interactions are not expected in the present case, an explicit account for the solvent molecules remains essential for a correct description of EES, where both ME and EE altered remarkably the in vacuo computed spectral line shape</p><p>As already observed in previous works, the present protocol is well suited also for the description of more complex systems,34,35,38,39 where the surrounding medium consists of advanced materials such as nano-particles, polymers, lipid bilayers and other types of molecular environments. Indeed, the extension of the present computational procedure to the study of the NfO-TEMPO dye in a polymeric matrix is currently in progress.</p><p>Finally, it may be important to note that, in the in silico design of such complex materials, the choice of a reliable combination of DFT functional and basis set is of fundamental importance. The present work suggests that the selection of the most suitable combination in the TDDFT calculations of the spectral properties should be driven by the comparison of the experimental spectra with broadened theoretical lines, rather than by single point vertical excitations. According to the more or less flexible structure of the target dye and/or to the strenght of the interactions with the surrounding medium, this can be done in a time independent fashion, by accounting for the vibronic structure,27 or through a time-dependent statistical approach, as done in the present work.</p>