Modeling electrostatic interactions in complex systems

DEPARTMENT OF CHEMISTRY, LIFE SCIENCES AND ENVIRONMENTAL SUSTAINABILITY Doctoral Programme in Material Science and Technology XXXIII Cycle MODELING SUPRAMOLECULAR ELECTROSTATIC INTERACTIONS IN COMPLEX SYSTEMS Coordinator: Prof. Enrico Dalcanale Tutor: Prof. Anna Painelli PhD student: Mattia Anzola 2017-2020 Contents Introduction 1 1 Resonant Energy Transfer: a dynamical approach 5 1.1 Introduction and plan of the work . . . . . . . . . . . . . . . . . . . . . . 5 1.2 Force Field validation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 1.2.1 NBD................................... 11 1.2.2 NR.................................... 13 1.2.3 Parametrization of the excited state FF: QM calculations . . . . . 15 1.3 Molecular Dynamics simulations . . . . . . . . . . . . . . . . . . . . . . . 18 1.3.1 The free DA pair in chloroform . . . . . . . . . . . . . . . . . . . . 18 1.3.2 The bound DA pair in chloroform . . . . . . . . . . . . . . . . . . 21 1.3.3 The bound DA pair in DMSO . . . . . . . . . . . . . . . . . . . . . 21 1.3.4 Conformational Analysis . . . . . . . . . . . . . . . . . . . . . . . . 24 1.3.5 Charateristics timescales . . . . . . . . . . . . . . . . . . . . . . . . 25 1.4 A dynamical model for excited state decay . . . . . . . . . . . . . . . . . . 28 1.4.1 DecayTimes .............................. 33 1.4.2 Simulations in DMSO . . . . . . . . . . . . . . . . . . . . . . . . . 38 1.4.3 Decay rates in static configurations . . . . . . . . . . . . . . . . . . 39 1.5 Conclusions................................... 43 2 Molecular Aggregates 45 2.1 Introduction................................... 45 2.2 Aggregates of polar and polarizable dyes . . . . . . . . . . . . . . . . . . 47 2.2.1 DA chromophores as polar and polarizable dyes . . . . . . . . . . . 47 I 2.2.2 The model for the aggregate . . . . . . . . . . . . . . . . . . . . . 50 2.2.3 Electronic Hamiltonian: the rotation on adiabatic states . . . . . . 51 2.2.4 Accounting for vibrations: the Lang-Firsov transformation . . . . 54 2.2.5 Computational Strategy . . . . . . . . . . . . . . . . . . . . . . . . 55 2.2.6 Results ................................. 57 2.2.7 Discussion................................ 63 2.3 Aggregates of non-polar molecules . . . . . . . . . . . . . . . . . . . . . . 65 2.3.1 The model Hamiltonian . . . . . . . . . . . . . . . . . . . . . . . . 66 2.3.2 Exact diagonalization approach: calculating absorption and fluorescencespectra ............................ 68 2.3.3 Exact results on finite size systems: validating the Heitler-London approximation ............................. 71 2.3.4 Testing approximation schemes . . . . . . . . . . . . . . . . . . . . 74 2.3.5 2Daggregates.............................. 89 2.4 Conclusions................................... 96 3 Chiral aggregates of αand β-dicyanostylbenes: chiroptical properties 99 3.1 Introduction................................... 99 3.2 DCSB: experimental results . . . . . . . . . . . . . . . . . . . . . . . . . . 101 3.3 MolecularDynamics ..............................103 3.3.1 Definition of the force field . . . . . . . . . . . . . . . . . . . . . . 103 3.3.2 MD simulations on aggregates . . . . . . . . . . . . . . . . . . . . 105 3.4 Absorption and CD spectra of DCSB aggregates . . . . . . . . . . . . . . 107 3.4.1 Themodel ...............................107 3.4.2 Calculated Spectra . . . . . . . . . . . . . . . . . . . . . . . . . . . 109 3.5 Conclusions...................................111 4 Chiral aggregates of squraine dyes 113 4.1 Introduction...................................113 4.2 The three state model for squaraine dyes . . . . . . . . . . . . . . . . . . . 115 4.3 Aggregate model: electrostatic interactions . . . . . . . . . . . . . . . . . 118 4.4 The role of intermolecular charge transfer interactions . . . . . . . . . . . 123 II 4.5 Calculation of CD spectra of aggregates with delocalized electrons . . . . 129 4.6 Molecular Dynamics modeling of chiral squaraine aggregates . . . . . . . . 133 4.6.1 Definition of a force field . . . . . . . . . . . . . . . . . . . . . . . . 133 4.6.2 Enhanced sampling MD on tetramers . . . . . . . . . . . . . . . . 133 4.6.3 Selection of statistically representative structures . . . . . . . . . . 134 4.7 Chiral aggregates of squaraine dyes: absorption and CD spectra . . . . . . 135 4.8 Conclusions...................................136 General conclusions 143 A Molecular Aggregates 145 A.1 Oscillator Strength: sum rule . . . . . . . . . . . . . . . . . . . . . . . . . 145 A.2 Aggregates of polar dyes: additional results . . . . . . . . . . . . . . . . . 147 B Molecular dynamics 151 B.1 Foundamental parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . 151 B.2 ForceFields...................................153 B.3 Ensembles....................................153 B.4 Enhanced Sampling Techniques . . . . . . . . . . . . . . . . . . . . . . . . 155 B.4.1 Umbrella Sampling . . . . . . . . . . . . . . . . . . . . . . . . . . . 155 B.4.2 Hamiltonian Replica Exchange . . . . . . . . . . . . . . . . . . . . 156 B.5 Simulation conditions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 157 B.5.1 Chapter1................................157 B.5.2 Chapter3................................158 B.5.3 Chapter4................................158 Acknowledgements 161 Bibliography 163 List of Publications 177 III Introduction Every known interaction between objects or particles can be ascribed to one of the so called four fundamental forces: strong, weak, gravitational, and electromagnetic. They are separated according to the average strength of the force, the kind of particles involved in the interactions, and the range of effectiveness. The eletromagnetic force, and specifically electrostatic force in chemistry, represents undoubtedly the most important. Almost the totality of transformations involving atoms and molecules, from the freezing of water to complex energy transfer mechanisms, obey its rules. Electrostatic interactions depend on attractive or repulsive forces between atoms or molecules generated by their electrical charges. This interaction is also known as Coulomb interaction, named after physicist Charles-Augustin de Coulomb, who first characterized it in 1785.[1]. Even though the majority of non-bonded forces are electrostatic, chemists usually discriminate between them on the basis of the types of charges involved. Charged particles are responsible for the strongest non-bonded forces, and include electron-nucleus, ion-ion and ion-dipole interactions. Conversely, weaker electrostatic forces manifest between non-charged objects, such as polar or polarizable molecules, and comprehend, among many others, dipole-dipole interactions, Van der Waals and London dispersion forces. In this work we will focus on intermolecular interactions that play a major role in the definition of properties in molecular materials. The study of novel molecular materials provides great opportunities to solve a large amount of needs. A careful arrangement of molecules at the nanoscale can lead to the formation of systems with unique optical, electrical and magnetic properties. Molecules arranged in supramolecular aggregates interact through weak non-covalent forces, giving rise to collective electronic states that radically change spectroscopic features of these systems. A large variety of fields exploit 1 INTRODUCTION this high versatility, from organic solar cell [2, 3, 4] to OLED technology.[5, 6, 7] Another noticeable example that leverages intermolecular electrostatic interactions is made up of the so called light harvesting complexes, acting in the very firsts step of photosynthesis, in which solar light is absorbed and efficiently carried as excitation energy through resonant energy transfer.[8, 9] These systems have been investigated for a long time since their artificial replications find numerous applications as power sources, sensor systems and nanofabricated devices.[10, 11, 12, 13] Of crucial importance in researching new materials, as well as in developing new features in existing one, is the theoretical modeling of these systems and the rationalization of their spectroscopic properties. Purpose of this thesis is the study and integration of novel methodologies to get a deep understanding of intermolecular electrostatic interactions in complex nanosized systems. In Chapter 1 we investigate F¨orster RET for a pair of dyes linked to a calixarene scaffold. In systems where the geometry of the chromophore pair is not known a priori, the rate of energy transfer is estimated considering two limiting regimes, static and dynamic, offering an oversimplified view of the mechanism. For this reason, we make extensive use of MD simulations in two different solvents to describe the effect of solvation and conformational motions on the rate of energy transfer. Moreover, we develop a fully dynamical model, based on Monte Carlo method, to analyze the characteristic timescales of such processes and compare them with the experimental picture. Subsequent chapters are all focused on molecular aggregates. In Chapter 2, we examine spectroscopic properties of molecular aggregates, testing new approaches and approximation schemes for polar and non-polar supramolecular assemblies. The first part is focused on the discussion of aggregates of polar and polarizable dyes, improving already existent models to account for vibrational coupling and hence for spectral band-shapes. We then turn attention to aggregates of non-polar chromophores, addressing the reliability of the Heitler-London approximation and presenting a model for two dimensional aggregates. Chapters 3 and 4 are focused on chiral aggregates. In Chapter 3 we investigate aggregates formed by dicyanostilbenes decorated with chiral pendants. Through the use of an hybrid approach, involving MD simulations and exciton modeling, we are able to get a deep understanding on both aggregation and spectroscopic features of these system, questioning the effectiveness of widely adopted 2 INTRODUCTION rules to assess the system chirality from Circular Dichroism (CD) spectra. In Chapter 4 we focus attention on non-symmetric squaraine aggregates. An extensive theoretical work is discussed, devoted to the study of spectroscopic features of squaraine assemblies in solution. While the possibility of these system to form stable Charge Transfer (CT) states has already been extensively discussed,[14, 15, 16] literature lacks of a methodology able to describe rotatory power induced by CT transitions. We here present a new model for the calculation of absorption and CD spectra of squaraine complexes using a delocalized electrons approach, taking into account for both intraand intermolecular charge transfer mechanisms. During my thesis period I spent 3 months at the International Centre for Theoretical Physics in Trieste under the guidance of Dr. Ali Hassanali and 3 months at the Ruder Boˇskovi´c Institute in Zagreb under the supervision of Dr. Luca Grisanti where I learnt and mastered all the Molecular Dynamics techniques implemented in this essay. 3 CHAPTER 1: RESONANT ENERGY TRANSFER: A DYNAMICAL APPROACH of the system were unsuccessful, with theoretically estimated timescales ranging from a few tenths to a few hundreds fs.[34] Here, while adopting the same approach to the VDA estimate, we are able to properly simulate RET timescales thanks to our fully dynamical approach to RET. We will demonstrate that the RET process is governed by a complex interplay between different competing dynamical processes that include not just the D radiative and non-radiative relaxation, but also the conformational and solvation degrees of freedom of the system that, modulating VDA on similar timescales as RET, cannot be neglected in the description of this dynamical phenomenon. We can summarize our roadmap as follows: •build a computationally reliable dynamical system for the D,A molecules in the presence of the calixarene (clx) scaffold •develop a fully dynamical model for the energy transfer process and how this compare with experimental picture •analyze the characteristics timescales for such processes to understand the role of different physical components 1.2 Force Field validation Because of the complexity of the system, we decided to study first the two free dyes in water (a very well studied solvent for MD). Then we proceeded according to the following scheme: (DA-c) Dand Aunbound pair in chloroform (clx-DA-c) Dand Aconnected to calix-[4]-arene in chloroform (clx-DA-d) Dand Aconnected to calix-[4]-arene in DMSO A reliable computational model of the above systems can be built by taking advantage of MD. In a MD simulations, atoms mutually interact, accounting for attractive and repulsive forces as well as for chemical bonds constraints, generating a potential field. 10 1.2 Force Field validation The numerical solution of the Newton’s equations of motion generates a trajectory, showing the dynamic evolution of the system. See Appendix B for more details. The first step in any MD investigation is the definition of an optimal Force Field (FF) for the system at hand. Typically, non-biological organic molecules of moderate size are modelled using force fields that includes all atoms and describe their interactions, starting from ab-initio optimized geometries. Commonly adopted FF are GROMOS, CHARMM and GAFF. GAFF will prove the most appropriate choice for our system. A specific requirement for our work is a good representation of the donor dye both in the ground state, Dand in the relaxed excited state D∗. Different charge distributions for the same molecule in the two different states can be readily converted into GAFF topologies, enabling a fine tuning of parameters. 1.2.1 NBD Figure 1.3: 4-amino-1-nitrobenzoxadiazole, NBD The donor molecule of our system, NBD, is a small molecule with a fairly rigid structure. However, some geometrical parameters must be monitored to ensure a realistic 11 CHAPTER 1: RESONANT ENERGY TRANSFER: A DYNAMICAL APPROACH behaviour during simulations. We specifically investigate the geometry of amino and nitro groups. After topology files were created for each FF, three different simulations at constant number of molecules, volute and temperature (NVT ensamble, see Appendix B) simulations were performed in a 5 ×5×5nm3water box. To test the goodness of the topology the Radial Distribution Function (RDF) of water oxygens around the nitro and amine nitrogen atoms is calculated (Fig. 1.4). RDF describes how the density of surrounding atoms varies as a function of distance from a point. It is usually determined as the number of atoms A (in our case water oxygens) at a distance rfrom a given atom B (nitro and amino nitrogen), normalized for the density of A atoms in the system: gAB(r) = 1 ρ(A)hX A δ(~rA−~rB)i(1.6) RDF results for GROMOS differ from those obtained by other FFs, giving an incorrect Figure 1.4: Left: radial distribution function of water oxygens around amino N; right: radial distribution function of water oxygens around nitro N. distribution of water oxygens around the amino nitrogen, while CHARMM and GAFF results are comparable (Fig. 1.4). We modified some parameters in GAFF topology file of the NBD molecule in order to improve the consistency of results obtained with different force fields. We expect the dyhedral angle distribution of both the nitro and amino groups to be peaked at 0◦, since 12 1.2 Force Field validation Figure 1.5: H-N-C-C dihedral angle distribution for different FFs. the strong donor/acceptor character of the NH2/NO2groups ensures conjugation and hence planarity. While the H-N-C-C dihedral (shown in Fig. 1.5) already satisfies our predictions, the O-N-C-C dihedral distribution presents a minimum at 0◦(Fig. 1.6 left). For this reason we increased the relevant force constant, from 2.51 to 4.4 Kcal·mol−1·˚ A−2 (i.e. setting it to the same value as for the H-N-C-C force constant). The resulting distribution obtained properly peaks at 0◦, as shown in the right side of Fig. 1.6. 1.2.2 NR For NR we run analogous simulations as described for NBD. As shown in the left panel of Fig. 1.8, the magnitude of the dipole moment shows different time profiles depending on the adopted FF. A correlation is observed between the oscillations of the dipole magnitude and the orientation of the methoxy group (see Fig. 1.8). In GROMOS, where the methoxy group is confined in the [−45◦,45◦] region, dipole module is fixed around 13 CHAPTER 1: RESONANT ENERGY TRANSFER: A DYNAMICAL APPROACH Figure 1.6: Distribution of dihedral angle O-N-C-C in different Force Fields. Left: before correction; right: after correction. Figure 1.7: methoxy-Nile Red, NR 6.7 D, while it strongly oscillates when the methoxy group can interchange between “open” and “closed” conformations (defined respectively as conformations where the O−CH3group is in “cis” or “trans” configuration with respect to the carbonyl group). We performed two different Quantum Mechanical (QM) calculation (b3lyp/6-31g(d,p), frozen geometry) on the “closed” and “open” conformation of the methoxy group to 14 1.2 Force Field validation Figure 1.8: Left: NR dipole moment magnitude as a function of time. Right: methoxy group dihedral distribution. assess the magnitude of the electric dipole in the two conformers. We obtain 6.2899 D in the closed confirmation and 7.8773 Din the open conformation, in good agreement with MD results. The good agreement with DFT supported our choice of GAFF as the reference force field, with the minor correction of the force constant relevant to the OMe group that we set to 9.17 kJmol−1nm−2, to be compared with the orginal value 3.77. With this new constraint the methoxy group is still able to rotate but the distribution is more pronounced in the 0◦and 180◦region (Fig. 1.9), consistently with chemical intuition. 1.2.3 Parametrization of the excited state FF: QM calculations For the D molecule we need also a FF for the excited D* state. Working with excited states in classical molecular dynamics would in principle require a full reparametrization of the FF. This is a non-trivial task and definitely beyond the scope of our work. Therefore we adopted an empirical approach using the same GAFF parametrization defined for the ground state, but replacing the ground state equilibrium geometry and charge distribution with those relative to the Kasha state (the lowest vibrationally relaxed excited state). Indeed very minor conformational changes, are observed, as expected for 15 CHAPTER 1: RESONANT ENERGY TRANSFER: A DYNAMICAL APPROACH Figure 1.9: Left: comparison between methoxy dihedral distribution in modified-GAFF with FFs in Fig. 1.8; right: dipole magnitude oscillation of the modified GAFF force field compared with “closed” (orange) and “open” (red) conformation dipole magnitude as obtained from DFT calculations. NBD, a planar and rigid molecule. For excited state calculations we tested few levels of theory. We used standard HF as well as TD-DFT with two different functionals, with 6-31G(d,p) basis set. Calculations were run both in vacuum (vac) and in chloroform (clf, PCM). Table 1.1 show relevant results. The diagram in Fig.1.10 schematizes the process to extract the proper topology from the result of the ab-initio calculation. Molecular dynamics simulations address slow Figure 1.10: Calculation performed to obtain excited state charges on isolated chromophores 16 1.2 Force Field validation NBD CIS-HF B3LYP CAM-B3LYP Exp. clf (vac) clf (vac) clf (vac) Trans. En. (Abs)a4.11 (4.35) 3.14 (3.32) 3.39 (3.60) 2.73 µt(Abs) 5.978 (4.272) 3.964(2.222) 4.894 (2.885) Trans. En. (Emi) 3.70 (3.96) 2.56 (2.75) 2.91 (3.12) 2.38 µt(Emi) 5.972 (4.186) 2.682 (1.501) 4.289 (2.401) Table 1.1: Calculated absorption and emission energies (eV) and transition dipole moments (atomic units, a.u.) relevant to the lowest excited state NBD. Last column shows experimental transition energies from Ref. [34]. degrees of freedom, such as conformational changes and diffusion within the solvent, much slower than molecular electronic and vibrational degrees of freedom. Therefore the excitation is simulated by simply switching the molecular topology file from the one relevant to the ground state to that relevant to the relaxed excited state. As shown in the diagram above, we start from the optimized ground state, focus on the first excited state, relax its geometry, calculate restrained electrostatic potentials (RESP charges) and finally transfer this information to GAFF topology. Since best results for ground state calculations were obtained with B3LYP functional in vacuum, the same functional is adopted in steps 2 and 3. For step 4 instead, we adopted the HF-CIS functional, the one most suitable with RESP charges. The basis set is maintained as 6-31 G(d,p) in all calculations. Table 1.2 report the difference between ground and excited state electric dipole (both in magnitude ∆|µ|and orientation θµ), calculated for different levels of theory. Along with them, we listed the variation in charge upon excitation for a few atoms. 17 CHAPTER 1: RESONANT ENERGY TRANSFER: A DYNAMICAL APPROACH NBD ∆|µ|D 1.15 θµdeg. 12.65 O−NO ∆q0.031 NO2∆q-0.074 NH2∆q0.053 Table 1.2: Ground and excited state RESP charge distribution. Magnitude and orientation changes in electric dipole moment are shown along with charge differences for most relevant atoms. 1.3 Molecular Dynamics simulations 1.3.1 The free DA pair in chloroform As a first step, we address MD simulations of the unbound DA pair in chloroform. We focus on two different systems: -D−A -D∗−A We first generated a 5 ×5×5nm3chloroform box with the two dyes randomly placed in it. We then proceeded with a 1 ns simulation at fixed molecule number, pressure and temperature (NPT simulation), followed by a preliminary NVT run of 10 ns and a final NVT run of 100 ns. The Potential of Mean Force (PMF) curve is defined as the free energy surface along a chosen coordinate. Specifically, PMF P, along a generic coordinate ξ, can be defined as[35]: P(ξ) = P(ξ0)−RTln D(ξ) D(ξ0)(1.7) where Ris the universal gas constant , Tthe temperature ξ0is an arbitrary reference point. The function D(ξ), the distribution function along the coordinate ξ, is obtained directly from the MD trajectory. 18 1.3 Molecular Dynamics simulations Figure 1.11: PMF extracted from a 100 ns NVT simulation in chloroform for DA, D*A and DA*. PMFs shown in Fig. 1.11 (and all the others) were elaborated from the distributions D(d) of the distance (d) between the D and A center of mass (COM), calculated over long dynamical runs (at least 100 ns, see Appendix B for further details). The PMF profiles in Fig. 1.11 show a very broad minimum at 2.58 nm, which is however an artifact due to periodic boundary conditions (the box dimension is 5.16 nm). The peak at 0.39 nm corresponds instead to a π−πconfiguration (Fig. 1.12 left). A second peak at 1 nm acorresponds (Fig. 1.12 right) to a configuration characterized by a strong H-bond between the amino hydrogen of NBD with the quinonoid oxygen of NR. Umbrella Sampling For our purpose, since RET is strongly affected by the D-A distance, we must ensure that all possible distances are properly sampled. Therefore the PMF curve obtained 19 CHAPTER 1: RESONANT ENERGY TRANSFER: A DYNAMICAL APPROACH Figure 1.19: Bidimensional histograms of θµ(left) and θπ(right) over chromophore center of mass distance in chloroform. Top: results for clx-DA. Bottom: results for clx-D*A. where h...iξindicates the integral over time on the entire dynamics. Since in MD simulations we only have access to discrete data points separated by time intervals ∆t, the integration is substituted by a sum as follows: ACf(j∆t) = 1 N−j N−1−j X i=0 f(i∆t)f((i+j)∆t) where iand jrun over all Nframes of the simulations. The variables of interest are the distance, the angle bewteen the vectors normal to the molecular planes and the angle between permanent dipoles (d,θπand θµ), as defined previously. For each degree of freedom, an autocorrelation curve is calculated and then each AC(t) is fitted with a 26 1.3 Molecular Dynamics simulations Figure 1.20: Bidimensional histograms of θµ(left) and θπ(right) over chromophore center of mass distance in DMSO. Top: results for clx-DA. Bottom: results for clx-D*A. double exponential (Table 1.3) to extract the relevant time-scale that will be compared to the excited state decay in the next section. CXX(t) = a0exp(b0t) + (1 −a0) exp(b1t) with X=d, θπ, θµ. The results in Fig. 1.22 confirm a strong similarity of behaviour in chloroform for the intermolecular distance and the orientational motions, with all variables dynamically active on a time scale of ∼10 ns. In contrast, in DMSO the orientational motion occurs on a faster time scale (∼1ns) than the variation of intermolecular distances (∼5ns). 27 CHAPTER 1: RESONANT ENERGY TRANSFER: A DYNAMICAL APPROACH Figure 1.21: MD simulation snapshot of two of the most visited conformations in chloroform: left, π−πstacking; right, dye-linker H-bond 1.4 A dynamical model for excited state decay Simulating the relaxation of D∗is quite tricky. In fact, once the donor has been excited, it can decay along different paths: •Non-radiative decay, with kinetic constant knr •Radiative decay, fluorescence, with kinetic constant krad •Energy Transfer, with kinetic constant kRET The first two processes are marginally affected by the system dynamics, while the RET rate, being strongly dependent on the distance and mutual orientation of the dyes, is highly variable. We are particularly interested to understand and model the concurrent dynamics of RET and conformational motion, well beyond standard treatments 28 1.4 A dynamical model for excited state decay clx −DA Chloroform clx −DA DMSO dfit range :100000 :100000 a00.722 0.664 b−1 011.63 5.47 b−1 11.95 0.883 θπfit range :100000 :100000 a00.678 0.748 b−1 010.87 0.747 b−1 10.37 0.052 θµfit range :100000 :100000 a00.550 0.877 b−1 011.76 1.02 b−1 10.68 0.067 Sampling frequency: 10 ps. Table 1.3: Parameters extracted from exponential fitting of autocorrelation functions of d,θπand θµ. All time constants are in ns. that typically address two limiting cases. Specifically, if the molecular motion is much slower than the transfer process, dyes orientation is approximated as frozen and the orientational factor κ2is considered constant and set to a value ranging from 0 to 4, depending on the mutual orientation of the dyes. In the opposite case, the donor and acceptor have enough time to explore all their variable-space before exchanging energy, and κ2=2 3. In realistic situations, and especially for bound DA pairs, the situation is actually intermediate and calls for more detailed model, as made possible by MD. In our approach we convert the rate of each decay path for D∗(radiative/nonradiative decay, RET) into a probability. The constant rates for the radiative and non-radiative decay are set to the the experimental values. Instead, the probability decay along the the RET channel is controlled by the instantaneous value of VDA. Plotting 29 CHAPTER 1: RESONANT ENERGY TRANSFER: A DYNAMICAL APPROACH Figure 1.22: . Autocorrelation functions of distance d,θµand θφ, between D and A for clx-DA in chloroform (left) and DMSO (right) the value of interaction, along with D-A separation, as a function of time (Fig. 1.23), it is clear how long inter-dyes distance is always related to a negligible magnitude of VDA while, as chromophores get closer, the interaction can assume large values, depending on relative orientation of transition dipole moments. Knowing, for a generic process, the characteristic decay time τ, defined as: τ(t) = 1 krad +knr +kRET (t),(1.9) we can assess an infinitesimal probability Pdt to every “instant of time”: Pdt =1 τdt →Zτ 0 Pdt = 1 (1.10) Since MD simulation uses numerical integration, we only have access to finite time-steps ∆t, and we need to translate Eqn. 1.10 into a discrete equation: P∆t=1 τ∆t→ τ X 0 P∆t= 1 (1.11) In practice, each timestep carries a fraction of the probability for each of the possible paths. 30 1.4 A dynamical model for excited state decay Figure 1.23: D-A interaction VDA (black lines) and Donor Acceptor distance d(red lines) plotted as a function of time for three sample excited state trajecotries. Top panels: full dynamics; bottom panels: magnified sections of the aforementioned trajectories. Energy transfer rate kRET is easily calculated from Eqn. 1.2, while krad and knr are extracted from lifetime (τD) and quantum yeld (Φ) experimentally measured on the isolated chromophore: krad =1 τD·Φknr =1 τD−krad (1.12) Once ∆tand τare determined relaxation time is calculated as follows: 1. Initialization A finite group of chromophores is selected; number of simulation steps is set to Ncycle =−1 and initial time is defined as t= ∆t·Ncycle 2. Time step Time is increased by ∆t(Ncycle =Ncycle + 1) 31 CHAPTER 1: RESONANT ENERGY TRANSFER: A DYNAMICAL APPROACH 3. Decays For each molecule, a random number 0 ≤Rn<1 is generated: - if Rn≤1 τ, a relaxation occur - else, nothing happens 4. Continuation Relaxed molecules are removed from the initial group, and their decay time (Ncycle· ∆t) is saved; remaining chromophores are sent back to Step 2 5. End When all molecules are relaxed, the simulation stops A graphical representation is shown in Fig. 1.24. In this example a generic step N for an ideal system is illustrated. Remembering the relation between lifetime and rate constants (Eq.1.9), P∆tcan be written as: P∆t= (krad +knr +kRET )·∆t= 3.5ns−1·0.1ns = 0.35 (1.13) Each time step is represented by a box with Nb= 200 balls, of different colours. The colour of each ball represents the ”action” pursued by the excited moiety for each timestep, as depicted in Fig. 1.24. Red, green and cyan balls all represent a de-excitation event, respectively RET, radiative decay and non-radiative decay, while the black balls accounts for timesteps in which the donor stays in the excited state. The portion of balls of each colour are directly related to the probability the excited donor has to take that path. The simulation consists in consecutively extracting a ball from each box, starting from box number 0, until a non-black ball is picked. The output of the simulation is a distribution of relaxation times, one for each molecule. Plotting the fraction of molecules in the excited state per time, an exponential profile is obtained. Moreover, being able to discriminate between the three different paths of relaxation, our model allows to determine the energy transfer efficiency for each system as the ratio between molecules relaxing through RET channel over all the decayed molecules: ΦRET =NRET Ntot (1.14) 32 1.4 A dynamical model for excited state decay Figure 1.24: A simple graphic representing a generic step in the simulation. 1.4.1 Decay Times To start with, we consider the isolated D∗, setting kRET = 0. The experimental lifetime of the free donor is τ=1 knr+krad = 7.0ns. We selected a population of 105molecules and set the time step for the simulation to ∆t= 10ps. In absence of RET the decay probability is P∆t= 1.429 ·10−3and we obtain quite naturally an exponential decay profile. The calculated distribution of decay-times is shown in Fig. 1.25. The cumulative difference over the histogram, generated subtracting to the total excited state population the number of molecule relaxing at each timestep, gives the characteristic exponential decay shown in the right panel. We now turn on the donor-acceptor interaction, and hence RET. Since the interaction between the chromophores depends on the separation and relative orientation between chromophores, that vary during dynamics, the new probability will be itself function of time (P∆t(t)). We run a long clx-DA simulation (1µs) to have a statistically relevant number of configurations, and from this set we randomly select 1000 configurations, making sure that the subset gives a good replica of the original distribution (Fig. 1.26). For each configuration a 10 ns NVT dynamics 33 CHAPTER 1: RESONANT ENERGY TRANSFER: A DYNAMICAL APPROACH Figure 1.25: Left: distribution of relaxation times for τ= 7.0ns; right: normalized excited population as function of time for the same system. Figure 1.26: PMF curves calculated from the 1000 frames used to RET analysis (violet) and distribution calculated for all frames in the trajectory (green). 34 1.4 A dynamical model for excited state decay Fragment µt xµt yµt z|µt|atom 1 atom 2 Free NBD -1.1920 0.2826 0.0000 1.23 N1 N4 Free NR -3.0375 0.5667 0.0783 3.09 C51 C62 Table 1.4: Transition dipole moment components for the two chromophores. In MD simulations, the direction is determined by the versor connecting highlighted atoms in Fig. 1.27. is calculated, simulating the donor excitation by switching the charge distribution with the one of clx-D*A. The decay times following the instantaneous excitation were estimated calculating VDA(t), using transition dipole moments from the relaxed excited state to the ground state and from the ground state to the bright excited state, for donor and acceptor respectively. We assume that molecular properties, including the transition dipole moments, are marginally affected by the dynamics, so their magnitudes and orientations relative to the corresponding chromophoric units are constant. For each chromophore, two atoms are selected whose distance is parallel to the transition dipole moment orientation evaluated with ab-initio calculations (Table 1.4 and Fig. 1.27). For each instantaneous excitation dynamics, we were able to calculate kRET (t) as follows: kRET =V2 DA ~2cJDA (1.15) where cis the speed of light (in m s−1) and JDA is the spectra overlap between NBD emission spectrum and NIR absorption spectrum, each normalized for unit area, considered constant in our model[34] (1.95 ·10−6m). Following the strategy outlined above, we finally obtain a decay time for each derived trajectory. To further improve statistical accuracy we repeated the calculation 100 times for each dynamics, ending up with 105 decays. Fig. 1.28 shows the calculated decay of the excited donor, and compares the results of the dynamical calculation for the RET pair with those obtained for the isolated D dye. To fit the exponential decay of D∗we used a multiple exponential fitting, with 35 CHAPTER 1: RESONANT ENERGY TRANSFER: A DYNAMICAL APPROACH channels. On the contrary, in a dynamical model, the population is continuously transferred to replenish these hot spots, so that fast channels stay active all along the process. Quite interestingly, the different results obtained in the static and dynamic calculations cannot be ascribed to different distributions of kRET . Indeed the histograms calculated along the static and dynamic trajectories shown in Figure 1.31 are very similar. As the static scenario only refers to a sampling through graound-state MD, while the dynamic one is the results of the excited state (non-equilibrium) dynamics, the similarity of the two distributions suggests that the sources of non-equilibrium effects are actually rather modest. Figure 1.31: Left panels: decay of the D∗population calculated in chloroform (top) and DMSO (bottom). Dashed lines refer to the exponential decay of the isolated D species; red lines show static results (obtained neglecting the conformational motion of the excited RET pair after excitation), and the black lines show the full dynamical result. Right panels show the (area normalized) distributions of the RET rates, resulting from static and dynamic calculations for both solvents 42 1.5 Conclusions 1.5 Conclusions We proposed an original computational protocol for the calculation of the RET dynamics and quantum yields for a donor-acceptor pair in two different solvents, highlighting the importance of conformational fluctuations. The proposed approach is validated against an extensive set of experimental results available for a RET-pair bound via flexible links to a calixarene scaffold. We obtained a very good comparison with the experiment and solved a theoretical problem associated with this system, where estimates of VDA from TD-DFT calculations on few preselected representative configurations of the D-clx-A system underestimated the RET lifetimes by several orders of magnitude.[34] Our approach combines equilibrium and non-equilibrium MD calculations with TDDFT results, leading to a detailed description of the concurrent processes: D∗decay, energy transfer and conformational dynamics. The conformational motion modulates VDA, the intermolecular interaction responsible for RET. In a fully dynamical picture the system after photoexcitation is allowed to explore conformational regions where VDA is large, then opening fast RET channels and leading to faster RET than in a static picture, where the conformational motion is frozen. RET in disordered systems is a delicate issue: the standard approach relying on the use of an average κ2value has been questioned in several ways. In the first place the factorization of the VDA interaction in a term κ2that only depends on the intermolecular orientation and in a term that only depends on the intermolecular distance, is incorrect - particularly if the system can explore regions where intermolecular distances are comparatively short. More generally those approaches represent a too crude approximation to describe DA interactions in systems, like the one investigated here, where a complex supramolecular structure poses serious constraints to the mutual arrangements of the Dand Amoieties. We show that these problems can be easily addressed by ground state MD calculations, that offer reliable information on the conformational heterogeneity of the RET pair. However, we also demonstrated that this is not sufficient in the case investigated here. In such flexible systems, the conformational motion modulates intermolecular interactions on a timescale relevant to RET, leading to important effects that cannot be accounted for through an orientational average in either in the static or ultrafast regime. 43 Chapter 2 Molecular Aggregates 2.1 Introduction Molecular materials are characterized by intermolecular forces much weaker than the chemical bonds inside each individual molecular unit. In spite of that, in molecular materials intermolecular interactions deeply affect optical spectra, that therefore cannot be calculated as the sum of molecular spectra. Intermolecular charge transfer (CT) was early recognized as a source of impressive spectroscopic phenomena in absorption spectra of molecular materials, both in the visible and near-IR spectral regions, where so-called CT absorption bands appear[36, 37]. Vibrational spectra are also affected by CT, with the appearance of strong features due to large charge fluxes driven by molecular vibrations or lattice modes (phonons).[38, 39, 40, 41] In molecular materials with intermolecular distances larger than the sum of Van der Waals radii, electrons are localized within each molecular unit and CT interactions are negligible. Even in these condition, electrostatic intermolecular interactions may have prominent effects, driving resonance energy transfer among different molecular species[17, 42, 43, 44] and energy delocalization among equivalent (or nearly so) molecules in molecular crystals and aggregates.[45, 46, 47] The physics of excitons and of optical spectra in molecular crystals was first addressed in the seminal works of Craig,[48] Davidov[49] and Agranovich[50]. The same physics also applies to molecular aggregates: the ground-breaking discovery 45 CHAPTER 2: MOLECULAR AGGREGATES of the anomalous spectra of cyanine dyes in poor solvents[51] opened the research field of molecular aggregates, with the seminal theoretical work of Kasha.[52] An enormous body of experimental and theoretical work can be found in the literature, as recently reviewed.[46, 47] Current understanding of optical spectra of molecular aggregates and crystals is mainly based on the so-called exciton model that, only accounting for electrostatic interactions among degenerate states, leads to a large reduction of the basis set. Analytical solutions have been obtained of the electronic problem that offer a reliable basis for understanding spectral properties of molecular aggregates and crystals. The approximations of the exciton model were discussed in the original papers[49, 50] and have been recently addressed in relation with aggregates of polarizable molecules with polar[53, 54, 55] or quadrupolar character.[56, 57, 58, 59] Molecular vibrations add another layer of complexity to the physics of molecular aggregates and crystals: the deformation of the molecular structure upon excitation is responsible for the Franck-Condon structure of absorption and fluorescence spectra of isolated molecules in solution, but the shape of absorption and fluorescence spectra of molecular aggregates often largely deviates from the Franck-Condon behavior,[60, 46] as first recognized in the narrow and structurless absorption and emission spectra of aggregates of cyanine dyes.[51] Delocalization energies and vibrational energies are often comparable in molecular aggregates and the adiabatic approximation must be abandoned. Treating the electronic and vibrational degrees of freedom on the same foot is a formidable task. Indeed analytical solutions of the coupled electronic and vibrational problem are available for an infinite one-dimensional array of molecules in the exciton approximation and accounting for a single coupled vibration.[60, 61] Most often, approximations schemes have been proposed to treat the problem, whose validity and applicability need a careful discussion. In this chapter we will address spectral properties of molecular aggregates, discussing first aggregates of polar and polarizable dyes, extending a previous work[53, 62] to account for vibrational coupling and hence for spectral band-shapes. We will then turn attention to aggregates of non-polar, yet polarizable, dyes, addressing the reliability of the Heitler-London approximation. Finally we will briefly address the behavior of two dimensional aggregates. This chapter is based on a published work on aggregates 46 2.2 Aggregates of polar and polarizable dyes of non polar dyes (Anzola M, Di Maiolo F, Painelli A. Optical spectra of molecular aggregates and crystals: testing approximation schems, PCCP 2019), and a second paper (in preparation) on aggregates of polar and polarizable dyes. 2.2 Aggregates of polar and polarizable dyes 2.2.1 DA chromophores as polar and polarizable dyes Conjugated dyes with an electron-donor (D) and acceptor (A) group represent a large family of dyes of interest for several applications, ranging from non-linear optics[63], molecular electronics[64, 65], OLED[5] etc. Interamolecular interactions in aggregates of these dyes considerably alter their spectral properties[54, 55, 66] One of the most widely used approximation to describe DA dyes is the so called Mulliken model[36, 39, 67], in which a dye is described in terms of two electronic states, a neutral (|DAi) and a zwitterionic (|D+A−i) state. The underlying hypothesis is that higher excited states are located at too large energies to be relevant. Indeed the dyes of interest typically have the lowest excited state in the visible region of the spectrum, while higher excited states are in the ultraviolet region. Defining 2z0as the energy difference between the two states and −τthe mixing matrix element, the Hamiltonian reads ˆ H= 0−τ −τ2z0!= 2z0ˆρ−τˆσx(2.1) In the following we will set τ= 1 as the energy unit. The Hamiltonian can be readily diagonalized, giving the eigenvalues Eg=z0− pz2 0+τ2and Ee=z0+pz2 0+τ2and eigenstates |gi=p1−ρ|DAi+√ρD+A− |ei=−√ρ|DAi+p1−ρD+A−(2.2) where ρ=hg|ˆρ|gimeasures the weight of the zwitterionic state in the ground state, and hence the molecular polarity. Its dependence on z0is (see also Fig 2.1a): ρ=1 2−z0 2pz2 0+ 1 (2.3) 47 CHAPTER 2: MOLECULAR AGGREGATES Figure 2.1: The isolated (gas-phase) dye. Top: the two resonating structures. (a) The ρ dependency from z0. (b) the transition energy Ω and the transition dipole moment as a function of ρ. (c) The potential energy surfaces for a system with εv= 0.4 and z0= 0.7. (d) The adiabatic PES calculated for the same system as in panel (c). 48 2.2 Aggregates of polar and polarizable dyes To address optical spectra we need a definition of the dipole moment operator (ˆµ). Following Mulliken [36], we neglect all matrix elements of the dipole moment operator in the chosen basis but µ0, the large dipole moment associated with the zwitterionic state. Accordingly, the dipole moment operator reads: ˆµ=µ0ˆρ= 0 0 0µ0!(2.4) All spectral properties can be expressed as function of ρ. Specifically, the transition energy and the transition dipole moments read: ~Ω = 1 ρ(1 −ρ)=1 β(2.5) hg|ˆµ|ei=µ0β(2.6) and the mesomeric dipole moment, the difference between the permanent dipole moments in the excited and ground states reads ∆µ=µ0(1 −2ρ) = µ0α(2.7) The dependence of transition energy and dipole moment on ρis shown in Fig. 2.1b. The ∆µ(ρ) dependence, a straight line, is not shown. To properly address the shape of optical spectra of polar aggregates, molecular vibration coupled to the electronic system must be taken into account. As sketched in Fig 2.1(c,d), we introduce a single effective vibrational coordinate Q, and assign to the two basis states two harmonic potential energy curves with the same curvature and displaced minima. The relevant Hamiltonian reads: ˆ H= 2z0ˆρ−τˆσx+~ωv(ˆa†ˆa+1 2)−g(ˆa†+ ˆa)ˆρ. (2.8) where ωvand gare the frequency and the electron-vibration coupling of the effective vibrational mode, respectively (see Fig. 2.1c). The coordinate and its conjugated momentum are expressed in second quantization as follows Q=r~ 2ωv (ˆa† i+ ˆai) P=ir~ωv 2(ˆa† i−ˆai) (2.9) 49 CHAPTER 2: MOLECULAR AGGREGATES The vibrational relaxation energy (see Fig. 2.1c) is εv=g2/ωv. In the adiabatic approximation the two-dimensional electronic Hamiltonian is diagonalized for each Qto get the Q-dependent adiabatic eigenstates whose energy as a function of Qis shown in Fig. 2.1d. If we are only interested in adiabatic results at equilibrium Q ¯ Q=r2ωv ~ g ω2ρ(2.10) we may solve the equilibrium adiabatic Hamiltonian that reduces to the two-state electronic Hamiltonian in Eq. 2.1 but with the energy gap between the N and Z states that self-consistently depends on the ρ: 2z= 2z0−2εvρ[68]. This self-consistent problem is easily solved to get the ρ(z0) for fixed v(see Fig. 2.1a). 2.2.2 The model for the aggregate To model an aggregate we consider a linear ordered array of NDA dyes, imposing periodic boundary conditions to minimize finite size effects, while preserving translational symmetry. In line with the standard model to describe aggregates, we assume that intermolecular distances are large enough to neglect the overlap between orbitals on nearby molecules, so that electrons are fully localized on the molecular units, and intermolecular interactions are only due to electrostatic forces. The aggregate Hamiltonian is the sum of the molecular Hamiltonians plus a term that accounts for electrostatic intermolecular interactions: ˆ H=X i (2z0ˆρi−τˆσx,i) + 1 2X i,j Vij ˆρiˆρj+~ωvX i (ˆa† iˆai+1 2)−gX i (a† i+ ˆai)ˆρi(2.11) where Vij is the electrostatic interaction between two molecules in their zwitterionic state, and i,jrun over all aggregate units. The adopted model for dipolar aggregates, despite its semplicity, implies a problem related to its dimension. Indeed, the dimension of the electronic basis grows as 2N, allowing to treat fairly large aggregates (up to more than 10 molecules). But vibrations must be accounted for in a non-adiabatic approach to properly deal with aggregates, and if just 5 vibrational quanta are considered per molecule, one ends up with a basis 50 2.2 Aggregates of polar and polarizable dyes increasing as 10N, making the problem intractable already at N∼3−4. In the following we will discuss strategies to deal with this problem, before showing selected results. 2.2.3 Electronic Hamiltonian: the rotation on adiabatic states We start our discussion focusing on the electronic Hamiltonian, corresponding to Eq. 2.8 with g=0: ˆ Hel =X i (2z0ˆρi−τˆσx,i) + 1 2X i,j Vij ˆρiˆρj(2.12) This Hamiltonian, written on the diabatic DA e D+A−basis, is fairly simple, but we need to account for all 2Nstates already to describe the ground state. In fact each molecular unit is described in the ground state by a linear combination of the DA and D+A−state (see equation 2.2). Therefore it is convenient to rotate the electronic basis and rewrite the Hamiltonian on the adiabatic basis, gand e. Of course the rotation does not affect the basis dimension. However, as it will be discussed below, a proper choice of the adiabatic basis will help us to better understand the physics of the system, also allowing for a controlled truncation of the basis, while maintaining the quality of the results. Following a previous work[69] we define two new operators, linear combinations of Pauli operators σx,z, whose projections are defined by the parameter ρ: ˆ Sx,i =−2βˆσz,i +αˆσx,i ˆ Sz,i =αˆσz,i + 2βˆσx,i (2.13) where iruns on the molecular sites and α= 1 −2ρand β=pρ(1 −ρ) are the same variables introduced previously. The above operators are then expressed in terms of creation and annihilation operators: ˆ Sx,i = 1 −2ˆ b† iˆ bi ˆ Sz,i = (ˆ b† i+ˆ bi) (2.14) The operator ˆ b† icreates an excitation on site i, by turning the molecule from state |gi to |ei, while ˆ bidestroys the excitation. As discussed by Agranovich,[70] these operators 51 CHAPTER 2: MOLECULAR AGGREGATES with low and intermediate polarity (left and central panels), but red-shifts in the case of a largerly polar dye (right panel). These apparently crazy results, possibly suggesting the failure of the exciton picture, are indeed related to a bad choice of the reference state. A large part of the shift in fact is not excitonic in origin, but is related to the effects that surrounding charges have on the energy of the states. This is easily calculated in the mf approximation, solving the problem for an isolated dye feeling the potential from the surrounding dyes. Repulsive intermolecular interactions leads to a reduction of the polarity of each dye in the aggregate (see fig. 2e) and hence to a variation of the frequency of the absorption band. The proper reference for the exciton model is indeed represented by the mf frequency. Specifically, for the dye in the left panels of Fig. 2.4, the ionicity decreases from 0.19 in the gas phase to a mf value of 0.17. Accordingly, the maximum of the absorption blueshifts, slightly reducing the exciton shift. Similar considerations apply to the dye in the middle panels, whose ionicity is reduced from 0.64 in the gas phase to 0.5 in the mf approach. For ρ= 0.5 (the so-called cianine limit) the Condon vibrational coupling (proportional to the squared mesomeric dipole moment) vanishes leading to the disappearance of the vibronic structure for absorption and fluorescence bands. More inspiring is the case of the zwitterionic dye in right column of Fig. 2.4. Here the decrease of the ionicity from 0.76 in the gas phase to 0.61 in the mf approximation is responsible for a large red-shift of the absorption band. Taking as proper reference the mf frequency, a blue-shift of the absorption band is observed for the aggregate, fully in line with its H character, as due to repulsive (V > 0) intermolecular interactions. Fluorescence in H-aggregates comes from electronic states at the border of the Brillouin zone and are only allowed due to the coupling to vibrational modes. As a result, very weak and largely red-shifted bands are observed, but what we notice here is that, since the dominating (Condon) term accounting for electron-vibration coupling vanishes in the cyanine limit, the fluorescence intensity is vanishingly small in this limit. A similar analysis can be done for the aggregates in Fig. 2.5, corresponding to the case of weak attractive intermolecular interactions (V=−1). Intense emission bands and vanishing Stokes shifts in the aggregate are fully in line with J-aggregate behavior. The redshift of absorption (and emission) bands observed for the dyes in the left and 58 2.2 Aggregates of polar and polarizable dyes Figure 2.4: H aggregate, V= 1, εv= 0.4, ωv= 0.17: top and bottom panel show calculated absorption and fluorescence spectra. Intensities per molecules are reported in arbitrary units. The weak fluorescence spectra of the aggregate are multiplied by the factor shown in the figure. Left panels refer to a system with z0= 0.8, corresponding to an ionicity for the isolated dye ρ= 0.19 that deacreses in the mf approximation to ρ= 0.17. Middle panels: z0=−0.3, gas phase ρ= 0.64, mf ρ= 0.50. Right panels: z0=−0.6, gas phase ρ= 0.76, mf ρ= 0.61. 59 CHAPTER 2: MOLECULAR AGGREGATES Figure 2.5: J aggregate, V=−1, εv= 0.4, ωv= 0.17: top and bottom panel show calculated absorption and fluorescence spectra. Intensities per molecules are reported in arbitrary units. Left panels refer to a system with z0= 1.0, corresponding to an ionicity for the isolated dye ρ= 0.15 that increase in the mf approximation to ρ= 0.21. Middle panels: z0= 0.7, gas phase ρ= 0.21, mf ρ= 0.50. Right panels: z0= 0.3, gas phase ρ= 0.36, mf ρ= 0.82. 60 2.2 Aggregates of polar and polarizable dyes middle panels of Fig. 2.4 are again in line with a J-aggregate behavior. The most striking results is however recognized again for the most polar molecule (ρ= 0.36 in the gas phase) in the right panels of Fig. 2.4: here in fact the exciton band moves to the blue with respect to the gas-phase molecule. But again this anomalous behavior is simply related to the choice of a wrong starting point. In the aggregate, the mf solution of the problem drives the molecule deep in the ionic regime with ρ= 0.82. This implies a large blue shift of the absorption and fluorescence bands, so that, when taking as reference the gas phase molecule, an apparent blue-shift of the exciton band is observed, that actually corresponds to a red-shift when the proper reference is considered, in line with the attractive nature of the interactions. We also notice that for the zwitterionic system, when the wrong reference state is considered, the intensity of the transitions (both absorption and fluorescence) decreases and the vibronic structure becomes more prominent, in striking contrast with the J-nature of the aggregate. This inconsistency is however quite naturally solved if the proper mf reference is considered: in all cases the spectral intensity increases when going from the mf dye to the aggregate, while the vibronic structure becomes less and less prominent. Quite interestingly, results in the central panel fo Fig. 2.4 refer to a dye with ionicity ρ= 0.16 in the gas phase that is driven to the cyanine limit, ρ= 0.50 when embedded in the aggregate. Once again, in the cianine limit the vibronic structure of absorption and fluorescence bands disappears. Medium and Strong Coupling We will now address the cases of medium and strong coupling. Fig. 2.4 show absorption spectra calculated for H-aggregates in the medium (V= 1.6) and strong-coupling (V= 2.0) regimes. It turns out that Ne= 4 is the minimum number of exciton states needed for convergence, Ne= 3 results are totally untenable, with the only exception of the systems that in the mf approaximation have ρ= 0.5. In these conditions in fact all terms in Eq.2.18 proportional to 1−2ρvanish. Accordingly, the vibronic structure disappears, as discussed above, as well as all terms related to the mesomeric dipole moment. For the electronic part, the Hamiltonian in the ρ= 0.5 limit reduces to that relevant to nonpolar aggregates and most of the anomalous effects associated with aggregates of polar and polarizable dyes are washed out. Once convergence is reached, finite size effects are 61 CHAPTER 2: MOLECULAR AGGREGATES marginal for largely neutral dyes, as well as for dyes in the cyanine limit, but become relevant for zwitterionic dyes. Figure 2.6: H aggregate absorption spectra. All results refer to a system with v= 0.4 and ωv= 0.17. Top panels show results for V=−1.6, from left to right: z0= 0.8, gas phase ρ= 0.21, mf ρ= 0.15; z0=−0.6, gas phase ρ= 0.84, mf ρ= 0.5; z0=−1.0, gas phase ρ= 0.9, mf ρ= 0.62. Bottom panels show results for V=−2.0, from left to right: z0= 0.8, gas phase ρ= 0.21, mf ρ= 0.14; z0=−0.8, gas phase ρ= 0.88, mf ρ= 0.5; z0=−1.2, gas phase ρ= 0.92, mf ρ= 0.61. For N= 4 the full electronic basis is considered, Ne= 4. More interesting is the case of J-aggregates, where electrostatic intermolecular interactions lead to intriguing phenomena.[53, 62] Fig. 2.7 shows absorption and fluorescence spectra calculated for a system with V=−1.6, corresponding to the curve in Fig. 2.2 that marks the boundary between the normal (weak coupling) and the bistable (strong coupling) regime. Much as in the weak case, the apparently anomalous behavior observed when comparing aggregate spectra with spectra calculated for the isolated dye 62 2.2 Aggregates of polar and polarizable dyes are relieved if the proper reference system is considered, corresponding to the mf solution. In all cases in fact the aggregate spectrum is red-shifted with respect to the relevant mf spectrum. The most important difference with respect to the weak coupling is the appearance of finite size effects, with N= 6 results differing from N= 4, pointing to excitons with large delocalization. Moreover, to get convergence for N= 6 at least Ne= 4 is needed (see Appendix A) in sharp contrast with the weak coupling case. Quite interestingly, finite size effects are marginal for the system described in middle column of fig. 2.7 where the mf ionicity is 0.5. As discussed above, in this limit, the vanishing of terms proportional to 1 −2ρnot only kills the main vibronic coupling term, but also reduces the electronic part of the Hamiltonian to that of aggregates of non-polar dyes. This is even more evident in the strong coupling limit in fig. 2.8, showing spectra calculated for V=−2.0. Similar considerations apply as in the medium-coupling regime, but in this case N= 6 results do not converge until the maximum number of excitons Ne= 6 is accounted for in the calculation, or in other terms, the complete electronic basis is considered (see Appendix A). This immediately tells us that the exciton-exciton interaction term (the first term in the last line of Eq.2.18, lowers the energy of multiexciton states that get mixed with the lowest excited states giving a sizable multiexcitonic character to the state, as extensively discussed in refs. [53, 62] Again, this term vanishes for a system with a mf ionicity ρ= 0.5, so that for this system (middle panel of fig. 2.8) the N= 6 results already converge at Ne= 3. 2.2.7 Discussion Extending a previous work[53] to account for molecular vibrations, as needed to properly address spectral bandshapes, a two-step approach is introduced for the description of optical spectra of aggregates of polar and polarizable molecules. The first step is the definition of the proper reference state as the mf solution of the problem. Basically, the ground state polarity of each dye is self-consistently defined by the polarity of the surrounding dyes, leading to increased polarity for attractive intermolecular interactions and reduced polarity for repulsive interactions. Of course all molecular properties (including transition frequencies and dipole moments) are affected by this variation. These mf states define the proper reference states for the exciton model. The molecular ge63 CHAPTER 2: MOLECULAR AGGREGATES ometry is also affected by the molecular polarity and the correct reference state for the vibrational problem is defined via a Lang-Firsov transformation that simply translates the origin of the vibrational coordinates to the equilibrium position relevant to the charge distribution of the molecule inside the aggregate. Since molecular vibrations in turn affect the molecular polarity, this leads to another self-consistent interaction. While this may look as a difficult problem, it boils down to a simple self-consistent diagonalization of a two by two Hamiltonian if the molecules are described in an essential state picture.[77, 55] The essential state model adopted here has been extensively validated against experiment and describes in a very effective way the low-energy spectral properties of push-pull dyes accounting for environmental effects in solution[68, 78, 79] and aggregates[80, 55], films[81] and crystals.[72, 77] In the context of this work, however, we underline that the model relies on similar approximation as the standard exciton model, accounting for a single electronic excitation and a single vibrational mode per molecule. At variance with the standard exciton model, however, the proposed model fully accounts for the molecular polarizability and for the dependence of the ground and excited state molecular geometry on the molecular polarity. Once the proper reference state is defined, several interaction terms are recognized in the Hamiltonian that can be classified as excitonic, when conserving the exciton number, and ultraexcitonic when mixing states with a different number of excitons.[53] The vibrational coupling lead to an excitonic term that corresponds to the Condon coupling in the exciton model, but turns out proportional to 1−2ρ, and therefore vanishes in system whose mf ionicity is close to 0.5: in these system the vibronic bandshape is washed out. The ultraexcitonic term exchanges vibrational quanta and excitons and has marginal spectroscopic effects in the weak coupling limit as shown in Fig. 2.9 and 2.10, that compare exact results obtained in the weak and strong coupling regimes for H and J aggregates with those obtained suppressing the non-Condon vibronic coupling in the Hamiltonian in Eq. 2.18. Non-Condon corrections give rise to sizable effects only in the strong coupling regime. As for excitonic terms originating from electrostatic interactions, we recognize terms ∝ρ(1−ρ), i.e. proportional to the squared transition dipole moment of the mf molecules: 64 2.3 Aggregates of non-polar molecules these terms are responsible for the exciton hopping. Other terms appear ∝(1 −2ρ), i.e. proportional to the mesomeric dipole moment, that account for exciton-exciton interactions. These last terms vanish when the mf molecular ionicity is close to 0.5, and the system reduces to an aggregate of non-polar dyes. The exciton approximation works reasonably well for weak coupling but becomes clearly untenable in the strong coupling regimes (see fig. 2.9 and 2.10). Indeed, when increasing the coupling, ultraexcitonic terms enter into play with particularly impressive effects in J-aggregates, where bistability regions are observed in the mf solution.[53, 62] Finite size effects become important in these conditions and the exciton basis cannot be reduced to account for just the first few exciton states (up to 3 excitons are enough to get converged results in the weak coupling limit). Indeed, the lowest excited state in these conditions cannot be described, not even approximately, as a state with a single exciton, rather it corresponds to a state where several excited molecules cluster together in a multiexciton state.[53, 62] 2.3 Aggregates of non-polar molecules Aggregates of non-polar dyes are typically described in terms of the exciton model. Specifically, the model assumes that each molecule in the aggregate can be either on the ground |gior excited |eistate, both states having a negligible permanent dipole moment. Electrostatic intermolecular interaction then only imply transition dipole moments, hg|ˆµ|ei=µt. Relevant interactions enter the aggregate Hamiltonian in two different terms: an Heitler-London (HL) term that is responsible for the exciton hopping, and a non-HL term that mixes states whose exciton number differs by two units. When intermolecular interactions are much smaller than the exciton energy, one can neglect all terms in the Hamiltonian that mix states with different energy. The non-HL terms are then neglected, leading to the standard exciton model. The HL approximation is very useful and widely adopted as it leads to an enormous reduction of the electronic basis. Indeed, while for the complete model one should account for 2Nstates (Nis the number of molecules in the aggregate) in the HL approximation, as long as one is interested to linear spectral properties, only the subspace with a single exciton is of relevance and the 65 CHAPTER 2: MOLECULAR AGGREGATES basis has dimension N. In physical terms, imposing the HL approximation amounts to fully neglect the molecular polarizability, imposing that the nature of the molecular ground and excited states is not affected by intermolecular interactions. As it will be discussed below, this approximation, while usually leading to acceptable spectra, leads to some fundamental problem. Specifically the sum rule for the oscillator strength and for the rotational strength in chiral aggregates[82] are broken. In the following, we will see how these issues can be solved relaxing the HL approximation. 2.3.1 The model Hamiltonian We consider a one dimensional array of Nnon-polar molecules assuming periodic boundary conditions. Each molecule is either in the ground, |gior excited |eistate, the two states being separated by an energy ~ω0. An internal vibrational coordinate ˆqiis introduced per molecule and the two electronic states are assigned harmonic potential energy surfaces with the same frequency, ωv, but displaced minima. The strength of the electron-vibration coupling is measured by the vibrational relaxation energy, λ, that, being related to the Huang-Rhys factor, S=λ/~ωv, can be extracted from the analysis of the absorption or fluorescence bandshape for the isolated dye in solution.[25] The electrostatic interaction between two molecules at sites iand jreads:[50] Ji,j =µ2 t 4πεd3 ij Dij,(2.21) where Dij is a geometrical factor that only depends on the relative orientation of the transition dipole moments on sites iand j, while dij is the distance between the two molecular sites. Dij can assume positive values (repulsive interactions) or negative values (attractive interactions). In the following, we will address one-dimensional molecular aggregates with one molecule per unit cell. We will impose periodic boundary conditions as to maintain translational symmetry. In these conditions, intermolecular interactions only depend on the relative distance between molecules and we define Jm=Ji,i±m. Two extreme cases will be considered: (a) nearest-neighbor interactions with J1=Jand Jm= 0 for m > 1; (b) unscreened long-range Coulomb interactions with J1=Jand Jm=J[sin(π/N)/sin(mπ/N)]3for m > 1, with Nbeing the number of molecular units. 66 2.3 Aggregates of non-polar molecules In either case a single parameter, J, measuring the strength of the nearest-neighbor interaction, fully defines the model. The magnitude of µtis experimentally accessible from the oscillator strength of the g→etransition measured for the isolated dye in solution: fge =2 3 me ~e2ω0µ2 t,(2.22) where meis the electron mass, ethe electron charge and ω0the frequency of the g→e transition. Alternatively, the transition dipole moment can be obtained from quantum chemical calculations, with the added value of getting information about the dipole moment orientation with respect to the molecular frame. With these definitions, the Hamiltonian for a linear array of Nmolecules reads: ˆ H=X ihE−λ(ˆa† i+ ˆai)iˆni+~ωvX iˆa† iˆai+1 2 +X m Jm(ˆ b† iˆ bi+m+h.c.) + X m Jm(ˆ b† iˆ b† i+m+h.c.),(2.23) where iand jrun on the Nmolecular sites. The operators ˆa† i, ˆaiare the boson operators associated with the harmonic oscillator on site i. The operator ˆ b† icreates an excitation on site i, by turning the molecule from state |gito |ei, while ˆ bidestroys the excitation (these operators obey Paulion algebra as discussed in Section 2.2). The first term in the above Hamiltonian describes the molecular problem, where E is the vertical excitation energy of the molecule in the aggregate. It may differ from ~ω0, the transition energy in the isolated molecule, due to local field effects, but we will neglect these corrections in the following, setting E=~ω0. The last two terms account for intermolecular interactions: the same interaction Jmis responsible for the hopping of the excitation from site ito i+m(and viceversa) and for the simultaneous creation (destruction) of two excitations on sites iand i+m. As mentioned above, the hopping term mixes states with the same number of excitons, i.e. states having the same diagonal energy, while the last term in the above Hamiltonian describes the interaction among states whose energy differs by 2Eand, in the HL approximation, it is neglected. Of course, the HL approximation works well for J2E. Before closing this Section, it is instructive to compare the Hamiltonian for non - polar dyes in Eq. 2.23 with the Hamiltonian describing polar and polarizable dyes in 67 CHAPTER 2: MOLECULAR AGGREGATES 2.3.4 Testing approximation schemes Exact diagonalization approaches are limited to small aggregates, up to 7 sites for the complete model and up to 10 sites in the HL approximation. The electronic basis is comparatively small, growing with Nin the HL approximation and as 2Nin the complete model. Indeed the basis blows up because of the vibrational states: accounting for just three vibrational quanta per site would multiply by a 3Nfactor the basis dimension. It is therefore very important to discuss approximation schemes to cut the basis dimension and particularly so for the vibrational states. Indeed the HL electronic basis is already very small, while we already discussed how the electronic basis in the complete model can be reduced by fixing a maximum number of excitons (Me= 3 seems to work pretty well in most cases of interest in this study, even if this approximation is untenable for clusters of polar and polarizable dyes in medium or strong coupling regimes. Recently,[55] discussing J-aggregates of polar dyes, we realized that for largely delocalized excitons only the vibrational modes in the close proximity of the center of the Brillouin zone are effectively coupled to the electronic degrees of freedom so that, instead of accounting for Nlocal harmonic oscillators, reasonable results are obtained accounting for the single oscillator with q= 0. This of course leads to an enormous reduction of the basis dimension. The Hamiltonian in Eq. 2.26 shows that accounting for just the q= 0 mode one obtains a similar coupling Hamiltonian as for the isolated molecule, but with the strength of the coupling reduced to λ/√N. As a result, in this approximation the same bandshape is calculated for J and H aggregates, as shown in the top panels of Fig. 2.13, where we show absorption spectra calculated in the HL approximation for the same model parameters as in Fig. 2.11. A full decoupling of the q= 0 vibrational mode is expected in the infinite chain limit and a washing out of the vibronic structure in either J or H aggregates of infinite size. Adding the two nearest modes to the q= 0 mode in the Brillouin zone (middle panels of Fig. 2.13) improves the agreement and adding 2 more modes (for a grand total of 5 delocalized vibrations, bottom panels) gives very good results for J aggregates and an acceptable agreement for H-aggregates. Cutting vibrational modes in the Fourier space works in principle for the complete as well as for the standard exciton model, but, apart from the simplest case where only the q= 0 mode is accounted for, the approach is 74 2.3 Aggregates of non-polar molecules difficult to implement in the complete model. Moreover this approximation is expected to work well for largely delocalized excitons. Of course for localized excitons or in the presence of disorder the approach could only work if many modes (possibly all) in the reciprocal space are introduced, making the approximation useless. A useful and widely adopted approach to reduce the vibrational space is the socalled few-particle approximation,[83, 84] that has been extensively applied by Spano[46, 47] in the 2-particle approximation (2PA) or 3-particle approximation (3PA) flavors. It works in the real space, so that it does not require a symmetric or ordered system, but only applies in the HL approximation where all relevant basis states have a single exciton. In the 2PA approximation, the basis is cut imposing a maximum number neof vibronic excitations on the electronically excited states (notice that these vibronic states refer to the displaced harmonic oscillator as relevant to the electronically excited state). Moreover a maximum number of vibrational quanta nvcan be present in just a single additional site, different from the site bearing the exciton. In the three particle approximation (3PA), one accounts for vibrational excitations occurring on up to two sites. Two different approximation schemes are possible for both 2PA and 3PA, a small-range (sr) scheme, where vibrational excitations are only allowed in the nearest sites of the site bearing the exciton (Fig. 2.15), or a long-range (lr) scheme, where vibrational sites can be spread all over the aggregate. Of course the 2PA-sr or the 3PA-sr only apply when the exciton model accounts for nearest neighbor interactions, while one must resort to lr-schemes when accounting for long-range electrostatic interactions. To be specific, the 2PA basis set is: |ψ2PAi=|n, ˜ν, νli, l 6=n, (2.32) where nmarks the site where the exciton resides and ˜νcounts the number of vibrational quanta in the displaced oscillator associated with the same site. The numbers νl6=ncount the vibrational quanta in the undisplaced harmonic oscillator on site l. In the sr flavor of 2PA, l=n±1, while in the lr flavor, lcan assume any value different from n. The diagonal energy of the 2PA states is easily calculated as ~ω0+ (˜ν+νl)~ωv. Off diagonal 75 CHAPTER 2: MOLECULAR AGGREGATES matrix elements, accounting for the interaction between different sites, are: hn, ˜ν, νl|ˆ H|m, ˜µ, µki=Jf˜ν,µnf˜µ,νmY i6=n,m δνi,µi,(2.33) where f˜ i,j is the Franck-Condon factor measuring the overlap between the vibrational level of excited state ˜ iand vibrational level of ground state j. Finally, the transition dipole moment is calculated as follows: µtrans i=hG|ˆµ|ψii=X n,˜ν cn,˜νhG|ˆµ|n, ˜ν, 0i=X n,˜ν cn,˜νµ0f˜ν,0.(2.34) Moving to the 3PA, the relevant basis set reads |ψ3P Ai=|n, ˜ν, νl, νl0i,l, l06=n. Accordingly, the diagonal energy is ~ω0+ (˜ν+νl+νl0)~ων, while the Hamiltonian offdiagonal matrix elements are: hn, ˜ν, νl, νl0|ˆ H|m, ˜µ, µk, µk0i=Jf˜ν,µnf˜µ,νmY i6=n,m δνi,µi.(2.35) Fig. 2.16 compares absorption spectra calculated via exact diagonalization and with the 2PA-sr and 3PA-sr for an aggregate of 10 molecules, described by the standard exciton model, with increasing strength of nearest-neighbor interactions, J. For J-aggregates the 2PA and 3PA approximations work pretty well up to medium-large interactions, but for H-aggregates, the approximation is poor already for interactions of medium strength. Similar results hold true for emission spectra in Fig. 2.17. Moving to 2PA-lr or 3PA-lr does not change the picture, as expected. The lr extensions of the 2PA and 3PA approaches has to be invoked for systems where long-range Coulomb interactions are accounted for. However the 3PA-lr basis is very large, making it impossible to deal with aggregates with more than 6 sites. Therefore Fig. 2.18 and Fig. 2.19 compare exact and 2PA-lr results for absorption and emission spectra, rescpectively, of 10 site aggregates with long-range intermolecular interactions. 76 2.3 Aggregates of non-polar molecules Figure 2.7: J aggregate, V=−1.6, v= 0.4, ωv= 0.17: top and bottom panel show calculated absorption and fluorescence spectra. Intensities per molecules are reported in arbitrary units. Left panels refer to a system with z0= 1.5, corresponding to an ionicity for the isolated dye ρ= 0.09 that increases in the mf approximation to ρ= 0.10. Middle panels: z0= 1.0, gas phase ρ= 0.16, mf ρ= 0.50. Right panels: z0= 0.5, gas phase ρ= 0.33, mf ρ= 0.90. For N= 4 the full electronic basis is considered, Ne= 4. For N= 6 only converged results are shown with Ne= 4 77 CHAPTER 2: MOLECULAR AGGREGATES Figure 2.8: J aggregate, V=−2.0, v= 0.4, ωv= 0.17: top and bottom panel show calculated absorption and fluorescence spectra. Intensities per molecules are reported in arbitrary units. Left panels refer to a system with z0= 1.5, corresponding to an ionicity for the isolated dye ρ= 0.09 that increases in the mf approximation to ρ= 0.11. Middle panels: z0= 1.2, gas phase ρ= 0.12, mf ρ= 0.50. Right panels: z0= 1.0, gas phase ρ= 0.16, mf ρ= 0.87.For N= 4 the full electronic basis is considered, Ne= 4. For N= 6 only converged results are shown with Ne= 4 78 2.3 Aggregates of non-polar molecules Figure 2.9: H aggregates with N= 6. Top panel show weak-coupling results, V= 1 for the same values of model parameters as in Fig. 2.4; bottom panels show results for strong coupling, V= 2, for the same parameters as in the bottom panels of Fig. 2.6. In all panels blue lines show converged results for the total Hamiltonian, dashed black curves show results obtained neglecting the non-Condon electron-vibration coupling term, continuous black lines show results for the exciton model, i.e. suppressing all ultraexcitonic terms in the Hamiltonian. 79 CHAPTER 2: MOLECULAR AGGREGATES Figure 2.10: J aggregates with N= 6. Top panel show weak-coupling results, V=−1 for the same values of model parameters as in Fig. 2.5; bottom panels show results for strong coupling, V=−2, for the same parameters as in Fig. 2.8. In all panels blue lines show converged results for the total Hamiltonian, dashed black curves show results obtained neglecting the non-Condon electron-vibration coupling term, continuous black lines show results for the exciton model, i.e. suppressing all ultraexcitonic terms in the Hamiltonian. 80 2.3 Aggregates of non-polar molecules Figure 2.11: Absorption spectra calculated for the complete Hamltonian in Eq.2.23 only accounting for nearest neighbor interactions. Results are obtained for N=6, ~ω0= 2.0 eV, λ=0.17 eV, ~ωv=0.17 eV and |J|=0.255 eV. The vibrational basis is truncated setting the maximum number of total vibrational quanta Mv=6. The electronic basis is truncated fixing the maximum number of excitons Meto a value ranging from 1 to 4. Results for Me=1 coincide with those obtained in the HL approximation. Top and bottom panels refer to H-aggregates (J > 0) and J-aggregates (J < 0), respectively. In both panels the dash-dotted line shows the spectrum relevant to non-interacting molecules (i.e., the monomer limit). 81 CHAPTER 2: MOLECULAR AGGREGATES Figure 2.12: Emission spectra calculated for the the same system as in Fig. 2.11 82 2.3 Aggregates of non-polar molecules Figure 2.13: Calculated absorption spectra for J and H aggregates (left and right panels, respectively) calculated in the HL approximation for the same model parameters as in Fig. 2.11. Results refer to aggregates of 10 molecules, black lines show numerically exact results, obtained with Mv= 4, magenta lines show results calculated only accounting for the q= 0 mode (top panels), q∈ {−π 5,0,π 5}(middle panel), q∈ {−2 5π, −π 5,0,π 5,2 5π} (bottom panels). 83 CHAPTER 2: MOLECULAR AGGREGATES Figure 2.21: Crystal structure of BF38: (a) layered structure; (b) highlights of the two directions of interactions inside a single layer. Structures of BF31 and GK08 are completely analogous. (kx, ky), where kxand kyare defined as: kx=2π Nx sxky=2π Ny sy(2.36) Nx(y)represents the number of molecules in the x(y) dimension, and sx(y)is anologous to its monodimensional counterpart. The Hamiltonian is written in the reciprocal space, as in Eq.2.26. The Paulion creation and annihilation operators are: ˆ bK=ˆ bkxky=1 pNxNyX m,n eimkxeinkyˆ bmn (2.37) ˆ b† K=ˆ b† kxky=1 pNxNyX m,n eimkxeinkyˆ b† mn (2.38) where i(j) runs on the Nx(Ny) molecules. The exciton Hamiltonian in the reciprocal space reads: H2D=X K ˆ b† KbKE+ 2cos(kx)Jx+ 2cos(ky)Jy+~ωvX Qˆa† QˆaQ+1 2 −λ √NX K,Q (a† Qb† KbK+Q+h.c.) (2.39) 90 2.3 Aggregates of non-polar molecules where b† Kcreates an exciton with wavevector ~ K, whose energy is represented by the quantity in square parenthesis, Nrepresent the total number of molecules forming the aggregate (N=Nx·Ny) and ~ Q= (qx, qy) represents the two dimensional counterpart of the reciprocal space vibrations described in Eq. 2.26. Selection rules for absorption and emission are the same as in one dimension. Absorption is only possible from the ground state to a total-symmetric state ( ~ K= (0,0)), while the symmetry of the fluorescent state (the lowest excited state) depends on the sign of the interactions. Specifically, there are four possibilities for the wavevector of the fluorescent state, as shown in Table 2.1. where JxJy~ K + + (0,0) + - (0, kmax y) - + (kmax x,0) - - (kmax x, kmax y) Table 2.1: Possible values of ~ Kfor the lowest excited state in bidimensional aggregates. kmax i, the highest value for the wavevector, depends on the number of molecules in the idimension (Ni): kmax i=   π, if Nieven Ni−1 Niπ, if Niodd. (2.40) Reliable values for the Jxand Jyinteractions are extracted in two steps. First we use TD-DFT (B3LYP functional, 6-31g(d,p) basis set) to calculate the magnitude and orientation of the transition dipole moment for the isolated molecules (see Table 2.2). Then, crystallographic data are exploited to calculate the interactions between transition dipole moments on different molecules, in the dipolar approximation: J12 =1 4πε0η2r3 12 [(~µ1·~µ2)−3 r2 12 (~µ1·~r12)(~µ2·~r12)]|2(2.41) Specifically, we calculated interactions Jbetween all nearest neighbor couples in the crystal. For all the three systems, two major interaction are found as reported in Table 2.3 for η2= 2. 91 CHAPTER 2: MOLECULAR AGGREGATES Molecule Bright Exc. State Energy (eV) Trans. Dipole Moment (D) BF31 3.27 5.056 BF38 2.69 5.259 GK08 3.42 4.578 Table 2.2: Excited state energy and relative transition dipole moment magnitude for the 3 molecules under investigation. Molecule Jx(eV) Jy(eV) BF31 0.0059 -0.030 BF38 0.053 0.028 GK08 0.061 0.021 Table 2.3: Interaction energies calculated for the three molecular crystals Figure 2.22: Normalized emission spectra calculated for the target molecules in THF solution. The single layers are then modeled using a 5 ×3 bidimensional aggregate for a grandtotal of 15 molecules. Vibrational energies ωvand electron-phonon couplings λ are extracted from experimental monomer emission spectra (shown in Fig. 2.22). The difference in energy between the 0 −0 and 0 −1 vibronic peaks gives the vibrational energy, while the ratio between 0 −0 and 0 −1 peak intensities is proportional to the 92 2.3 Aggregates of non-polar molecules coupling λ. In fact, recalling the expression for the Franck-Condon coefficients: h0|ni2=Sne−S n!with S=εv ωv ,(2.42) and knowing εv=λ2/ωv, the electron-photon coupling is obtained as: λ=√Sωvwith S=h0|1i2 h0|0i2.(2.43) In Fig. 2.23, 2.24 and 2.25, experimental absorption and emission spectra are compared to calculated spectra. The monomer absorption energies (Ein Eq. 2.39) are selected Figure 2.23: Comparison between experimental measured spectra, in green, and calculated spectra from Hamiltonian in Eq. 2.39, in violet, for the bf31 crystal. Left: absorption spectra; right: emission spectra. in order to match the position of experimental crystal absorption peaks. For both BF38 and GK08 the comparison between calculated and experimental spectra show a remarkably good match, both in absorbtion and emission. The bandshapes of the transitions, as well as the Stokes shift is satisfactorily reproduced. The results acquire additional value remembering that the parametrization of the Hamiltonian is obtained, except for the exciton energies, using electron-phonon couplings directly extracted from monomer spectra in solution and evaluating the Jxand Jyinteractions from ab-initio calculations and experimental crystallographic data. On the contrary the modelization of BF31, 93 CHAPTER 2: MOLECULAR AGGREGATES Figure 2.24: Comparison between experimental measured spectra, in green, and calculated spectra from Hamiltonian in Eq. 2.39, in violet, for the BF38 crystal. Left: absorption spectra; right: emission spectra. Figure 2.25: Comparison between experimental measured spectra, in green, and calculated spectra from Hamiltonian in Eq. 2.39, in violet, for the GK08 crystal. Left: absorption spectra; right: emission spectra. while leading to an acceptable match in absorption spectrum, fail in the description of fluorescence, regarding both position and shape of the spectrum. The problem with this system can be due to the presence, in the two direction, of a positive and a negative 94 2.3 Aggregates of non-polar molecules interaction, which leads to a more cumbersome landscape that require bigger dimensions of the simulated aggregate. 95 CHAPTER 2: MOLECULAR AGGREGATES 2.4 Conclusions In this chapter we have discussed effective models for molecular aggregates. First we have addressed linear aggregates of polar and polarizable molecules, extending a previous work[69, 62] to account for electron-vibration coupling. A clever choice of the basis, defined on the adiabatic electronic states in the mf approximation and accounting for the vibrational displacement via a Lang-Firsov transformation, allows to largely reduce the basis dimension, so that, also exploiting translational symmetry, we are able to treat fairly large systems, with up to 6 molecules. Apart from this important technical result, the proper choice of the basis amounts to a proper choice of the reference state. The apparently erratic behavior observed for aggregates of polar and polarizable molecules, where red-shifts are sometimes observed for H-aggregates and blue-shifts are sometimes observed for J-aggregates,[66] actually results from a poor choice of the reference state. Excitonic (and ultraexcitonic) effects, including shifts, must be evaluated against the transition frequencies obtained in the mf approximation, i.e. taking into account the large variation of the nature of the polarizable dye when inserted in a lattice of polar dyes. Once the proper reference is selected, J and H aggregates always give rise to red and blue-shifted absorption band. Similarly, anomalous effects on bandshapes and intensities are relieved when the proper reference state is considered. Our analysis demonstrates that, for aggregates of polar and polarizable dyes, when the proper reference is taken, the excitonic approximation works reasonably at least in the weak coupling regime. Ultraexcitonic effects are important in the strong couplig regime and particularly so in J-aggregates, where multiexcitonic states become prominent.[69, 62] For aggregates of non-polar dyes we were able, exploiting symmetry and a clever computational implementation, to treat fairly large systems, as needed to validate several approximation schemes. More importantly we were able to address a fundamental problem of the exciton model that generates spectra that do not obey the sum rule for the oscillator strength. This is safely ascribed to the HL approximation. We demonstrated that the HL approximation, the main approximation adopted in the exciton model for aggregates of non-polar dyes, leads to good estimates of transition energies because of a cancellation of errors, while already in the weak interaction regime, the HL approximation leads to overestimated absorption (and fluorescence) intensities for 96 2.4 Conclusions H-aggregates and underestimated intensities for J aggregates. Of course for large interactions, ultraexcitonic effects are also recognized in the frequencies of the bands. Having demonstrated that the exciton model leads to reasonable results for weak interactions, we finally applied the model to describe optical spectra of three crystals of organic dyes, with good results, particularly in view of the minimal number of adjustable model parameters. 97 Chapter 3 Chiral aggregates of αand β-dicyanostylbenes: chiroptical properties 3.1 Introduction Objects that cannot be superimposed to their mirror image, technically objects whose symmetry group does not contain any improper axis (including S1the mirror plane and S2the inversion center), are called chiral (from the greek word for hand). Chiral systems show special properties when interacting with light. At the linear order (weak electromagnetic fields) the Optical Rotatory Dispersion (ORD)[86, 87] measures the (frequency-dependent) difference between the refractive index for the left and right circularly polarized light and Circular Dichroism (CD)[88, 89] measures the difference in the corresponding extinction coefficients.[90, 91, 48] Of course ORD and CD spectra are related through integral expressions similar to the Kramers Kr¨onig relations.[92] Chirality and chiroptical activity are observed in intrinsically chiral molecules, including molecules with chiral centers as well as chiral structures like helicenes. However, nonchiral molecules may give rise to chiral responses is some conditions as (a) under the 99 CHAPTER 3: CHIRAL AGGREGATES OF αAND β-DICYANOSTYLBENES and β-DCSB molecules with non-chiral substituents[99, 100] and start our MD simulation from a configuration obtained arranging the eight molecules of the αand β-DCSB derivatives mimiking stacks extracted from the relevant crystal structures, as shown in Fig. 3.6. We notice that, since the crystal structures are obtained for non-chiral systems, the starting configuration is non-chiral. Figure 3.6: Initial configurations of αand β-DCSB. Molecular dynamics simulations are then run with a first minimization step, followed by a 30 ns NPT (see Appendix B) equilibration (needed to allow molecules to arrange in the thermodynamically stable structure) and a final 200 ns NPT main run with a 1 ps timestep. In the poor solvent (water) the molecules reorganize to minimize the interaction with the solvent. The central chromophoric units interact via π-πstacking while the side chiral chains, with asymmetric interactions, impart an overall torsion to the whole supramolecular aggregates. The high flexibility of the molecules together with the large number of units, result in a somewhat messy conformational landscape. Nevertheless, a preferential packing is clearly detected for the 4 investigated systems. In Fig. 3.7 the most representative structures arising for each of the long trajectories are presented. All systems are characterized by a well-defined helical structure, shaped by interactions of side pendants. The RR enantiomers (either αor β-DCSB) show a left-handed 106 3.4 Absorption and CD spectra of DCSB aggregates helix organization, while SS enantiomers form a right handed helix. The MD results, Figure 3.7: Most relevant structures of the 4 systems under investigation. showing that αand β-DCSB substituted with the same chiral pendant aggregate in supramolecular structures with the same handedness but showing opposite CD signal is a first confirmation of the validity of the new chirality rule, stating that the sign of CD spectra depends not only on the chirality of the system, but also on the sign of intermolecular interactions. However, more reliable results and a better understanding of the phenomenon can be obtained via a detailed calculation of absorption and CD spectra of the aggregates. 3.4 Absorption and CD spectra of DCSB aggregates 3.4.1 The model To model our aggregate we adopt the standard exciton model for non-polar dyes (i.e. we impose the HL approximation, see Chapter 2) and neglect the coupling to molecular vibration. As discussed previously, the exciton model gives acceptable results for optical spectra of non-polar dyes and, for largely disordered system, the information about the vibronic structure is smeared out by inhomogeneous broadening effects. 107 CHAPTER 3: CHIRAL AGGREGATES OF αAND β-DICYANOSTYLBENES With these approximations, the excitonic Hamiltonian for a generic aggregate of N identical molecules reads: ˆ H=X n E0|nihn|+X nm Vnm|nihm|(3.1) where |niidentifies the basis state in which the nth molecule is in its excited state (|g1, g2. . . en. . . gNi), with energy E0, and Vnm is the dipole-dipole interaction term: Vnm =1 4πε0n2|rn,m|2~µ n 0·~µ m 0−3(~µ n 0·~rnm)(~µ m 0·~rnm) r5 n,m (3.2) where ~ µn 0is the transition dipole moment for the excitation on the n-th molecule and ~rnm is the vector distance between sites nand m. In the following we will use the molecular transition dipole moments obtained by TD-DFT (B3LYP/6-31G**) as discussed in Section 3.3, anchored to the molecular site. The diagonalization of the excitonic Hamiltonian on the Nbasis gives the excitonic eigenstates accounting for the linear combination of singly excited states. Transition dipole moment, from the ground state to the keigenstate, ~µ k t, is obtained as linear combinations of molecular transition dipole moments: ~µ k t=X n ~µ kn t=X nhk|ni~µ n 0(3.3) where kruns on the eigenstates and ~µ n 0is the dipole moment associated to the |gi → |eitransition in the nth molecule. The absorption spectrum is calculated assigning a Gaussian lineshape to each transition: A(ω) = ~ωX k~µ k t 2e (~ω−Ek)2 2σ2(3.4) where Ekis the energy of the transition from the ground to the keigenstate. In the proposed model, electrons are localized on each molecular unit. In this case, as first demonstrated in a classical work by Condon[91], CD spectra can be calculated from the knowledge of transition dipole moments. In fact, in this limit, transition magnetic dipole moments, needed to calculate CD spectra (see Chapter 4 for further details), can be expressed in terms of transition dipole moments. Following Condon we define the 108 3.4 Absorption and CD spectra of DCSB aggregates rotational strength for each transition from the ground to the k-eigenstate as Rk=−EkX nm ~rnm ·(~µ kn t×~µ km t) (3.5) where ~µ kn tis the dipole moment of the |0i → |kitransition relative to molecule n. Replacing Eq. 3.3 in Eq. 3.5, we obtain the explicit expression for Rk: Rk=−EkX nm ~rnmhn|kihm|ki·(~µ n 0×~µ m 0) (3.6) CD spectra, measuring the difference between exctinction coefficient of left and right circularly polarized light, is simply obtained assining a Gaussian lineshape with width σto each transition: ∆(ω)∝X k Rke (~ω−Ek)2 2σ2(3.7) 3.4.2 Calculated Spectra The output of MD simulations is a collection of configurations. Assuming that the transition dipole vectors calculated for the isolated molecules stay constant during the dynamics with respect to the molecular frame of each monomer, we can calculate, for each configuration of the aggregate, the electrostatic interactions entering the Hamiltonain (see Eq. 3.2). The values of E0, corresponding to the transition energy for the isolated monomer, are extracted from experimental values in Fig. 3.2 and are set to 3.37 eV for αderiatives and 3.10 eV for βderivatives. The Hamiltonian is then diagonalized and transition energies, transition dipole moments and rotational strengths are finally obtained. With these information, absorption and CD spectra are calculated (Eq. 3.7, setting σ=0.07 eV). The resulting spectra are then averaged over a large portion of dynamics (60 ns). Results for the absorption and and CD spectra of the 4 systems under investigation are shown in Fig. 3.8 and 3.9. Calculated spectra and experimental data presented in Fig. 3.2 are generally in good agreement. Rotational strenghts extracted from trajectories mimic accurately experimental CD spectra. Specifically, we obtain the correct signs for CD spectra for the 4 simulations. This result, together with the most relevant structures presented in Fig. 109 CHAPTER 3: CHIRAL AGGREGATES OF αAND β-DICYANOSTYLBENES Figure 3.8: Top panels: absorption spectra of both enantiomers of α-DCSB compared to the monomer (black curve). Bottom panels: corresponding CD spectra. 3.7, demonstrate the reliabilty of our hybrid method in the modelization of these complex systems. Results for β-DCSB are very good: absorption spectra of both enantiomers are very similar, comfirming that the equilibrium geometry is obtained for both species, and show a very clear blue-shift with respect to the monomer. Results for α-DCSB are less satisfactory. While the spectrum for the R enatiomer has its main peak slightly red shifted with respect to the monomer, both R and S systems present a not negligible band at higher energy, denoting the formations of H type clusters. More critically, the absorption spectra of the two α-enantiomers differ considerably, suggesting that the MD simulation has not yet fully converged to the equilibrium. For this reason, while results for all systems are promising, calculations on αderivatives are still on going. 110 3.5 Conclusions Figure 3.9: Top panels: absorption spectra of both enantiomers of α-DCSB compared to the monomer (black curve). Bottom panels: corresponding CD spectra. 3.5 Conclusions In this Chapter we suggested a novel approach to the spectroscopic characterization of complex aggregates in solution. We proposed a hybrid method, based on the combination of MD simulations and exciton models to calculate absorption and circular dichroism spectra of a family of dicyanostilbenes (DCSB) derivatives. In particular, we focused our investigation on 4 systems: α-DCSB (R and S enantiomers) and β-DCSB (R and S enantiomers) that show different spectroscopic behavior depending on the nature of the chromophoric unit and the chirality of the side chains. With the use of MD we studied the aggregation of these systems in solution. Starting from completely achiral arrangements of monomers, extracted from crystallographic data on related systems, well defined clusters are obtained, with an overall helicity that only depends on the chirality of side chains. Geometric informations were extracted from classic trajectories and then 111 CHAPTER 3: CHIRAL AGGREGATES OF αAND β-DICYANOSTYLBENES processed in order to parametrize an exciton model. Absorption and CD spectra are then calculated for the 4 systems of interest, averaging over thousands of conformations. We obtain CD spectra in good agreement with experimental data. The hybrid technique proposed in this Chapter represents a reliable method for the replication of chiroptical properties in unusual chiral assemblies, where the fine interplay between chromophores nature and peripheral ligand chirality contributes to complex outcomes. 112 Chapter 4 Chiral aggregates of squraine dyes 4.1 Introduction Squaraines are a widely investigated family of organic dyes, characterized by a resonance stabilized structure where the central squaryl ring, a strong electron-acceptor unit, is conjugated to two electron-donating groups, as in the examples shown in Fig. 4.1. Squaraines have interesting spectroscopic properties, with intense and sharp absorption and emission bands in the red and near-IR spectral regions and large two-photon cross-sections. They are investigated for applications in dye-sensitized solar cells[3, 4], in colorimetric sensors [109] and non-linear optics.[110] Moreover, their rigid conjugated structure allows the formation of stable aggregates, [111, 112] whose unsusual spectroscopic properties further widen the possible range of applications. The study of selforganization and aggregation behaviour of squaraine dyes is a fairly hot topic in recent years. [113, 16, 114, 14] Chiral aggregates of squaraine dyes were obtained ∼15 years ago via the supramolecular arrangment of squaraine dyes decorated at the two sides by chiral groups[115, 112]. The chiral supramolecular arrangment is demonstrated by the observation of characteristic CD spectra for the aggregate structure, in the region of the squaraine absorption, CD-silent for the non-aggregated dye. New attention on chiral squaraine aggregates was called by a recent paper from the group of Manuela Schieck [114], showing how 113 CHAPTER 4: CHIRAL AGGREGATES OF SQURAINE DYES Figure 4.1: a): example of two squaraine dyes; b): schematic representation of the squaraine D-A-D structure. aggregates deposited on films give rise to very large CD signals. In this chapted we describe an extensive theoretical work that, devoted to the analysis of absorption and CD spectra of chiral squaraine aggregates in solution, represents a first fundamental step to understand the behavior of aggregates deposited on films. This work is done in collaboration with experimentalists and specifically with researchers in the group of Professors Manuela Schiek (Johannes Kepler University of Linz, Austria) and Arne L¨utzen (University of Bonn, Bonn). Fig. 4.2 (courtesy of PhD student Jennifer Zablocki) shows representative absorption and CD spectra for aggregates in solution. Indeed a very large number of results are available on aggregates of squaraine dyes with different substituents (Fig. 4.3). For all aggregates absorption spectra show aggregations features both to the blue and to the red of the monomer band, located at 650 nm. We will dub the two features as H and J bands, respectively, even if their nature (as it will be demonstrated below) is different. Initially, the two bands were ascribed to two features associated with the exciton splitting of the monomer band. This however would require a fairly large intermolecular interactions. Moreover, this interpretation contrasts sharply with the observed CD spectra. Indeed if the H and J bands were the two features due to the exciton splitting, one would have observed in CD spectra a single bisignated signal with opposite sign at the 114 4.2 The three state model for squaraine dyes Figure 4.2: Absorption (left) and CD (right) spectra measured for the ProSQ-SS-C10 molecule in methanol upon addition of water. H an J band locations. Instead, the CD spectrum shows two distinct bisignated signals, one in the region of the H and and one in the region of the J band. The possibility that the H and J bands are due to the simultaneous presence of different kinds of aggregates in the same sample is not very likely, due to the consistent observation of the two features in samples obtained in different experimental conditions and using dyes with different substituents. Absorption and CD spectra of chiral squaraine aggregates then call for a careful modelization that will be the topic of this chapter. 4.2 The three state model for squaraine dyes The spectroscopic behaviour of quadrupolar DAD dyes, and specifically of squaraines, can be rationalized adopting an essential state model, based on 3 electronic basis states [116]. The three orthogonal states are a neutral state |Ni, and two degenerate states, |Z1iand |Z2i, corresponding to the two zwitterionic structures D+A−D and DA−D+. On this basis, the electronic Hamiltonian is: ˆ Hel = 2z0ˆρ−τˆσ(4.1) 115 CHAPTER 4: CHIRAL AGGREGATES OF SQURAINE DYES Figure 4.6: Absorption (top-left) and CD (bottom-left) spectra for the same system in Fig. 4.5, with a finite tilt angle α= 25◦. 122 4.4 The role of intermolecular charge transfer interactions Figure 4.7: The model dimer studied by Collison and Spano with highlighted intramolecular and intermolecular charge transfers. 4.4 The role of intermolecular charge transfer interactions In a recent paper Spano and Collison studied non-chiral aggregates of squaraine dyes characterized by an absorption spectrum showing prominent absorption features both to the red and to the blue of the monomer absorption spectrum.[14]. They ascribed this observation to the presence of a charge transfer (CT) interaction between adjacent molecules. This work inspired us to investigate intermolecular CT as a source of the anomalous spectral features observed in our chiral aggregates. Specifically, Spano and Collison, considered a dimer structure as in Fig. 4.7, and accounted for a CT interaction between the adjacent D and A sites located in the two different molecules. Indeed other structures can also be considered, including e.g. those drawn in Fig. 4.5a. Of course, when accounting for intermolecular CT interaction the basis must be enlarged to account for states with charge separation between the two molecules, so that in addition to the 9 localized states in Eq. 4.9, other states enter into play. Specifically, Spano and Collison proposed to add the states obtained considering that each dye could 123 CHAPTER 4: CHIRAL AGGREGATES OF SQURAINE DYES be in any of the following states: D+AD |C1i DAD+|C2i D+A−D+|C3i(4.10) DA−D|Ai with Cnstates bearing a positive charge (cations) and Abearing a negative charge (anion). Of course dimer states must be overall neutral so that one adds 6 states to the 9 localized states in Eq. 4.9, for a grand total of 15 states. While Spano and Collison definitely had the good intuition about the possible role of intermolecular CT in squaraine aggregates, the adopted basis is far from complete. Indeed they fully disregard the electronic spin degrees of freedom, basically building a model that implicitely treats electrons as spinless fermions. This is easily understood adopting the valence-bond basis where electrons are paired in singlet states.[123, 124] Fig. 4.8 shows a very simple example referring to a state with two zwitterionic molecules, indeed a state that is already present in the localized basis as state |Z2Z1i. If CT interactions are disregarded, electrons can only be paired within each single molecular unit (the basis states defined as direct product of the three basis state for each dye constitute a complete basis set), but, if CT interactions are accounted for, a second state must be considered to account for the two different possibilities to pair electrons in singlet states.[123, 124] The valence bond basis, constructed accounting for the subspace with total spin zero, would correspond to the smallest basis for our system. However a problem arises since the valence bond basis is non-orthonormal, making the procedure to calculate observables and transition dipole fairly cumbersone. We therefore adopt a simpler approach, using the so-called real-space basis, where electrons are accomodated in the different D/A site, selecting only states with Sz= 0, where Szmeasures the total spin component along the z direction. The real space basis is larger than the valence bond-basis, but, in view of the fairly complex calculations required for CD spectra (see Section 4.7) we stick on it for our calculations. Specifically, the real space basis for a squaraine dimer comprises 53 basis states, represented in the bit-representation by integer numbers, as shown in Fig. 124 4.4 The role of intermolecular charge transfer interactions Figure 4.8: a) Example of bi-zwitterionic state as proposed by Spano and coworkers; b) When intermolecular charge delocalization is taken into account, two singlet states are needed to represent the same charge distribution. 4.9 To proceed, we must define the Hamiltonian accounting for intra and intermolecular CT and for electrostatic intermolecular interactions. We propose a Hubbard-like Hamiltonian as follows: ˆ H=X mX i εiˆnm i+U 2X mX i ˆnm i(ˆnm i−1) + 1 2X mn X ij Vmn ij ˆqm iˆqn j(1 −δij) −tX mX i6=jˆ bm ij +ˆ b†m ij −βX m,m+1 X ij ˆ bm ij +ˆ b†m ij (4.11) where ˆnm icounts the electrons on site i(= 1,2,3) of molecule m, ˆqm imeasures the onsite charge with ˆqm i= 2 −ˆnm ifor sites D (i=1 or 3) and ˆqm i=−ˆnm ifor sites A (i=2). Moreover, nm i=Pσˆc†m iσ ˆcm iσ counts the electrons on site iof molecule mand ˆc†m iσ and ˆcm iσ create and annihilate and electron with spin σon site i, m. Intramolecular and intermolecular CT are defined by the hopping integrals tand βwith the hopping operator (an anti-Hermitian operator) defined as ˆ bmn ij =Pσˆc†m iσ ˆcn jσ. At variance with the Hamiltonian in Eq. 4.6, in this Hamiltonian both inter and intramolecular electrostatic 125 CHAPTER 4: CHIRAL AGGREGATES OF SQURAINE DYES Figure 4.9: The 53 diagrams composing the real state basis for a squaraine dimer. The first column numbers the diagrams, the second and third columns are the integer numbers representing each basis state in decimal and binary representation, respectively. Each of the 6 sites is represented by two digits in the binary number, 00 stands for a vaccum site, 01 and 10 for a spin up and a spin down electron, respectively, 11 for a doubly occupied site. In other terms, for each site, the first digit counts the number of spin down, the second digit counts the number of spin up. 126 4.4 The role of intermolecular charge transfer interactions interactions are accounted for, with Vmn ij designed as in Eq. 4.8. To ensure that the Hubbard Hamiltonian in the above equation reduces to the Hamiltonian in Eq. 4.6 for β= 0 we set t=τ/√2 and fix on site energies so that 2z0= 2εi−U−Vwhere V=Vmm 12 =Vmm 23 is the modulus of the electrostatic interaction between charges on D and A sites in the same molecule (of course fully defined by the D-A distance). Finally, we set U, the repulsion between two electrons residing on the same site, to a very large value as to make its value irrelevant, and to ensure that states with doubly charged site D2+ and A2−have very large energies and can be safely disregarded. The above Hamiltonain is written and diagonalized on the real-space basis to get relevant eigenstates and eigenvalues. The calculation of absorption spectra only requires the definition of the dipole moment operator, that reads: ˆ ~µ =−X mX i ~r m iˆnm i(4.12) where ~r m iis the position vector for site iof molecule m. On the chosen basis the dipole moment operator is diagonal. To calculate absorption spectra we calculate transition dipole moments from the ground to the excited states and calculate absorption spectra assigning a Gaussian lineshape with width σ= 0.08 to each transition: A(ω)∝~ωX E|hE|ˆµ|Gi|2e −(~ω−~ωEG)2 2σ2(4.13) Fig. 4.10 shows absorption spectra calculated for two different dimers with aligned molecules with exactly the same model parameters as in Fig4.5 but introducing a finite β= 0.6. It is clear that for both dimers, accounting for intermolecular CT interactions leads to a splitting of the band with the appearance of two transitions, one to the blue and one to the red of the monomer absorption band. We observe that the absorption spectrum arising from the configuration B is very similar to experimental data presented in Fig. 4.2. To better understand the nature of the double picked absorption spectrum we give a closer look to the transion dipole moments. Specifically, we calculate absorption spectra of the dimer in Fig. 4.10B accounting for polarized radiation. We oriented the dimer in order to have the molecules aligned along the x-axis and packed along z. The dipole moment of transitions lead by intramolecular CT (term τ) will be parallel to x-axis 127 CHAPTER 4: CHIRAL AGGREGATES OF SQURAINE DYES Figure 4.10: Absorption spectra (left) for a model squaraine dimer in the same configurations introduced in Fig. 4.5, compared to the monomer. The CT interactions taken into account for each configuration are depicted with green dashed lines. 128 4.5 Calculation of CD spectra of aggregates with delocalized electrons (same direction of the D-A arms), while intermolecular CT (term β) will be concordant with the packing direction, the z-axis. In Fig. 4.11 the total absorption spetrum of the dimer (top panel) is plotted together with absorption spectra obtained with xand z-polarized light (middle and bottom panels). It is evident how the intensity of both transitions is governed by intramolecular, rather than intermolecular, charge migration. The two characteristic transitions seen in CT dimers arise from the mixing of CT and exciton states. However the CT weight in the ground state is negligible so that CT states do not contribute to the intensity of transitions from the ground state. In other terms the intensity of the localized transition redistributes in two states, the red and the blue band. 4.5 Calculation of CD spectra of aggregates with delocalized electrons Results in the previous section suggests that accounting for intermolecular CT interaction is probably the key to understand the strange spectroscopic behavior of chiral squaraine aggregates. However, to proceed we must attack a non-trivial problem, i.e. the calculation of CD spectra in a molecular aggregate where, due to intermolecular CT interactions, electrons are not localized on the molecular units. This makes it impossible to use the approach described in Section 4.3 that, based on the classical Condon work,[91] relies on the definition of the aggregate dipole moment as the sum of molecular dipole moments. This approach of course is inadequate for delocalized electrons. According to Condon[91], CD spectra can be calculated from the rotational strength associated to each transition Ri. The rotational strength plays for CD spectra the same role as the oscillator strength in absorption spectra, so that CD spectra can be calculated, once the Riare known for each transition assigning a spectral bandshape (typically Gaussian or Lorentzian) to each transition. If |giis the ground state the rotational strength for the |gi→|fitransition is Rf= Imnhg|ˆµ|fi·hf|ˆ M|gio(4.14) where ˆµand ˆ Mare the electric and magnetic dipole operators, respectively. While 129 CHAPTER 4: CHIRAL AGGREGATES OF SQURAINE DYES Figure 4.11: Top panel: absorption spectrum of dimer in Fig. 4.10B. Middle panel: absorption calculated for the same system with radiation polarized along the x-axis. Bottom panel: absorption calculated for the same system with radiation polarized along the z-axis. the definition of electric dipole moment operator is straightforward (see Eq, 4.12), the definition of the magnetic dipole moment operator is more delicate. The magnetic dipole 130 4.5 Calculation of CD spectra of aggregates with delocalized electrons is proportional to the total angular momentum: ˆ ~ M∝ −ˆ ~ L=−X k ˆ ~ lk(4.15) where kruns on all electrons. The angular momentum of k-electron is related to its linear momentum ˆ ~pkby ˆ ~ lk=ˆ ~rk׈ ~pk(4.16) Finally the linear momentum can be obtained as ˆ ~pk=−i ~hˆ ~µk, Hi(4.17) The problem is that in a real-space description we cannot address the properties of a single electron. To overcome the problem we start evaluating the total linear momentum. For the sake of clarity we explicitly write the expression for only one of the three components: ˆ Px=−i ~[ˆµx, H] (4.18) where ˆµxis the xcomponent of the total dipole moment in Eq. 4.12. The dipole moment operator commutes with the first three terms of the Hamiltonian in Eq. 4.11, while it does not commute with the hopping terms. Overall ˆ Px=i ~X mn X ij δi,j±1(xn j−xm i)ˆvmn ij τδmn +β(1 −δmn) (4.19) where the bond velocity is an anti-Hermitian operator ˆvmn ij =ˆ bmn ij −ˆ b†mn ij =ˆ bmn ij −ˆ bnm ji (4.20) If we calculate ˆ ~ R׈ ~ P, where ˆ ~ R, the global position operator, is identical to the total dipole moment operator in Eq. 4.12, we do not get the total angular momentum in Eq. 4.15, as ~ R×~ Pdoes contain both one-electron and two-electron terms. We then define ˆ ~ Lselecting out of the ˆ ~ R׈ ~ Poperator only the one-electron terms. Accordingly: ˆ Lx=−i ~X mn X ij δi,j±1(zn j−zm i)(ym iˆ bmn ij −yn jˆ bnm ji ) −(yn j−ym i)(zm iˆ bmn ij −zn jˆ bnm ji )τδmn +β(1 −δmn) (4.21) 131 CHAPTER 4: CHIRAL AGGREGATES OF SQURAINE DYES Figure 4.14: Modelization of Prosq-SS into a schematic representation of DAD quadrupolar dye (same modelizazion applies to ProSQ-RR). Figure 4.15: The three parameters used to define relative orientation of squaraine pairs. From lefto to right: plane distance d, shift angle θSand tilt angle θT 138 4.8 Conclusions Figure 4.16: Histograms representing the distribution of the three geometric parameters. From top to bottom: d,θSand θT. In each panels the distribuions relative to the 3 nearest neighbors pairs of the tetramer are shown. 139 CHAPTER 4: CHIRAL AGGREGATES OF SQURAINE DYES Figure 4.17: Dimer absorption (left) and CD spectra (right) calculated for different value of intermolecular Charge Transfer integral β. 140 4.8 Conclusions Figure 4.18: Trimer absorption (left) and CD spectra (right) calculated for different value of intermolecular Charge Transfer integral β. 141 CHAPTER 4: CHIRAL AGGREGATES OF SQURAINE DYES Figure 4.19: Tetramer absorption (left) and CD spectra (right) calculated for different value of intermolecular Charge Transfer integral β. 142 Conclusions In this Thesis we extensively discussed intermolecular electrostatic interactions with special emphasis on molecular aggregates and FRET. The presented theoretical work provides a series of tools with the purpose of improving the description of complex supramolecular systems. Specifically, we demonstrated how the combined use of MD, exciton and essential state models can provide interesting results and a reliable description of optical properties of a variety of systems. In the first part of the dissertation we investigated, for a selected pair of chromophores, the relations between energy transfer efficiencies and slow degrees of freedom connected to conformational motion and solvation. To this aim we proposed a new computational method for the calculation of the dynamics and RET quantum yields for a RET pair in two different solvents, highlighting the importance of conformational fluctuations. Our approach, combining equilibrium and non-equilibrium MD simulations with ab-initio results, offers a detailed description of all the processes that simultaneously affect the decay rates of our system. We validated the reliability of the protocol against an extensive set of available experimental results, with very good results. This work solved a major issue related to the characterization of RET rates in systems with high flexibility, emphasizing the importance of a dynamical approach to the description of complex decay mechanisms and opening new opportunities for their theoretical treatment. The second part of the Thesis, focused on molecular aggregates, aims to extend the portfolio of relevant models and computational techniques. We conducted a detailed study on aggregates of polar and polarizable molecules, testing different approximation approaches and their limitations. Moreover, we proposed a simple yet accountable method for the calculation of optical spectra of bidimensional aggregates, whose descrip143 CONCLUSIONS tion is usually very complex. With respect to chiral molecular aggregates, we extensively investigated spectroscopic and supramolecular features of two systems, where the subtle interplay between intermolecular electrostatic interactions and chiral supramolecular arrangements results in non-obvious outcomes. Specifically, we focused attention on aggregates of non-polar and quadrupolar chromophores. We suggested a novel approach for the characterization of a family of dicyanostilbenes derivatives (non-polar). The hybrid method proposed, based on the combination of MD simulations and exciton models, simulates well the chiroptical properties observed experimentally. As for aggregates of quadrupolar dyes, a group of squaraines chiral derivatives, showing unique absorption and circular dichroism spectra upon aggregation, have been investigated. Through accurate modelization, we were able to mimic the supramolecular arrangements of the molecular units by the use of molecular dynamics simulations, We then presented a new model, able to account for CT intermolecular interactions, that offers an accurate tool for the description of spectroscopic properties of squaraine aggregates in solution. To conclude, this Thesis offers new theoretical approaches for the description of intermolecular electrostatic interactions that are of outstanding importance in the field of innovative molecular materials. Complex mechanisms, such as resonant energy transfer and charge transfer interactions, have been investigated in a completely different light, leading to a more profound understanding of their properties and opportunities. 144 Appendix A Molecular Aggregates A.1 Oscillator Strength: sum rule To relate the total oscillator strength of a system to a ground-state expectation value, we define the velocity dipole operator, ˆv: i~ˆv= [ˆµ, ˆ H],(A.1) where, to simplify the notation, we have suppressed the vector notation on both ˆvand ˆµoperators. The oscillator strength associated with the G→Etransition is: fEG =2 3 me ~e2ωEGhG|ˆµ|EihE|ˆµ|Gi.(A.2) To eliminate the transition frequency from the above expression we use: ihG|ˆv|Ei=ωEGhG|ˆµ|Ei(A.3) and its complex conjugate, thus getting: F=X E fEG =−ime 3~e2hG|[ˆµ, ˆv]|Gi,(A.4) thus proving Eq.2.28 in the main text. Furthermore, using Eqs.2.26 and 2.27 (main text) and remembering the Paulions algebra (see Eq. 2.15 in the main text), the commutator 145 APPENDIX A [ˆµ, ˆv] can be easily calculated. In particular, we have: [ˆµ, ˆv] = 1 i~hˆµ, [ˆµ, ˆ H]i, =1 i~µ2 0X i [~ω0−λ(ˆa† i+ ˆai)](4ˆni−2),(A.5) that is Eq.2.30 in the main text. 146 APPENDIX A A.2 Aggregates of polar dyes: additional results Figure A.1: The same esults as in Fig. 2.7, main text (N= 6), accounting for different value of Ne. 147