Modelling exciton dynamics inlight-harvesting moleculesbyLeonard RuoccoM.Phys, The University of Sussex, 2012A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THEREQUIREMENTS FOR THE DEGREE OFDOCTOR OF PHILOSOPHYinTHE FACULTY OF GRADUATE AND POSTDOCTORAL STUDIES(Physics)The University of British Columbia(Vancouver)February 2019c© Leonard Ruocco, 2019The following individuals certify that they have read, and recommend to the Faculty of Graduate and Postdoctoral Studies for acceptance, the dissertation entitled: Modelling exciton dynamics in light-harvesting molecules submitted by Leonard Ruocco in partial fulfillment of the requirements for the degree of Doctor of Philosophy in Physics Examining Committee: Prof. Philip Stamp Supervisor Prof. Ian Affleck Supervisory Committee Member Prof. Gordon Semenoff Supervisory Committee Member University Examiner University Examiner Additional Supervisory Committee Members: Prof. Robert Kiefl Supervisory Committee Member Supervisory Committee Member Prof. Roman KremsProf. George SawatzkyAbstractI investigate the dynamics of multi-state central systems coupled bilinearly to an externaloscillator bath within the noninteracting-blip approximation. I focus on both a 3-siteconfiguration, as well as a 2-site model for the central systems of interest. The 2-sitemodel, dubbed the dual-coupling spin-boson (DCSB) model, includes both diagonal andnon-diagonal system-bath couplings, whereas the 3-site model considers only diagonalcouplings. The bath spectral densities considered in this work include both Ohmic andsuper-Ohmic forms, as well as single optical phonon peaks. This work is motivated by therecent observance of long-lived quantum coherence effects in the photosynthetic organismknown as the Fenna-Matthews-Olson (FMO) complex. The models investigated in thisthesis are applied to this system in an attempt to explain its remarkably efficient excitontransfer mechanism, as well as to shed light on the functionality of coherence. The DCSBmodel is shown to reproduce the rapid exciton transfer times as well as the relatively longcoherence times observed in the FMO complex. The non-diagonal system-bath couplingis shown to play a crucial role in this process. This can be attributed to the inelasticphonon-assisted tunnelling (IPAT) mechanism arising from the presence of significantnon-diagonal system-bath interactions. Conversely, the 3-site model predicts rapid butincoherent exciton transfer. This can be attributed to the presence of a resonant state inthe 3-site architecture, resulting in a relatively slow exciton transfer mode in the system.Therefore efficient exciton transfer requires a careful configuration of the chromophoreenergy landscape to avoid a resonant 3-site-V configuration. Furthermore, I concludethat coherence effects arising from excitons delocalised across multiple chromophores,promotes IPAT processes arising from non-diagonal system-bath couplings, producingrapid exciton transfer between chromophores. This offers a potential explanation as to thefunctional role that coherence plays in the energy transfer mechanism of photosynthesis.iiiLay summaryA number of landmark experiments in the last decade have suggested that quantummechanics may be responsible for the remarkably efficient energy transfer in photosyn-thesis. Theoretical research has ensued hoping to explain the role of quantum mechanicsin this process; however, the exact mechanism responsible for these observations remainsunexplained. In this thesis I investigate certain mathematical models that could poten-tially explain this mechanism. I compare the results of these calculations with thosedetermined experimentally on certain photosynthetic organisms. In doing so I manageto closely reproduce the observed coherence times with a model that incorporates criticalphysical features that, as of yet, have not been applied to the photosynthesis theoreticalmodelling process. Nature has managed to produce remarkably efficient light harvestingorganisms. A better understanding of the underlying mechanisms responsible for this mayenable us to harness this knowledge towards the improvement of our own light-capturingtechnologies.ivPrefaceThis thesis includes both original work, as well as introductory material based upon asummary of the relevant literature. The work presented in Chapters 3, 4 and 5 is origi-nal, unpublished work, carried out by the author, Leonard Ruocco, under the guidanceof supervisor Prof. Philip Stamp. The work presented in Chapters 4 and 5 is underpreparation for publication in the near future.vContentsAbstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiiLay Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ivPreface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vTable of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viList of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xList of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiGlossary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xv1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.1 Background and motivation . . . . . . . . . . . . . . . . . . . . . . . . . 11.2 Quantum phenomena in photosynthesis . . . . . . . . . . . . . . . . . . . 71.2.1 Light-harvesting molecules . . . . . . . . . . . . . . . . . . . . . . 81.2.2 The Fenna-Matthews-Olson complex . . . . . . . . . . . . . . . . 101.3 Microscopic origin of chromophore-environment interactions . . . . . . . 131.3.1 Diagonal and non-diagonal couplings in chromophores . . . . . . . 141.4 Exciton dynamics in photosynthetic systems: traditional theories and theirlimitations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 161.4.1 Resonance Energy Transfer . . . . . . . . . . . . . . . . . . . . . 171.4.2 Forster theory . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 181.4.3 Quantum master equations in molecular systems: the Redfieldequation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 191.5 3-site configurations in the Fenna-Matthew-Olson complex . . . . . . . . 222 A path-integral approach to exciton transfer . . . . . . . . . . . . . . . 27vi2.1 Truncation procedure for the 3-site-boson model from an extended system 282.2 Modelling the environment as a macroscopic harmonic oscillator bath withlinear system-bath coupling . . . . . . . . . . . . . . . . . . . . . . . . . 292.3 The oscillator bath spectral density . . . . . . . . . . . . . . . . . . . . . 312.3.1 Ohmic damping: s = 1 . . . . . . . . . . . . . . . . . . . . . . . . 332.3.2 Super-Ohmic damping (acoustic phonons): s = 3 . . . . . . . . . 332.3.3 Optical phonon damping with spectral broadening . . . . . . . . . 352.3.4 Spectral density functions for light-harvesting molecules . . . . . . 352.4 The path-integral approach to modelling open-quantum mechanical sys-tems and the Feynman-Vernon influence functional . . . . . . . . . . . . 392.4.1 The path integral formalism of quantum mechanics . . . . . . . . 392.4.2 The influence functional . . . . . . . . . . . . . . . . . . . . . . . 412.5 The noninteracting-blip approximation (NIBA) . . . . . . . . . . . . . . 462.6 Validity of NIBA: a quantitative measure . . . . . . . . . . . . . . . . . . 503 Limiting and perturbative analysis of the 3-site-boson model . . . . . 523.1 The central 3-site system . . . . . . . . . . . . . . . . . . . . . . . . . . . 523.1.1 The 3-site-V configuration and population trapping . . . . . . . . 523.1.2 A path-integral formalism for the bare 3-site system . . . . . . . . 553.1.3 Population trapping in a 3-site-V system . . . . . . . . . . . . . . 603.2 Weak system-bath coupling regime for the 3-site-V system . . . . . . . . 663.2.1 Linear response and the fluctuation-dissipation theorem . . . . . . 713.2.2 Results of the 3-site-boson model in the perturbative regime . . . 733.2.3 Validity of NIBA in the perturbative regime for Ohmic and super-Ohmic damping . . . . . . . . . . . . . . . . . . . . . . . . . . . . 743.2.4 Coherent phase space in the perturbative regime . . . . . . . . . . 753.2.5 Decay times in the perturbative regime for Ohmic spectral densities 763.2.6 Decay times in the perturbative regime for super-Ohmic spectraldensities . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 783.3 Perturbative analysis: optical phonon bath with Lorentzian lineshape . . 793.4 The 3-site-independent-boson model and dephasing . . . . . . . . . . . . 813.4.1 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 824 Non-perturbative analysis of the 3-site system coupled diagonally toan Ohmic oscillator bath . . . . . . . . . . . . . . . . . . . . . . . . . . . . 83vii4.1 NIBA in the 3-site-boson model . . . . . . . . . . . . . . . . . . . . . . . 834.1.1 Ground-state propagator for the unbiased 3-site-boson system . . 854.1.2 Excited-state propagator for the unbiased 3-site-boson system . . 864.1.3 The influence functional for a biased 3-site system . . . . . . . . . 864.2 Non-perturbative analysis for Ohmic spectral densities . . . . . . . . . . 874.2.1 Validity of NIBA for the 3-site-Ohmic-boson model . . . . . . . . 904.2.2 Mid-high temperatures in the non-perturbative regime . . . . . . 924.2.3 Coherent phase space and relaxation times: Ohmic regime for mid-high-temperatures . . . . . . . . . . . . . . . . . . . . . . . . . . . 944.2.4 3-site-Ohmic bath model in FMO . . . . . . . . . . . . . . . . . . 954.2.5 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 965 A dual-coupling-spin-boson model for light-harvesting molecules . . . 985.1 The dual-coupling polaron transformation . . . . . . . . . . . . . . . . . 995.2 Dynamics of the dual-coupling-spin-boson model . . . . . . . . . . . . . . 1005.3 Super-Ohmic spectral densities in the DCSB model . . . . . . . . . . . . 1055.4 High-temperature limit in the super-Ohmic NDSB model . . . . . . . . . 1075.4.1 DCSB model for the FMO system with acoustic phonons . . . . . 1115.5 The DCSB model coupled to optical phonons . . . . . . . . . . . . . . . 1125.5.1 DCSB model for the FMO system with optical phonons . . . . . . 1145.5.2 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1166 Summary, conclusions and future work . . . . . . . . . . . . . . . . . . . 1176.1 Further Work . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 122Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 125Appendices . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 137A Calculation of the electron-phonon correlation function for DCSBmodel . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 137B Connection between continuous and discrete electron-phonon corre-lation function . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 141C Full eigenvalues of detuned-3-site-V-system . . . . . . . . . . . . . . . . 143viiiD Fluctuation-dissipation theorem . . . . . . . . . . . . . . . . . . . . . . . 145E Cubic polynomials . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 146F Saddle point integration of optical phonon correlation functions . . . 149G Pure dephasing dynamics in the 3-site-V system . . . . . . . . . . . . . 151ixList of Tables1.1 Tuned V-system parameters. From Adolphs et al [46] . . . . . . . . . . . 262.1 FMO parameters for optical phonons characterised by their peak heightλ, full-width-half-maximum ξ and dimensionless coupling parameter ν =λξ/piω2o in units of meV . Taken from Olbrich et al. [114]. . . . . . . . . . 38xList of Figures1.1 Single BChla chlorophyll molecule. Adapted and reprinted with permis-sion from [63]. c© 2018 DOAJ . . . . . . . . . . . . . . . . . . . . . . . . 81.2 Composition of light harvesting systems within Chlorobaculum tepidum . 111.3 Structural depiction of the Fenna-Matthew-Olson complex from Chlorobac-ulum tepidum. (a) Full view of FMO trimer including 3 monomer subunits.(b) Single view of monomer sub unit plus surrounding proteins (grey lines).(c) Single monomer of 8 BChla pigments excluding surrounding proteinscaffold. Reprinted with permission from [71]. . . . . . . . . . . . . . . . 111.4 Quantum beating signatures for a 77 K FMO 2D spectrum. Axes ωτand ωt represent the Fourier transformed excitation and detection pulsesrespectively, of a pump and probe laser directed on the FMO sample. Thethird, vertical axis, ωT , represents the Fourier transform of the signalsevolving over time. The presence of cross-peak frequencies, evolving overtime indicates exciton delocalisation across pigments, and thus coherence.Reprinted with permission from [74]. . . . . . . . . . . . . . . . . . . . . 121.5 Pictorial representation of Hamiltonian for FMO complex, taken fromspectroscopic studies [46], in units of meV. Coloured lines indicate excita-tion energies of different chromophores (diagonal Hamiltonian entries).Arrows indicate inter-chromophore couplings (off-diagonal Hamiltonianentries). Numbers included next to lines and arrows indicate the asso-ciated energy with that excitation energy of coupling respectively. Onlydominant energy transfer pathways included; couplings > 4meV . . . . . 231.6 Effective 3-site-V configuration for FMO complex. . . . . . . . . . . . . . 24xi1.7 (A) Arrangement of BChla pigments within the FMO units. BChla sitenumbering according to Fenna is in black Roman numerals. Schematicrepresentation of spatial extent of the excitons according to Adolphs etal. [46] is shown by shaded areas. Exciton numbering (red numbers) isgiven in order of increasing energy. (B) Linear absorption spectrum ofFMO at 77 K with excitonic transitions (vide infra) represented by verti-cal bars. (C) Normalized absorptive 2D spectra at increasing populationdelays with dashed lines indicating excitonic transition energies. All spec-tra were recorded in 1:2 aqueous buffer:glycerol mixture at 77 K. Reprintedwith permission from [91]. . . . . . . . . . . . . . . . . . . . . . . . . . . 252.1 Spectral densities for the FMO trimer used in studies by Ishizaki andFleming [21], Cho et al. [103], Adolphs and Renger [46], as well as Nalbachet al. [26]. Spectral densities fitted to experimental data by Olbrich [114]also included, indicated by the lines ‘present - BChl 1-6’, ‘present - BChl7’, ‘present - BChl 8’, corresponding to the average over BChl molecules1-6, 7 and 8 respectively. The inset shows an enlarged energy and spectraldensity range. Reprinted with permission from [114]. . . . . . . . . . . . 372.2 An example of possible paths x, x′ ∈ {1, 0,−1} for the general 3-site systemand their corresponding symmetric/anti-symmetric paths χ ∈ {2, 1, 0,−1,−2},ξ ∈ {2, 1, 0,−1,−2}. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 442.3 Diagrammatic depiction of possible coherent states (blips) permitted in 3-site-V model under NIBA. ξ10 represents a wavefunction overlapping withstates |1〉 and |0〉, while ξ10 represents an overlap with |2〉 and |0〉 . . . . 472.4 Diagrammatic representation of the interaction terms contributing to theblip-blip interaction propagators Λjk. The contribution depicted is a nearest-neighbor interaction of blips. . . . . . . . . . . . . . . . . . . . . . . . . . 493.1 Eigenvalues λ0(red), λ+(blue), λ−(orange), as a function of (a) bias energy for the tuned system where 1 = − δ/2, 2 = + δ/2, (b) tunnellingmatrix element separation ∆ where ∆10 = −∆/2,∆20 = ∆/2, (c) detuningδ where ∆10 = 0.5,∆20 = 0.5, (d) detuning δ where ∆10 = 0.5,∆20 = 1 .In units of . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 543.2 Diagrammatric representation of 3-site-V model coupled to a sink via anirreversible decay channel from the ground state . . . . . . . . . . . . . . 60xii3.3 Return probability to site 1: P11(t), as a function of time for tuned (red)and detuned (blue) cases. In units of 1 = 1. 2 = 1.51(blue), 2 = 1 =1(red). (a): Γ = 0.1, ∆10 = 0.5, ∆20 = 0.3, (b): Γ = 0.1, ∆10 = 1,∆20 =0.8, (c): Γ = 0.1, ∆10 = 0.5,∆20 = 0.5, (d): Γ = 0.5, ∆10 = 0.2,∆20 = 0.1,(e): Γ = 0.5, ∆10 = 0.1,∆20 = 0.1, (f): Γ = 0.5, ∆10 = 1,∆20 = 0.1 . . . . 613.4 Return probability to ground state |0〉: P00(t), as a function of timefor tuned (red) and detuned (blue) cases. In units of 1 = 1. 2 =1.51(blue),2 = 1 = 1(red). (a): Γ = 0.1, ∆10 = 0.5, ∆20 = 0.3, (b):Γ = 0.1, ∆10 = 1,∆20 = 0.8, (c): Γ = 0.1, ∆10 = 0.5,∆20 = 0.5, (d):Γ = 0.5, ∆10 = 0.2,∆20 = 0.1, (e): Γ = 0.5, ∆10 = 0.1,∆20 = 0.1, (f):Γ = 0.5, ∆10 = 1,∆20 = 0.1 . . . . . . . . . . . . . . . . . . . . . . . . . . 633.5 Return probability to state |1〉: P11(t), as a function of time for tunedbiases (red) and detuned biases (blue), as well as tuned tunnelling (green)and detuned tunnelling (purple), in the small bias regime. 2 = 0.2, 1 =0.1 (blue), 2 = 1 = 0.1(red). (a): Γ = 0.1, ∆10 = 1, ∆20 = 0.8, (b):Γ = 0.01, ∆10 = 1,∆20 = 0.8, (c): Γ = 0.5, ∆10 = 1,∆20 = 0.8, (d):Γ = 0.5, ∆10 = 0.2,∆20 = 0.1, (e): Γ = 0.1, ∆10 = 1,∆20 = 0.8 . . . . . . . 643.6 Return probability to site 1 for unbiased system. ∆10 = 1, ∆20 = 0.5(purple), ∆10 = 1, ∆20 = 1(green). (a): Γ = 0.1, (b): Γ = 0.5 . . . . . . . 653.7 Return probability to site 1 at t→∞: P11(∞), as a function of detuningδ = 1 − 2. ∆10 = ∆20 = 0.1, = 1 (a) Γ = 0 (Bath off) (b) Dark-statepeak for various decay rates . . . . . . . . . . . . . . . . . . . . . . . . . 653.8 Long time return probability P11(∞) to site 1 vs ∆10: ∆20 = ∆20 = 0.1, = 1; (a) δ = 1 (b) δ = 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . 663.9 Coherent-incoherent phase space of tuned-perturbative 3-site-bath model.Parameters: ∆10 = 0.1, ∆20 = 0.12 γ1 = 0.01, γ2 = 0.05, in units of .(blue) Ohmic, (orange) super-Ohmic. . . . . . . . . . . . . . . . . . . . . 763.10 (a) Pure relaxation times vs kBT ; τ±r = 3/Φ(E) (orange/green), τDr =3/Φ() (blue), (b) Full relaxation times vs kBT ; τ± = 1/Reλ˜± (orange/green),τD = 1/Reλ˜D (blue). Parameters: γ1 = 0.01, γ2 = 0.05, = 1, ∆10 = 0.1,∆20 = 0.2, . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 773.11 (a) Pure relaxation times vs kBT ; τ±r = 3/Φ(E) (orange/green), τDr =3/Φ() (blue), (b) Full relaxation times vs kBT ; τ± = 1/Reλ˜± (orange/green),τD = 1/Reλ˜D (blue). Parameters: γ1 = 0.01, γ2 = 0.05, = 1, ∆10 = 0.1,∆20 = 0.2, . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79xiii3.12 Decay times for 3-site model and optical phonons. Parameters: ∆10 = 0.1,∆20 = 0.2, Ξ1 = 0.01, Ξ2 = 0.01, Γ = 0.001, ω0 = 10. . . . . . . . . . . . . 803.13 Phase space for 3-site model and optical phonons. τ± (orange/green),τD (blue). Parameters: ∆10 = 0.1, ∆20 = 0.2, Ξ1 = 0.01, Ξ2 = 0.01,Γ = 0.001, ω0 = 10. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 814.1 F1 vs γ1 figures (a) and (b), F1 vs kBT figures (c) and (d). Parameters:γ0 = 0, γ2 = 0.2, 1 = 1, 2 = 1.5. (a) ∆10 = 0.1,∆20 = 0.2, kBT = 10 (b)∆10 = 0.1,∆20 = 0.2, kBT = 0.1, (c) ∆10 = 0.1,∆20 = 0.2, γ1 = 0.2, (d)∆10 = 0.1,∆20 = 0.2, γ1 = 0.8 . . . . . . . . . . . . . . . . . . . . . . . . 914.2 F1 vs γ2. Parameters: γ0 = 0, γ1 = 0.2, 1 = 2 = 1, ∆1 = ∆2 = ∆ = 1(a) γ2 = 0.2, (b) kBT = 10 . . . . . . . . . . . . . . . . . . . . . . . . . . 924.3 Phase space, decay rates, relaxation times and oscillation frequency vs γ(coherent regime). ΓD, τD (blue), Γ−, τ+ (green), Γ+, τ− (orange). Param-eters: ∆1 = 0.1, ∆2 = 0.05, = 1, kBT = 10 . . . . . . . . . . . . . . . . 954.4 Phase space and relaxation times of FMO complex at T=300K in units offemtoseconds with Ohmic bath. FMO parameters: ∆˜10 = 2meV , ∆˜20 =3meV , = 20meV . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 965.1 Phase space, decay rates and oscillation frequency for Non-diagonal-SBmodel vs ζx. ζz = 1. Parameters: = 1, ∆ = 0.3, kBT = 5, ωc = 0.1. Inunits of . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1105.2 Exciton transfer time τr, coherence time τφ, and dimer oscillation fre-quency ω as a function of non-diagonal coupling strength ζx. Parameters: = 20meV, ζz + ζx = 1, kBT = 27meV, ωc = 8.7meV. Bath spectraldensity of super-Ohmic form: J(x,z)(ω) = pi~ ζ(x,z)(ω3/ω2c )exp(−ω/ωc). . . 1115.3 (a) ∆˜ = 0.3, = 1, kBT = 1, ωo = 0.4, ξ = 0.01, λz = 5 in units of (meV). 1145.4 Exciton transfer time τr (blue line), coherence time τφ (orange line), anddimer oscillation frequency ω (purple line) as a function of non-diagonalcoupling strength νx. Parameters: = 20meV, ∆ = 6meV, ζz=1, kBT =27meV, ω0 = 6meV. Bath spectral density of Gaussian-optical form: J(x,z)(ω) =λ(x,z)exp(−(ω − ω0)2/2ξ2). . . . . . . . . . . . . . . . . . . . . . . . . . . 115xivGlossaryAntenna complex : A system containing many individual chromophores which actsas the first point of entry for photon absorption in photosynthesis.Bacteriochlorophyll A (BChla): The principle light harvesting molecule present ingreen sulfur bacteria. Contains a magnesium atom at its centre, surrounded by 4nitrogen atoms. A delocalised pi-electron spans the entire molecule.Bath correlation function: Determines how the environment fluctuations affect thecentral system through the system-environment coupling λ. This tells us how aperturbation of the environment caused by the central system through λ, affectsthe system at a later time through λ. Fluctuations of the bath coupled to thecentral system at some time are correlated with fluctuations at some later time.Boson: A particle with integer spin that obeys Bose-Einstein statistics, so there isno restriction of the number of them occupying the same quantum state.Coherence: Describes a degree of correlation between physical quantities of a sin-gle wave, or between several waves. In in the context of wave mechanics, andby extension the wave-like aspect to quantum mechanics, coherence describes theconstructive interference of waves (wavefunctions). Therefore coherence dependson the relative phase of the two waves. Temporal/spatial coherence describes thecorrelation beween waves observed at different moments in time/points in space.Therefore coherence between energy states in a quantum mechanical system de-scribes the overlap of the wavefunctions corresponding to each state. A particletunnelling coherently between states is said to be ’delocalised’ across the statesinvolved in the tunnelling process.Central system: In the language of open quantum systems, the central system de-scribes the system of interest. It usually contains only a few degrees of freedomxvcompared to its environment which contains many. We are usually interested inobserving physical quantities pertaining to the central system under the influenceof its environment. Therefore the variables pertaining to the environment are ’in-tegrated out’, leaving our dynamical equations in terms of variables correspondingto the central system only.Chromophore: A light-absorbing-conjugated molecule that contains extended pi-orbitals, the energy difference of which falls within the visible spectrum. Ofteninterchanged with ’pigment’.Cut-off frequency ωc: A high frequency mode of the bath that regularises the bathspectral density to negligibly small values for ω ωc. This fixes the ultra-violetdivergences present in power-law spectral densities and represents a physical limitto the frequency of bath fluctuations.Dark state: A term coined in quantum optics, a dark state, represents an atom ormolecule that cannot absorb or emit photons. In the molecular context, a dark-state represents an eigenstate of a molecular system (or interacting molecules), thatdoesn’t couple to any external bath. Therefore an exciton that occupies the dark-state will remain ’trapped’. The external bath does not ’see’ the dark-state, as itdoes not couple to it, and therefore offers the exciton no means of leaving that statewithin the system.Decoherence: The loss of quantum coherence from a central system of interestdue to interactions with an external environment. In the language of quantummeasurement theory, decoherence represents the process of wavefunction collapsewithin the central system due to measurements made by an external environment.Density matrix : A quantum mechanical operator, or matrix, that describes the sta-tistical state of a system. The density matrix is capable of describing a statisticalmixture of states whereas the quantum mechanical state vector can only describea ’pure state’ with no mixture. Mixed states arise in situations where the experi-menter does not know which particular states are being manipulated. Therefore thedensity matrix contains information regarding superpositions of states and there-fore quantum mechanical coherences. The diagonal elements of the density matrixrepresent the probability of finding the system in a particular state, and are there-fore dubbed the ’populations’. The off-diagonal elements on the other hand containxviinformation regarding interference effects and are correspondingly dubbed ’coher-ences’. A particularly illuminating form for the density matrix comes by way ofthe path-integral formalism of quantum mechanics. In this case the density matrixrepresents two independent paths including possible interference effects betweenthem. Therefore the diagonal elements represent both paths starting and endingin the same state, whereas the off-diagonal elements represent different endpointsof the paths and therefore superpositions between states.Diagonal coupling : A system-environment interaction that couples the environmentto the diagonal elements of the Hamiltonian. Therefore the diagonal couplingcorresponds to fluctuations in the on-site (or kinetic energies) of the central system.Often referred to as ’local’ coupling due to the environment coupling to quantitiespertaining to individual states in the system..Dimer : A molecular aggregate containing two constituent molecules, one acting asthe ’donor’ and the other the ’acceptor’.Direct-Coulomb interaction: A pair interaction term in the Hamiltonian. It de-scribes the classical Coulomb interaction between charge densities due to electronsoccupying orbitals on different molecules. For the long range Coulomb interaction,the point-dipole approximation is invoked which results in the dipole-dipole ap-proximation of pair interactions. In molecular transport theory, retention of onlythis potential term in the Hamiltonian results in Fo¨rster theory.Donor/acceptor molecules : A donor molecule represents the point of entry of anexciton to the system. The acceptor molecule, coupled to the donor, represents theadditional state in the two-level system that the exciton can tunnel to.Dual-coupling-spin-boson model (DCSB): An extension to the well known Spin-boson (SB) model that includes both diagonal and non-diagonal system-bath cou-plings. The SB model only includes diagonal couplings.Environment (bath): The environment surrounding the central system of interest,usually composed of a very large number of degrees (d.o.f) of freedom with respectto the central system d.o.f. When coupled to the central system, the bath d.o.f in-fluence the dynamics, and quantum properties like coherence, of the central systemd.o.f. Usually one is concerned with the time-evolution of the central system d.o.fonly, so the bath d.o.f are ’integrated out’.xviiExciton: A bound state of an electron-hole pair. The electron and hole interactvia the Coulmomb interaction forming a charge-neutral quasiparticle that has lessenergy than the unbound electron and hole. Exciton’s can be highly delocalised dueto the screening force of surrounding electrons. However, in amorphous biosystems,exciton delocalisation is suppressed due to disorder by way of Anderson localisation.Exchange interaction: A pair interaction term in the Hamiltonian that includesan appreciable overlap of neighbouring molecular wavefunctions. Inclusion of thisterm permits inelastic tunnelling processes involving particle exchange. Unlike thedirect-Coulomb interaction, this term has no classical analog. With the inclusionof an external bath coupled to the system, this term permits tunnelling processesinvolving particle exchange with the bath.Fenna-Matthews-Olson complex : A water-soluble light-harvesting complex that ex-ists in green sulfur bacteria. It is composed of BChla molecules and mediates theenergy transfer from light-harvesting chlorosomes in the antenna complex to thereaction center.Fermion: A half-integer-spin particle that is subject to the Pauli exclusion principleand therefore obeys Fermi-Dirac statistics.Fermi’s Golden Rule (FGR): A formula describing the transition rate from oneenergy eigenstate of a quantum system to a continuum of energy eigenstates. Theformula is derived using time-dependent perturbation theory in the tunnelling ma-trix element between initial and final states. Therefore it is assumed that thetunnelling time is much longer than the transition time to the continuum whichputs the theory in the regime of very strong system-bath coupling.Fo¨rster theory : A molecular theory of energy transport between two chromophores.A donor molecule may transfer an exciton to an acceptor molecule via the direct-Coulomb interaction (dipole-dipole coupling). With the inclusion of spectral broad-ening to the spectral lineshapes of the donor and acceptor molecule, the equationsof motion are derived with use of a FGR approach, treating the exciton tunnellingterm as a perturbation. Fo¨rster theory therefore describes incoherent hopping ofan exciton between chromophores.Green sulfur bacteria: Photosynthetic organisms that contain the FMO complexas part of their light-harvesting mechanism. They can be found living in the low-xviiiest light intensity regions of any known photosynthetic organism and exhibit aremarkably efficient energy transfer process, from the point of photon capture tophotochemsitry.Green’s function: See Propagator.Hamiltonian: A quantum mechanical operator, or matrix, corresponding to thetotal energy of the system.Highest-occupied-molecular orbital (HOMO): A molecular orbital representing thehighest energy state of the molecule occupied by an electron.Ideal-dipole approximation: An approximation valid when the wavelength asso-ciated with the energy of the dipole is much larger than the size of the atom.Therefore the dipole is considered to be a point in space.Inelastic scattering process : A scattering process in which the kinetic energy of aparticle is not conserved. The energy can be lost or increased by exchange with itsenvironment.Instanton: A semi-classical solution to the equations of motion regarding tunnellingbetween potential wells. In the path integral formalism, the instanton solutiondescribes instantaneous jumps between potential wells, hence the name.Light-harvesting complex/molecule (LHC): A photosynthetic complex of subunitproteins that may be part of a larger supercomplex. A LHC therefore forms part ofthe functional process of photosynthesis. They contain a number of chromophores,used to capture a photon of light and/or shuttle an exciton through the complexto the next part of the photosynthesis process.Localisation: The process by which an atomic or molecular wavefunction is confinedto the atom or molecule due to wavefunction collapse. This is the opposite of thedelocalised limit where a wavefunction spans several or more atoms or moleculesand retains coherence effects. Localisation therefore reflects the incoherent limit ofparticle hopping.Lowest-unoccupied-molecular orbital (LUMO): A molecular orbital representing thelowest energy state of the molecule’s excited states. The LUMO state is unoccupiedin the ground state of the molecule, and occupied by an electron in the 1st-excitedstate.xixMarkovian approximation: An approximation used to simplify the equations ofmotion for an open quantum system. It assumes that the correlation time of thebath dynamics is much smaller than the characteristic timescale of central systemdynamics. Therefore ’memory effects’ induced by the bath are neglected. In otherwords a bath mode excited by the central system immediately relaxes to its groundstate, from the perspective of the central system, and does not affect the centralsystem dynamics at a later time.Monomer : A single molecule that is classified in this way as it can undergo ’poly-merisation’ in certain circumstances leading to polymer formation involving multi-ple monomers.Non-diagonal coupling : A system-environment interaction that couples the envi-ronment to the off-diagonal elements of the central system Hamiltonian. Thereforefluctuations in the environment modulate the tunnelling matrix elements (potentialenergy terms/transfer integrals) of the central system. Often referred to as ’nonlo-cal’ coupling due to the environment coupling to quantities pertaining to transitionsbetween states in the system.Non-interacting-blip approximation (NIBA): An approximation used to simplify theequations of motion (e.o.m) of a central system interacting with its environment(bath). A general e.o.m will include bath correlations that extend to t → ∞. Ina diagrammatic language this means that all higher order system-bath interactiondiagrams as well as bath-bath interaction diagrams are considered. NIBA usesthe assumption that the system spends vastly more time in diagonal states ofthe reduced density matrix for the system, than off-diagonal states. In a pathintegral language, it permits the system to spend one ’time block’ in an off-diagonalstate before returning to a diagonal state. In a diagrammatric language NIBAretains only the self-interaction term of the bath correlation function and excludesany higher order bath-bath interaction diagrams. The approximation also retainsthe system-bath interaction diagram pertaining to neighbouring interactions only.Since these diagrams are exponentiated in the propagator describing the centralsystem dynamics, they are still considered to all orders in the system-bath couplingenergy. Therefore NIBA represents a non-perturbative, quasi-coherent limit to thedynamics of the system.Ohmic/super-Ohmic/Optical phonon spectral density : Ohmic and super-Ohmic de-xxscribe algebraic frequency dependence of the spectral density with linear orderand powers higher than one respectively. A regularisation function with a high-frequency cut-off is required to avoid the ultraviolet divergences. An optical spectraldensity represents a peak in the profile. The peak is centered around a particularfrequency and therefore represents a collection of bath excitations with just onefrequency. Hence the application to optical phonon spectral densities.On-site (bias) energy : The kinetic energy term in the Hamiltonian. For interactingmolecules in the two-state limit, the on-site energy represents the energy differencebetween the two molecules LUMOs, which form the basis for the interacting system.Open quantum system: A global quantum system that involves some ’central sys-tem’ of ’interest’, and a classical or quantum environment that may be coupled tothe central system.Oscillator bath: An environment (bath) surrounding a central system, comprisedof a large number of simple-harmonic oscillators.Path-integral : An alternative to the classical single trajectory of a system. Insteadthe trajectory involves a sum, or functional integral, over an infinite number ofpossible trajectories permitted by quantum mechanics.Perturbation theory : An approximate solution to a system where the full, com-plex, system has a small quantity relative to other comparable quantities. Thesmall quantity is considered separately first, with the remaining part of the systemcontaining a known solution. Then the small quantity is reintroduced, and the ap-proximate solution to the full problem is calculated using the fact that correctionis a ’perturbation’ from the known solution.pi-bond/pi-electron: A covalent bond where two ’lobes’ of an orbital on one atomis shared with two lobes on another atom. pi-bonds tend to be much weaker thanother covalent bonds and can result in a delocalised pi-electron across the molecule.Pigment : See ChromophorePhonon: acoustic: An ’in-phase’ collective motion of atoms out of their equilib-rium position. Acoustic phonons exhibit a linear relationship between frequencyand phonon wavevector for long wavelengths which tends to zero for the longestwavelengths; the limit relevant to amorphous materials.xxiPhonon: optical : An ’out-of-phase’ collective movement of atoms where one atommoves in the opposite direction to its neighbour. Optical phonons show no dis-persion, and therefore constant energy, in the long-wavelength limit which is thelimit relevant to amorphous materials. In molecular theory, optical phonons inter-acting with a molecule, usually come from intramolecular bonds. For e.g. BChlamolecules contain a Carbon-Carbon bond that oscillates within the optical rangeand couples strongly to delocalised pi-electron over the molecule.Phonon-assisted tunnelling (PAT): A tunnelling process that is assisted by inter-actions with phonons in the environment. PAT is usually characterised as elasticor inelastic phonon-assisted tunnelling. Elastic PAT results from only diagonalsystem-bath couplings in a certain parameter space where the tunnelling rate is fa-cilitated by the bath. Inelastic PAT arises from non-diagonal system-bath couplingswhich permits a tunnelling process via phonon exchange.Polaron: A quasi-particle formed from an electron and phonon ’cloud’. As anelectron moves through a solid, it displaces the atoms around it from their equi-librium positions, dragging them along with it. This results in a bound state ofthe electron-phonon cloud, lowering the energy compared to the non-interactingsystem.Propagator : A function that specifies the probability amplitude of a particle to gofrom one state to another. State in this case is purely general and could describeanything from position to abstract energy eigenstates. Propagators in the contextof quantum field theory are often referred to as Green’s functions.Quantum beating signal : Time-dependent oscillations in the cross-peak correla-tions measured using 2D-spectroscopy. Quantum beating signals are interpreted asmeasuring coherence in a system.Quantum master equations (QME): A completely general QME simply representsthe dynamics of the density matrix as opposed to just a quantum state vector. Thisway coherences, as well as site populations, can be tracked over time. The termQME are often however used interchangeably with Lindblad, or Redfield, equationsand therefore represents a restricted class of equations subject to the secular andMarkovian approximations.xxiiReaction center (RC): A sub unit of the photosynthesis mechanism of Green sulfurbacteria. It is situated after the FMO complex, harvesting the exciton to createa charge separation. This provides the energy required for the photochemistrypertaining to photosynthesis. The RC therefore represents the target destinationof the exciton.Redfield equation: An approximate QME that utilises the Markovian and secularapproximations.Relaxation time: The timescale associated with the inverse of a decay rate. Thepure relaxation time refers to the inverse of the decay rate out of the high-energystate to the lower-energy state. Dephasing time is the timescale associated withthe decay rate of oscillatory terms in the exciton dynamics.Reorganisation energy : The energy associated with the bath reconfiguration back toequilibrium. This is a measurable quantity in many cases, and is the energy emittedby the bath upon relaxing to its ground state after excitation due to interactionwith an exciton.Resonance/exciton-energy transfer (RET/EET): A process describing energy, orexciton, transfer from one (donor) molecule to an (acceptor) molecule in biologicallight-harvestin complexes.Resonant states : A condition when two energy eigenstates are equal.Secular approximation: An approximation to the dynamics of a central systemin which rapidly oscillating terms in the Markovian quantum master equation aredisregarded.Solvent : The liquid solution that permeates the surrounding regions of the chro-mophores in photosynthetic organisms. Usually constitutes water and electrolytesgiving the solvent dialectic properties.Spectral broadening : The broadening of an unphysical delta function peak in aspectrum to a physical peak with some finite linewidth. The linewidth comes frominteraction processes with external systems, representing an exchange of energywith the environment. A delta function peak is said to be unphysical in principlebecause a quantum system can never completely isolated from its environment.xxiiiAt the very least, there will also exist some broadening due to the Heisenberguncertainty principle.Spectral density : A function that maps out the system-environment interactionsacross the range of bath frequencies. It therefore depends on the distribution ofbath frequencies, the density of states, and the system-bath interaction energy toeach bath mode.Spin-boson model (SB): A mathematical model that includes a spin-1/2 particle(two-state system) coupled non-perturbatively to an oscillator bath. NIBA is in-voked to close the system of equations of motion and achieve results in various lim-its, and for various bath spectral densities. The spin-boson mode is valid in largeregions of the parameter space due to the non-perturbative path-integral techniquesemployed. However, for arbitrary system-bath couplings, it can generally be saidthat the SB model is valid for large bias energies relative to tunnelling energies.Three-site-V system: A quantum mechanical central system that contains threestates, two of which permit tunnelling between them. Configured in a V shape, theupper two sites are non-interacting and no tunnelling can occur between them. Thelower site is coupled to a continuum of states, or in the context of quantum optics,a laser field. In either case the lower site can be depleted due to its interaction withan external source.Transition dipole moment : An electric dipole moment associated with the transitionbetween two states. Molecular orbitals have different charge distributions and thustwo orbitals, or states, will have a transition dipole moment associated with them.Trimer : A molecular aggregate containing three units, each unit containing one ormore molecules.Tunnelling matrix element ∆: A term in the central system Hamiltonian represent-ing a coupling between different states. Often used interchangeably with transferintegral or inter-state coupling.Two-dimensional spectroscopy : Correlates excitation and emission energies of asample as a function of delay time between events. Spectra are plotted as a func-tion of absorption and emission where diagonal peaks represent on-site energiesand cross-peaks represent transitions. This essentially measures the Hamiltonian.xxivThe time-dependence of cross-peak rephasing can be observed to measure coher-ence. This essentially measures the density matrix and the off-diagonal componentsidentified as ’quantum beating signals’.xxvChapter 1Introduction1.1 Background and motivationUnderstanding the role of quantum coherence effects in molecular systems has been thesubject of intensive research for a number of decades. The obvious applications to quan-tum computing [1] and energy-transport mechanisms at the nanoscale [2] have been aprimary motivating factor, not to mention the general interest that physicists have inunderstanding the quantum-to-classical transition. Quantum coherence effects are wellunderstood at the atomic scale, and it is generally accepted that they can manifest atthe molecular scale under optimum conditions. However for most systems comprised of alarge number of molecules, in contact with a densely packed local environment at physi-ological temperatures, we expect quantum coherence effects to be absent. The ‘hot’ and‘messy’ environments that characterise many molecular systems involve a large numberof rapidly fluctuating degrees of freedom which overall lead to a relatively strong system-environment interaction. These are far from optimal conditions for quantum coherenceeffects and one would expect molecular wavefunctions to be fairly well localised on eachmolecule, representing the incoherent regime of energy transport. In other words, thetimescale for coherence effects should be relatively short compared to the timescale ofdynamics of energy transport and thus coherence effects are not expected to play a sig-nificant role in the process. This is indeed the case for systems such as betaine dyemolecules [3] and charge transport in trans-polyacetylene [4], which see coherence timesof the order of femtoseconds, in line with inter-molecular transfer times. As a result, aclassical description of the energy transfer mechanisms usually suffices in such systemsand the role of quantum coherence effects can largely be ignored.1A landmark experiment performed by Engel et al. in 2007 [5] reported the discovery ofrelatively ‘long-lived quantum-beating signals’ in the prototypical photosynthetic systemChlorobaculum tepidum (green sulfur bacteria) at cryogenic temperatures. Other exper-iments found similar coherence effects in the system Rhodobacter sphaeroides [6], also atcryogenic temperatures, and furthermore, subsequent experiments [7, 8] found coherenceeffects persisted even at physiological temperatures. To this date, the functional role ofquantum mechanics in the energy transport process of these particular photosyntheticsystems remains unclear [9], and considering the delicate nature of quantum coherenceeffects at such scales and temperatures, it is a surprising discovery to find the presenceof these so-called ‘quantum-beating signals’ in biomolecules.It has long been known that certain photosynthetic organisms can harvest sunlightwith near perfect efficiency; however, the application of classical energy-transport modelshas proven unsuccessful at reproducing the observed efficiencies [10, 11, 12, 13, 14]. Thequestion remains, to what degree do quantum effects play a functional role in optimisingthe photosynthesis process?Since the publication of the experiments mentioned above, a plethora of theoreticalwork has ensued aiming to understand these findings. However, the degree to whichquantum coherence effects actually facilitate the energy transfer process in photosyn-thesis remains controversial. Recently a number of theoretical studies have suggestedthat quantum effects might actually play a key role in optimising the efficiency of energytransfer in biomolecules [15, 16, 17, 18, 19, 20, 21, 22, 23]. Most theoretical work todate has relied upon a number of assumptions that are arguably invalid when it comesto the specific systems of interest: the mathematical models employed have either beenperturbative, or assumed incoherent energy transfer dynamics from the outset. Neitherare able to predict the relatively long coherence times observed in experiment.The task of including the full breadth of coherence effects in such a complex systemis formidable, and most studies so far have resorted to fairly extreme approximations inorder to make the problem tractable. These approximations tend to exclude any possibil-ity of coherent dynamics from the outset, producing models that are ill-equipped to findcoherent behaviour. Attempts to go beyond these limiting models have had promising re-sults [24, 25, 26, 27], and the predicted coherence times are now comparable to those seenin experiments. However, much of this progress has so far been facilitated to a greateror lesser extent by numerics, and a rigorous analytical framework of understanding islacking. While the numerical results are in accordance with observation, the underlyingphysical mechanisms remain elusive and a deeper analytical understanding is required.2One mathematical model that treats the interaction of a central spin-half system withits environment non-perturbatively is the spin-boson (SB) model. Originally introducedby Emery and Luther [28] to shed light on the Kondo problem, they demonstrated thatthe phase space of the low-temperature Kondo model can be understood with the useof the equivalent SB model. Caldeira and Leggett [29, 30, 31, 32], as well as a numberof contemporaries [33, 34, 35, 36], used the model to study the effects of dissipation inthe dynamics of a particle hopping between two states. The model used the influencefunctional technique developed by Feynman and Vernon [37] to analyse the effects ofthe environment non-perturbatively on the central two-state system. Functional integralmethods were also employed to study a particle tunnelling in a periodic potential inter-acting with its environment [38].Quantum coherence effects are retained up to a certain degree in the SB model, andit has been a powerful tool for modelling 2-state systems interacting with their environ-ments. The model has had many applications ranging from quantum computation [39]to condensed matter and solid-state physics [40]. It has also seen some cursory appli-cations to the photosynthesis mechanism [41, 42, 43, 27]. Here, certain 2-level systemswith appreciable wavefunction overlap and the appropriate parameters were selected fromspecific photosynthetic systems. It was also found recently that in a slight extension tothe SB model, where certain higher order coherence effects were taken in to account,decoherence times calculated in the system matched those observed experimentally [27].Such a direct application of the results of the SB model to photosynthesis is of courserestricted to just two states interacting with a shared environment. The energy-transportstructure in photosynthesis generally involves more than just two sites and the larger,more representative structure, has generally been approached numerically due to its com-plexity. It has however been argued that for the photosynthetic systems of interest, thefull multi-level structure effectively reduces to a 3-site model [44, 45]. Furthermore, ex-citon wavefunctions in the FMO complex have been shown to be delocalised across twoto three molecules at a time [46]. Therefore in the interest of accurately modelling theenergy transport mechanism in photosynthesis analytically, we propose the extension ofthe SB model to three sites for one model of study in this thesis. This will not onlymodel coherent effects analytically in photosynthesis, but will give a more accurate rep-resentation of the structure of the photosynthesis mechanism.Another limiting aspect of the SB model is the fact that it only considers ‘diagonal’system-environment couplings. This kind of system-environment interaction couples theenvironment to the diagonal elements of the central system Hamiltonian only. In other3words, the diagonal interaction induces fluctuations in the potential energy terms of thetwo-state system, but not the kinetic energy terms, which constitute the non-diagonalcomponents of the Hamiltonian. In the original theoretical treatment of the SB model,it was argued that one can ignore couplings to the non-diagonal-Hamiltonian elements,provided the tunnelling matrix elements–the kinetic energy terms–are small compared tothe on-site energies–the potential energies. These non-diagonal couplings were first intro-duced by Holstein in 1959 [47] in the context of polaron motion within a one-dimensionalcrystal lattice. Later they would be discussed in the context of the Peierls transition[48] and subsequently in the Su-Schrieffer-Heeger model of polyacetylene [49]. Thesenon-diagonal couplings correspond to transitions in the central system that involve anabsorption or emission of a boson from the environment. As the non-diagonal system-environment coupling energy should be proportional to the tunnelling energy [50], oneexpects these couplings to be relatively small in the limit of small tunnelling energy.Indeed, the limit of large on-site energy to tunnelling energy is also the limit in whichthe noninteracting-blip approximation—one of the key approximations used to solve theSB model—is valid.The problem with this heuristic justification for the exclusion of the non-diagonalcouplings in the original SB model is that, while the tunnelling energy is indeed smallrelative to the bias energy, it is not negligibly small. In fact, the SB model is consid-ered such a powerful model because it retains coherence effects (up to a degree), andcontributing to this, is a small but significant tunnelling energy. Numerical studies mod-elling polaron formation on a lattice, with the inclusion of both non-diagonal as well asdiagonal couplings, demonstrated how even small values of the non-diagonal couplingstrength greatly alter the polaron properties [50]. Experimental evidence for significantnon-diagonal exciton-phonon couplings in photosynthetic systems has also surfaced [51].Despite this, the inclusion of non-diagonal couplings in theoretical modelling of photo-synthesis has largely been ignored. As Mahan [52] remarked, the system-environmentcouplings of this type are usually excluded due to the difficulty of obtaining reliable so-lutions with them [53].Some preliminary theoretical work has been done exploring the role that non-diagonalcouplings play in photosynthesis [18, 53]. However, in these cases, some fairly limitingassumptions were employed. In [18] each state in the system was assumed to interactonly with its own local, and independent, environment–meaning that the possibility ofenvironment induced transitions between states was ignored. For the effects of non-diagonal couplings on coherent dynamics to be truly resolved, one needs to consider a4shared environment between the states in the system, as is done in this work. In [53],a shared environment was included, however a quantum-master-equation approach wasused, assuming that the system-environment interaction could be treated perturbatively.A perturbative approach has been argued to be valid for strong diagonal system-bathcouplings modelling the residual system-bath interaction for a system transformed intothe polaron frame [54]. However, the inclusion of non-diagonal system-bath couplings–which can be much smaller than diagonal couplings–can render this approach invalid.Therefore a non-perturbative approach, even in the polaron transformed frame, is neces-sary. The importance of non-diagonal couplings in the photosynthesis modelling processis apparent now, however a proper treatment of the system-environment interaction in-cluding these couplings is still lacking.Aside from the general interest that one might have in the inclusion of non-diagonalcouplings to the SB model and photosynthesis, the need for such an extension to themodel has found applications in other areas too, for example in the modelling of chargetransfer in polymeric solar cells (PSC) [55]. Non-diagonal couplings have also been shownto play an important role in charge trasnfer modelling of organic semiconductors [56].Recent experiments performed on PSC’s have found ‘ultrafast quantum beating signals‘with surprisingly short periods [57] not unlike those found in photosynthesis. It has beensuggested that phonon-assisted transfer may be responsible for this, as the tunnellingenergies present in the system cannot account for the ultrafast oscillation frequenciesalone [58]. Perhaps the same mechanism is responsible for the long-lived coherence ef-fects observed in photosynthesis too?One major question that arises from all these considerations is: to what degree doescoherence play a functional role in exciton transport in biomolecules? While it is de-sirable to explain the emergence of coherence in photosynthesis, the coherent nature oftransport in these systems could be just a coincidence unless a specific functional roleis identified. Coherence-enhanced functionality in biological systems has been discussedpreviously involving constructive interference of exciton pathways through molecules,leading to improved transfer efficiency [59]. However, how these phenomena apply tophotosynthesis remains unclear. In this work, the importance of non-diagonal couplingsin modelling exciton transfer is demonstrated. Since couplings of this nature arise whenthere is appreciable exciton wavefunction overlap across neighbouring molecules, the func-tional role of coherent exciton dynamics is possibly due to this effect. In other words, thepropagation of coherent excitons through the system, facilitates the process of phonon-assisted-transport between molecules, greatly improving transport efficiency.5In summary, the original work presented in this thesis primarily constitutes the fol-lowing: a new analytical model for the photosynthesis mechanism involving a dimer non-perturbatively coupled both diagonally and non-diagonally to an oscillator bath, and a3-site model non-perturbatively coupled diagonally to an oscillator bath. Both modelsshed light on important processes present in the photosynthesis mechanism. The dualcoupling dimer model, where both diagonal and non-diagonal system-bath couplings areincluded, reveals the importance of non-diagonal couplings in not only assisting excitontransfer through the system, but also retaining coherence effects at the same time. The3-site model demonstrates the potential importance of the dark-state in the photosyn-thesis mechanism. The dark-state is a state intrinsic to the 3-site-V configuration, andthe application of this model to the FMO complex reveals this state to have the longestrelaxation times in the system.Relaxation times in both models are found to be comparable to those observed ex-perimentally. The inclusion of non-diagonal couplings is shown to have a particularlystrong influence on both the rate of exciton transfer through the system, as well as thecoherence of the exciton. Three different spectral densities, that characterise the envi-ronment, are used and their effects on exciton dynamics explored. These constitute twolow-frequency distributions of the environment oscillator modes: known as Ohmic andsuper-Ohmic, and one discrete optical mode with spectral broadening. The inclusion ofoptical phonons in the environment is shown to be especially important for the FMOcomplex, which exhibits an experimentally determined spectrum with multiple discretepeaks. Both optical and acoustic phonons are shown to aid in exciton transfer througha photosynthetic dimer system. The non-diagonal system-bath coupling is shown to notonly decrease exciton transfer times for relatively large coupling strengths, but also in-troduce a further decoherence mechanism in to the system. However, while the excitoncoherence is shown to decrease with increasing non-diagonal coupling, the coherence timeremains longer than the exciton transfer time. This indicates that the exciton motionthrough the dimer remains coherent for its full duration.An additional novel feature of the present work involves the retention of the on-siteenergies in modelling the photosynthesis mechanism with a dual coupling approach. Thespin-boson model, with diagonal couplings only, is often simplified to the resonant casewhere the energy difference between the two central system states is zero [33, 60, 61].Resonant two-level systems, coupled diagonally to the bath, produce a simple, quadraticpole structure to the propagator which is easily dealt with analytically. Dual couplingmodels for photosynthesis have also looked at the resonant case [53], however the validity6of such a simplification is dubious. The non-diagonal coupling serves to renormalise theon-site energies in the polaron frame, and setting these energies to zero after the polarontransformation trivialises the non-diagonal interaction. In this work, the on-site energyis retained, and the cubic pole structure of the propagator evaluated analytically. Thecubic formalism also serves as the solution to the pole structure of the 3-site configurationin this work.This thesis is organised as follows: the rest of the introduction chapter will introducethe biological system of interest, the FMO complex, and review some of the ways ithas been modelled so far. I discuss why they are lacking in their application to pho-tosynthesis and motivate a more advanced approach. Finally I motivate the study ofa 3-site system to model the photosynthesis mechanism and discuss the experimentalsignatures for a dominant 3-site system operating within the larger seven-site system.In Chapter 2 I introduce the mathematical models to be used in later chapters, and thepath integral techniques used to analyse their dynamics. In Chapter 3 I analyse first the‘bare’ 3-site model i.e. just 3-sites and their respective tunnelling terms, isolated fromany environment; and I then incorporate the environment perturbatively. In Chapter4, I treat the system-environment interaction non-perturbatively, producing my generalresults for the 3-site-boson model. In Chapter 5 I look at the inclusion of non-diagonalsystem-bath couplings to a 2-site model. At the end of each results chapter, I apply myresults to the case of photosynthesis. This involves identifying the parameters relevantto the respective model and calculating the relevant relaxation times, decay rates andcoherent-incoherent phase spaces of the system. In Chapter 6 I summarise my findings,draw some conclusions and suggest avenues of future work. Appended at the end of thisthesis is a section on nomenclature. Due to the interdisciplinary nature of this thesisand the many esoteric concepts introduced and referenced later on, an alphabetical listwhere each term is summarised is included to assist the reader.1.2 Quantum phenomena in photosynthesisOne of the key motivations of this thesis is to understand the emergence of quantumeffects in photosynthetic biological systems. While there have been a number of photo-synthetic organisms shown to exhibit quantum effects at the molecular scale, I focus mystudy on one such organism, namely that of green-sulfur bacteria Chlorobaculum tepidum.Green-sulfur bacteria has been extensively studied using state-of-the-art spectroscopymeasurements and much is known about its internal energetic structure therefore. It is7also one of the key photosynthetic organisms known to exhibit very high energy-transferefficiency as well as long-lived quantum coherence effects in certain parts of its excitontransfer pathway.1.2.1 Light-harvesting moleculesEvery photosynthetic organism contains light-harvesting complexes (LHC) comprised ofchromophores: light-absorbing molecules attached to a protein structure holding them inplace. These molecules typically consist of a small collection of atoms with separationsbetween molecular-orbital energy levels that fall within the visible spectrum [46, 62, 11,12]. There are various photosynthetic organisms with different structural arrangements;however, the general process starts with the absorption of light by a chromophore inan antenna-like structure consisting of many individual chromophores. In Figure 1.1 oneimportant chromophore—the bacteriochlorophyll A (BChla) molecule—and its moleculararrangement are depicted.Figure 1.1: Single BChla chlorophyll molecule.Adapted and reprinted with permission from [63]. c©2018 DOAJBChla is the principle chlorophyll-type pigment in the majority of photosynthetic bac-teria [12]. Spanning around 10A˚ on its side, it contains a magnesium atom in the centerand is surrounded by four nitrogen atoms. An extensive delocalised pi-electron systemextends over the molecule. The chlorophyll contains two major absorption bands, one inthe near ultra-violet region and one in the near infra-red region. These absorption bandscome from the molecular transitions within the molecule. Each transition arises due to8different electron density distributions across the molecule, resulting in different dipolemoments. In BChla these charge density distributions are primarily governed by theacetyl group at the C-3 position in the molecule, and the single bond in ring B betweenC-7 and C-8 [12]. This reduces the degree of conjugation as well as the symmetry of themolecule.The spectra are understood theoretically using a 4-molecular-orbital model. Theseare made up from the two highest occupied molecular orbitals (HOMOs) and the twolowest unoccupied molecular orbital (LUMOs). While this picture is an approximationof a complex relationship between electronic states and orbital energies, it is generallyunderstood that two dominant transitions arise from these states [12]. Associated witheach transition is a transition dipole moment with different strengths. Spectroscopic dataof the BChla molecule resolves two fairly dominant peaks for intramolecular transitionscorresponding to the HOMO and LUMOs of the molecule and hence form the basis forthe 2-state model considered for each chromophore [12, 64].The photons absorbed in the antennae chromophores create a molecular excitationknown as an exciton [65]. The exciton is a bound state created from the transition fromone molecular orbital to another; an electron found in the LUMO and an electron holefound in the HOMO, bound together via the Coloumb interaction [66, 67]. After itscreation by photon absorption in the antenna complex, the exciton travels through anumber of structures, eventually reaching its target destination where it creates a chargeseparation used for the biochemistry of photosynthesis [11]. This final destination is usu-ally known as the Reaction Center (RC) and the biological structures preceding it serveto ‘funnel’ the captured photon energy to the RC in as efficient a manner as possible.Molecular excitons typically have characteristic lifetimes on the order of nanoseconds,after which relaxation to the ground state of the molecule occurs with photon or phononemission. Associated with each molecular state is a specific electronic density and thusdifferent charge distributions across the molecule. This gives rise to a transition dipolemoment that interacts with either the photon of light incident on a LHC or with the dipolemoment of a neighbouring chromophore. The chromophore can therefore be modelledas a two-level system and for example, two interacting chromophores would comprise a4-D Hilbert space with states |LUMO〉1 , |HOMO〉1 , |LUMO〉2 , |HOMO〉2 correspondingto the ground (HOMO) and excited (LUMO) states of chromophores 1 and 2. Typi-cal relaxation times of chromophore excited states (radiative lifetime) are of the orderτrad ∼ 10ns whereas inter-chromophore excitation transfer τ ∼ 10ps happens on the orderof picoseconds for BChla molecules in the light harvesting complex LH-2 [68]. Therefore9the ratio τrad/τ 1 tells us that we can exclude the ground states and consider only theinteractions between excited states, effectively truncating the Hilbert space to includeonly the LUMOs of each chromophore.It turns out that almost every exciton created in the antennae complex by photonabsorption eventually makes it to the RC [12]. This is the ’quantum efficiency’ usuallyquoted in the literature pertaining to photosynthesis; the ratio of excitons that reach theRC over the number created by photons. This is simply reflected in the fact that theradiative lifetime (or recombination time) of the exciton τrad is very large compared tothe exciton transfer time. In this sense, we can see that the quantum efficiency of photo-synthesis is facilitated by the processes that lead to making the ratio τrad/τ so large. Wewill see in this thesis that the presence of non-diagonal couplings is one of the drivingmechanisms in reducing the exciton transfer time τ . Perhaps a more meaningful measureof the efficiency of photosynthesis is the percentage of energy absorbed by a photosyn-thetic organism that is eventually stored in the RC. This number is actually around 27%[12], which is comparatively low compared to the ’near unity quantum efficiencies’ quotedin the literature. Therefore the results of this thesis make progress towards explainingthe quantum efficiency of photosynthesis more so than the energy efficiency.1.2.2 The Fenna-Matthews-Olson complexOne such example of an intermediate sub-system that exhibits remarkably efficient ex-citon transfer in photosynthesis is the Fenna-Matthew-Olson (FMO) complex, namedafter the researchers who discovered and determined its structure [69]. Forming part ofthe full light-harvesting infrastructure of Chlorobaculum tepidum, it essentially serves asa ‘molecular wire’, transferring an exciton created by photon absorption in the antennacomplex to the RC (see Figure 1.2). The FMO complex has long been known to ex-hibit near unity exciton transfer efficiency1 [12] to the RC; however our understandingof exactly how this high efficiency is achieved is lacking.1The ‘near unity’ quantum efficiency usually quoted in the literature refers to the ratio of photonsabsorbed to excitons that reach the RC. The amount of energy that reaches the RC compared to thatabsorbed is actually much lower, reflecting the dissipative effects that happen along the way. This canbe anywhere from 30%-50% in the literature [12], however this is still remarkably high and unexplainedby the exciton transfer modeling techniques to date.10Chlorosome antennaBChl IncidentphotonBChl BChl BChl BChl BChl BChl DecayFMO complexReaction Centre(RC)InitialexcitationFigure 1.2: Composition of light harvesting systemswithin Chlorobaculum tepidumThe FMO complex consists of a trimer, formed of three identical monomers, each con-taining seven BChla molecules as depicted in Figure 1.3. An eighth site has recently beenidentified with energy ∼ 500cm−1, which sits close to the antenna complex and possiblyacts as the entry point to the system [70].Figure 1.3: Structural depiction of the Fenna-Matthew-Olson complex from Chlorobaculum tepidum. (a) Fullview of FMO trimer including 3 monomer subunits. (b)Single view of monomer sub unit plus surrounding pro-teins (grey lines). (c) Single monomer of 8 BChla pig-ments excluding surrounding protein scaffold. Reprintedwith permission from [71].As the couplings between BChla molecules in neighbouring monomers are negligiblecompared to couplings between chromophores in the same monomer, exciton transfer isonly modelled through an individual monomer [72, 46]. The trimer therefore is consideredto contain 3 independent paths of energy transfer to the RC complex.We see in Figure 1.3 how the eight chromophores in the FMO complex are situated11inside a protein scaffold indicated by the grey lines surrounding the chromophores. Thisprotein structure surrounding the chromophores serves to hold the chromophores at theright distances and orientations for efficient energy transfer [46]. In addition to thislocal environment of proteins is a surrounding aqueous solvent that permeates the wholestructure [73].It is the FMO complex that was investigated in the landmark 2007 paper by Engel andFleming [5], which reported on experimental observations of long-lived coherence effectsin the exciton transfer pathways. Here, using 2D femtosecond nonlinear spectroscopy,they observed quantum coherence effects in the energy transfer process to persist for∼ 660fs at T = 77K [5] and 300fs at T = 300K [7, 8].Figure 1.4: Quantum beating signatures for a 77 KFMO 2D spectrum. Axes ωτ and ωt represent the Fouriertransformed excitation and detection pulses respectively,of a pump and probe laser directed on the FMO sample.The third, vertical axis, ωT , represents the Fourier trans-form of the signals evolving over time. The presence ofcross-peak frequencies, evolving over time indicates exci-ton delocalisation across pigments, and thus coherence.Reprinted with permission from [74].In Figure 1.4 we see an example of 2D spectroscopic FMO data that is interpreted asevidence of coherence in the system. The ωτ and ωt axes represent the Fourier trans-formed excitation and detection pulses respectively, of a pump and probe laser directedon the FMO sample. Frequencies where ωτ = ωt are associated with chromophore on-siteenergies and the cross-peaks ωτ 6= ωt represent tunnelling energies in the system. The12plane ωτ , ωt therefore maps out the Hamiltonian of the FMO complex. The third axis,ωT , represents the Fourier transform of the signals evolving over time. Time evolution ofthe cross-peaks is interpreted as a ’quantum beating’ signal and thus evidence of excitondelocalisation across chromophores.1.3 Microscopic origin of chromophore-environmentinteractionsThe BChla molecules of the FMO complex exist not only in an aqueous solution butare also embedded in a protein scaffold. The effect of this is to shift the excitation en-ergies of the BChla molecules as well as modulate the distance between the molecules.For the solution surrounding the molecule, it is the solvent’s dialectric that changes itsenergy, and for the protein scaffold it is predominantly charged amino acids that interactwith the molecule [46]. The interaction of the chromophores with the various parts oftheir local environment leads to fluctuations of not only the excitation energies, but thecouplings between chromophores too. On the one hand, environmental changes inducefluctuating local electric fields that induce a Stark shift in the excitation energies andthus the transition dipoles of the chromophores to fluctuate. This results in the diago-nal system-environment coupling. On the other hand, changes in the environment willaffect the relative spacing of the chromophores, introducing a fluctuating component tothe inter-chromophore couplings. This results in the non-diagonal system-environmentcoupling.While most studies to date have focused on the diagonal couplings, a number of recentstudies have indicated the importance of including both types of coupling in photosynthe-sis [75, 76]. While it should be noted that some experiments have found the non-diagonalcouplings to be small compared to the diagonal ones [77], it has nonetheless been demon-strated that even small non-diagonal couplings can have a considerable impact on theexciton dynamics [78, 53]. Furthermore, in the study of polaron dynamics, Marchand etal found similarly strong effects arising from non-diagonal couplings to phonons [79, 50].It is worth remarking at this point that determining the exact size of the system-environment coupling parameters from experiments of the FMO complex is difficult.Analysing spectroscopic data requires complex fitting algorithms based upon models ofsystem-environment interactions [46, 14]. The size of the parameters obtained are highlydependent on the models employed. Furthermore, non-diagonal couplings are usually13excluded from these models meaning that the magnitudes of these couplings are undeter-mined. This also highlights an issue with the sizes of the diagonal coupling parametersdetermined to date. If non-diagonal couplings should be present in the models, and thetwo types of couplings distinguished, then the models that assume only diagonal cou-plings have presumably determined the size of a composite diagonal and non-diagonalcoupling parameter instead.1.3.1 Diagonal and non-diagonal couplings in chromophoresTwo interacting chromophores will have associated with them: excitation energies per-taining to each molecule and tunnelling energies between the molecules. The effect of theenvironment induces fluctuations in the chromophore’s energies, which leads to a dynam-ical modulation of these parameters. Therefore we can think about two distinct types ofsystem-environment interactions. The first type of interaction is one that modulates theexcitation energies of the chromophores. This is known as a ‘diagonal’ interaction as theenvironment is coupling to the diagonal elements of the central system’s Hamiltonian.The second type of interaction is the ‘non-diagonal’ interaction, where the environmentmodulates the tunnelling energies between chromophores. In this section I discuss themicroscopic origin of these two types of interaction in the context of biomolecules.As I discussed in Section 1.2.1, chromophores contain two dominant molecular states:the HOMO and LUMO, which form the basis for their ground and excited states re-spectively. We recall that the radiative lifetime of the excited (LUMO) state of eachchromophore relaxing to their ground (HOMO) state is orders of magnitude smaller thanthe interchromophore transfer time. Therefore we only consider the interaction of twochromophores’ LUMO states. Associated with each chromophores’ LUMO states is adipole moment. If the chromophore’s surrounding environment, comprised of a solventand protein scaffolding, is sufficiently far away, we can treat the chromophore withinthe dipole approximation. This leads to a central assumption in the Fo¨rster picture ofenergy transfer between chromophores, which we discuss below. Charges present in thechromophore’s environment, coming from the surrounding solvent or amino acid residuesin proteins, produce an electric field. This leads to a fluctuating energy term for theexciton of the form∆Eα =1∑jqj~µi · ~xijR3ij(1.1)14where ~µi is the dipole moment of an exciton on the ith chromophore, is the optical di-aletric constant, Rij the distance from the center of the chromophore to the environmen-tal charge, and xij is the fluctuation about the equilibrium position of the chromophore.Since the fluctuating energy shift depends on the relative distance between chromophoresand amino acid residues, we see how phonon vibrations present in the surrounding proteinstructure couple to the system as they modulate this distance. In a second quantisedlanguage the site displacement operator for the phonon xq =√~/2Mωq(bq + b†q) (seeMahan [52]) gives∆Eα =∑qλα,q(bq + b†q) (1.2)where α here labels a chromophore and q a phonon mode of the bath. The chromophore-bath coupling λα,q, is a function of the charges present in the surrounding medium andis therefore thought of as a backreaction resulting from the polarisation of the environ-ment by the chromophore (see section 2.5.1.2 in [64]). This is the origin of the diagonalcouplings.We now turn to the microscopic origin of the non-diagonal chromophore-environmentcouplings. Provided the interchromophore spacing is large compared to the width of achromophore i.e. more than several Angstroms [12], then the exciton transfer betweenchromophores can be considered a long range dipole-dipole interaction. This leads to theideal-dipole approximation of molecular interactions often employed in the literature,where non-diagonal couplings can be excluded. If however the interchromophore spacingis such that there exists appreciable wavefunction overlap between molecular orbitalsthen this approximation breaks down.Calculations performed using the ideal-dipole approximation of chromophore-chromophoreinteraction, compared to experiment, demonstrate the failure of the ideal-dipole approx-imation for BChla molecules [80]. The sizes of the chromophores in this case are ∼ 9A˚,which is comparable to the inter-chromophore distance determined to be ∼ 15A˚, whichputs the system out of the range of applicability of the ideal dipole approximation.Wavefunction overlap of neighbouring chromophores contains an exponential depen-dence on separation distance. Vibrations in the local environment can cause the twochromophores to move closer or further away from each other, thereby modulating theseparation distance. For an environment with quantised vibrations, this introduces aprocess by which an exciton can hop between chromophores by way of an inelastic scat-tering process [52]. Phonons are not conserved in this case and an exciton can tunnel15between chromophores by absorption or emission of phonons. Since the tunnelling termdepends exponentially on separation now due to the molecular wavefunction overlap ofthe two interacting chromophores, it follows that the tunnelling frequency between donord and acceptor a molecules ∆ad is modulated by∆˜ad = ∆ad exp[−∑qFad,q~ωqxq](1.3)where Fad,q is the force between the two chromophores due to the displacement xq, and qlabels the phonon mode. In a second quantised language the site displacement operatorfor the phonon position xq =√~/2Mωq(bq + b†q) gives us∆˜ad = ∆adexp[−∑qVq~ωq(bq + b†q)](1.4)with Vq representing the product of the force applied by the phonon mode q and thedistance between the two chromophores, with ωq the frequency of phonon mode q, andbq, b†q the destruction/creation operators for the phonons. We can expand this to linearorder in Vq/~ωq, provided the phonons modulate the distance between chromophores onlyslightly. In this case we have∆˜ad = ∆ad −∑qλad,q(bq + b†q) (1.5)where the coupling strength λad,q is now a function of the tunneling matrix elementbetween the donor and acceptor molecules and the phonon frequency. So we see that thephonons modulate the distance between chromophores, thereby modulating the tunnelingamplitude. This is the origin of the non-diagonal coupling between chromophores.1.4 Exciton dynamics in photosynthetic systems: tra-ditional theories and their limitationsI now discuss the conventional theories of excitation energy transfer (EET) that use themolecular parameters described above to obtain rate constants determining dynamicswithin FMO and associated relaxation effects. There is of course no general theory (yet)that applies across the whole parameter space, but the choice of theory depends on therelative strength of pigment-pigment and pigment-environment couplings.16Traditionally, EET in FMO has been modelled in two limits. One limit, knownas Fo¨rster theory, considers small chromophore-chromophore couplings ∆ compared tochromophore-environment couplings Λ (I define the coupling here to be Λ in anticipationof its relationship to the experimentally determined reorganisation energy of the bath).In this case the exciton is localised on each chromophore and ∆ can be treated perturba-tively. In the opposite limit, the exciton-bath coupling Λ is treated perturbatively, leadingto Redfield theory. Associated with each of these couplings are therefore two separatetimescales. One is the inverse of the bath coupling 1/Λ, which describes the timescaleover which energy is transferred to the environment. The other is the inverse of theexciton-exciton coupling 1/∆, which describes the timescale of exciton transfer betweenchromophores. Therefore in the Forster limit 1/Λ 1/∆ we have a rapid energy lossto the environment occurring over a timescale much shorter than the inter-chromophoretransfer. Conversely in the opposite limit, exciton transfer occurs at a much faster ratethan dissipation to the environment.In the case of FMO the tunnelling energies are of order ∆ ∼ 1meV corresponding toa transfer time τ ∼ 4ps, while reorganisation energies are found to be around Λ ∼ 4meVcorresponding to a dissipation time τΛ ∼ 1ps[62]. The reorganisation energy is an ex-perimentally determined value related to the system-bath coupling strength (see Section2.3) So we see that in FMO the size of these parameters are actually similar and onemust therefore go beyond the limiting cases and utilise a non-perturbative theory. Itwas shown by Sharp et al. [81] as well as Ishizaki et al. [24] that a calculation of the2D optical spectrum assuming weak coupling to the environment fails to reproduce theoscillations seen by Engel et al. [5]. Therefore the application of these models to date hasbeen arguably inappropriate to the FMO complex, and a non-peturbative analysis suchas the one I have employed here is warranted. Nonetheless I briefly summarise some ofthese traditional methods before discussing the results of my non-perturbative methodsand their application to FMO.1.4.1 Resonance Energy TransferFirst I will describe the process of exciton couplings across different chromophores. This isknown as resonance energy transfer (RET) and describes the transfer of energy from onemolecule (the Donor) to another (the Acceptor) in biological chromophores. A molecule(or chromophore) is initially excited by a laser pulse in a lab setting or a monochromaticlight source in a natural setting, creating an exciton on the first molecule. Splitting the17chromophore coupling into its long and short range components, i.e. dipole-dipole Vc andexchange terms Vex respectively, readsV = Vc + Vex. (1.6)It is assumed in RET that the exciton is transferred between molecules via a dipole-dipoleinteraction, therefore V ≈ Vc whereVc =14pi0κµDµAR3(1.7)where κ = µD · µA − 3(µD ·R)(µA ·R) is an orientation factor that reflects the fact thatno interaction is observed between perpendicularly oriented chromophores. ~R is the sep-aration vector between the molecule’s centres and ~µd and ~µa are the donor and acceptortransition dipole moments respictively. The transition probability is proportional to |∆|2and therefore goes as 1/R6. The dipole-dipole interaction expression is obtained fromthe expansion of the Coulomb interaction in the acceptor-donor distance parameter toget a multipole series, while retaining only the dominant term.1.4.2 Forster theoryForster theory presented a key step in the understanding of EET at the molecular scalewhen it was introduced in the last century [82]. Its applications have had success notonly in understanding light-harvesting mechanisms in biochemistry but also in the de-velopment of artificial organic-based light-emitting diodes [83]. Its application to theexperimental sphere facilitated the development of fluorescence-resonance-energy trans-fer (FRET) [67], a technique that has been widely used to detail the structure of biologicalsystems. Central to the theory, however, is the approximation of weak inter-chromophorescoupling ∆da (between donor and accepter chromophores) compared to environmentalcouplings γ. Forster therefore derived his EET rate expression with the use of the Fermi-Golden-rule approach, treating ∆ as the perturbation in time-dependent perturbationtheory. The rate of exciton transfer between donor and acceptor states was derived byForster [82] to bekd→a = ∆2daˆdω2piRe{Ad}Re{Fa}(1.8)where Ad and Fa are the absorption and fluorescence lineshapes of the donor and acceptorchromophores. The resulting rate expression is expressed as the overlap integral between18the fluorescence spectrum of a donor and the absorption spectrum of an acceptor [84].The limitations of this approach to modeling EET rest upon some key, rather limitingassumptions. The donor and acceptor molecule must be coupled sufficiently weaklycompared to chromophore-bath coupling such that Fermi’s Golden Rule applies. Thisleads to incoherent hopping between donor-acceptor molecules such that the dissipationto the bath is one-way. No feedback to the system from the bath is permitted and thetheory excludes possible coherence effects in principle due to the assumption of verystrong system-bath coupling relative to the tunnelling energy.1.4.3 Quantum master equations in molecular systems: theRedfield equationIn the opposite limit to Forster theory, quantum master equations are used to describeexciton dynamics. The fundamental assumptions here are that the exciton-environmentcoupling is small compared to inter-molecular interactions and so the system-environmentcoupling is treated perturbatively. The most commonly used version of this approach inthe context of biomolecules is the Redfield equation [64, 85]. Master equations have alsoseen widespread application in the field quantum optics [86].Derivations of quantum master equations for specific systems usually start from theLiouville equation [87]∂tρ(t) = − i~ [H, ρ] (1.9)where ρ(t) is the time-dependent density matrix of the entire system: central system Splus bath B such that the Hamiltonian is H = HS + HB + HI where HI is the system-bath interaction. The Liouville equation is essentially the density matrix form of theSchrodinger equation and so allows the time-evolution of a statistical mixture of states.One is generally interested in the dynamics of the central system so we trace over thebath degrees of freedomρS(t) = TrB[ρ(t)]. (1.10)In deriving the Redfield equation one starts from the assumption that the initial systemand bath states factoriseρ(0) = ρS(0)⊗ ρB (1.11)19This assumption is valid in the context of exciton dynamics as the creation of the excitonin the system at t = 0 is incredibly fast compared to the inverse bath couplings. Werecall that in light-harvesting complexes, an exciton is created either by the suddenabsorption of light, or the injection of an exciton from the antenna complex. Anotherkey assumption in Redfield theory takes the system-environment coupling to be weak andsecond-order perturbation theory is used accordingly. In addition to these assumptionsis the use of the Markov approximation, which makes the master equation local in time.This latter simplification, it should be noted, removes any possibility for feedback in tothe central system from the bath i.e. any energy transfer to the bath is one way andcannot re-enter the central system at a later time. The relevant physical condition forthe Markov approximation is when the memory timescale τB of the bath is very smallcompared to the timescale t of the dynamics of the central system. This means that thebath correlation function decays sufficiently during an appreciable timescale over whichthe central system evolves, removing the possibility for bath memory effects [88]. Theresulting Redfield equation is [24]∂tραβ(t) = −iωαβρ(t) +∑α′,β′Rαβ,α′β′ρα′β′(t) (1.12)where ~ωαβ = α − β describes the energy gap between chromophores and the indicesα, β run over chromophore sites. The second term describes the relaxation dynamicswhere the Redfield tensor is given byRαβ,α′β′ = Γβ′β,αα′ + Γ∗α′α,ββ′ − δββ′∑µΓαµ,µα′ − δαα′∑µΓ∗βµ,µβ′ (1.13)The damping terms are given byΓαβ,α′β′ =1~2∑qq′λαβ,qλα′β′,q′Cqq′(ωi′j′) (1.14)where λαβ,q = 〈α|λq |β〉 are the chromophore-bath couplings. The Fourier transformedbath correlation function is defined asCqq′(ω) =ˆ ∞0dteiωtCqq′(t) (1.15)where the time-domain correlation function is20Cqq′(t) =~piˆ ∞∞dωJqq′(ω)[n(ω) + 1]e−iωt (1.16)and n(ω) = 1/(eβ~ω−1) is the Bose-Einstein distribution function and Jqq′(ω) the spectraldensity of the bath (see Section 2.3). An additional approximation applied to equation1.12, which casts it in the quantum master equation form commonly used in the liter-ature [88, 64], is known as the secular approximation. This amounts to neglecting therapidly oscillating terms that contain ei(ωi−ωi′ ), which effectively amounts to a RotatingWave approximation. This is generally justified under the Riemann-Lebesgue-Lemmathat states that the integral of a rapidly oscillating function goes to zero in the limitwhen in this case ωi − ωi′ →∞.Along with the fundamental assumption that the system-bath interaction is perturba-tive to second order, we see that the Redfield equation describes a fairly restricted classof systems. Despite the popularity of quantum master equation approaches to modellingphotosynthesis [15, 23, 16, 17], a number of studies have shown explicitly the failure ofthis approach [24, 84, 89] in the context of photosynthesis.There have been some attempts to go beyond the traditional theories described how-ever, often these changes have involved small deviations from the master equations models[21, 84]. In this case non-Markovian effects are reintroduced, however the system-bathcoupling is still treated perturbatively. Some analytical studies on photosynthetic sys-tems have utilised non-perturbative methods. The spin-boson model, with just diagonalsystem-bath couplings, was applied to the FMO complex and found persistent oscillationsin the exciton dynamics in accordance with the experimentally determined values [27].However, this model assumes a relatively large tunnelling energy between chromophoreswhich can account for the predicted rapid exciton transfer rates. Furthermore, the modelassumes just an Ohmic form for the bath spectral density in order to utilise the resultsfrom the spin-boson literature. Despite these limitations the initial success of utilisingnon-perturbative methods in studying the FMO system is promising.In the next chapter, I set up the path-integral approach to analysing exciton dynamics.In subsequent chapters, this will be used to model the FMO system in photosynthesiswith the intention of accurately predicting relaxation rates compared to experiment aswell as providing a better understanding of the underlying physics.211.5 3-site configurations in the Fenna-Matthew-OlsoncomplexI now turn to the motivation behind studying a 3-site model in the context of photosyn-thesis in more detail. I begin with some preliminary, qualitative considerations: lookingat the energy structure of the FMO system and how the branched nature of the systemleads to an intuitive 3-site reduced model. I also look at some spectroscopic data, whichhas hinted at the presence of a dominant 3-site structure in the system. While therehave been a number of studies investigating a 3-site-ladder model for the FMO complex[44, 90], where the site energies are therefore different, here we focus on the 3-site-Vconfigurations [45].To see how a 3-site-V configuration can be effective at modelling the dynamics ofthe full FMO system, we first look at the on-site energies and transition matrix ele-ments between sites (dipole-dipole couplings between chromophores), established fromspectroscopic data [46].HFMO =26.7 −12.9 0.6 −0.5 4.7 −0.58 −1.0−13.0 27.3 4.0 0.9 0.7 1.0 0.10.6 4.0 0 −5.8 0.12 −1.0 0.63−0.5 0.9 −5.8 15.5 −8.8 −1.8 −7.60.58 0.7 0.12 −8.8 55.8 89.7 −2.5−15.1 8.3 −8.1 −14.7 11.1 41 4.1−1.0 0.1 0.63 −7.6 −0.31 4.1 34.7(meV) (1.17)In Figure 1.5 we see a pictorial representation of the FMO Hamiltonian equation 1.17,including only dominant couplings greater than 3.7meV.22Figure 1.5: Pictorial representation of Hamiltonian forFMO complex, taken from spectroscopic studies [46], inunits of meV. Coloured lines indicate excitation ener-gies of different chromophores (diagonal Hamiltonian en-tries). Arrows indicate inter-chromophore couplings (off-diagonal Hamiltonian entries). Numbers included nextto lines and arrows indicate the associated energy withthat excitation energy of coupling respectively. Onlydominant energy transfer pathways included; couplings> 4meVHaving only included dominant energy transfer pathways above 4meV we hope to demon-strate the branched nature of the system. It’s clear that there are two dominant branchesin the system with relatively small tunnelling energies between them. One connects states|1〉 ↔ |2〉 ↔ |3〉 and the other connects states |3〉 ↔ |4〉 ↔ |5〉 ↔ |6〉 ↔ |7〉. If the statesthat comprise the excited state portion of each branch can be described by a single ‘ef-fective’ state, with some correponding ‘effective’ tunneling energy with the ground stateat |3〉, then the remaining system would be of the 3-site-V form depicted in Figure 1.6.23Figure 1.6: Effective 3-site-V configuration for FMOcomplex.Schmidt and Renger, in fact, demonstrated that there is no qualitative difference inthe site-population dyamics between the full FMO Hamiltonian and an ‘effective’ 3-sitemodel [90].After some qualitative considerations of the possible 3-site nature of the FMO system,I look now at some spectroscopic data suggestive of a dominant 3-site structure in FMO.The energy spacing of excitonic levels is relatively small compared to the inhomogenousbroadening present in the optical transitions, so much so that the resulting linear ab-sorption spectrum is almost entirely featureless at physiological temperatures [91]. Atcryogenic temperatures of ∼ 77K however, 3 distinct peaks emerge in the spectra ofFMO. Simulations performed in the literature [46, 91, 71] suggest exciton delocalisationacross combinations of 2-3 chromophores in FMO. If coherence effects are to play a role inefficient energy transfer in the FMO complex, then it is across at most three chromophoresthat this effect will be present. Beyond that the system presumably evolves accordingto the Forster limit corresponding to incoherent hopping between chromophores. Thevarious 3-chromophore delocalisation configurations are shown in Figure 1.7. Here wesee the results of experimental spectroscopic data on FMO showing the emergence ofthree distinct peaks at cryogenic temperatures. We also see the extent of wavefunctiondelocalisation across the various chromophores.24Figure 1.7: (A) Arrangement of BChla pigments withinthe FMO units. BChla site numbering according to Fennais in black Roman numerals. Schematic representation ofspatial extent of the excitons according to Adolphs et al.[46] is shown by shaded areas. Exciton numbering (rednumbers) is given in order of increasing energy. (B) Lin-ear absorption spectrum of FMO at 77 K with excitonictransitions (vide infra) represented by vertical bars. (C)Normalized absorptive 2D spectra at increasing popula-tion delays with dashed lines indicating excitonic transi-tion energies. All spectra were recorded in 1:2 aqueousbuffer:glycerol mixture at 77 K. Reprinted with permis-sion from [91].The central 3-site-V system parameters to be used in our analysis we take from [46],where the electronic energies and couplings between them were determined from a time-dependent density functional theory (DFT) on the chromophores with an electrostaticPoisson-Boltzmann type calculation. Similar studies have also used a DFT approach todetermine the chlorophyll energies immersed in organic solvents to great effect, comparingtheir results to experimentally observed spectra [92]. The calculated diagonal site energiesand tunneling terms in [46] are shown in Table 1.1.25BChl site On-site energy (meV) BChl-BChl transitionTunneling(meV)Ground state: |0〉 0 |0〉 → |1〉 1.8Site 1: |1〉 20 |0〉 → |2〉 3.2Site 2: |2〉 20Table 1.1: Tuned V-system parameters. From Adolphset al [46]We see in Table 1.1 how the on-site energies are large compared to the tunnelling energies,putting us within the strong-localisation region. Just based upon these numbers onewould expect the wavefunctions of each chromophore to be fairly-well localised aroundtheir respective BChla molecules, with possibly some delocalisation spreading out toneighbouring chromophores. This is somewhat in line with the data presented in Figure1.7, however there is clearly still some delocalisation across chromophores.26Chapter 2A path-integral approach to excitontransferIn this chapter, I set up the 3-site model in a path-integral formalism. While I willprimarily focus on the 3-site application here, the resulting formalism will be applicableto the 2-site model used in Chapter 5 with the appropriate modifications. Indeed, the2-site path integral formalism has been described extensively in the literature and so inthis chapter I will focus on the extensions of the 2-site formalism to the 3-site case.While my final results will be non-perturbative in the system-bath coupling, the auxil-iary approximations made along the way will build upon a number of techniques employedin the famous spin-boson model. The spin-boson model [31], considered a spin-1/2 systemcoupled to an oscillator bath, developed in the path integral formalism along with thenoninteracting-blip approximation (NIBA). While this approach is fully non-perturbativein the system-bath coupling, and permits memory effects in the system, the underlyingapproximation restricts the longevity of coherences in the central system. My 3-site-boson model will involve an extension to this model, by way of the addition of not only athird site to the system but an extra tunnelling matrix element. Although this extendedmathematical model will also share many of the varied applications that the spin-bosonmodel has, it is the recent surge of theoretical physics research in the area of quantumtransport phenomena in photosynthesis that primarily motivates the development of thismodel.There have been a few preliminary studies, looking at the role that 3-level systemsmight potentially play in photosynthesis [45, 44]. However they have so far relied uponnumerical techniques that assume one-way dissipation to the bath, which are insufficient27at modelling coherence effects. Therefore in order to model the photosynthesis mechanismmore accurately, such that the possibility of relatively long-lived coherences can be estab-lished, one must adopt a more advanced approach to modelling the system-environmentinteraction—one that retains the desired coherence effects. A number of numerical ap-proaches are suited to this; however, in this thesis we are interested in the possibility ofobtaining analytic results and ultimately applying them to the case of photosynthesis.2.1 Truncation procedure for the 3-site-boson modelfrom an extended systemPrior to performing any analysis I must set up the problem, and like most physical systemsof interest, I am considering an approximation of some more general system. In my casethis amounts to the truncation of a more general potential energy landscape, with threedistinct potential wells that permit tunnelling between them, down to three distinct,interacting ground states of each potential well. The problem that I wish to investigatetherefore involves a quantum-mechanical system that occupies a 3-dimensional Hilbertspace.In the more general case, the system has a continuous degree of freedom q with acorresponding potential energy function V (q), that forms three separate potential wells.Each well supports multiple quantum mechanical states, but if we assume that the barrierheights between wells are large compared to the energy separation between the groundand first excited state of each well, then the states should localise within each well.Moreover if I assume that the temperature of the system kBT , as well as bias energiesseparating the ground states of the potential minima, are small compared to the excitedstate energies, then states in each well should still be restricted to just their groundstates. In other words, the bias energies are small enough such that the system does notlocalise in any one well, and the temperature kBT is not too large so as to thermallyexcite any higher energy states. In effect, the system is now described by a 3D Hilbertspace spanned by these three ground states. In this regime, the dynamics of the systemis governed by tunnelling between the wells associated with a tunnelling matrix element~∆, which must also be small compared to the energy level spacings within wells so asnot to mix the ground and excited states.Isolated systems like this are known to exhibit quantum mechanical interference effectsbetween wells. Realistically though any physical realisation of such a system is not28isolated, and exists in the presence of some environment. So far I have only considereda 3-level system isolated from any environment. I refer to this as the ‘central system’ forwhich the truncated Hamiltonian isH0 = 1 |1〉〈1|+ 2 |2〉〈2|+ 0 |0〉〈0|+ ∆10(|1〉〈0|+ h.c.) + ∆20(|2〉〈0|+ h.c.) (2.1)where the lower site energy 0 = 0.The complete 3-state central system would of course include a tunnelling term betweenthe upper two levels as well. We exclude this term here due to our interest in this partic-ular configuration of the 3-site system, often referred to as a V-system in the quantumoptics literature where it sees widespread application [86]. Our interest in this partic-ular 3-site-V-system configuration is motivated by its application to the photosynthesismechanism discussed in Chapter 1.2.2 Modelling the environment as a macroscopic har-monic oscillator bath with linear system-bathcouplingI now discuss the interaction of the central system with an external environment. Here Iam considering a central system that contains only a few degrees of freedom with respectto the environment, which is assumed to be large with many modes. As the interactionbetween the central system and each environmental mode is inversely proportional to thesize of the bath, for a macroscopic environment we can assume weak coupling to eachindividual mode. However, the total coupling to the bath can be strong since the fullcoupling is due to a summation over the couplings to each bath mode. Therefore we canassume that the coupling is linear in bath coordinates, since each individual coupling is29small, but the total coupling can be arbitrarily large [61].Any linear system composed of a distribution of interacting oscillators can be rep-resented by an equivalent set of independent oscillators [37]. In other words, a lineartransformation of the coordinates of the interacting oscillators, recasts the system interms of its eigenmodes. Since we’re assuming that the bath is comprised of a set ofsimple harmonic oscillators, perturbed only slightly from their equilibrium position, thenit is permissible to consider the equivalent set of independent oscillators. Since the bathcoordinates will ultimately be summed over, and therefore appear as dummy variablesas part of the path-integral measure, it is inconsequential whether the coordinates repre-sent positions or otherwise. The linear transformation coefficients will of course introducesome new normalisation term due to the path-integration measure, however, this is takenin to account in the path integral formalism [37]. I therefore write the Hamiltonian asH = HS +HB +HI (2.2)where HS is the Hamiltonian of the central system, the 3-site system, HB is the Hamil-tonian describing the environment and HI the coupling between the two. In the statespace the central system is a function of the state label α and for a bath of N HarmonicoscillatorsHB =N∑q=1(p2q2mq+12mqω2qx2q)(2.3)andHI = −∑αβN∑q=1Fαβ,q |α〉〈β|xq(t) (2.4)where the interaction is linear in bath coordinate xq(t) and we have excluded a counter-term in the interaction that offsets the renormalisation of the potential in HS. Fαβ,q isthe state-dependent force between state αβ and an oscillator mode q. Therefore we haveassociated with this force some state-dependent fluctuating energy Eαβ(t). We have keptHI general for now, permitting both diagonal α = β and non-diagonal α 6= β couplings.Second quantising the bath and interaction terms givesHB +HI =N∑q=1ωqb†qbq +∑αβ|α〉〈β|N∑q=1λαβ,q(bq + b†q) (2.5)30where Fαβ,q =√2mqωqλαβ,q. For the case of just diagonal couplings α = β in the 3-sitemodel the state-dependent couplings would be λα,q ∈ {λ1,q, λ2,q, λ0,q}. Physically, we havedescribed here a bath that is sensitive to the degrees of freedom of the central system.In the polaron literature for example, this usually describes phonons in the environmentthat are sensitive to the spatial coordinate of the central system particle of interest; i.e.as the particle moves through the environment of phonons, the environment and theparticle experience a force between them as a function of their spatial separation [52].Here, since we’re considering a central system in state space, the environment couplesindividually to each state in the system and there exists a force between each state andeach oscillator mode corresponding to a coupling constant λαβ,q.2.3 The oscillator bath spectral densityI have so far only stipulated that the environment is comprised of a large set of har-monic oscillators and restricted the system-bath coupling to a linear form. However theenvironment, being the large, complicated structure that it often is, requires further spec-ification as to the distribution of its modes across frequency space. For the case wherea thermal equilibrium statistical average is taken over the initial and final states of theenvironment, then the system-environment coupling can be completely characterised bythe spectral density [31, 61]Jα(ω) =pi2∑qλ2α,qωqδ(ω − ωq) (2.6)which contains the density of states of the bath weighted by the state-dependent system-bath coupling λα,q. So we have a spectral density corresponding to each state in thecentral system with the bath mode part of the coupling summed over. When consideringboth diagonal and non-diagonal bath couplings we will also have to differentiate betweendiagonal and non-diagonal coupling strengths. If one is considering some vibronic de-grees of freedom in the environment, such as phonons, then to linear order, the bathmodes couple to some coordinate in the central system. In the standard literature this isusually derived as a particle, moving with a geometrical coordinate interacting with itsenvironment by displacing the positions of the atoms around it. Therefore the interac-tion potential depends on the position of the particle from the atoms in the environmentaround it. To linear order in displacement one arrives at the Hamiltonian above. Howeversince we’re considering a state space for the central system, which does not necessarily31involve spatial coordinates, then we don’t describe the coupling as being dependent onthe central system coordinate, but instead dependent on the state α. In the originalspin-boson model for the TSS, the system is said to suddenly jump between positions±q0/2, where the bath is coupled to these relative positions. This arises as a couplingto q0σz in the interaction term of the Hamiltonian. In fact, the original motivator forLeggett et al. to work on the SB model was to apply it to the SQUID system wherethe central system coordinate was the flux ϕ and not at all a geometrical coordinate[93]. Nevertheless they used the coordinate q0 throughout, along with that caveat. Hereinstead I consider a coupling λα,q which contains both the state dependence and bathmode dependence.In the limit of a large number of bath oscillators, where the frequencies ωq are suffi-ciently dense so as to form a continuous spectrum, it is appropriate to replace the discretesum in Jα(ω) by a continuous integral. In the low frequency limit of the spectral densitythe spectral density takes the general formJα(ω) = Aα,sωsω1−sph e−ω/ωc (2.7)For s = 1 we have the Ohmic case, for s > 1 we are in the super-Ohmic regime. In thelatter case we see that the bath frequency is defined in relation to some characteristicbath frequency ωph that fixes the dimensionality of the spectral density. The state de-pendence is now contained within the parameter Aα,s, which will change depending onthe value of s and which I will define explicitly in the next sections.A quantitative measure of the system-bath coupling strength usually used in the lit-erature is the reorganisation energy Λ [84, 27, 94]. This describes the energy released bythe bath when relaxing to its equilibrium ground state [64].Λ =1piˆ ∞0dωJ(ω)ω(2.8)This parameter can be determined experimentally using spectroscopy [95] and is foundto be around Λ = 4meV [24, 26]. The justification for quantifying the coupling strengthwith the measurable reoganisation energy is usually understood in the literature, by firstconsidering the case of zero tunnelling matrix element [56]. Then the chromophore re-duces to the well known Independent-boson model [52], an exactly diagonalisable model,which results in site energies shifted by the distortion of the local environment. Asso-ciated with this deformation is the reorganisation energy, and the system-environmentcoupling in this limit is clearly defined. Reintroducing small, but appreciable, tunnelling32matrix elements is then assumed to leave the reoganisation energy picture of system-bathcouplings unaltered.2.3.1 Ohmic damping: s = 1For the special case of Ohmic spectral densities, which will be one form of J(ω) con-sidered in this work, the damping in the system is independent of frequency [34, 29].This situation is known to describe the environmental effects of conduction electrons insolids [40, 96] and has been argued to model dissipative energy transfer mechanisms inphotosynthesis as well [42, 24, 27] (see also Section 2.3.4). In the Langevin equation thiscorresponds to a frictional force proportional to the velocity of the path (hence the termOhmic) [61]. In the high temperature limit, the classical Langevin equation is recovered,which describes a heavy Brownian particle immersed in a fluid of light particles [61]. Inthis limit the power spectrum of the thermal fluctuations can be described by Johnson-Nyquist noise [97], and thus is of Ohmic form. that In this case the spectral density takesthe formJα(ω) = ηαωe−ω/ωc (2.9)which is valid for all frequencies much less than some cutoff frequency ωc. As the spectraldensity is a function of both the density of states of the bath modes as well as the state-dependent system-bath coupling, the parameter ηα has an index corresponding to eachchromophore. This parameter therefore has units of Joule · seconds and is often referredto as the viscosity coefficient [61]. At this point I can define a dimensionless couplingconstantγα =ηαpi~(2.10)which will be used in proceeding calculations.The measured value of the dimensionless coupling constant is measured to be γ = 0.22for an Ohmic bath in the FMO complex [98]. The corresponding cut-off frequency wasfound to be ~ωc = 21.1meV.2.3.2 Super-Ohmic damping (acoustic phonons): s = 3For the case of s = 3 the spectral density is known as super-Ohmic. The physicalbasis for a spectral density of this form can be realised in the coupling of the system to33a 3-dimensional acoustic phonon bath [61, 31]. One can arrive at this cubic frequencydependence with a system-bath coupling λq ∝ ω−1/2q for elastic waves and a Debye densityof state∑q δ(ω − ωq) ∝ ω2. Recall that, for low enough frequencies, it is assumed thatphonons still have a well defined wavevector for amorphous solids. Further summationover the one longitudinal, and two transverse branches of the dispersion, results in aspectral density to the 3rd order in ω [99]. Therefore, for the low frequency distributionof the bath, the spectral density takes the formJα(ω) = ραω3ω2phe−ω/ωc (2.11)Once again the coefficient ρα has units of Joule · seconds. I can therefore define adimensionless coupling parameter that this time takes the formζα =ραpi~(2.12)which is independent of the characteristic bath frequency ωph. Instead, we’ll see thatduring the calculations in subsequent chapters, that I instead end up with a ratio of cut-offfrequency to characteristic phonon bath frequency in our expressions. The dimensionlesscoupling constant and cut-off frequency were found to be ζ = 0.31 and ~ωc = 8.7meVrespectively [98].It is worth briefly discussing the significance of phonons and super-Ohmic spectraldensities in the context of amorphous solids. Phonons are normally discussed in thecontext of ordered lattices where there exists a well-defined translational symmetry. Inthis case we can describe the distribution of phonon modes with a dispersion relationaccording to a well-defined plane-wave excitation momentum. In disordered (amorphous)solids, despite the lack of translational symmetry, low-frequency vibrational modes similarto those of a crystalline lattice [100, 101] are present. This can be conceptually understoodby first considering an ordered lattice with a phonon dispersion. Introducing disordershould affect the high-frequency modes first as the low-frequency modes do not ’see’ thedisorder due to their long wavelength. For higher, and higher disorders, we expect thehigh frequency phonons to be scattered and only very low frequency modes to survive.The environment of a chromophore is of course not an ordered lattice but an amorphousone with some degree of disorder. Nevertheless it is clear that vibrational excitationsexist in such amorphous systems arising from the protein structure that surrounds eachchromophore [46, 102].342.3.3 Optical phonon damping with spectral broadeningFor the case of optical phonons the spectral density contains just one frequency ω0 andthe zero-temperature spectral density includes just a delta function. However for realisticapplications we would like to add a small linewidth to simulate some damping by theenvironment. In this case I choose a Gaussian lineshapeJα(ω) = λαexp[−(ω − ω0)22ξ2](2.13)with fullwidth-half-maximum ξ and coupling strength λα that determines the height ofthe peak and has units of Joules. In this case the coupling strength λα correspondsdirectly to the coupling strength in the system-bath interaction part of the Hamiltonian(see equation 2.5). The dimensionless coupling parameter in this case isνα =λαξpi~ω20(2.14)where ξ and ωo have units of frequency.Optical phonon peaks appearing in the structured spectra of the FMO complex, tendto be relatively narrow and oscillate at very high frequency relative to the FMO systemtunnelling frequencies. Linewdiths are of the order ~ξ ∼ 1 − 10meV, and oscillationfrequencies ~ω0 ∼ 50− 250meV. This range of values describes narrow peaks centered atrelatively high energies.2.3.4 Spectral density functions for light-harvesting moleculesHaving introduced various general forms for spectral density functions, I now discuss theirapplication to the case of photosynthesis. Previous analytical approaches to modellingchromophore dynamics have made various assumptions as to the nature of chromophore-environment interactions. Approaches by Gilmore [41, 42, 43], Ishizaki [21] and Cho[103] have assumed a Debye-solvent approximation to the susceptibility function of theenvironment. The resulting spectral density is of the Ohmic form with a Drude-Lorentzregularisation [87], which has the same low-frequency behaviour as the overdamped Brow-nian oscillator model [104, 105]. The Ohmic form for the spectral density describes aheavy Brownian particle immersed in a fluid of relatively light particles [61]. The modelis linear in bath-mode-frequencies with a long-tail high frequency cut-off. For the caseof a Debye-solvent, the model assumes that the environment is a homogeneous dielectric[41, 42] based upon the Onsager model of solvation [106, 107]. This combined with the35Debye-solvent approximation models the dielectric solvent well, provided it is highly di-lute [108].While this model is useful in considering a specific microscopic model, of given so-lute dimensions as well as solvent dielectrics, it is generally restricted to the class ofsolvents and excludes any rigid protein structure in the environment. In either case onewould expect the influence of both acoustic and optical phonons to play a significantrole in chromophore dynamics. Nonetheless, light-harvesting molecules are suspended inan aqueous solution and we expect the solvent to affect the solute dynamics as well aspossibly explaining why the Ohmic form for the spectral density has had some successin modelling light-harvesting molecules [109, 110, 111].Light-harvesting molecules are surrounded by a protein structure that contains a largenumber of degrees of freedom. We are often interested in the response of the system tothe low-frequency portion of the spectrum and thus a super-Ohmic spectral density canbe employed. A number of studies on these systems have used a super-Ohmic spectraldensity [112, 46, 113]. One can of course turn to experimentally determined spectra foran understanding of the nature of system-environment interactions. An accurate de-termination of the energetic structure and system-environment interactions present in acomplex system such as FMO is difficult. As a result, the general method of determiningthese parameters is to perform complex fitting algorithms on experimentally determinedoptical spectra, and various approximations have been proposed and implemented onFMO.One of the canonical works on FMO line-spectra calculations is that by Adolphs andRenger [46]. They estimated the exciton-environment coupling based on the fluorescenceline-narrowing spectrum. The result was a spectral density of cubic form with a peak ataround 0.02eV and an exponential cut-off for high frequencies. Nalbach et al. extendedthis to a ω5 frequency dependence as well as the inclusion of a single optical transitionpeak. They found the inclusion of the peak to alter the exciton dynamics only slightly[26]. Flemming’s group [21] calculated the exciton dynamics using a low-frequency Ohmicapproximation for the spectral function that relies on the Drude approximation for thesolvent environment [41]. Olbrich et al. [114] have calculated the spectral function basedon a combination of exponential and damped oscillations for the correlation function. Theresulting function is combination of a temperature-dependent function and Lorentzians.The number of peaks corresponds to the number of terms included in their fitting func-tions to experimental data. Depending on their algorithm, they find around 13 peakswith various spectral weights.36Figure 2.1: Spectral densities for the FMO trimer usedin studies by Ishizaki and Fleming [21], Cho et al. [103],Adolphs and Renger [46], as well as Nalbach et al. [26].Spectral densities fitted to experimental data by Olbrich[114] also included, indicated by the lines ‘present - BChl1-6’, ‘present - BChl 7’, ‘present - BChl 8’, correspondingto the average over BChl molecules 1-6, 7 and 8 respec-tively. The inset shows an enlarged energy and spectraldensity range. Reprinted with permission from [114].We see from Figure 2.1, where the results of the various approximation methods are com-pared, that they differ not only in qualitative features, but also in magnitude. The studiesof Cho [103], Adolphs [46] and Ishizaki [21] assumed only low-frequency vibrational modesof either Ohmic or super-Ohmic form. Nalbach [26] included a single internal vibrationalmode, evident by the single peak on top of the low-frequency distribution. However, themost up-to-date study by Olbrich et al. [114] not only found many more optical phonontransitions in the spectra, but found a much larger amplitude to the spectral densityacross frequencies, especially for the outlying BChla molecules 7 and 8 that sit on theperipheries of the FMO structure and are thus weakly bound. The core of the FMOstructure, comprised of BChla 1-6, which were the focus of the other studies, are morerelevant to compare. Olbrich et al. still finds an elevated amplitude across frequencies37but less significant optical phonon transitions. This is with the exception, however, of astrong peak at around 0.22eV. This peak has been attributed to several carbon-oxygenand carbon-carbon bonds, present in the BChla molecules internal structure that vibratewith similar frequencies [115, 114]. For low-energies below 0.01eV, Olbrich find system-bath interactions approximately 2-3 times larger than other studies. In Table 2.1, twooptical phonon peaks are selected from the FMO spectra calculated by Olbrich [114],corresponding to the 0.22eV peak and a lower frequency peak at 0.075eV. We see howthe dimensionless coupling parameter associated with these peaks is much less than onewhich suggests a perturbative approach to modelling the system-bath dynamics withthese optical phonon peaks should be sufficient. This is done in Section 3.3. In theinterest of investigating the non-perturbative effect of optical phonon peaks on the FMOdynamics, I also select peaks at 6, 8 & 10meVOptical phononsν λ(meV ) ξ(meV ) ω0(meV )0.003 35 12.5 2200.004 12.5 5 750.1 12.5 0.5 5Table 2.1: FMO parameters for optical phonons charac-terised by their peak height λ, full-width-half-maximumξ and dimensionless coupling parameter ν = λξ/piω2o inunits of meV . Taken from Olbrich et al. [114].In the interest of clarity I’ve converted all the parameters to units of meV (N.B. Planck’sconstant in the above units is ~ = 658meV · fs).For the acoustic phonons I take the FMO values determined from Jang and Mennucci[98], which provides an up-to-date, comprehensive review of fitted spectral functions forthe system-environment interactions in FMO. Here the spectral density function used isof the formJα = ραω3ω2ce−ω/ωc (2.15)where the corresponding dimensionless coupling is ζα = 0.31 and cut-off frequency ωc =8.7meV.For an Ohmic bath, the spectral density function used is of the form38Jα = ηαωe−ω/ωc (2.16)where the corresponding dimensionless coupling is γα = 0.22 and cut-off frequency ωc =21.1meV.2.4 The path-integral approach to modelling open-quantum mechanical systems and the Feynman-Vernon influence functionalWhen investigating how a quantum mechanical system interacts with its environment,we are usually interested in tracking the degrees of freedom within the ‘central’ quantummechanical system under the influence of the degrees of freedom within the environment.This ultimately means that the mathematical formalism I wish to work with must cast theeffects of the external system only in terms of the coordinates of the central system. Thisway I can explore the free parameters of the central system of interest under the influenceof its environment without having to track every degree of freedom in the problem,system-plus-environment at once. One of the most successful formalisms that achievesthis is that of the Feynman-Vernon influence functional approach [37]. The resultingmathematical object that facilitates these calculations is known as the Feynman-Vernoninfluence functional. This method allows us to consider an arbitrarily strong coupling tothe environment and retain ‘memory effects’, such that the system can feed energy tothe bath and experience feedback as well.2.4.1 The path integral formalism of quantum mechanicsCentral to this method of modelling open quantum systems is the path-integral approachto quantum mechanics devised by Feynman [116]. When the formalism was introducedin the mid-20th century it provided physicists with the understanding that a quantummechanical system can be thought of as exploring every possible path available to it,with each path weighted by a phase factor. The overall probability amplitude for a givenprocess therefore amounts to integrating (or summing) over all possible paths the systemmight take. The resulting functional integral can be thought of as a sum over possibleconfigurations of the system at each consecutive time. I briefly describe this conceptbelow due to its relevance to our calculations in subsequent chapters.39Consider a time-independent Hamiltonian for a particle of massM in a one-dimensionalpotential V (x), with coordinate xH = T + V, T =p22M(2.17)where p is the particle’s momentum. The solution of the Schro¨dinger equation can bewritten in the form|Ψ(t)〉 = e−iHt/~ |Ψ(0)〉 〈x|Ψ(t)〉 =ˆdx′G(x, t;x′, t′ = 0) 〈x′|Ψ(0)〉 , t > 0 (2.18)with the introduction of the propagator (or Greens function)G(x, t;x′, 0) = 〈x| e−iHt/~ |x′〉 (2.19)Dividing the time interval into infinitely small pieces and utilising the Trotter product (seeSection 3.2 of [117]), which allows one to ignore the non-commutivity of the kinetic andpotential operators in the Hamiltonian, leads to the path integral form of the propagatorG(x, t;x′, 0) =ˆ xx′DxeiScl (2.20)whereScl =ˆ t0dt′Mx˙22− V (2.21)is the classical action andDx = limN→∞ˆdx1...dxN−1(MN2piit)N2(2.22)is the Feynman path integral measure. We see from equation 2.20 how the problem ofcalculating the Greens function for a system is reduced to ‘summing over all possiblepaths’. In the case of our relatively simple system, where there are only a few degrees offreedom in the central system, this functional integral amounts to the product of time‘slices’, where we sum over all possible configurations of the system at each time slice.Such a formalism for the transition amplitude involves a single path integral. In theinterest of modelling decoherence processes, one must consider a double path integralformalism for the density matrix. For the isolated case, where there is no environment40and just a central system, the two paths are independent. In realistic physical systemswe have an environment interacting with our ‘central system’ and this couples the twopaths. Therefore I must include some system-bath interaction term in the Hamiltonianand then deal with all the degrees of freedom introduced to the problem with the bath.The method I employ here to do so is the influence functional method.2.4.2 The influence functionalI now include the environment by way of the influence functional method developed byFeynman and Vernon [37]. This method would later see application to the spin-bosonmodel analysed by Leggett et al. [31]. Here I briefly summarise the technique.Suppose that at time t = 0 the system and environment are uncoupled and the en-vironment is in thermal equilibrium. In this case the density matrix of the compositesystem is in a product stateρ(0) = ρS(0)⊗ ρE(0) (2.23)where ρS is the reduced density matrix for the central system and ρE that of the envi-ronment. It should be noted that this is of course a somewhat artificial situation butnonetheless one that can bear direct experimental relevance to a system if it is initiallyprepared in this state with sufficiently strong bias forces, which are subsequently switchedoff at t = 0. It’s also worth noting, within the context of the applied model to photosyn-thesis, that this initial condition has been argued to be valid [21]. Since the electronicexcitation process in EET corresponds to an excited (or ground) state prepared by pho-toexcitation in accordance with the Franck-Condon transition, this factorised initial stateshould be applicable. Incidentally the product initial state is the simplest starting pointfor these kinds of calculations and so is utilised in this work. I therefore assume thatthe system-bath coupling is suddenly switched on at t = 0 and I consider the dynamicsof ρ(t),∀ t > 0. Feynman and Vernon found the resulting form for the reduced densitymatrix to beρS(xf , x′f ; t) =∑xi,x′iKFV (xf , x′f ;xi, x′i; t)ρ(xi, x′i; 0)KFV (xf , x′f ;xi, x′i; t) =ˆ xfxDxˆ x′fx′Dx′exp{i~(SS[x]− SS[x′])}FFV [x, x′] (2.24)41where KFV is the propagating function determining the time evolution of the centralsystem under the influence of the bath. The object FFV [x, x′] is the Feynman-Vernoninfluence functional, which contains all the information about the bath as a function ofthe central system coordinates x, x′. Therefore the coordinates pertaining to the bathvariables are integrated out within the influence functional itself and we are left withonly the system coordinates as desired. One of the main results of the work of Feynman-Vernon was to calculate the exact form of the influence functional for a bath of harmonicoscillators coupled linearly to the central system. The result isFαβγδ[x, x′]= exp{− 1pi~ˆ tt0dτˆ τt0ds[− iL′(τ − s)(xα(τ) + x′β(τ))(xγ(s)− x′δ(s))(2.25)+ L′′(τ − s)(xα(τ)− x′β(τ))(xγ(s)− x′δ(s))]} (2.26)where xα, x′α represent the two paths of the density matrix that can visit each state inthe system [37, 61]. The full influence functional involves a sum over all possible pathsF =∑αβγδFαβγδ (2.27)The bath correlatorsL′(τ − s) =ˆ ∞0dωJ(ω) sinω(τ − s) (2.28)L′′(τ − s) =ˆ ∞0dωJ(ω) cosω(τ − s) coth(β~ω/2) (2.29)describe the time evolution of the bath, and contain the bath spectral density J(ω). L′(t)is related to the damping kernel in the classical Langevin equation [61] and physicallydescribes the coherent exchange of energy with the environment. L′(t) determines the lossof phase coherence in the system due to fluctuations in the environment. It is convenientto transform the double path integral over x and x′, into a single path integral that visitsall the possible states. For this I define ξ, χ to be the anti-symmetric and symmetricpaths respectivelyξαβ(t) ≡ xα(t)− x′β(t), χαβ(t) ≡ xα(t) + x′β(t) (2.30)42The state ξ(t) represents the off-diagonal terms of the rotated density matrix–the quan-tum coherences of the state between different positions. The state χ(t) indexes how fardown the diagonal of the density matrix we move and so tracks the incoherent hoppingof the system. So ξ(t) measures the difference in coordinates x−x′ and is zero wheneverthe system is in a diagonal (onsite) state and vice versa for χ. For a 3-site system, whichallows the values xα ∈ {x1, x0, x2}, we have for the density matrix thereforeρ =ρx1,x1 ρx1,x0 ρx1,x2ρx0,x1 ρx0,x0 ρx0,x2ρx2,x1 ρx2,x0 ρx2,x2 (2.31)Which is reparameterised in terms of the functions χ, ξ asρ = ρχ11,0 ρχ10,ξ10 ρχ12,ξ12ρχ01,ξ01 ρχ00,0 ρχ02,ξ02ρχ21,ξ21 ρχ20,ξ20 ρχ22,0 (2.32)10−1x(t)10−1x′(t)43210−1−2χ(t)210−1−2ξ(t)Figure 2.2: An example of possible paths x, x′ ∈{1, 0,−1} for the general 3-site system and theircorresponding symmetric/anti-symmetric paths χ ∈{2, 1, 0,−1,−2}, ξ ∈ {2, 1, 0,−1,−2}.In Figure 2.2 we see an example of possible paths taken by x and x′ for a general 3-sitesystem and the corresponding symmetric and anti-symmetric paths χ, ξ.Substituting the symmetric and antisymmetric paths into the influence functionalgivesFαβγδ[χ, ξ]= exp{ 1pi~ˆ tt0dτˆ τt0ds[2iL′(τ − s)χαβ(τ)ξγδ(s)− L′′(τ − s)ξαβ(τ)ξγδ(s)]}(2.33)I assume here that the paths take the form of an instantaneous flip whereby the pathjumps between the available locations instantly. Assuming the system begins in a sojournstate, and flips between successive sojourn states via blips in either transition |1〉 ↔ |0〉or |2〉 ↔ |0〉, I can parameterise the symmetric/anti-symmetric paths according to thesudden-flip approximation [31]44χαβ(t) ≡∑i=0χαβj[Θ(t− t2i)−Θ(t− t2i+1)], ξαβ(t) ≡∑j=1ξαβj[Θ(t− t2j−1)−Θ(t− t2j)](2.34)This particular parameterisation in terms of the above combination of step functionsproduces sojourns that live for even-times t2i+1 − t2i and blips for odd-times t2j − t2j−1N.b. I assume we start in a sojourn at t = 0. This way, the system moves from sojourn-blip-sojourn and so forth. We’ll see below that the construction of the influence functionalis such that the effect of the bath is to suppress time spent in the blip states, i.e. thecoherent/off-diagonal states. And so in the classical limit, where the bath has fullysuppressed the blips, the system hops incoherently between sojourns i.e. diagonal stateswithin the system. Substituting into the influence functional and performing the integralsgivesFαβγδ[χ, ξ]= exp{ 1pi~[i∑i=0∑j=1χαβi ξγδj Xij −∑j(ξαβj)2Q′′2j,2j−1 −∑i=1∑j=2ξαβi ξγδj Λij]}(2.35)where I have defined the bath correlations Λij between blip-blip pairs {i, j} and blip-sojourn correlations Xij between a sojourn at i and a blip at jXij ≡ Q′2i,2j+1 +Q′2i−1,2k −Q′2i,2k −Q′2i−1,2k+1Λij ≡ Q′′2i,2j−1 +Q′′2i−1,2k −Q′′2i,2k −Q′′2i−1,2k−1 (2.36)where I have used the compact notationQ′i,j ≡ Q′(t2i − tj) Q′′i,j ≡ Q′′(ti − tj) (2.37)and whereQ′(t− t′) =ˆ ∞0dωj(ω)ω2sinω(t− t′) (2.38)Q′′(t− t′) =ˆ ∞0dωj(ω)ω2(1− cosω(t− t′)) coth(β~ω/2) (2.39)45So a sojourn at time interval i can interact with a blip at j and a blip can interact withitself as well as other blips.It’s worth noting that the prescription to sum over all paths in the path integralformalism involves an integration over each ‘time-slice’ as well as a sum over all thepossible blip and sojourn states that are available within each transition. The first termin the influence functional represents interferences between sojourn-blip-sojourn paths,i.e. paths that start in some sojourn, undergo a blip, and return to the same sojourn.The second term is the self-interaction of the blips, so the same blip interacting withitself. The third term is the interaction between different blips, therefore a blip at someinterval tj− tk in one transition can interact at a later time with another blip from eithertransition.The sojourns themselves represent the diagonal elements of the density matrix whilethe blips represent the off-diagonal elements. Therefore one should think of the blipsas the coherences between states in the system. We see that the second term in theinfluence functional serves to suppress the weight of these coherent paths. Physically,this represents the bath ‘measuring’ the system and leading to the destruction of phasecoherence between paths through the system. The effect of the bath therefore is to inducedecoherence in the system, such that the paths eventually are restricted to the classicalcase of hopping between discrete sites.2.5 The noninteracting-blip approximation (NIBA)The full influence functional as it stands includes all possible pairings of blips and so-journs as well as blip-blip interactions. This means that all time-non-local interactionsare included such that a blip-sojourn/sojourn-sojourn interaction can occur for blips andsojourns separated across the entire time domain. Not only does this include an inordi-nate number of pairings for large number of flips n, but the exponentiation of these, inthe path integral formalism, creates all possible diagrams. Therefore in its general form,the influence functional presents a formidable mathematical object to evaluate and someapproximations are required. In the interest of retaining a non-perturbative system-bathcoupling, one approach to truncating the number of processes is to consider a systemthat dwells mostly in diagonal states of the density matrix. The system is still permittedto occupy the off-diagonal elements, and therefore we retain those quantum coherenceeffects; however, these excursions are infrequent and short-lived. This means that time-non-local interactions between these off-diagonal states with each other are suppressed46due to the fact that they are separated by large time-intervals, and the memory effectsintroduced via the bath are lost. This amounts to the well-known noninteracting-blipapproximation (NIBA) [31].The underlying assumption here is that the average time spent by the system in adiagonal element (sojourn) of the density matrix 〈s〉 is much larger than that spent inan off-diagonal element (blip) 〈b〉. This ultimately amounts to neglecting interactionsbetween different blips except for the self-interaction of blips. This leads to a strongsuppression of time spent in the off-diagonal terms of the density matrix as the bath israpidly measuring the state of the central system and forcing it back to diagonal states.For the 3-site-system of interest, this mathematically corresponds to considering only theinteraction between diagonal states χ11, χ22, χ00 and off-diagonal blips ξ10, ξ20.Figure 2.3: Diagrammatic depiction of possible co-herent states (blips) permitted in 3-site-V model underNIBA. ξ10 represents a wavefunction overlapping withstates |1〉 and |0〉, while ξ10 represents an overlap with|2〉 and |0〉To see this we first inspect the term involving the bath correlator Q′(t− t′) in equation2.39. In the influence functional, this term appears in the part involving sojourn-blipinteractions. Upon inspection we see that the full function Xij that contains all the bathcorrelators Q′, reduces to just the single term Q′(t2j − t2j−1), in the limit of 〈b〉/〈s〉 → 0as the length of sojourns dominates. Terms that involve interactions between blips andsojourns and are not nearest-neighbours in time, contain arguments that are very largeand they become rapidly oscillating functions. By the Riemann-Lebesgue lemma theseterms go to zero with the integral over ωˆ ∞0dωJ(ω)ω2sinω(τ)→ 0 as |τ | → ∞ (2.40)47The blip-blip interactions are handled in a less direct way within NIBA. The third termin equation 2.35 represents the interaction of non-time-local blips, and for this to beneglected, the blip-blip-propagators contained in Λj,k ∀j 6= k must all be much less thanthe blip-sojourn propagators contained in Xj,j−1 as well as the blip-self-interaction termcontaining Q′′j,j−1. This amounts to minimising the ratioQ′′〈s〉+ 〈b〉)Q′(〈s〉) (2.41)where 〈s〉, 〈b〉 are the average sojourn and blip times respectively. Q′′〈s〉+ 〈b〉) containsas its argument the time between nearest-neighbour blips, separated by the sojourninterval. In order to better understand this part of NIBA, I pre-emptively quote theforms of the bath correlators Q′(t), Q′′(t) for an Ohmic spectral density so that I canjustify the approximation in the appropriate regime. The Ohmic form of the correlatorsare [31]Q′(t) = arctan(ωct) (2.42)Q′′(t) =12ln(1 + ω2c t2) + ln[~βpitsinh(pit~β)](2.43)The length of a sojourn is of the order 1/∆ and if one considers timescales in the problemωc/∆ 1, then the blip-sojourn propagator Q′(〈s〉) = arctan(ωc/∆) ∼ pi/2. In the limitof 〈s〉 〈b〉, the argument of the blip-blip-propagator term, which contains non-localinteractions between blips, vanishes, i.e.lim〈s〉〈b〉Λjk → 0 (2.44)This is because the propagators in Λij involve non-local blip-blip interactions whosetime arguments are of the order of 〈s〉 + 〈b〉 ∼ 〈s〉 and therefore cancel. In Figure2.4 we see a graphical representation of the nearest-neighbor contribution to the term∑i=1∑j=i+1 ξiξjΛij, where the four propagators lines correspond to the four terms inthe equationΛ12 ≡ Q′′(t4 − t1) +Q′′(t3 − t2)−Q′′(t4 − t2)−Q′′(t3 − t1) (2.45)and all the terms in this function cancel when the resolution of the blips 〈b〉 vanishes.48t1 t2 t3 t4〈b〉 〈s〉Figure 2.4: Diagrammatic representation of the in-teraction terms contributing to the blip-blip interactionpropagators Λjk. The contribution depicted is a nearest-neighbor interaction of blips.In summary, NIBA amounts to setting the factors Q′i,j to zero for j 6= i − 1 such thatwe’re just left with Q′2j,2j−1 = Q′(t2j − t2j−1). I also set all Q′′ij equal to zero as theyrepresent blip-blip interactions.As it stands, the full density matrix in equation 2.32 presents a fairly formidablecombinatorics problem in terms of the possible paths that can be tracked through it.Since in NIBA we are permitted to spend only one time interval in an off-diagonal stateof the density matrix, this greatly reduces the number of possible paths. In this case thestates ξ12, ξ21 become inaccessible and the intermediate states χ10, χ20, χ12 are neglectedas well.ρ =ρχ11 ρξ10 0ρξ01 ρχ00 ρξ020 ρξ20 ρχ22 (2.46)So the blips and sojourns for the 3-site-V system within NIBA can take on the valuesξαβ ∈ {ξ10, ξ20}, χαβ ∈ {χ11, χ00, χ22} (2.47)Having established that NIBA corresponds to the situation where 〈s〉 〈b〉, it is prudentto ask in what regions of the parameter space of the system is it valid. A general centralsystem of interest can be parametrised by bias energy , tunnelling energy ∆, temperature49kBT and system-environment coupling λ. The first intuitive limit to consider, which ismost applicable to this thesis, is where ∆. In this limit, one intuitively expects thesystem to spend more time in the diagonal elements of the density matrix as opposedto the off-diagonal elements. This can be appreciated as the limit where the centralsystem is best described in the site basis as opposed to the eigenstate basis. In the lattercase, where tunnelling energies dominate, states are hybridised and one would expectthe system to spend relatively long time intervals in off-diagonal elements of the densitymatrix. However, NIBA is known to be applicable in other areas of the parameter spaceand is worth commenting on there.When the bias energy is not necessarily large, long blips are known to be suppressedat long times in conjunction with large damping and/or high temperatures for Ohmicbaths. Incidentally both of these criteria are actually met in this thesis when Ohmicbaths are considered, indeed they are necessary assumptions to achieve analytic solutions.For super-Ohmic spectral densities, blip-blip interactions tend to be suppressed relativeto intrablip interactions for large temperatures, and the NIBA condition is met [61].Finally, for zero bias, NIBA is also known to be valid in the limit of weak system-bathcoupling such that only one-phonon processes need be considered [31]. In fact it is exactin this limit, and can be used to evaluate the bath correlation functions when applyingperturbation theory in the system-environment interaction.2.6 Validity of NIBA: a quantitative measureI established in the preceding section that the condition that must be met in order toignore the time-non-local blip-sojourn and sojourn-sojourn interactions in the influencefunctional is 〈s〉 〈b〉. How this condition is met depends upon the various parametersin the model, as discussed above, and will change depending on the limits I inspect.In order to quantify the validity of NIBA in each limit more precisely, I outline themathematical condition that applies. As the influence functional contains all the bathparameters, including the system-bath coupling, it dictates the nature of the blips andsojourns in the system. I am interested specifically in the average length of a blip 〈b〉, andI can extract this quantity from a consideration of the various moments to the probabilitydistribution that the influence functional represents. To see this I first expand the generalinfluence functional in a power series about λ = 0F (λ) = F0 + F1λ+ ... (2.48)50whereF0 = limλ→0F (λ), F1 = limλ→0∂∂λF (λ) (2.49)The first term, F0, in the expansion represents the purely incoherent decay rate from thesystem and the ratio, F1/F0 = 〈b〉, is the first moment of the system. To see this, recallthat the expansion of the moment-generating function of a random variable X isMX(t) = 1 +m1t+m2t22!+ ..., mn = limt→0dnMXdtn(2.50)where mn is the nth moment. We further recall that the term in the influence functionalthat’s a function of Q′′(t), controls the width of the blips. Therefore, aside from theadditional cosine terms coming from the tunneling matrix element renormalisation, thefirst moment of the influence functional tells us the average blip width 〈b〉. NIBA is validwhen 〈b〉/〈s〉 = F1 1 since 〈s〉 = 1/F0 [31].I have outlined above a prescription for calculating the mathematical condition forwhich NIBA is valid. This means that for a given central system, with certain spectraldensity function for the environment, one can calculate the quantity F1 and inspect theregimes in which it is minimised for the various parameters in the model. In subsequentchapters I will perform this calculation explicitly each time NIBA is invoked. For thecase of photosynthesis, I have already identified the values of all the relevant parameters,which means a calculation of F1 will result in a final number for inspection.51Chapter 3Limiting and perturbative analysisof the 3-site-boson modelBefore entering into the full, non-perturbative calculations of the 3-site system plus os-cillator bath, I first inspect some of the limiting cases of the model so as to betterunderstand the system. To the same end I also then look at the inclusion of certainparameters perturbatively.3.1 The central 3-site systemThe first case that I look at in detail is that of the ‘bare’ 3-site system which comprisesthe central system of interest. It is described as ‘bare’ as I exclude the oscillator bath fornow.3.1.1 The 3-site-V configuration and population trappingIt is well understood that the addition of a perturbation to the degenerate two statesystem (with no coupling between states) , gives rise to an avoided crossing in the energydispersion diagram. This is often discussed in terms of lowered energy eigenstates andthus more stable configurations in chemical physics [118], level repulsion and tunnelling.The fact that the eigenenergy dispersion lines avoid any crossing reflects the fact thatan excitation in the system can tunnel between the available eigenstates. This is due tothe introduction of a potential energy term that is present even if the states are broughtinto resonance at zero energy. In certain configurations of 3-level systems; however, thesituation is more interesting. Let us consider the Hamiltonian52H3 =+ δ ∆10 0∆10 0 ∆200 ∆20 − δ (3.1)which describes a 3-site system with the zero energy set at the second site, two differentcouplings and some detuning between the upper levels. I’ll refer to this configuration asa V-system [86, 119]. Finding the eigenvalue expressions for the general case of δ 6= 0requires solving a cubic equation and produces long, unilluminating expressions. If weset δ = 0 however, one can easily find the eigenenergies to beλ± =12(±√2 + 4(∆210 + ∆220)), λ0 = (3.2)With δ 6= 0, a cubic polynomial provides the energies, the exact expressions of which I’vechosen not to include here due to their lengthy nature; however, I’ve included them inAppendix C.53-4 -2 0 2 4-4-2024ϵλ(a)-1.0 -0.5 0.0 0.5 1.0-10123Δλ(b)-2 -1 0 1 2-2-1012δλ(c)-2 -1 0 1 2-2-1012δλ(d)Figure 3.1: Eigenvalues λ0(red), λ+(blue), λ−(orange),as a function of (a) bias energy for the tuned systemwhere 1 = − δ/2, 2 = + δ/2, (b) tunnelling matrixelement separation ∆ where ∆10 = −∆/2,∆20 = ∆/2,(c) detuning δ where ∆10 = 0.5,∆20 = 0.5, (d) detuningδ where ∆10 = 0.5,∆20 = 1 . In units of .We see that the avoided crossing between the two eigenergies λ1,2 is equal to 4(∆210 +∆220)but there also exists an eigenenergy λ0 with a linear dispersion. Considering that onecould always make a shift in the zero energy, let us refer to this state as the zero energy, aconcept often discussed within the field of quantum optics in terms of 2-photon resonanceof dressed states [120]. This represents an anti-symmetric coherent superposition of theupper two states with eigenvectorδ = 0 : |D〉 = 1√∆210 + ∆220(∆20 |1〉 −∆10 |2〉)(3.3)and for completeness I will also include the symmetric combination, often referred to asthe ‘bright state’54δ = 0 : |B〉 = 1√∆210 + ∆220(∆10 |1〉+ ∆20 |2〉)(3.4)The reason this has been ascribed the population trapping state (or dark state) monikeris because even if one considers some decay mode out of the ground state |0〉, the state|D〉 never actually ‘sees’ this decay mode. In other words, as equation 3.3 contains nooverlap with the lower-decaying site, any excitation propagating through the system thatenters in to this state will remain there. From equation 3.3 we see that, if we start thesystem of in one of the upper states, for e.g. state |1〉, then the final long-time populationleft in the system after |±〉 have decayed, is the overlap of state |1〉 with the dark-state| 〈1|D〉 |2 = 11 +(∆10∆20)2 (3.5)We see that the ‘trapped’ population is maximised for ∆10 = ∆20.As is the case with the 2-level system, the ratio of the tunnelling elements to biasenergies ∆/ describes the general nature of the dynamics. If this ratio is small, then theeigenstates are small deviations from the basis states |α〉. In the opposite limit, wherethe bias energies go to zero, the eigenstates are mixtures of the basis states in the system.The resulting eigenstates include the dark-state |D〉 included above as well as the states|±〉 = 1√∆210 + ∆220 +N2+(∆10 |1〉 −N+ |0〉+ ∆20 |2〉)(3.6)In equation 3.6 one notices the link between |±〉 and the bright-state |B〉. While thebright-state is not an eigenstate of the 3-site-V-system, it is an eigenstate of the unbiased2-level system. As |D〉 remains unchanged upon the introduction of a 3rd state (|0〉) anda bias energy , it’s the bright state that must split into |±〉 becoming a mixture of thebasis states of the 3-basis-states. While in the limit /∆ 1, it would be appropriateto cast the 3-site-boson Hamiltonian in the eigenstate basis {|D〉 , |+〉 , |−〉} and solvethe system from there. However because we’re mostly concerned with the opposite limit,/∆ 1, in the application to photosynthesis, we therefore remain in the state basis.3.1.2 A path-integral formalism for the bare 3-site systemI now develop the path-integral formalism for the bare 3-site system–the isolated central3-site system without any environment–by first calculating the transition amplitude, or55Greens function [121]. This will form the basis of my analysis using the path integral for-malism to include the system-bath coupling both perturbatively and non-perturbatively.As discussed in Section 2.4.1, the bare amplitudes in the Hamiltonian formalism canbe related to the transition amplitudes for the bare system in the path integral formalism[61] byA[σf , σi] =ˆ σfσiDσeiS0[σ] = 〈σf | e−iH0t/~ |σi〉 (3.7)and then for an unbiased system I rewrite the exponential as an infinite sum〈σf | e−iHˆ0t |σi〉 = 〈σf |( ∞∑n=0(−iHˆ0t/~)nn!)|σi〉 (3.8)where I’ve been careful to identify the fact that Hˆ0 is in fact an operator. Therefore I’veavoided the power series expansion in terms of individual terms within the Hamiltoniandue to the non-commutativity of the constituent operators eA+B 6= eAeB for [A,B] 6= 0.Instead we inspect the even and odd power contributions to equation 3.8,A[σf , σi; t] = 〈σf |( ∞∑n=0(−iHˆ0t/~)2n(2n)!+∞∑n=0(−iHˆ0t/~)2n+1(2n+ 1)!)|σi〉 (3.9)I findA[σf , σi; t] = 〈σf |(1+∞∑n=1(∆210+∆220)n−1 (−it)2n(2n)!(Hˆ0)2+∞∑n=0(∆210+∆220)n (−it)2n+1(2n+ 1)!Hˆ0)|σi〉(3.10)where for the unbiased systemH0 = 0 ∆10 0∆10 0 ∆200 ∆20 0 , (H0)2 = ∆210 0 ∆10∆200 ∆210 + ∆220 0∆10∆20 0 ∆220 (3.11)The bare amplitude can then be rewritten as an infinite series of products of time-intervals. For the A11(t) element this readsA11(t) = 1 + ∆210∞∑n=1(−i)2n(∆210 + ∆220)n−1ˆ t0dt2nˆ t2n0dt2n−1...ˆ t20dt1 (3.12)Once I consider the full system, including the effects of the bath by way of the influence56functional, I’ll proceed with the analysis by way of taking Laplace transforms. Here Ipreemptively take the Laplace transform of 3.12 and also check the method against asimple matrix inversion to calculate the Greens functionA11(λ) =1λ+ ∆210∞∑n=1(−i)2n(∆210 + ∆220)n−1L[ˆ t0dt2nˆ t2n0dt2n−1...ˆ t20dt1](3.13)which yieldsA11(λ) =1λ+ ∆210∞∑n=1(−i)2n (∆210 + ∆220)n−1 1λ2n+1 (3.14)and the summation can be performed to giveA11(λ) =λ2 + ∆220λ3 + λ(∆210 + ∆220)(3.15)The inverse-Laplace transform yields the dynamicsA11(t) =∆220 + ∆210 cos(√∆210 + ∆220t)∆210 + ∆220(3.16)It’s worth remarking at this point that the (∆210+∆220)n−1 term arises due to the ambiguityas to whether the system flips in either branch. This is where we pick up the abovesuperposition of ∆210 and ∆220 terms due to the summation over both outcomes, n − 1times, until we reach the final state.To check our result here we can compare to the calculation of the bare 3-site Greensfunction by simple matrix inversionAσ,σ′(λ) = 〈σ|(1λ1+ iHˆ0)|σ′〉 (3.17)where the A11(λ) element is given byA11(λ) =∣∣∣∣∣∣∣λ −i∆10 0−i∆10 λ −i∆200 −i∆20 λ∣∣∣∣∣∣∣−1 ∣∣∣∣∣ λ −i∆20−i∆20 λ∣∣∣∣∣ (3.18)which does indeed give equation 3.15 confirming the above formulation of the 3-site path57integral.The simplest case, which is useful to note here for future calculations, is the ground-state propagatorA00(λ) =1λ+∞∑n=1(−i)2n(∆210 + ∆220)nL[ˆ t0dt2nˆ t2n0dt2n−1...ˆ t20dt1](3.19)which givesA00(λ) =λλ2 + (∆210 + ∆220)(3.20)andA00(t) = cos(√∆210 + ∆220t)(3.21)Now I would like to complete the bare-system by adding the biases 1 and 2. Theintroduction of the kinetic energy term is using the Suzuki-Trotter decomposition of thepropagator [117]. The simplest way to proceed is to add the appropriate bias terms toequation 3.13 and explain the physical basis for each termA11(λ) =1λ+ i1+ ∆210∞∑n=1(−i)2nL[ˆ t0dt2ne−i1(t−t2n)ˆ t2n0dt2n−1ˆ t2n−10dt2n−2(∆210e−i1(t2n−1−t2n−2) + ∆220e−i2(t2n−1−t2n−2)) ˆ t2n−20dt2n−3...ˆ t20dt1e−i1(t2−t1)](3.22)The first term represents the free-particle Greens function: the transition amplitude forthe system to remain in state |1〉 stipulated by the initial conditions. The second term,which contains all higher order processes, begins with a transition matrix element ∆210,coming from the initial and final conditions. The system starts and ends in state |1〉,therefore with one or more flips in the system, there will at least be two coming fromthe left branch so as to bring the system back again. So we see that the n = 1 termsatisfies this condition, and only contains the propagator for state |1〉 integrated over oneintermediate time interval that lets the system sit in state |0〉 for some time (weightedby unity since 0 = 0) and then return to |1〉. All higher order terms contain n − 1 ofthe superposition (∆210exp[−1(t2n−1 − t2n−2)] + ∆220exp[−2(t2n−1 − t2n−2)]). This comes58from the fact that whenever the system enters the lower site |0〉, it has the option to eithersubsequently flip in either branch. We also see that there will be n contributions fromthe number of times the system visits site |3〉, which will show up as a 1/λn term in theLaplace transform. A sum over all possible paths will include any combination of theseprocesses, i.e. the system should be allowed to propagate as many times as it likes ineither branch and with any possible combination of flips in both branches. Computingthe Laplace transform I findA11(λ) =1λ+ i1+ ∆210∞∑n=1(−i)2n(1λ+ i1)2(∆210λ+ i1+∆220λ+ i2)n−11λn(3.23)and computing the series summation then givesA11(λ) =λ(λ+ i2) + ∆220λ(λ+ i1)(λ+ i2) + (λ+ i1)∆220 + (λ+ i2)∆210(3.24)For the case where the upper two levels are ‘tuned’ such that 1 = 2, the transitionamplitude can be inverse Laplace transformed to giveA11(t) =∆220eit∆210 + ∆220+2∆210eit2− 12it√4∆210+4∆220+2√4∆210 + 4∆220 + 2(√4∆210 + 4∆220 + 2 + )+2∆210e12it√4∆210+4∆220+2+ it2√4∆210 + 4∆220 + 2(√4∆210 + 4∆220 + 2 − ) (3.25)Once again we can check this result against a simple matrix inversion which does indeedproduce equation 3.24, confirming our procedure for the biased case as well.I similarly calculate the ground state propagatorA00(λ) =1λ+∞∑n=1(−i)2n(∆210λ+ i1+∆220λ+ i2)n1λn+1(3.26)which yieldsA00(λ) =(λ+ i1)(λ+ i2)λ(λ+ i1)(λ+ i2) + (λ+ i1)∆220 + (λ+ i2)∆210(3.27)593.1.3 Population trapping in a 3-site-V systemIn this section I would like to demonstrate the phenomenon of population trapping inthe dynamics of the 3-site-V system. The simplest way to do this, is by considering adecay channel coupled only to the ground state |3〉 that manifests from an environmentwith a white-noise spectrum [61]. This represents a form of irreversible decay [120] fromthe system. For this, the Greens function isA11(λ) =(λ+ Γ)(λ+ i2) + ∆220(λ+ Γ)(λ+ i1)(λ+ i2) + (λ+ i1)∆220 + (λ+ i2)∆210(3.28)and the configuration with the introduction of the sink is depicted in Figure 3.2Figure 3.2: Diagrammatric representation of 3-site-Vmodel coupled to a sink via an irreversible decay channelfrom the ground stateAs before I calculate the inverse-Laplace transform to find the time-domain amplitude.Then the probability to start and return to site one can be calculated fromP11(t) = K11;11(t)ρ11(0) (3.29)where the density matrix propagator isK11;11(t) = A∗11(t)A11(t) (3.30)60(a) (b)(c) (d)(e) (f)Figure 3.3: Return probability to site 1: P11(t), as afunction of time for tuned (red) and detuned (blue) cases.In units of 1 = 1. 2 = 1.51(blue), 2 = 1 = 1(red).(a): Γ = 0.1, ∆10 = 0.5, ∆20 = 0.3, (b): Γ = 0.1, ∆10 =1,∆20 = 0.8, (c): Γ = 0.1, ∆10 = 0.5,∆20 = 0.5, (d):Γ = 0.5, ∆10 = 0.2,∆20 = 0.1, (e): Γ = 0.5, ∆10 =0.1,∆20 = 0.1, (f): Γ = 0.5, ∆10 = 1,∆20 = 0.1In Figure 3.3 we see the return probabilities to site 1 for various areas of the parameterspace. The blue line corresponds to the detuned system with a δ = 0.25 and the red linesthe tuned case with δ = 0. The irreversible decay from state 0 is quantified in the rate Γ.61For the underdamped case, where Γ is sufficiently small to permit the system to executevarious cycles, we see damped coherent oscillations in the dynamics. The beat signalscharacteristic of the two interfering frequencies in the 3-site system are present in theunderdamped case and visible for a few cycles (see Figures 3.3b and 3.3c). However thebeats are evidently suppressed after a few cycles. Population trapping is visible in theunderdamped cases with the tuned (red) lines maintaining some asymptotic occupationprobability for long times, while the detuned case (blue), decays eventually.The efficacy of the population trapping effect is evident even in the overdamped case(see Figures 3.3d and 3.3e), where the system is unable to even complete one cycle in thedetuned case, while the tuned case makes it through about one cycle before settling into the ‘dark state’. Furthermore we see the effect of tuned tunnelling matrix elements onthe dynamics. In the underdamped case, Figure 3.3b, demonstrates the complimentaryeffect that the tuning of the tunneling energy has on population trapping. We see that itleaves the tuned system with increased occupation probability at long times. This effectis also clearly visible in the overdamped regime too (see Figure 3.3d). We can also seethe destruction of the dark-state properties of the system when the second tunnellingmatrix element ∆20 of the system is too small relative to the decay Γ out of site-3, andthe system is unable to enter into dark-state. Here we see that the long-time populationdecays to zero in a similar fashion to the 2-state system where population trapping isabsent.62(a) (b)(c) (d)(e) (f)Figure 3.4: Return probability to ground state |0〉:P00(t), as a function of time for tuned (red) and detuned(blue) cases. In units of 1 = 1. 2 = 1.51(blue),2 =1 = 1(red). (a): Γ = 0.1, ∆10 = 0.5, ∆20 = 0.3,(b): Γ = 0.1, ∆10 = 1,∆20 = 0.8, (c): Γ = 0.1,∆10 = 0.5,∆20 = 0.5, (d): Γ = 0.5, ∆10 = 0.2,∆20 = 0.1,(e): Γ = 0.5, ∆10 = 0.1,∆20 = 0.1, (f): Γ = 0.5,∆10 = 1,∆20 = 0.1In Figure 3.4 I calculate the ground-state dynamics P00(t) for the same parameters asthe P11(t) dynamics of Figure 3.3. We immediately notice how the ground state is63unaffected by the dark-state occupancy, considering that it has no overlap with |D〉.The tuned biases do however introduce a phase-shift in the dynamics as can be seen inFigures 3.4a,3.4b, 3.4c. In accordance with this, the dynamics decay faster, and coherentoscillations persist for shorter times. Coherent oscillations are entirely washed out for astrong ground-state-decay rate relative to the tunnelling matrix elements, and we are inthe fully incoherent regime.The other limit of interest within the biased regime is that of small bias relative totunnelling matrix element. In this regime we expect the tunnelling terms to dominate,such that we observe strong coherent oscillations persisting as well as beat frequenciesarising from the interference of the two branches.(a) (b)(c) (d)Figure 3.5: Return probability to state |1〉: P11(t), asa function of time for tuned biases (red) and detunedbiases (blue), as well as tuned tunnelling (green) anddetuned tunnelling (purple), in the small bias regime.2 = 0.2, 1 = 0.1 (blue), 2 = 1 = 0.1(red). (a): Γ = 0.1,∆10 = 1, ∆20 = 0.8, (b): Γ = 0.01, ∆10 = 1,∆20 = 0.8,(c): Γ = 0.5, ∆10 = 1,∆20 = 0.8, (d): Γ = 0.5,∆10 = 0.2,∆20 = 0.1, (e): Γ = 0.1, ∆10 = 1,∆20 = 0.864(a) (b)Figure 3.6: Return probability to site 1 for unbiasedsystem. ∆10 = 1, ∆20 = 0.5 (purple), ∆10 = 1, ∆20 =1(green). (a): Γ = 0.1, (b): Γ = 0.5We see from the unbiased plots in Figure 3.6 the presence of the zero energy ‘dark-state’eigenvalue λ = 0 in the long-time population trapping. The effect of tuning the transitionmatrix elements to aid in population trapping is also evident even in the case of strongdissipation.From the perspective of population trapping, it is instructive to investigate the longtime dynamics P11(t→∞).(a) (b)Figure 3.7: Return probability to site 1 at t → ∞:P11(∞), as a function of detuning δ = 1 − 2. ∆10 =∆20 = 0.1, = 1 (a) Γ = 0 (Bath off) (b) Dark-statepeak for various decay ratesFor the case of zero-decay rate to the reservoir 3.7a, we see how the detuning parameterδ introduces a phase shift in the dynamics. In the incoherent regime, with non-zerodecay to the reservoir, we see the effect that continuously varying the detuning has on65P11 for long times. For large detuning, again we see how for large level asymmetry (largedetuning), the system fails to leave |1〉. This is symptomatic of a large level-spacingleading to localisation and not population trapping. For small detuning, the system iswell within the parameter regime for the 3-site model. For small, but non-zero δ, thesystem is depleted due to Γ and there is no state 1 population for long times. For δ → 0however, the system enters in to the dark state and population remains for long times,indicated by the peak in Figure 3.7a.(a) (b)Figure 3.8: Long time return probability P11(∞) to site1 vs ∆10: ∆20 = ∆20 = 0.1, = 1; (a) δ = 1 (b) δ = 03.2 Weak system-bath coupling regime for the 3-site-V systemIf the system bath couplings are small then I can perform perturbation theory in theratio uαq = λα/ωq. Here I apply the Lang-Firsov polaron transformation [52] to thebare 3-site Hamiltonian coupled diagonally to an oscillator bath, which decouples thecentral system from the bath. To do this I define the unitary operator U = eS, whereS = −∑αq uαq(bq − b†q) |α〉〈α| and uαq = λαq/ωq such that H → UHU † = eSHe−S.Physically, this represents the shifting of the boson cloud to its new equilibrium position.So S can be thought of as a shift operator, for the diagonal system-bath interaction.Applying this transformation, and using the Baker-Campbell-Hausdorff formula, yields66H˜ = 1 |1〉〈1|+ 2 |2〉〈2|+ 0 |0〉〈0|+ ∆10(|1〉〈0|B10 + |0〉〈1|B01)+ ∆20(|2〉〈0|B20 + |0〉〈2|B02)+HB (3.31)where Bαβ = exp (φα − φβ), φα =∑q uαq(bq − b†q) and a constant (−1)∑q uαqλαq |α〉〈α|has also been dropped in H˜ such that the zero energy has been shifted accordingly1.Since the bath is decoupled from the central system via the polaron transformation, Ican calculate the ground-state density matrix propagator of H˜ viaK˜00;00(t) = 〈A˜00(t)A˜∗00(t)〉B (3.32)where the statistical average is over the bath coordinates and I assume that the initialdensity matrix components for the bath and central system factorise at t = 0 such thatρ(0) = ρS(0)⊗ ρB(0). Since the transformation has decoupled the bath from the centralsystem, I can infer the transition amplitude A˜(t) for the transformed Hamiltonian H˜,from the form determined in Section 3.1.2. Now each tunneling matrix element alsocontains a corresponding term Bαβ that shifts the boson cloud as well. I present theform for the unbiased transformed Greens function here for brevityA00(t) = 1 +∞∑n=1(−i)2nˆ t0dt2nˆ t2n0dt2n−1[∆210B10(t2n)B01(t2n−1) + ∆220B20(t2n)B02(t2n−1)]× ...ˆ t20dt1[∆210B10(t2)B01(t1) + ∆220B20(t2)B02(t1)](3.33)Therefore the ground-state-density-matrix propagator isK00;00(t) =〈1−ˆ t0dt2ˆ t20dt1[k10(t1, t2) + k20(t1, t2)]+ˆ t0dt4ˆ t40dt3ˆ t30dt2ˆ t20dt1[k10(t1, t2)k10(t3, t4) + k10(t1, t2)k20(t3, t4)+ k20(t1, t2)k10(t3, t4) + k20(t1, t2)k20(t3, t4)]+ ...〉B(3.34)1I can rewrite the shift operators in terms of the momentum operator φα =∑q(xαcq/mqω2q )pˆα wherewe can see more explicitly how the transformation serves to displace the oscillators about the distancexα − xβ , corresponding to the separation between wells67where the kernelskαβ(t, t′) =∆2αβ2[Bαβ(t)Bβα(t′) +Bβα(t)Bαβ(t′)](3.35)describe the ‘blips’ in the system, in the language of the spin-boson model. The termsBαβ(t) are the shift operators for the bath modes. They move the bosons between thepotential minima in the system and are dependent on the system-bath coupling, and thetemperature of the bath.The statistical average over the bath modes leads to the coupled-bath correlationfunctions containing the coupled-boson propagator, which I introduce shortly. We seethat in the general case of equation 3.34 above, where blip interactions can be long rangein time, that we have all higher order combinations of blip-blip interactions in the fullsummation. Since we’re operating in the weak-coupling limit, and the shift operatorsBαβ(t) are functions of the dimensionless coupling parameter uαq, we can partially sumthis series by only including the zeroth order blip interaction, i.e. the self-interactionof the blip. Blip-blip interactions contain terms that are higher order in uαq and aretherefore neglected. All of this is suggestive of a physical description of the situationthat involves the bath rapidly measuring the state of the system and thereby suppressingcoherent states within it. A blip–as discussed in the context of the spin-boson model–represents an off-diagonal excursion within the reduced density matrix of the centralsystem. In the current formalism it takes the equivalent form of a tunnelling process(back and forth between wells) that interacts with the bath along the way. With this inmind the partial re-summation isK˜00;00(t) = 1 +∞∑n=1(−i)2nˆ t0dt2nˆ t2n0dt2n−1(Σ10(t2n, t2n−1) + Σ20(t2n, t2n−1)) ˆ t2n−10dt2n−2× ...ˆ t30dt2ˆ t20dt1(Σ10(t2, t1) + Σ20(t2, t1))(3.36)whereΣαβ(t, t′) =∆2αβ2[e−i(α−β)(t−t′)/~〈Bαβ(t)Bβα(t′)〉B+ ei(α−β)(t−t′)/~〈Bβα(t)Bαβ(t′)〉B](3.37)and the bias terms have been reintroduced. Since I have decoupled the bath from the68central system, what I have here is a similar situation to the bare-3-site model of Section3.1.2. The system can flip in either transition and return to the ground state each time,but now each time the system flips, it displaces the bath and drags a ‘cloud’ of bosonswith it. This is accounted for by the averages 〈Bαβ(t)Bαβ(t′)〉 which displace the boson-cloud to site α at time t and then back to site β at time t′. The bath correlator Σ(t, t′)includes both the forwards and backwards paths due to the density matrix formalism.Implicit in the analysis so far is the assumption that bath correlations die off for longtimes, which is why I have only included nearest-neighbour-bath correlations. This isequivalent to the application of NIBA [99]. I will also investigate the regions of validityof this approximation later on in this chapter. Computing the Laplace transforms tosolve the system yieldsK˜00;00(λ) =1λ+∞∑n=1(−i)2n(Σ10(λ) + Σ20(λ))n 1λn+1(3.38)where the self-energies areΣαβ(t, t′) = ∆2αβ〈Bαβ〉〈Bβα〉 cos[(α − β)(t− t′)/~]eϕαβ(t−t′) (3.39)and 〈Bαβ〉 give rise to the Debye-Waller factor that renormalises the tunnelling energy([99]) ∆˜αβ = ∆αβ√〈Bαβ〉〈Bβα〉. This term represents the adiabatic renormalisation ofthe tunnelling energy due to high-frequency bath modes much greater than the tunnellingfrequency. The bath correlations are contained in the second, exponential term, whoseexpansion contains all possible diagrams of electron-boson interaction. The phase factoris given by the well-known coupled phonon propagator [52]ϕαβ(t) =∑q(uαq − uβq)2[nq(1− eiωqt)+ (1 + nq) (1− e−iωqt) ] (3.40)with the Bose occupation numbers nq = (eβωq − 1)−1, and the dimensionless couplingparameter uαq = λαq/ωq. In Appendix B I show how in the continuum limit, the phaseis equivalent toϕαβ(t) =1~ˆ ∞0dωJαβ(ω)ω2[i sin(ωt)− (1− cos(ωt)) coth(~βω/2)](3.41)69For weak system-bath interaction uαq 1, the 1st-order expansion of the self-energy interms of the phase ϕ, givesΣαβ(t) ≈ ∆˜2αβ(cos[(α − β)t/~](1 + ϕαβ(t))(3.42)Absorbing the bias term into the phase, such thatΣαβ(t) ≈ ∆˜2αβ(cos[(α − β)t/~]+ ϕ˜αβ(t, α − β))(3.43)whereϕ˜αβ(t, α − β) = cos[(α − β)t/~]ϕαβ(t) (3.44)The Laplace-transformed self-energy is thereforeΣαβ(λ) ≈ ∆˜2αβ(2λλ2 + [(α − β)/~]2 + ϕαβ(λ, α − β))(3.45)and the biased ground-state propagator becomesK˜00;00(λ) =1λ+ (Σ10(λ) + Σ20(λ))(3.46)Equation 3.46 has a formidable pole structure in the arbitrary bias regime 1 6= 2. In theinterest of producing tractable analytic results here, I inspect the case of tuned biases1 = 2 = , with the ground state 0 = 0. In this case the pole structure produces a poleat λ = 0 and solutions to the cubic polynomialλ3 + λ2Φ(λ) + λE2 + 2Φ(λ) = 0 (3.47)whereΦ(λ) = ∆˜210ϕ10(λ, ) + ∆˜220ϕ20(λ, )E =√2 + Ω˜2, Ω˜2 = ∆˜210 + ∆˜220 (3.48)N.B. when the bath couplings are zero (or equal), the bath decouples from the centralsystem and the phase factors go to zero. In this limit I recover the Laplace-transformeddynamics for the bare-site system of equation 3.15, as one would expect. In accordance70with the perturbative approach I have taken here, I assume that the phase factors aresmall fluctuations about the bare poles, and so I Taylor expand these functions to 1st-order and evaluate them at the pointsλD = 0, λ± = ±iE/2 (3.49)which are the solutions to the pole structure of equation 3.47 in the absence of dampingϕ = 0. Recall that in the bare system, is the energy required to excite the dark state|D〉 and ( ± E)/2 for the states |λ±〉. The term involving the bath correlators Φ(λ)(evaluated at the bare poles) can be evaluated with use of the Fluctuation-Dissipationtheorem [122, 31] (see Appendix D).3.2.1 Linear response and the fluctuation-dissipation theoremI now make a quick digression here to provide a heuristic derivation of the bath correlatorsabove in the linear response regime and show how this is a form of the well knownfluctuation-dissipation theorem.The fluctuation-dissipation theorem is a central feature of linear response theory andis applied here in the 3-site perturbative model. The theorem relates the relaxation ofa weakly perturbed system to the thermal fluctuations in the environment. The mainresult of the theory relates the power spectrum S(ω) of the fluctuations, to the Fouriertransform of the susceptibility χ(ω) (the linear response function). We begin with thebath correlation functionϕ(t) = cos(t/~)ˆ ∞0dωJ(ω)~ω2(i sinωt− (1− cosωt) coth(~βω/2))(3.50)the exponentiation of which produces all orders of possible bath interactions with thecentral system. In the perturbative limit the exponential is expanded to linear orderand the Laplace transform of the function reduces to just the transform of the bathcorrelation function. The phase therefore describes the corrections to the system intro-duced by the bath. The time-independent term can be factorised in the exponential asa renormalisation of the tunneling energy. So the remaining fluctuating part is71ϕ(t) = cos(t/~)ˆ ∞0dωJ(ω)~ω2(i sinωt+ cosωt coth(~βω/2))(3.51)In the perturbative limit we expect the solutions to the system to be small deviationsfrom the ’bare-system’ eigenvalues and thus I Taylor expand the phase about these points.In this case the complex frequencies of the transform take on the bare-eigenvalues λ0 andwe haveϕ(λ0) =1~ˆ ∞0dte−λ0t cos(t/~)ˆ ∞0dωJ(ω)~ω2(i sinωt+ cosωt coth(~βω/2))(3.52)For the bare-system, uncoupled from the bath, the complex frequencies contain no realpart and are related to the Fourier-transform real frequencies by λ0 ∼ iω0. Performingthe time integration first, and using the definition of the delta-function, we findϕ(ω0) =14~ˆ ∞0dωJ(ω)ω2{[δ(ω − ω0 + /~) + δ(ω − ω0 − /~)− δ(ω + ω0 − /~)− δ(ω + ω0 + /~)]+[δ(ω − ω0 + /~) + δ(ω − ω0 − /~)+ δ(ω + ω0 − /~) + δ(ω + ω0 + /~)]coth(~βω/2)}(3.53)and performing the frequency integralsϕ(ω0) =14~[J(ω0 − /~)(ω0 − /~)2 +J(ω0 + /~)(ω0 + /~)2− J(/~− ω0)(ω0 − /~)2 −J(−ω0 − /~)(ω0 + /~)2+J(ω0 − /~)(ω0 − /~)2 coth(~β(ω0 − /~)/2) +J(ω0 + /~)(ω0 + /~)2coth(~β(ω0 + /~)/2)− J(/~− ω0)(ω0 − /~)2 coth(~β(ω0 − /~)/2)−J(−ω0 − /~)(ω0 + /~)2coth(~β(ω0 + /~)/2)](3.54)A comparison with the formal result of the fluctuation-dissipation theorem providedin Appendix D shows the similarity between the two results. What I have done, iseffectively derive the power spectrum of the fluctuations within linear response theory72for the density matrix propagator. The linear response function quoted in the formaltheory usually applies to the 2-point propagator (Greens function). We’ve left the aboveresult general for any form of the spectral density, and we see below how it simplifies forthe various specific forms of interest.3.2.2 Results of the 3-site-boson model in the perturbative regimeThe solutions to equation 3.47 are obtained using the method outlined in AppendixE for solving cubic polynomials [123, 124, 125]. When the discriminant D of the cubicpolynomial is greater than zero: D > 0, we have 1-real and 2-complex roots. I can ascribethe real pole to the Dark-state, as with tuned upper energy levels we have no oscillationfrequency between them, so the resultant energy is purely real. The polynomial is solvedin the most illuminating form asλ˜D = Γγ(0)− Γr(0)λ˜± = −Γγ(±E)2− Γr(±E)± iω(±E) (3.55)whereΓr(z) =Φ(z)3, Γγ(z) = u(z)− v(z), ω(z) =√32(u(z) + v(z))u =3√√D − (2Φ3 − 9ΦE2 + 272Φ)54, ∀D > 0v =3√√D +(2Φ3 − 9ΦE2 + 272Φ)54, ∀D > 0 (3.56)So I have one entirely real pole λD that contains two decay rates Γr and Γγ. The remain-ing complex-conjugate poles contain a combination of the decay rates with a differentargument and the frequency of oscillation ω. The discriminant D tells us the nature ofthe polesD =(2Φ3 − 9ΦE2 + 272Φ2916)2+(3E2 − Φ281)3(3.57)For Ohmic spectral densities the bath correlation function takes the form73ΦD(0) =1[∆˜210γ10 + ∆˜220γ20] coth(ωc)coth(~β2)Φ±(±E) = [∆˜210γ10 + ∆˜220γ20]{coth(E−ωc)E − [1 + coth(~β(E − )2)]+coth(E+ωc)E + [1 + coth(~β(E + )2)]}(3.58)For super-Ohmic spectral densities the bath correlation function takes the formΦD(0) = [∆˜210γ10 + ∆˜220γ20] coth(ωc)coth(~β2)Φ±(±E) = [∆˜210γ10 + ∆˜220γ20]{coth(E − ωc)(E − )[1 + coth(~β(E − )2)]+ coth(E + ωc)(E + )[1 + coth(~β(E + )2)]}(3.59)3.2.3 Validity of NIBA in the perturbative regime for Ohmicand super-Ohmic dampingIn Section 2.5 I introduced NIBA and discussed the mathematical condition for its va-lidity. This amounted to the calculation of the 1st-moment of the self-energy (influencefunctional) which tells us the ratio of the blip/sojourn times F1. In order for the coher-ences in the system to be short-lived, this 1st-moment must be small, such that F1 1,and NIBA applies. Here we quantify it for the various spectral densities of interest withinthe perturbative regime. For Ohmic spectral densities the bath correlators take the formϕαβ(λ) =1pi~(ηα − ηβ)2λ coth(βλ2)(3.60)and for the super-Ohmic caseϕαβ(λ) =1pi~(ρα − ρβ)2λ3 coth(βλ2)(3.61)where λ are the complex-frequencies of the bare-3-site system. The total self-energy isΣ(λ) = Σ10(λ) + Σ20(λ), where74Σαβ(λ) ≈ ∆˜2αβ(λλ2 + (α − β)2 + ϕαβ(λ− i(α − β)) + ϕαβ(λ+ i(α − β)))(3.62)In order to calculate F1, I compute the derivative of the self-energy, which for both Ohmicand super-Ohmic spectral densities, I findF1 = limλ→0∂∂λΣ(λ) =∆21021+∆22022(3.63)The condition for F1 here in the perturbative regime, demonstrates the dependency ofNIBA on the ratio ∆/. We see that as the bias energy goes to zero, the 1st-momentdiverges, and thus we invalidate NIBA in this regime. Therefore NIBA only appliesto biased systems in the perturbative regime, specifically small values of the tunnellingfrequency relative to the bias energy. Qualitatively this is intuitive; as for large tunnellingfrequency, we expect the system to spend more time in the off-diagonal states, andcorrespondingly for small bias energies, the system becomes less localised in each well.The coherences become more pronounced in this limit, and thus interactions betweencoherences to higher orders must be considered.3.2.4 Coherent phase space in the perturbative regimeI now investigate the coherent-incoherent phase space of the model here in the perturba-tive limit. The sign of the discriminant D (see equation 3.57) determines what region ofthe phase space we’re in.75(a) (b)Figure 3.9: Coherent-incoherent phase space of tuned-perturbative 3-site-bath model. Parameters: ∆10 = 0.1,∆20 = 0.12 γ1 = 0.01, γ2 = 0.05, in units of . (blue)Ohmic, (orange) super-Ohmic.In Figure 3.9 I investigate only the region where NIBA is valid, i.e. for large bias relativeto tunnelling terms. We see that for all temperatures, the discriminant remains positive,and therefore the system is always in the underdamped regime. Therefore there are alwaystwo complex conjugate solutions and one real solution to the poles of the propagator.3.2.5 Decay times in the perturbative regime for Ohmic spec-tral densitiesHere I investigate the relationship of the various decay times in the system versus temper-ature and system-bath coupling for Ohmic spectral densities. In the high-temperaturelimit I can expand the temperature dependent function in terms of the dimensionlessparameter z/KBT 1coth(z2kBT)=2kBTz+2z3kBT+ ... (3.64)where z ∈ {, ± E}. If I take just the 1st order term I can investigate this for variousspectral densities. The relaxation time τr is the inverse of the pure relaxation rateΦ ≡ 1/τr and for an Ohmic spectral density, the times are76τ±r =1(± E)[(γ1 − γ3)2 + (γ2 − γ3)2] , T = 0K, (Ohmic)τDr =1[(γ1 − γ3)2 + (γ2 − γ3)2] , T = 0K, (Ohmic)τ±r = TDr =12kBT [(γ1 − γ3)2 + (γ2 − γ3)2] , kBT ± E, , (Ohmic) (3.65)The relaxation time τγ is far more complicated in analytic form so instead I graph thefull decay times associated with the frequencies λ˜D, λ˜±(a) (b)Figure 3.10: (a) Pure relaxation times vs kBT ; τ±r =3/Φ(E) (orange/green), τDr = 3/Φ() (blue), (b) Full re-laxation times vs kBT ; τ± = 1/Reλ˜± (orange/green),τD = 1/Reλ˜D (blue). Parameters: γ1 = 0.01, γ2 = 0.05, = 1, ∆10 = 0.1, ∆20 = 0.2,We see that the pure relaxation times corresponding to λ±, λD are different for lowtemperatures. As the temperature is increased however, the two times converge andfollow an inverse function of kBT , independent of , E as indicated in equation 3.65. Thesystem rapidly relaxes to the bath in this limit.For all kBT the dark-state dominates the dynamics as its relaxation time is much lowerthan the other relaxation time in the system. Furthermore, the dark-state relaxation timetends towards zero very rapidly with increasing temperature, whereas τ± tends towardsa constant asymptotic value and the two relaxation times τ± become unresolved in thislimit. This suggests that for high enough temperatures, the dark state is hardly populatedand it is |±〉 that dominates the dynamics. Additionally, in this limit we expect only twodistinct frequencies to be present as Γ± have become unresolved. Of course as we take77this limit further /∆ → ∞, the relaxation times should tend to infinity as we recoverthe Independent-boson-model discussed above, where we have only pure dephasing andno energy relaxation.3.2.6 Decay times in the perturbative regime for super-Ohmicspectral densitiesHere I investigate the relationship of the various decay times in the system versus tem-perature and system-bath coupling for super-Ohmic spectral densities. Once again Iinvestigate the zero T = 0K and high-T limits of the pure decay times but in the super-Ohmic regime now. For the decay times I findτ±r =1(± E)3[(γ1 − γ3)2 + (γ2 − γ3)2] , T = 0K, (super-Ohmic)τDr =13[(γ1 − γ3)2 + (γ2 − γ3)2] , T = 0K, (super-Ohmic)τ±r =12kBT (± E)2[(γ1 − γ3)2 + (γ2 − γ3)2] , kBT ± E, (super-Ohmic)τDr =12kBT2[(γ1 − γ3)2 + (γ2 − γ3)2] , kBT , (super-Ohmic) (3.66)Immediately we notice that the full decay times do not converge and τDr , τ±r remaindistinct for all temperatures in the super-Ohmic regime. We also see that the decaytimes fall off as a power-law in the bare-energies , ± E.78(a) (b)Figure 3.11: (a) Pure relaxation times vs kBT ; τ±r =3/Φ(E) (orange/green), τDr = 3/Φ() (blue), (b) Full re-laxation times vs kBT ; τ± = 1/Reλ˜± (orange/green),τD = 1/Reλ˜D (blue). Parameters: γ1 = 0.01, γ2 = 0.05, = 1, ∆10 = 0.1, ∆20 = 0.2,Once again, the dark-state decay channel is evidently the most sensitive to the bath forall but the lowest temperatures here in the super-Ohmic regime.3.3 Perturbative analysis: optical phonon bath withLorentzian lineshapeI now investigate the 3-site model coupled perturbatively to an optical phonon bath. Thespectral density for a realistic optical phonon bath will contain some line-broadening andhere I choose a Lorentzian lineshape. Normally Lorentzians are difficult to deal withdue to their long tails. However the result of equation 3.54 shows that the phase isdependent only on the spectral density evaluated at the bare poles. Therefore, thedivergent frequency integral of the Lorentzian is avoided in this case and I indulge in theuse of a Lorentzian spectral density here.The spectral function for a Lorentzian lineshape centered around a frequency ω0 andwidth ξ, including our dimensionless parameter is1pi~Jα(ω) = Ξαω20ξξ2 + (ω − ω0)2 (3.67)where the dimensionless coupling parameter for Lorentzian-optical phonons is defined asΞα =ναξpi~ω20(3.68)79Applying my perturbative analysis above, the bath correlators for the 3-poles are calcu-lated to beΦD(0) =8ω30ξ(∆˜210Ξ10 + ∆˜220Ξ20)coth(~β/2)(ξ2 + 2 − 2ω0 + ω20)(ξ2 + 2 + 2ω0 + ω20)Φ±(±E) = ξω20(∆˜210Ξ10 + ∆˜220Ξ20)(coth(~β(±E − )/2)( 1(∓E+ω0+)2+ξ2 − 1(±E+ω0−)2+ξ2)(±E − )2+coth(~β(±E + )/2)(1(±E−ω0+)2+ξ2 − 1(±E+ω0+)2+ξ2)(±E + )2+1(±E + )2 ((±E − ω0 + ) 2 + ξ2) −1(±E − )2 ((±E + ω0 − ) 2 + ξ2)+1(±E − )2 ((∓E + ω0 + ) 2 + ξ2) −1(±E + )2 ((±E + ω0 + ) 2 + ξ2))(3.69)In Figure 3.12 I calculate the three decay times for the optical phonon spectral densityin the perturbative limit. I split the two temperature scales 0 < kBT < 10 and 10 <kBT < 100 to show how the decay times τ± are indistinguishable for low temperatures,and distinct for high temperatures. We see how for low-T, the dark-state decay channeldisplays the longest relaxation time, but falls below τ± for high-T.(a) (b)Figure 3.12: Decay times for 3-site model and opticalphonons. Parameters: ∆10 = 0.1, ∆20 = 0.2, Ξ1 = 0.01,Ξ2 = 0.01, Γ = 0.001, ω0 = 10.In Figure 3.13 I plot the phase-space of the 3-site-optical-phonon model showing how the80system remains in the coherent phase for all temperatures in this perturbative limit.Figure 3.13: Phase space for 3-site model and opti-cal phonons. τ± (orange/green), τD (blue). Parameters:∆10 = 0.1, ∆20 = 0.2, Ξ1 = 0.01, Ξ2 = 0.01, Γ = 0.001,ω0 = 10.For the optical-phonon spectral density, the dark-state is the most robust state for rela-tively low-temperatures. We see that as the temperature increases however the 3 decayrates converge.3.4 The 3-site-independent-boson model and dephas-ingI now make a short comment on the limiting case of zero tunnelling ∆10 = ∆20 = ∆ = 0.In this limit we only have bias energies α as well as the bath HB including the system-bath coupling. This happens to be an exactly solvable model, given the application ofthe well known Independent-Boson model [52]. This amounts to applying the polarontransformation to each site independently. What we see in the dynamics amounts tothe case of ‘pure dephasing’, where only the off-diagonal density matrix elements aredepleted, while the diagonal elements remain unchanged. The Hamiltonian in this caseisH(∆ = 0) = 1 |1〉〈1|+ 2 |2〉〈2|+ 3 |3〉〈3|+∑q3∑α|α〉〈α|λα,q(bq + b†q) +∑qωqb†qbq(3.70)81Since the bath is sensitive to the position of the particle, i.e. which well it occupies, andthe tunnelling between wells has been turned off, the particles position does not changeand there is no energy dissipated to the bath. We can see that the time-evolution of thediagonal populations is accordingly zero, such that ∂t〈Pαα〉 = 0 and∂∂tPαα = −i[Pαα, H(∆ = 0)] = 0, ∂∂tPαβ = −i[Pαβ, H(∆ = 0)] 6= 0 (3.71)whereas the evolution of the off-diagonal populations is non-zero. In the latter case wehave what is known as ‘pure dephasing’. As this is an exactly solvable model, with justthree copies of the well known Independent-boson model, I will only comment on theresults here in the text. I have included a calculation of the dynamics for the Ohmicand super-Ohmic spectral density in Appendix G for reference. The dynamics of thesystem are qualitatively the same for both spectral densities. In the zero-temperature tolow-temperature limit, the off-diagonal density matrix elements decay only algebraicallyin time. For mid-to-high temperatures however, we see an exponential decay.3.4.1 SummaryIn this Chapter I have analysed the 3-site-V system in a number of limiting cases. Firstof all I considered the case of no environmental coupling such that the 3-site-V systemis isolated but for a decaying ground state. This served to demonstrate the phenomenonof population trapping due to the formation of the dark-state. I then included an oscil-lator bath coupled weakly to the diagonal elements of the 3-site-V Hamiltonian. Thiscorresponds to the diagonal system-bath interaction and was treated using perturbationtheory due to the weak nature of the coupling. I investigated this system for three differ-ent spectral densities: Ohmic, super-Ohmic and optical. The resulting relaxation timesin the system were compared and the dark-state relaxation time, associated with theinverse of the pure decay rate, was shown to be the shortest in the system for all spectraldensities. This suggests that for a perturbative system-bath coupling, the dark-state isactually the most sensitive to the effects of the bath and experiences the fastest decay.82Chapter 4Non-perturbative analysis of the3-site system coupled diagonally toan Ohmic oscillator bathI now proceed to treat the oscillator bath coupling to the 3-site-V system non-perturbatively.In Chapter 3, I analysed the 3-site model perturbatively in the system-bath coupling pa-rameter. Here I relax that condition such that the system-bath coupling parameter canbe arbitrarily strong. To do this I return to the analysis framework of the 3-site modeloutlined in Section 3.2, but this time evaluate the influence function non-perturbatively.In this chapter I consider just an Ohmic form for the bath spectral density. However inthe next chapter, where I include non-diagonal couplings as well as diagonal system-bathcouplings for a 2-state-spin-boson model, I will consider both a super-Ohmic form for thespectral density and a line-broadened optical phonon spectral density.4.1 NIBA in the 3-site-boson modelWe learned from the calculation of the bare transition amplitude for the system, that thepath integral formalism for the bare-3-site-system involves a path that returns sequen-tially to the ground state, and for each subsequent flip, has the option of entering eitherthe left or right branch before once again returning to the ground state. The symmetryabout the ground state |0〉 means that we can be certain that with every flip the systemmakes from the ground state, it can either go to state |1〉 or state |2〉 and return. In thedensity matrix formalism, my return state is therefore the sojourn χ00. This is possible83to deduce because I operate within NIBA, which favours occupation of diagonal statesand only allows an excursion in the off-diagonal portion of the density matrix over onetime interval in the path-integral formalism (see section 2.5). Otherwise, the possibleblip-sojourn paths would involve a nearest-neighbour hopping around a three-by-threesquare with all possible combinations of hops. This is a substantially harder combi-natorics problem and essentially comes down to the fact that there would exist statesχ10, χ20, χ12 which would form intermediate symmetric paths in the off-diagonal portionof the density matrix. As I am working in NIBA, which suppresses the time spent by thesystem in the off-diagonal terms of the density matrix, these states are also suppressed.Mathematically they would manifest in the influence functional as terms that go likecos(χ10ξ12) and would physically represent subsequent coherent-superposition-states.With this in mind the influence functional describing this system isF ξ10;ξ20χ11;χ22(t) = fξ10χ11(t) + fξ20χ22(t) + fξ10χ00(t) + fξ20χ00(t),= 4∆210[cos(χ11ξ10Q′j)+ 1]e−(ξ10)2Q′′j + 4∆220[cos(χ22ξ20Q′j)+ 1]e−(ξ20)2Q′′j= 8∆210 cos2(χ11ξ10Q′j/2)e−(ξ10)2Q′′j + 8∆220 cos2(χ22ξ20Q′j/2)e−(ξ20)2Q′′j (4.1)whereQ′(t) =ˆ ∞0dωj(ω)ω2sinωtQ′′(t) =ˆ ∞0dωj(ω)ω2(1− cosωt) coth(β~ω/2) (4.2)The factor 8 in equation 4.1 comes from the summation of various blips and sojourncombinations: χx,x′ and χx′,x, as well as ξx,x′ and ξx′,x. Remember that we’re effectivelysumming the probabilities for the system to be in all possible states. The function j(ω)in this case includes only the density of states of the bath and the system-bath couplingparameter has been absorbed in to the blip and sojourn terms in equation 4.1.The influence functional describes the correlations between the outbound and returnpath, in this case, to return to state |0〉 where the weight to remain in this state is unity.The two terms represent both paths x, x′ propagating through either branch, via the blips84ξ10 and ξ20 respectively1. I’ve simplified the influence functional in this case by settingγ0 = 0, and thus the sojourn χ00 = 0. One can think of this as setting the zero spatialposition at the ground state.At this point it is worth checking that the influence functional reduces to the Caldeira-Leggett SB influence functional in the limit that ∆10 = ∆,∆20 = 0 and γ2 = γ0. In thislimit the influence functional becomesfSB(t2j − t2j−1) = ∆2 cos(γ2Q′j)e−γ2Q′′j (4.3)which is indeed the SB model influence functional. Here I have made the same substi-tutions γ1 = γ/2, γ0 = γ2 = −γ/2 in order to recover the Caldeira-Leggett influencefunctional (see equation 4.33 in [31].4.1.1 Ground-state propagator for the unbiased 3-site-bosonsystemThe simplest propagator to calculate is K00;00(t) which describes the return probabilityto the ground state. In the interest of calculating relaxation times in the system, all therelevant information can be extracted from this propagator as the denominator containsthe pole structure of the system. In the path-integral-Influence-Functional formalism itreads in the time domainK00;00(t) = 1 +∞∑n=1(−i)2nˆ t0dt2nˆ t2n0dt2n−1F ξ10;ξ20χ11;χ22(t2n − t2n−1)× ...ˆ t30dt2ˆ t20dt1Fξ10;ξ20χ11;χ22(t2 − t1) (4.4)where we’ve collected the state couplings ∆10,∆20into the relevant influence functionals.The Laplace transform is1I’ve used the shorthand notation Q′j = Q′2j,2j−1, Q′′j = Q′′2j,2j−185K00;00(λ) =1λ+∞∑n=1(−i)2n 1λn+1(F ξ10;ξ20χ11;χ22(λ))n=1λ+ F ξ10;ξ20χ11;χ22(λ)(4.5)where F ξ10;ξ20χ11;χ22(λ) is the full influence functional that allows flips in both branches.4.1.2 Excited-state propagator for the unbiased 3-site-bosonsystemWhile this thesis will focus on calculating the ground-state propagator, I include theexcited-state propagator here in the interest of potential further work that requires suchpropagators. One must take care in the construction of the excited-state propagator asone must take in to the account the final return ’flip’ to the excited state, as we did insection 3.2, for the bare-3-site Greens function. Here we look at the return probabilityto site |1〉, defined by the propagator K11;11(t). A further complication arises in the3-site-boson propagator calculation as the influence functional for the system containstime-non-local interacting blips. This forces us to separate terms in the propagatorK11;11(λ) =λ+ f ξ20χ22(λ)λ2 + λF ξ10;ξ20χ11;χ22(λ)K22;22(λ) =λ+ f ξ10χ11(λ)λ2 + λF ξ10;ξ20χ11;χ22(λ)(4.6)The inverse Laplace transform of equations 4.6 provides the dynamics of the excitedstates in the system.4.1.3 The influence functional for a biased 3-site systemI previously demonstrated the path integral decomposition for the calculation of the bare-3-site system in Section 3.1.2. For the full system, including the bath, I calculate thedensity matrix; which compares two paths coupled via the influence functional. Whenthe paths x, x′ visit the same state, we see that the bias contribution from the baretransition amplitudes A[x]A∗[x] produces a factor 1. Whereas when x(x′) are, for e.g.86in states |1〉 (|2〉), we pick up factors in the transition amplitudes (ignoring factors dt2from the Suzuki-Trotter decomposition) of exp(−i1t)[exp(i2t)]. I therefore define thebias factors B{ξ} acting on a blip between times t2j − t2j−1 to beBξ10(t2j − t2j−1) = exp[− i1ξ10(t2j − t2j−1)],Bξ20(t2j − t2j−1) = exp[− i2ξ20(t2j − t2j−1)](4.7)The K00;00 propagator starts and ends in the sojourn χ00, where I set the zero pointenergy 0 = 0. Therefore, for the 3-site system in question, the system is constrainedto evolve from that point to either χ11 or χ22 via ξ10 or ξ20. For the final time-intervalt − t2n, the system is constrained to be in state |0〉, so both propagators in the densitymatrix combine to give just 1. I’ve collected the bias factors in to the influence functionalterms according toF˜ ξ10;ξ20χ11;χ22(τj) = 8∆210 cos(i1τj)cos2(χ11ξ10Q′j/2)e−ξ210Q′′j+8∆220 cos(i2τj)cos2(χ22ξ20Q′j/2)e−ξ220Q′′j (4.8)where the time interval is defined as τj = t2j − t2j−1. In order to proceed with theanalysis I must define the form of the spectral densities which then allows me to evaluatethe correlation functions.4.2 Non-perturbative analysis for Ohmic spectral den-sitiesFor Ohmic dissipation we have the following form for the spectral density (discussed inSection 2.3.1)Jα(ω) = ηαωe−ω/ωc (4.9)The dimensionless coupling parameter was defined asγα =ηαpi~(4.10)87however in our formulation of the influence functional above, we have for simplicity, setthe zero-position to x0 = 0.The bath correlators here contain just the bath density of states portion j(ω) of thespectral density. Performing the integrals in Q′(t) and Q′′(t) one gets [31]Q′(t) = arctan(ωct) (4.11)Q′′(t) =12ln(1 + ω2c t2) + ln[~βpitsinh(pit~β)](4.12)The Laplace-transformed influence functional readsfα(λ) = 8ˆ ∞0dt e−λt cos2[γα tan−1(ωct)](1 + ω2c t2)−γα[~βpitsinh(pit~β)]−2γα(4.13)As we’re considering times much greater than ω−1c , i.e. times greater than the character-istic timescale of the bath, the quantity ωct is much greater than one. Therefore applyingthe approximation ωct >> 1 to zeroth order provides the approximations1 + ω2c t2 ≈ ω2c t2, cos2[a arctan(ωct)] ≈ cos2(api) (4.14)such thatfα(λ) ≈ 8(~βpiωc)−2γαcos2(piγα) ˆ ∞0dt e−λt sinh(pit~β)−2γα(4.15)At this point one perform the Laplace transform here using the known integral [126]ˆ ∞0dxe−at sinhb cx =12b+1cB(a2c− b2, b+ 1)(4.16)where B(x, y) is the Beta function, producingfα(λ) = 8(~β2pi)1−2γαω−2γαc cos2(piγα)B(~βλ2pi+ γα, 1− 2γα)(4.17)Turning the Beta function in to Gamma functions usingB(x, y) =Γ(x)Γ(y)Γ(x+ y)(4.18)88finally yieldsfα(λ) = 8ω−2γαc cos2(piγα)Γ[1− 2γα](~β2pi)1−2γα Γ(γα + ~βλ/2pi)Γ(1− γα + ~βλ/2pi) (4.19)I can collect constant terms in 4.19 that contribute only to the renormalisation of thetunneling matrix elements in toCχ11ξ10 = 8ω−2γ1c cos2(piγ1)Γ(1− 2γ1), Cχ22ξ20 = 8ω−2γ2c cos2(piγ2)Γ(1− 2γ2) (4.20)These renormalisation constants renormalise the tunneling matrix elements according∆˜1 =√8ω−2γ1c cos2(piγ1)Γ(1− 2γ1)∆1, ∆˜2 =√8ω−2γ2c cos2(piγ2)Γ(1− 2γ2)∆2 (4.21)Extrapolating to the full biased influence functional F˜ with all terms included we findF˜ ξ10;ξ20χ11;χ22(λ) = ∆˜21µ1−2γ1 Γ(γ1 + ~β(λ+ i)/2pi)Γ(1− γ1 + ~β(λ+ i)/2pi) + ∆˜22µ1−2γ2 Γ(γ2 + ~β(λ+ i2)/2pi)Γ(1− γ2 + ~β(λ+ i2)/2pi)(4.22)where I’ve definedµ =~β2pi(4.23)The next steps in the analysis of the 3-site system coupled to an Ohmic bath requiresapproximations of the Gamma functions, which at this point create polynomials in thepole structure of the propagators up to infinite order as they stand. As the expansionin the coupling constant has already been treated in the perturbative analysis, we’llof course focus here on the different temperature regimes so as to expand the Gammafunctions. However as it turns out, the fully analytic form in the low-T (and zero-T byextension) regime is unfortunately unavailable to us in the case of the 3-site model. Thisis because the pole structure of the propagator takes on an arbitrary order due to thepresence of terms that go like λγ. The work on the original 2-site-spin-boson model wasable to get around this and reduce the pole structure to a tractable form because there89was only one exponent γ which could be factored out. However in our case, there aretwo, coming from the two paths in the system and the pole structure remains intractableanalytically. Of course one can easily solve these polynomials numerically, however thatis beyond the scope of this thesis which is primarily interested in analytic solutions.4.2.1 Validity of NIBA for the 3-site-Ohmic-boson modelI first inspect the regions of validity of NIBA so as to guide my exploration of theparameter space of results. In order to determine the validity of NIBA in this regime Ionce again calculate the ratio F1F1 =〈b〉〈s〉 = limλ→0∂∂λF ξ10;ξ20χ11;χ22(λ) (4.24)to beF1 = − 2µ2−γ1∆˜210Γ(1− γ1/2− iµ1)(γ1 − 2iµ1)2(Γ(γ1/2− iµ1)(2 + ψ0(1− γ1/2− iµ1)− ψ0(γ1/2− iµ1))γ1− 2iµ1(2 + ψ0(1− γ1/2− iµ1)− ψ0(γ1/2− iµ1)))− 2µ2−γ1∆˜210Γ(1− γ1/2 + iµ1)(γ1 + 2iµ1)2(Γ(γ1/2 + iµ1)(2 + ψ0(1− γ1/2 + iµ1)− ψ0(γ1/2 + iµ1))γ1+ 2iµ1(2 + ψ0(1− γ1/2 + iµ1)− ψ0(γ1/2 + iµ1)))− 2µ2−γ2∆˜220Γ(1− γ2/2− iµ2)(γ2 − 2iµ2)2(Γ(γ2/2− iµ2)(2 + ψ0(1− γ2/2− iµ2)− ψ0(γ2/2− iµ2))γ2− 2iµ2(2 + ψ0(1− γ2/2− iµ2)− ψ0(γ2/2− iµ2)))− 2µ2−γ2∆˜220Γ(1− γ2/2 + iµ2)(γ2 + 2iµ2)2(Γ(γ2/2 + iµ2)(2 + ψ0(1− γ2/2 + iµ2)− ψ0(γ2/2 + iµ2))γ2+ 2iµ2(2 + ψ0(1− γ2/2 + iµ2)− ψ0(γ2/2 + iµ2)))(4.25)We see from inspection of equation 4.25 how temperature, tunnelling and couplingstrength affect the NIBA condition F1 1. The parameter µ = ~/2pikBT contains90the temperature dependence. The function for F1 contains a fairly complicated depen-dence on the free parameters T, , γ,∆, due to the presence of the polygamma functions.F1 does however display a simple quadratic dependence in the the renormalised tunnel-ing terms ∆˜ where we see that by increasing the tunnelling strength brings us out of theregime of validity of NIBA. Intuitively this makes sense, as more dominant tunnellingleads to an increase coherence effects and thus counteracts the effect of the environmentto affect decoherence in the system. It’s for this reason that NIBA is indeed valid forsmall values of ∆/ωc. We see more clearly here how NIBA remains a valid approximationfor kBT ∆ as in this limit I minimise F1. Upon inspection of the function govern-ing the renormalisation of the tunneling matrix element (see equation 4.21), we see howNIBA becomes exact for γ = 1/4, verifying a result previously identified by Leggett et al.[31] to be the so called Toulouse limit. This happens because the cos(2piγ) term comingfrom the blip-sojourn interaction part of the influence functional vanishes in this limit.(a) (b)(c) (d)Figure 4.1: F1 vs γ1 figures (a) and (b), F1 vs kBTfigures (c) and (d). Parameters: γ0 = 0, γ2 = 0.2, 1 = 1,2 = 1.5. (a) ∆10 = 0.1,∆20 = 0.2, kBT = 10 (b) ∆10 =0.1,∆20 = 0.2, kBT = 0.1, (c) ∆10 = 0.1,∆20 = 0.2,γ1 = 0.2, (d) ∆10 = 0.1,∆20 = 0.2, γ1 = 0.891We see from Figure 4.1 how F1 is minimised for a range of γ1, predominantly in themid-strong coupling regime. For small values of γ1, as one would expect, F1 is large,and therefore outside the region of validity of NIBA. Here the system exhibits stronglycoherent oscillations, and one must go beyond NIBA in order to account for this. As γ1approaches 1 we also see a rapid divergence in F1 coming from the Γ[1− γ] factor in thetunnelling renormalisation term.Figure 4.2: F1 vs γ2. Parameters: γ0 = 0, γ1 = 0.2,1 = 2 = 1, ∆1 = ∆2 = ∆ = 1 (a) γ2 = 0.2, (b)kBT = 10In Figure 4.2 I inspect the regime where the tunnelling elements are comparable to thebias energies. We see that NIBA is essentially invalid here, as F1 remains much largerthan 1 throughout the range of γ even for high-temperature.4.2.2 Mid-high temperatures in the non-perturbative regimeIn the high temperature regime, kBT ~∆, I can expand the Gamma functions in theinfluence functional in terms of the parameter µλ, where µ = 2pi/kBT . First I expressthe influence functional in the formf ξχ(λ) =(µ1−2γγ + µλ)Γ(1 + γ + µλ)Γ(1− γ + µλ) (4.26)having used the functional relation of the gamma function Γ(1 + z) = zΓ(z) N.B. thisrelationship for the Gamma function takes it from its holomorphic form to its meromor-phic form which therefore allows us to access the pole structure of the propagator [127].Taylor expanding about µλ = 0 up to 2nd order92f ξχ(λ) ≈(µ1−2γγ + µλ)Γ(1 + γ)Γ(1− γ)[1 + µλ(ψ0(1 + γ)− ψ0(1− γ))+(µλ)22(ψ20(1− γ) + ψ20(1 + γ)− 2ψ0(1− γ)ψ0(1 + γ)−ψ1(1− γ) + ψ1(1 + γ))](4.27)where ψ0(z) is the digamma function defined asψ0(z) =Γ′(z)Γ(z), ψm =dmdzmψ0(z) (4.28)Explicitly the full influence functional in this regime is found to beF ξ10;ξ20χ11;χ22(λ) = ∆˜21(µ1−2γ1γ1 + µλ)ν1(1 + µΛ1λ+µ2Θ12λ2)+ ∆˜22(µ1−2γ2γ2 + µλ)ν2(1 + µΛ2λ+µ2Θ22λ2)(4.29)with the biases introduced viaF˜ ξ10;ξ20χ11;χ22(λ) = f˜ξ10χ11(λ+ i) + f˜ ξ10χ11(λ− i) + f˜ ξ20χ22(λ+ i) + f˜ ξ20χ22(λ− i) (4.30)whereν1 =Γ(1 + γ1)Γ(1− γ1) , ν2 =Γ(1 + γ2)Γ(1− γ2) ,Λ1 = ψ0(γ1)− ψ0(1− γ1), Λ2 = ψ0(γ2)− ψ0(1− γ2),Θ1 = ψ20(1− γ1) + ψ20(1 + γ1)− 2ψ0(1− γ1)ψ0(1 + γ1)− ψ1(1− γ1) + ψ1(1 + γ1),Θ2 = ψ20(1− γ2) + ψ20(1 + γ2)− 2ψ0(1− γ2)ψ0(1 + γ2)− ψ1(1− γ2) + ψ1(1 + γ2)(4.31)The pole structure of the propagator takes a cubic form if I restrict myself to the limit1 = 2 = and γ1 = −γ2. So I have for the pole structure of the propagator93λ+ Ω2(µ1−2γγ + µ(λ− i))ν1(1 + µΛ(λ− i) + µ2Θ2(λ− i)2)+Ω2(µ1−2γγ + µ(λ+ i))ν(1 + µΛ(λ+ i) +µ2Θ2(λ+ i)2)= 0 (4.32)where the coefficients of the cubic equationaλ3 + bλ2 + cλ+ d = 0 (4.33)area =γθνΩ2µ3−2γ + 2ΛνΩ2µ3−2γ + 2γµθνΩ2µ4−2γ + µ2(4.34)b =γ2 + 2γΛνΩ2µ2−2γ + 2νΩ2µ2−2γ + θνΩ22µ4−2γ + µ22θνΩ2µ4−2γ + µ2d =2γνΩ2µ1−2γ − γθνΩ22µ3−2γ + 2ΛνΩ22µ3−2γθνΩ2µ4−2γ + µ2(4.35)One can also verify, upon turning off the bath couplings γ1 = γ2 = 0, that we recoverthe bare-3-site frequencies of the density matrix propagator, which is an important con-sistency check.4.2.3 Coherent phase space and relaxation times: Ohmic regimefor mid-high-temperaturesNext I inspect the relaxation times in the system for an Ohmic spectral density in thehigh-T regime. For the case of small tunnelling relative to bias energies, the system is inthe coherent phase. This is characterised by a positive discriminant for all system-bathcoupling values up to a critical value, and so I solve the cubic equation using the solutionsoutlined in Appendix E.94(a) (b)(c) (d)Figure 4.3: Phase space, decay rates, relaxation timesand oscillation frequency vs γ (coherent regime). ΓD, τD(blue), Γ−, τ+ (green), Γ+, τ− (orange). Parameters:∆1 = 0.1, ∆2 = 0.05, = 1, kBT = 10In Figure 4.3 we see the results for the symmetric, high temperature, 3-site-boson modelwith Ohmic spectral density. I find for the parameter range considered here that thesystem remains in a persistent coherent phase across the range of the system-bath di-mensionless coupling parameter. As such I calculate the three oscillation frequencies inthe system and find them to be relatively constant for small system-bath couplings. Forlarge enough couplings they diverge however with the dark state having the smallestoscillation frequency. We also find that the dark state relaxation time dominates sub-stantially over small coupling strengths and for intermediate to high coupling strengthsit remains dominant but to a lesser degree.4.2.4 3-site-Ohmic bath model in FMOI now look at the results of the 3-site-boson model within the context of photosyn-thesis, for an Ohmic spectral density. In Section 4.2 I was able to solve the system95non-perturbatively in the Ohmic regime of the bath in the limit of γ1 = γ2 = γ and1 = 2 = . Therefore I apply these results here to the case of the FMO complex inphotosynthesis. The parameters used for the FMO system were outlined in Chapters 1and 2. It was determined that for an effective 3-site-model in the FMO complex the twotunnelling energies are ∆10 = 2meV and ∆20 = 3meV, the bias energy is = 20meV.Feeding these parameters in to the non-perturbative analysis of Section 4.2 and leavingthe system-bath coupling as a free parameter for now, I present the coherent phase andrelaxation times in Figure 4.4.(a) (b)Figure 4.4: Phase space and relaxation times of FMOcomplex at T=300K in units of femtoseconds with Ohmicbath. FMO parameters: ∆˜10 = 2meV , ∆˜20 = 3meV , = 20meVWe see from figure 4.4 how the 3-level system coupled to an Ohmic bath remains in thecoherent phase even at physiological temperatures for system-bath couplings below thecritical value γ = 1. Incidentally, this is the same critical value for the diagonal couplingwhere the coherent-incoherent phase transition was also observed in the original spin-boson model. For the FMO parameter γ = 0.22, we are therefore in the coherent phase.We see that the dark-state relaxation time dominates and in this case I find relaxationtimes τD ∼ 800fs, τ± ∼ 100fs.4.2.5 SummaryIn this Chapter I have treated the coupling of a 3-state-V system to an external oscillatorbath non-perturbatively. The oscillator bath spectral density was assumed to be of Ohmicform and the system was analysed using the influence functional techniques outlinedin Chapter 2. In this case it was found that the relaxation time associated with the96dark-state, the pure relaxation time with no oscillatory component, was found to bethe longest-lived of the three times in the system. This suggests that for the 3-site-Vsystem coupled non-perturbatively to an Ohmic oscillator bath, the dark-state is actuallythe most robust state and the population trapping mechanism described in Chapter 3,prevalent. The presence of the dark-state serves as a potential hindrance to efficientexciton transfer in the FMO complex. The 3-site-V system was shown in Chapter 3 tohave its applications to the FMO complex either as a configuration of three chromophores,or an effective 3-site model for the full FMO system. These results suggest that the FMOcomplex, with strong system-environment couplings, must avoid these possible trappingstates in order to facilitate efficient exciton transfer. Tuning of chromophore energy levelsleading to a 3-site-V configuration as used in this analysis has been shown to lead to theformation of a relatively long-lived dark-state in the system. Therefore, one can concludefrom this, that the avoidance of tuned chromophore energy levels in the FMO complexis a crucial aspect to facilitating efficient exciton transfer through the system.97Chapter 5A dual-coupling-spin-boson modelfor light-harvesting moleculesIn this chapter I investigate the effect of non-diagonal couplings on the 2-state-spin-boson(SB) model for both super-ohmic and optical spectral densities. The new model is there-fore dubbed the Dual-coupling-spin-boson (DCSB) model. This means that I include notonly the coupling of the bath to the diagonal elements of the central system Hamilto-nian (σz), as is present in the original spin-boson model, but also include couplings tothe off-diagonal elements (σx) too. The heuristic justification for excluding non-diagonalcouplings in the SB model normally relies on the assumption that small tunnelling raterelative to the on-site energy: ∆ , renders the non-diagonal coupling negligible [31].Since any interaction to σx should be proportional to the wavefunction overlap of thetwo wells–and is therefore proportional to ∆–the non-diagonal coupling is assumed to becomparatively small.The problem with this assumption in the context of the SB model is that, while thetunnelling energy is indeed required to be small relative to the bias energy for NIBA toapply, it is not negligibly small. In fact, the SB model is considered such a powerful modelbecause it retains coherence effects up to a degree, and this is because of a small butappreciable tunnelling energy. A perturbative tunnelling energy would actually justify aFermi’s-Golden-Rule approach as we saw in Section 1.4.2, which permits only a one waytransfer of energy to the bath and removes any coherence effects. All of this suggeststhat the SB model is somewhat incomplete without non-diagonal couplings.Indeed there are a number of problems where non-diagonal couplings are warranted.For example, tunnelling defects in solids interact with their surrounding atoms not only98by a fluctuation of their on-site energies (diagonal coupling) but can emit and absorbphonons (non-diagonal coupling) corresponding to atomic vibrations that modulate thedistance between defects. This latter effect should lead to an increase in the well knownphenomenon of phonon-assisted tunnelling [128, 99].Non-diagonal couplings were first discussed by Holstein (1959) [47] in the context of his‘small’ polaron model. The first application of these couplings was in the SSH model forpolyacetylene [129]. Silbey and Munn [130, 131, 132] later used a canonical transforma-tion method along with numerical approximations to investigate non-diagonal couplingsin molecular crystals. In their three part series on electron-phonon coupling, they foundthat the diffusion coefficient of charge carriers in molecular crystals is dominated by thecompetition of the two forms of coupling. Since then non-diagonal couplings have seenapplication to the SB model treated with a two bath approach [55] as well as numericalapproaches [133, 134]. Recently we’ve seen that the inclusion of both diagonal and non-diagonal couplings in the SSH model leads to a sharp transition between the behaviourat weak and strong couplings [50, 79].5.1 The dual-coupling polaron transformationFor the general case of a bath coupled to the diagonal and non-diagonal elements of the2-state tunnelling system we haveH = ({bq})σz2+ ∆({bq})σx2+HB (5.1)where HB =∑q ωqb†qbq is the kinetic energy of the bosonic bath. The Pauli spin matriceshere are defined as σz = |1〉〈1|− |2〉〈2| and σx = |1〉〈2|+ |2〉〈1|. Both the on-site energiesand transition matrix elements are modulated by the bath in this case. To linear orderin the bath couplings [50] (see Section 1.3) we then haveH =2σz +∆2σx +σx2∑qλx,q(bq + b†q) +σz2∑qλz,q(bq + b†q) +HB (5.2)We can decouple the bath from the central system with the appropriate unitary trans-formation, which in this case involves a product of the shift operators for the diagonaland non-diagonal couplings respectivelyU = e−S, S =σx2∑qux,q(bq − b†q) +σz2∑quz,q(bq − b†q) (5.3)99where ux,q = λx,q/ωq and uz,q = λz,q/ωq. Performing the transformation H → UHUT ≡H˜ using the Baker-Campbell-Hausdorff formulaH˜ = H − [S,H] + 12![S, [S,H]]− 13![S, [S, [S,H]]] +14![S, [S, [S, [S,H]]]] + ... (5.4)I find for the the transformed HamiltonianH˜ = ˆσz + Kˆ−σ+ + Kˆ+σ− +HB (5.5)where I have introduced the operator notation ˆ to reflect the fact that the on-site energy,in the transformed Hamiltonian picture, is fluctuating due to the non-diagonal coupling.The on-site energy is nowˆ =∆2+2cosh(φˆx)− ∆2cosh(φˆz)(5.6)The kernels are defined asKˆ± =2+∆2B±z −2B±x (5.7)whereφˆx =∑qux,q(bq − b†q), φˆz =∑quz,q(bq − b†q) (5.8)I identify the boson shift operators [52, 61, 99]Bx± = exp(±φˆx), Bz± = exp(±φˆz)(5.9)One can also verify that the transformed Hamiltonian reduces to the 2-state tunnellingHamiltonian with the system-bath coupling turned off.5.2 Dynamics of the dual-coupling-spin-boson modelThe time evolution of the system is governed by the reduced density matrix ρs(t). Theprobabilities for each state are therefore 〈σf | ρs(t) |σ′f〉 ≡ ρs(σf , σ′f ; t), and can be de-scribed in the Feynman-Vernon path integral formalism as100ρs(σf , σ′f ; t) =∑σi,σ′iˆ σfσiDσ(τ)ˆ σ′fσ′iDσ′(τ ′)eiS0[σ(τ)]−iS0[σ′(τ ′)]F [σ(τ), σ′(τ ′)]ρs(σi, σ′i; 0)(5.10)where σf = σ(τ = t), σ′f = σ′(τ ′ = t), σi = σ(τ = 0), σ′i = σ′(τ ′ = 0) are paths thatvisit the two states in the central system σ, σ′ ∈ {|1〉 , |2〉}. The paths σ, σ′ in the centralsystem are coupled via the influence functional F [σ, σ′], which incorporates all of theeffects of the bath. Since the bath has been decoupled from the central system via thedual-coupling polaron transformation, the transformed Hamiltonian H˜ of equation 5.5corresponds to a transformed action for the central system S˜0. Therefore the reduceddensity matrix for the transformed actionρs(σf , σ′f ; t) =〈∑σi,σ′iˆ σfσiDσ(τ)ˆ σ′fσ′iDσ′(τ ′)eiS˜0[σ(τ)]−iS˜0[σ′(τ ′)]ρs(σi, σ′i; 0)〉B(5.11)contains the effects of the bath in the transformed probability amplitudes, and the traceover the bath degrees of freedom leaves the dynamics in terms of central system variablesonly. The transformed probability amplitudes areA˜σf ,σi(t) ≡ˆ σfσiDσ(τ)eiS˜0[σ(τ)] = 〈σf | e−iH˜t/~ |σi〉 (5.12)The effects of the bath enter via the transformed on-site energies±/2 |±〉〈±| → ±ˆ |±〉〈±|,and transition matrix elements ∆σ± → K∓σ±, which are now operators containing thebath shift operators Bˆ±. In this formalism, the probability amplitudes can be calculatedindependently, and the trace over the bath degrees of freedom performed to obtain aclosed expression for the central system dynamics. Expanding each transition amplitudefor the density matrix, following the analysis outlined in Section 3.1.2 givesρs(1, 1; t) =〈(1−ˆ t0dt2ˆ t20dt1e−2i˜(t2−t1)K+(t2)K−(t1) + ...)×(1−ˆ t0dt2ˆ t20dt1e2i˜(t2−t1)K−(t2)K+(t1) + ...)〉B(5.13)where the density matrix element corresponding to the probability to start and end instate |1〉 was chosen. For clarity the ground state energy has been shifted to state |1〉,101such that 1 = 0, and the amplitude for no transitions A˜(0)11 (t)=1. Therefore the energyof state |2〉 is 2 = −2ˆ. Collecting termsρs(1, 1; t) = 1−ˆ t0dt2ˆ t20dt1Υ˜(t1, t2)+ˆ t0dt4ˆ t40dt3ˆ t30dt2ˆ t20dt1Υ˜(t1, t2)Υ˜(t3, t4) + ... (5.14)The kernel acting at times t, t′ is thereforeΥ˜(t, t′) =〈ei2ˆ(t−t′)/~K+(t)K−(t′)〉B+〈e−i2ˆ(t−t′)/~K−(t)K+(t′)〉B(5.15)where NIBA has been invoked such that〈K+(t1)K−(t2)K+(t3)K−(t4)〉B≡〈K+(t1)K−(t2)〉B〈K+(t3)K−(t4)〉B(5.16)and only correlations between successive bath fluctuations are permitted. This approxi-mation is valid when the system spends vastly more time in diagonal states of the densitymatrix compared to off-diagonal states (see Section 2.5). Separating the transformed on-site energy in to the sum of its average value ˜ ≡ 〈ˆ〉, and fluctuations about its meanvalue δˆ, leads to ˆ ≡ ˜+ δˆ. Now the Kernel readsΥ˜(t, t′) = ei2˜(t−t′)/~〈ei2δˆ(t−t′)/~K+(t)K−(t′)〉B+ e−i2˜(t−t′)/~〈e−i2δˆ(t−t′)/~K−(t)K+(t′)〉B(5.17)The system is solved by taking the Laplace transform of the density matrixρs(1, 1;λ) =1λ−L[ˆ t0dt2ˆ t20dt1Υ˜(t1, t2)]+L[ ˆ t0dt4ˆ t40dt3ˆ t30dt2ˆ t20dt1Υ˜(t1, t2)Υ˜(t3, t4)]+ ... (5.18)and the Laplace-transform of the Kernel is therefore102Υ˜(λ) =ˆ ∞0dte−λt[ei2˜t/~〈ei2δˆt/~K+(t)K−(t)〉B+ e−i2˜t/~〈e−i2δˆt/~K−(t)K+(t)〉B](5.19)The fluctuations of the on-site energy δˆ come from the shift operators δB±. Theyintroduce terms in the integrand that go like e−iδB±t/~. These terms become rapidlyoscillating functions, and by the stationary phase approximation, integrate to zero forlarge arguments of the phase. Therefore I can approximate this integral asΥ˜(λ) =ˆ ∞0dte−λt[ei2˜t/~〈K+(t)K−(t)〉B+ e−i2˜t/~〈K−(t)K+(t)〉B](5.20)The functions〈K+(t)K−(t′)〉Bdescribe the interactions of successive flips within thesystem - dragging their respective boson clouds with them, and are multiplied by a pre-factor corresponding to the kinetic energy term of the Hamiltonian. So I retain only theaverage of the fluctuating on-site energy˜ =∆2+2〈cosh(φx)〉B− ∆2〈cosh(φz)〉B=∆2+2〈Bx〉B − ∆2〈Bz〉B (5.21)Laplace transforming ρS(1, 1; t) leads to the geometric series summationρS(1, 1;λ) =∞∑n=0(−1)2n Υ˜(λ)nλn+1=1λ+ Υ˜(λ)(5.22)The Kernel is expanded toΥ˜(λ) =12ˆ ∞0dte−λt{cos(2˜t)[2(1− 2Bx)+ 2∆Bz]+ 2B2xfxx(t) + ∆2B2zfzz(t)− 2∆BxBzfxz(t)}(5.23)Paying close attention to the distinction here between diagonal and non-diagonal cou-plings, the bath correlation functions are calculated following Mahan [52] and are foundto be103〈Bα±(t)Bβ∓(0)〉B = 〈Bα±〉〈Bβ∓〉eϕαβ(t) (5.24)where the indices α, β ∈ {x, z}, and 〈Bα±(t)Bβ∓(0)〉B = 〈Bα∓(t)Bβ±(0)〉†B. The phaseϕαβ(t) = iQ′αβ(t) +Q′′αβ(t) (5.25)is comprised of the functionsQ′αβ(t) =ˆ ∞0dω√Jα(ω)Jβ(ω)~ω2sin(ωt)Q′′αβ(t) =ˆ ∞0dω√Jα(ω)Jβ(ω)~ω2cos(ωt) coth(~βω/2) (5.26)The biased influence functionals, within NIBA, are thereforefαβ(t) = BαBβ cos[2˜t+Q′αβ(t)]eQ′′αβ(t) (5.27)The spectral density function Jα(ω) = (pi/2)∑q(λ2α,q/ωq)δ(ω − ωq) describes the distri-bution of bath modes weighted by their coupling to the exciton. The terms Bx ≡ 〈Bx±〉B,Bz ≡ 〈Bz±〉B are the Debye-Waller factors, which can be cast in the continuum limit as[61]〈Bα±〉B = exp[−ˆ ∞0dωJα(ω)ω2coth(βω/2)](5.28)We see that the Debye-Waller factor has an infra-red divergence for Ohmic spectraldensities leading to an orthogonality catastrophe [61], however for super-Ohmic spectraldensities the integral is infrared convergent, reflecting the probability for the system totunnel between wells without exciting the bath.It’s prudent to check that when the non-diagonal coupling ζx → 0 is turned off, werecover the spin-boson-model influence function result with just diagonal coupling ζz.One can verify that this indeed the case. One can also verify that by turning off the bathentirely–both diagonal and non-diagonal couplings, we recover the eigenvalues of the 2-site tunnelling Hamiltonian. In this case the pole structure for the propagator reducesto104λ2 + 2 + ∆2 = 0,→ λ = ±i√2 + ∆2 (5.29)which is the pole structure for the propagator of the isolated 2-site system.5.3 Super-Ohmic spectral densities in the DCSB modelI now turn to the case of a bath characterised by a super-Ohmic spectral density. In thiscase the spectral density takes the formJx(ω) = ρxω3ω2phe−ω/ωc , Jz(ω) = ρzω3ω2phe−ω/ωc (5.30)where in section 2.3.2 I defined a dimensionless coupling constant for acoustic phononsasζα =ραpi~(5.31)The first step in evaluating the super-Ohmic influence functional is to do the frequencyintegrals in Q′(t), Q′′(t) given in their general form in equation 4.2. For the super-Ohmicphonon density of states they take the formQ′ph(t) =ˆ ∞0dω ωe−ω/ωc sin(ωt)Q′′ph(t) =ˆ ∞0dω ωe−ω/ωc(1− cos(ωt)) coth(~βω/2) (5.32)which appear in the individual influence functionalsf˜α(t) =12cos(2˜t+ ζαQ′ph(t)/ω2ph)exp(− ζαQ′′ph(t)/ω2ph) (5.33)and the indices form the set α ∈ {z, x, xz}. Next I split up the second bath correlator into its time-independent and time-dependent partsQ′′ph(t) = Q′′0 − Q˜′′ph(t) (5.34)The time-independent term is therefore105Q′′0 =ˆ ∞0dω ωe−ω/ωc coth(~βω/2) (5.35)which leads to the adiabatic phonon dressing of the tunnelling matrix element [34]. Theexponentiation of this time-independent term can be factorised out from the full influ-ence functional and the corresponding factor that is absorbed into the tunnelling matrixelement is known as the Debye-Waller factor [99]. This term serves to renormalise pa-rameters due to thermal excitations of the polaron cloud [34]. The integral I calculate tobeQ′′0 = −ω2c + 2(kBT~)2ψ′(kBT~ωc)(5.36)Expanding to 2nd-order in kBT ~ωc 1 producesQ′′0 ∼ ω2c +13(pikBT~)2+ 2ψ′′(1)ωc(kBT~)3+pi4(kBT )415~4ω2c(5.37)I absorb these terms into the renormalisation factors defined earlierBz = e−θzz , Bx = e−θxx , Bxz = e−θxz (5.38)whereθzz = ζzQ′′0/ω2ph, θxx = ζxQ′′0/ω2ph, θxz =√ζxζzQ′′0/ω2ph (5.39)The time-dependent part of the influence phase left over isϕ(t) = iQ′ph(t) + Q˜′′ph(t) (5.40)Defining a complex time-variableτ = t− i~β/2 (5.41)the phase can be cast in the convenient form for the phonon-phase in the complex time-domainf˜αβ(τ) = Re{exp[ϕαβ(τ)/ω2c + 2i˜(τ + i~β/2)]}(5.42)Computing the frequency integral using known integrals (Section 4.13 in [126]) yields106ϕ(τ) =1(~β)2[ψ′(12+1− iωcτ~ωcβ)+ ψ′(12+1 + iωcτ~ωcβ)](5.43)where the polygamma functions are definedψn(z) =dndznψ0(z), ψ0(z) =dndznln[Γ(z)](5.44)Analytic continuation of the phase back to the real axis yieldsϕαβ(t) =√ζαζβ(~β)2[ψ′(1− iωct~ωcβ)+ ψ′(1 + ~ωcβ + iωct~ωcβ)](5.45)As it stands, the Laplace transforms of the individual influence functionals in Equa-tion 5.23, are intractable. In the interest of exploring physiological temperatures, wherekBT = 27meV (T = 300K), and phonon cut-off frequencies ωc = 8.7meV [98], we canexpand in the parameter ~ωcβ 1. Asymptotic expansion of Equation 5.45 to thirdorder yieldsϕαβ(t) ≈ Θαβ + iφαβt− Φαβt2Θαβ =√ζαζβ(2ωcβ+ω3cβ3), φαβ = 2√ζαζβω3c ,Φαβ =2√ζαζβω3cβ(5.46)It now remains to evaluate the Laplace transform of the propagator.5.4 High-temperature limit in the super-Ohmic NDSBmodelAs they stand, the influence functionals produce a pole structure in the propagator ofEquation 5.22 of infinite order in λ, representing all possible excited states of the systemacross the whole temperature range. In the interest of keeping the system-bath couplingsto all orders, I can investigate certain temperature limits. In the high-temperature limitI can expand the influence phase in the small quantity λβ 1. the biased influencefunctional can be Taylor expanded to 1st-order as107f˜αβ = f˜0αβ − f˜ 1αβλ (5.47)where the coefficients aref˜ 0αβ = eΘαβˆ ∞0dt cos(2˜t+ φαβt)e−Φαβt2f˜ 1αβ = eΘαβˆ ∞0dt t cos(2˜t+ φαβt)e−Φαβt2 (5.48)The time interval in f˜ 0αβ can be evaluated with the known integral BI (263)(5) in [126]to bef˜ 0αβ =eΘαβ2√piΦαβexp[− (2˜+ φαβ)24Φαβ](5.49)The first order coefficient f˜ 1αβ can be evaluated using the known integral BI (362)(2) in[126] and the Maclaurin series of the Dawson function [135]D(x) =∑∞n=0(−2)nx2n+1/(2n+1)!!, to getf˜ 1αβ =eΘαβ2Φαβ− eΘαβ2Φ3/2αβ(2˜+ φαβ)D(2˜+ φαβ2√Φαβ)(5.50)The pole structure of the propagator now takes the formλ+(2(1− 2Bx)+ 2∆Bz)λ2(λ2 + 4˜2)+ Υ0 −Υ1λ = 0 (5.51)whereΥ0 = ∆2B2zf˜ 0z + 2B2xf˜ 0x − 2∆BxBzf˜ 0xzΥ1 = ∆2B2zf˜ 1z + 2B2xf˜ 1x − 2∆BxBzf˜ 1xz (5.52)Rearranging the pole structure to1082(1−Υ1)λ3 + 2λ2Υ0 + λ(22(1− 2Bx)+ 4∆Bz + 8˜2 − 4˜2Υ1)+ 8˜2Υ0 = 0 (5.53)The solutions of the cubic equation are obtained using the analysis framework of Ap-pendix E. For the poles I findλr = Γγ − Γθλφ = −Γγ2− Γθ ± iω (5.54)whereΓθ =Υ03(1−Υ1) , Γγ = u− v, ω =√32(u+ v)u = 3√√D − q2, ∀D > 0v = 3√√D +q2, ∀D > 0q =2(Υ0/(1−Υ1))3 − 9(Υ0Ξ)/(4(1−Υ1)2) + 108˜2Υ0/(4(1−Υ1))27p =3Ξ/(4(1−Υ1)−Υ20/(1−Υ1)23Ξ = 2(1− 2Bx)+ 2∆Bz + 4˜2 − 4˜2Υ1 (5.55)The entirely real pole λr ≡ Γr represents the pure relaxation rate. This pole leads to aterm describing exponential relaxation in the time regime upon Laplace inversion. Theinverse of Γr represents the inverse time scale for relaxation of the system to its groundstate, and is therefore interpreted as the exciton transfer time through the sysytem. Thereal part of the other two poles Γφ ≡ Re[λφ] is the decoherence rate describing the inversetimescale for the loss of phase coherence in the system. It manifests as the damping rateof the oscillatory terms in the dynamics upon Laplace-inversion [61].D =(q2)2+(p3)3(5.56)109(a) (b)(c) (d)Figure 5.1: Phase space, decay rates and oscillationfrequency for Non-diagonal-SB model vs ζx. ζz = 1. Pa-rameters: = 1, ∆ = 0.3, kBT = 5, ωc = 0.1. In units of.In Figure 5.1 the results for the DCSB model are presented as a function of non-diagonalcoupling strength ζx. The discriminant D remains positive, therefore the system remainsin the coherent phase for 0 < ζx ζz. However this just means that the exciton hassome oscillatory component to its dynamics. The duration of these oscillations can beseen in Figure c. The timescale associated with the decay of the exciton oscillations isτφ ≡ τ±, and we see that it remains longer than the exciton transfer time τr for smallvalues of ζx. In the intermediate regime of ζx, for the range considered here, the coherencetime falls below the exciton transfer time. For larger values of ζx we see the coherencetime is once again longer than τr. In Figure d the effect of phonon-assisted transport(normally associated with the oscillation frequency) can be seen by continuously varyingthe non-diagonal coupling parameter γx, and the dimer oscillation frequency rises. Forhigh enough values of γx, however, we see ω reduce rapidly due to an ‘over dressing’ of theparticle’s phonon cloud [52] i.e it becomes more difficult for the particle to tunnel betweenwells as it gets heavier. Therefore, while the non-diagonal coupling serves as an additional110decoherence mechanism in the system, it suppresses coherence at a slower rate than theexciton transfer time. This means that the exciton moves rapidly through the dimersystem while remaining coherent for the duration of its transfer between chromophores.5.4.1 DCSB model for the FMO system with acoustic phononsI now apply the DCSB model with super-Ohmic spectral densities to a dimer of twochromophores in the FMO complex. I use the excitation energies of the chromophoresand the inter-chromophore tunnelling energies determined in Section 1.5. The strengthof the diagonal system-environment coupling was discussed in Section 2.3.2.Figure 5.2: Exciton transfer time τr, coherence timeτφ, and dimer oscillation frequency ω as a function ofnon-diagonal coupling strength ζx. Parameters: =20meV, ζz + ζx = 1, kBT = 27meV, ωc = 8.7meV.Bath spectral density of super-Ohmic form: J(x,z)(ω) =pi~ ζ(x,z)(ω3/ω2c )exp(−ω/ωc).In Figure 5.2 we see how varying the system-bath coupling strength affects the excitontransfer time, coherence time and oscillation frequency of the dimer with both diagonaland non-diagonal couplings. We observe the increase in dimer oscillation frequency asζx is varied for small ζx. Along with this increase in oscillation frequency, the transfertime is seen to also increase. This can be interpreted as coherent tunnelling back and111forth within the dimer, as the density of states available to the exciton increases withphonon-assisted tunnelling. For large enough non-diagonal coupling strength, the trans-fer time is seen to rapidly decrease as the exciton tunnelling becomes one-way. At thesame time the exciton coherence time is seen to decrease with the additional decoherencemechanism provided by the non-diagonal coupling. However, the coherence time is seento remain persistently longer than the transfer time, indicating the exciton is coherentfor the duration of its motion.The sum of diagonal and non-diagonal coupling strengths used for these calculationswas ζz + ζx = 1. This number can be determined from the experimentally determinedreorganisation energy of a super-Ohmic spectral density bath model for the FMO com-plex. This was discussed in Section 2.3. The non-diagonal coupling strength for the FMOcomplex has not been experimentally measured yet and is therefore kept a free parameterin the results.5.5 The DCSB model coupled to optical phononsIn this section I model the environment as an oscillator bath comprised of optical phononswith frequency ω0. We saw in Section 2.3.4, how the experimentally determined spec-tra for biomolecules actually contain structure beyond a continuous distribution of low-frequency modes. Sharp optical transition peaks were observed revealing the presenceof discrete oscillation modes present in the environment of the biomolecule. Therefore,I would like to investigate the effect of including optical phonon modes in the spectraldensity of the DCSB model. Considering an optical phonon spectral density with spectralbroadening of, the spectral density function takes the formJα(ω) = λαe−(ω−ω0)2/2ξ2 (5.57)where ω0 is the optical phonon frequency, λx, λz the non-diagonal and diagonal couplingenergies respectively and ξ the peak width. The dimensionless coupling constant foroptical phonons is να = λαξ/piω20. We see that the optical phonon coupling scales linearlywith the spectral weight λξ and with the inverse of the optical phonon frequency squared.Therefore optical phonon peaks at high frequency with relatively small spectral weightcouple weakly to the exciton. In the interest of studying non-perturbative effects wetherefore consider relatively low-frequency optical phonons.The influence phase in complex time is now112ϕαβ(τ) =√νανβˆ ∞0dωω2e−(ω−ω0)2/2ξ2 cos(ωτ)csch(~βω/2) (5.58)The frequency integral can be evaluated with the use of the saddle-point approxima-tion since for light-harvesting molecules, the optical phonon peaks have a very narrowlinewidth ξ. This means that the prefactor to the argument of the Gaussian exponentialspectral function x = 1/2ξ2 is very large, and the peak therefore very narrow around thepoint ω0. We find for the Gaussian-optical phonon influence phaseϕαβ(τ) =ξ√2piλαλβpi~ω20cos(ω0τ)csch(~βω0/2) (5.59)and for the Debye-Waller factorsBα = exp(−να√2pi coth ~βω0/2)(5.60)Expanding the influence phase in the limit ω0τ 1 produces a Gaussian integralϕαβ(τ) ≈ξ√2piλαλβpi~ω20(1− ω20τ 2)csch(~βω0/2) (5.61)Analytically continuing ϕ(τ) to the real time axis and evaluating the first two terms ofϕ(λ), again expanded in λβ 1, yieldsf˜ 0αβ =eΛαβ2√piΛαβω20exp[−(2˜+ Λαβω20β)24Λαβω20]f˜ 1αβ =eΛαβ2(Λαβω20)3/2[√Λαβω20− (Λαβω20~β + 2˜)D(2˜+ Λαβω20~β2√Λαβω20)](5.62)where Λαβ =√2piνανβcsch(~βω0/2). The poles of the density matrix propagator withoptical phonons is then given by substituting Equations 5.62 and 5.60 in to the densitymatrix of Equation 5.22.113(a) (b)(c)Figure 5.3: (a) ∆˜ = 0.3, = 1, kBT = 1, ωo = 0.4,ξ = 0.01, λz = 5 in units of (meV).In Figure 5.3, the results for the DCSB model are presented in the low temperature regimekBT = 0.1. The discriminant D remains positive and the system is therefore persistentlyin the coherent phase. As the non-diagonal coupling energy λx is turned, the oscilla-tion frequency of the dimer is seen to increase in accordance with the phonon-assistedtransport mechanism. As λx is increased further however, the oscillation frequency dropsoff rapidly. The dimer relaxation rate Γr is seen to increase slightly for increasing λx,representing the inverse timescale for the dimer relaxation to its ground state.5.5.1 DCSB model for the FMO system with optical phononsIn this section I apply the DCSB model with an optical phonon bath to the FMO com-plex. It was determined in Section 2.3.4, for certain peaks selected from FMO spectraldata, that the corresponding optical phonon dimensionless coupling parameter was smallenough to permit a perturbative calculation of the system-bath interaction. This wasdone in Section 3.3. Here I use various optical phonon peaks at 6,8 and 10 meV, with114linewidth 0.5meV, and spectral weight 15meV to investigate the non-perturbative effectof optical phonons on the FMO system. This optical phonon frequency falls within thebandwidth of the FMO dimer and is expected to greatly affect the dynamics. This is incontrast to the perturbative results for very large optical phonon frequency relative toFMO dimer parameters.Figure 5.4: Exciton transfer time τr (blue line), co-herence time τφ (orange line), and dimer oscillation fre-quency ω (purple line) as a function of non-diagonal cou-pling strength νx. Parameters: = 20meV, ∆ = 6meV,ζz=1, kBT = 27meV, ω0 = 6meV. Bath spectral densityof Gaussian-optical form: J(x,z)(ω) = λ(x,z)exp(−(ω −ω0)2/2ξ2).In Figure 5.4 we see the effect of turning on the non-diagonal coupling on a dimer withFMO parameters. We see qualitatively similar results to the FMO system coupled to anacoustic phonon bath. The dimer oscillation frequency is once again seen to increase inaccordance with the phonon-assisted transport mechanism and drop off for large enoughvalues of the non-diagonal coupling energy. For small values of νx the exciton time is seento increase, alongside a an increase in the dimer oscillation frequency. As the non-diagonalcoupling energy is increased further, both the transfer time and coherence time, decreaserapidly as the additional coupling mechanism begins to decrease coherence effects and aidin transfer rate of the exciton through the dimer system. However, despite the additionaldecoherence observed with non-diagonal coupling, once again the coherence time remainspersistently longer than the transfer. Therefore we can conclude that the exciton remainscoherent for the duration of its motion.1155.5.2 SummaryIn this chapter I have analysed a two-state central system coupled both diagonally andnon-diagonally to an oscillator bath. I dubbed this model the Dual-coupling spin-boson(DCSB) model. I used both a super-Ohmic and optical spectral density for the bath,with the optical spectral density having a Gaussian form for the spectral broadening.The application of the DCSB model to photosynthesis was motivated in Chapters 2 and3 by the presence of delocalised excitons across two chromophores in the FMO complex.Oscillator baths consisting of acoustic phonons motivated the use of a super-Ohmic spec-tral density and optical phonons motivated the Gaussian-optical spectral density.For a super-Ohmic bath, the DCSB model found the exciton to remain in the coherentphase as the non-diagonal coupling was turned on and varied. The exciton transfer timethrough the dimer was found to decrease rapidly for larger non-diagonal couplings as didthe exciton coherence. Oscillations were found to persist in the exciton dynamics for therange of non-diagonal coupling strengths, with the coherence time found to be longerthan the exciton transfer time for the full range of non-diagonal couplings. Therefore thenon-diagonal coupling was found to not only facilitate exciton transfer through the dimersystem, but produce coherent exciton transfer for relatively large non-diagonal couplingstrengths. Similar results were found for a Gaussian-optical spectral density. Both theexciton transfer time and coherence time were found to decrease rapidly for significantnon-diagonal coupling strengths, with the coherence time again remaining longer thanthe transfer time.Both of these results can be interpreted within the context of the inelastic-phonon-assisted tunnelling mechanism. By turning on the non-diagonal coupling to the bath,inelastic processes are introduced to the dynamics. Excitons can tunnel between chro-mophores by emission and absorption of phonons from the bath. This produces anincreased density of final states available to the exciton, and therefore facilitates exci-ton tunnelling. As the coupling strength becomes large enough, the tunnelling becomesone-way as the exciton becomes heavily dressed with its phonon cloud.116Chapter 6Summary, conclusions and futureworkI now summarise the work presented in this thesis, and draw conclusions pertaining tothe results of the 3-site-boson model, and the dual-coupling-spin-boson (DCSB) model.I also conclude the investigation of applying both models to the photosynthesis mecha-nism. Finally I suggest a number of avenues of further research.In Chapter 3 I introduced some of the key aspects of the 3-site model isolated fromany environment, as well as some of the limiting and perturbative treatments of the 3-sitesystem interacting with a harmonic oscillator bath. First, I introduced the concept ofpopulation trapping in a 3-site-V configuration. I demonstrated that the addition of anarbitrarily strong decay mode out of the lower ground state still leads to a long-timepopulation present in the dynamics. This ‘trapping’ effect can be attributed to the pres-ence of an eigenstate, known as the ‘dark state’, that overlaps only with the upper twolevels despite there being no tunnelling term between these levels in the 3-site-V config-uration. I demonstrated how, by tuning the on-site energies of the upper two levels, thisdark-state is accessed and contains no overlap with the decaying ground-state, leading tothe long-time population trapping. To do this I established a path-integral formalism forthe ‘bare’ 3-site system, not only to calculate the Greens functions for the bare system,but also to set up the analysis for subsequent chapters where I treat the system-bathcoupling non-perturbatively. In the interest of exploring the 3-site system coupled tosome environment, I first introduced a perturbative system-bath coupling to an externaloscillator bath.I then solved the 3-site-V system coupled perturbatively to a harmonic oscillator bath117for 3-different spectral densities: Ohmic, super-Ohmic and optical phonons. I appliedthe noninteracting-blip-approximation (NIBA) to truncate the number of coherent-bathinteraction processes corresponding to higher-order coherence effects. The region of va-lidity of this approximation was not just discussed qualitatively—in the introductionchapter—but quantified. The quantity F1, describing the ratio of average time spent inan off-diagonal state to a diagonal state of the density matrix, was shown to be minimisedin the perturbative regime provided the tunnelling energies were small compared to theon-site energies. With this condition met, NIBA was shown to be valid in this regime.My final results for this model involved calculations of the three frequencies present inthe system and an exploration of the coherent-incoherent phase space. We saw thatthe system operates only in the coherent phase across the full temperature range (N.B.I don’t vary the system-bath coupling parameter here as it represents a perturbativequantity in this model). I found that the relaxation time associated with the dark-statewas the shortest of the three in both the Ohmic and super-Ohmic regime. Therefore, inthe perturbative limit at least, I show that the dark-state is the most sensitive state inthe system to the effects of the bath, and coherence effects attributed to this state arestrongly suppressed. However, the opposite is found to be true for the case of an opticalphonon bath for low-temperatures. Here, the dark-state dominates the relaxation timesand therefore is the most robust state in terms of exhibiting long-time coherence effects.For mid-high temperatures the decay rates are shown to converge and the correspondingrelaxation times in this limit are comparable.In Chapter 4 I moved beyond the perturbative analysis of the 3-site-boson model,incorporating the effects of the bath to all orders in the system-bath coupling parameter.Building on the path integral formalism for the 3-site propagator set up in Chapter 3, Iused an influence functional approach to model the system-bath interactions. NIBA wasonce again employed to simplify the number of coherent system-bath processes, howeverthe influence functionals still required expansions in small parameters in order to obtaintractable pole structures for the propagators. Therefore I chose to inspect the high-temperature limit in the interest of modelling physiological temperatures T = 300K inthe FMO complex. This provided me once again with a cubic pole structure for the prop-agator, this time to all orders in the system-bath coupling parameter. Incidentally, forthe case of Ohmic-bath spectral densities, I was unable to explore the low-temperatureregime for three distinct frequencies, as the pole structure in this limit is an arbitraryorder polynomial in the system-bath coupling parameter. The original 2-site spin-bosonmodel was able to explore this limit by utilising the special function: the Mittag-Leffler118function [31, 61], however the 3-site problem remains intractable analytically. A numeri-cal solution to the propagators pole structure is still permitted of course, but beyond thescope of this thesis. Nevertheless, the high-temperature regime was still accessible forthe Ohmic case and it is in this regime that we see the coherent-incoherent phase tran-sition. I was once again able to explore the phase space of the model by inspecting thediscriminant of the cubic polynomial, and in a similar fashion to the original spin-bosonmodel, we saw a phase transition beyond a critical value of the coupling parameter.I also quantified the validity of NIBA in the non-perturbative regime with a calcula-tion of F1 for various areas of the parameter space. I demonstrated how for the Ohmiccase, NIBA is primarily valid in the semi-classical regime of small tunnelling energy rel-ative to on-site energy. The complicated dependence of F1 on the system-bath couplingparameter as well as temperature required a graphical analysis. In general, we saw thatNIBA is valid in the mid-high temperature regime as well as the strong coupling regime.This is intuitive as NIBA is a strong-decoherence approximation and one would expectthis to apply in the case of strong system-bath coupling and/or high-temperature.The primary results of the non-perturbative analysis are the relaxation times. In theOhmic regime I calculated the three distinct frequencies coming from the cubic polestructure in the high-temperature limit. For the parameter space explored in this sec-tion: /∆ ∼ 10, kBT ωc, γ < 1, I found the system to exhibit a persistent coherentphase. The frequency associated with the dark-state of the system was found to have theslowest decay rate, and correspondingly the longest relaxation time, especially for smallto intermediate values of the coupling parameter. Therefore, it is this state in the systemthat dominates the long-time coherence in the system across the range 0 < γ < 1.The 3-site boson model, for an Ohmic bath spectral density, is ultimately applied tothe FMO complex. For the FMO parameters determined for an effective 3-site systemcoupled to an Ohmic bath, I found relaxation times τD ∼ 800fs, τ± ∼ 100fs at physio-logical temperatures. These results put the exciton transfer well within the incoherentregime, as the exciton coherence dies out faster than the time taken for the exciton topropagate through the system.In Chapter 5 I formulated the dual-coupling spin-boson model (DCSB) and investi-gated the results in both a general parameter regime–in units of the bias energy , andfor the FMO-complex parameters. This model describes a dimer system coupled bothdiagonally and non-diagonally to an oscillator bath, and I utilised an environment char-acterised by both super-Ohmic and optical phonon spectral densities. I once again usedNIBA to evaluate the influence functionals, meaning that this model applies in the regime119of small tunnelling energy to on-site energy.For the case of a super-Ohmic oscillator bath, I was able to explore results in thehigh temperature limit. We saw in this limit how the system remained in a persistentcoherent phase as the dimensionless non-diagonal coupling parameter is varied. For smallcouplings the dimer oscillation frequency rapidly increases along with a slight reduction inthe pure relaxation rate. The physical basis for this is due to the inelastic phonon-assistedtransport mechanism, whereby the transfer rate between dimer states is increased dueto the allowance of phonon emission/absorption processes. As the coupling is increasedbeyond a certain point, the oscillation frequency dies off as the particle’s phonon cloudbecomes heavier, making it harder for the particle to tunnel between wells. For theDCSB model with acoustic phonons applied to the FMO complex, similar effects werefound. The exciton transfer time–the inverse of the pure relaxation rate–was found toincrease slightly for small enough values of the non-diagonal coupling parameter as well.Correspondingly, the coherent relaxation time–the time associated with the relaxation ofthe oscillatory terms in the dynamics–was found to also increase for small non-diagonalcouplings. However, as the non-diagonal coupling was increased further, both the trans-fer and coherence times of exciton were seen to rapidly decrease. This happens alongsidea reduction in dimer oscillation frequency and can be interpreted within the context ofphonon-assisted tunnelling. For small enough couplings the exciton is able to tunnel backand forth coherently between the two states in the system. This is due to an increaseddensity of final states available to the exciton by way of phonon emission/absorption pro-cesses. However, as the coupling becomes sufficiently strong, the exciton is over ’dressed’by its phonon cloud and tunnelling becomes one-way.The addition of non-diagonal system-bath coupling serves to not only increase the rateat which the exciton is transferred through the dimer but also introduce an additionaldecoherence mechanism. However, the coherence time remains longer than the excitontransfer within a certain range of non-diagonal couplings. Therefore, one can concludethat, provided the non-diagonal coupling strength is significant, the exciton remains co-herent as it travels through the FMO dimer according to the DCSB model.The DCSB model has successfully managed to produce both rapid and coherent exci-ton transfer in the FMO complex. Traditional models for FRET fail to reproduce theseeffects, which suggests the inclusion of non-diagonal couplings is crucial to the excitontransfer mechanism in the FMO complex. Fo¨rster theory can produce rapid exciton trans-fer times but excludes any possibility of coherence. In the opposite limit, Redfield theoryfails to produce the rapid exciton transfer times observed in the FMO complex. Fur-120thermore, in both cases, the model parameters do not accurately reflect the real physicalsystems they model. Neither the assumption of perturbative tunnelling energy relative tosystem-bath coupling in Fo¨rster theory, nor the Markov approximation and perturbativesystem-bath coupling assumed in Redfield theory are valid. A non-perturbative theory,in any of these parameters, falls within the purview of the spin-boson model. Whenapplied to the FMO complex, the spin-boson model predicts incoherent exciton transferfor Ohmic spectral densities, and coherent but slow exciton transfer for super-Ohmicbaths. Therefore one can conclude that the inclusion of non-diagonal couplings is crucialto the maintenance of rapid and coherent exciton transfer in the FMO complex with asuper-Ohmic bath.Motivated by the evidence of structured spectral densities in the FMO complex, I alsoinvestigated the DCSB model for an optical phonon bath. A Gaussian lineshape wasused to model the optical phonon peaks, and the bath correlators were evaluated usingthe saddle-point approximation. This technique is valid for spectral lineshapes with verysmall width which is indeed the case for the FMO complex. We saw in Section 2.3.4,that experimentally determined spectra for the FMO complex found optical phononpeaks at relatively high frequencies with small linewidths. Therefore the saddle-pointapproximation was justified in this case. I also demonstrated in Section 2.3.4, how thedimensionless coupling parameters pertaining to the dominant peaks in the FMO spectrawere very small. This permitted a perturbative approach to modelling the system-bathcouplings. In the interest of investigating the non-perturbative effects of optical phononson the FMO system, a fictitious peak with frequency within the bandwidth of the FMOparameters was chosen for the DCSB model. The strong coupling to this optical phononpeak is hoped to reflect the overall coupling of the FMO system to the many weak opticalphonon vibrations.The results of the DCSB model with an optical phonon bath were similar to the acous-tic phonon (super-Ohmic) results. The dimer oscillation frequency increased for smallenough values of the non-diagonal system-bath coupling energy, and the coherence timeremained longer than the pure relaxation time for the range of couplings. Once again,while both the exciton transfer time and coherence time were greatly reduced for signifi-cant non-diagonal couplines, the coherence time remained persistently longer. Thereforethe optical phonon bath also aids in exciton transfer through the dimer in a similarlyefficacious way to the acoustic phonon bath, while the exciton remains coherent for theduration of its motion.It is prudent to ask the question at this point: what is the functionality of coherence121in the FMO exciton transfer mechanism? The experimental evidence for coherence in thesystem is clear, and the success of the DCSB model here has demonstrated the impor-tance of including non-diagonal couplings in the modelling of exciton transfer. However,what do these results say about the functionality of coherence in the system? Coherent(delocalised) excitons are characterised by overlapping exciton wavefunctions on nearbychromophores. The resulting exciton dynamics behaves in a wave-like manner, and asso-ciated with this will be the oscillatory ‘quantum beating’ signals observed in experiment.When there exists significant wavefunction overlap between chromophores, non-diagonalcouplings should be present. The DCSB model has demonstrated that significant non-diagonal couplings promotes exciton transfer through the dimer system in accordancewith the experimentally observed values. While the presence of non-diagonal couplingsserves to introduce an additional decoherence mechanism, the coherence time of excitonremains longer than the exciton transfer time. This is perhaps unsurprising, since thepresence of non-diagonal couplings is predicated on the coherent nature of excitons inthe first place. Therefore, the functionality of coherence in photosynthesis, if any, couldbe to introduce non-diagonal system-bath couplings between chromophores. This facil-itates the inelastic phonon-assisted transport mechanism and provides the exciton withadditional tunnelling pathways between chromophores.6.1 Further WorkI now discuss a number of avenues of further research that could build on the work pre-sented in this thesis. Neither the 3-site boson model presented in this thesis nor theDCSB model are exactly solvable due to the complicated nature of the system-bath cou-pling. The noninteracting-blip approximation (NIBA) employed here is a first step alongthe way of treating the system-bath interaction processes in a non-perturbative couplingapproach. A number of analytic studies have been done that go beyond NIBA. Thesetake the form of summing a larger class of diagrams, including higher order system-bathprocesses. The nearest-neighbour blip-approximation (NNBA) for example includes con-secutive blip-blip interaction terms, permitting the system to spend more time in theoff-diagonal state of the density matrix and therefore incorporate longer-lived coherenceeffects. This would allow one to consider stronger tunnelling matrix elements with re-spect to the on-site energies in both the 3-site and 2-site models. This is because thesystem is permitted to spend more time in off-diagonal states of the density matrix. Thispreserves longer coherence effects and one would expect to see an increase in relaxation122times calculated for these systems. For the FMO complex, there are a number of three-and two-site molecular configurations that contain a larger ratio of ∆/ than consideredin this work. Modelling these molecular configurations would require a model that goesbeyond NIBA, and NNBA could be a candidate for this. Of course, as one considersincreasingly higher order diagrams in this regard, the analytical complexity of the modelgreatly increases. However for the 2-site spin-boson model this has nevertheless beenshown to remain a tractable analytical problem yielding illuminating results [27, 99].For the 3-site boson model it remains unclear whether such a relaxed approxima-tion (NNBA) can be applied successfully analytically. The problem one faces in the3-dimensional Hilbert space of the model is substantially harder than that for the 2-D Hilbert space of the 2-site model. Excursions into the off-diagonal elements of the3-by-3 density matrix need to take into account all possible consecutive blip-blip andblip-sojourn interactions which make up the space of the 3-by-3 density matrix. Sincethe system is permitted to spend additional time in the off-diagonal space of the densitymatrix, the combinatorics problem grows substantially. Nevertheless, the solution to thisproblem would be a valuable one not just in the interest of relaxation NIBA and ob-serving the additional coherence effects in the model but also because it would allow thesystem to occupy the coherent state in the system that overlaps with both the upper twolevels. This possibility was excluded in our model because of the application of NIBAwhich quenches the off-diagonal excursions too rapidly before they can enter into thisfragile coherent state.One of the famous applications of the spin-boson model was to the Kondo problem.It was shown that the peculiar singular behaviour of charged interstitials in conductingmetals could be explained by the high density of electron-hole excitations around theFermi surface. This leads to a self-trapping phase transition at zero temperature above acritical value of the coupling strength, as well as an anomalous temperature dependenceon diffusion. The spin-boson model was appropriate in this case because the impuritiescan be represented by a spin-1/2 central system, and the surrounding Fermi gas can becharacterised by electron-hole excitations around the Fermi surface at such low tempera-tures. These electron-hole pairs are bosonic in character and obey an Ohmic form for thespectral density. Therefore, in the context of this research, one could ask the question:how does the Kondo effect change when the impurity is actually a 3-level system suchas the one studied in this system? As discussed in this thesis, the 3-level system con-tains intrinsic properties such as population trapping that are very different to a 2-levelsystem. This could have dramatic effects on the Kondo effect when applied to metals.123The 3-level model coupled to an Ohmic oscillator bath would need to be explored in thezero temperature regime, as only the high-T regime was considered in this thesis. Onewould also need to inspect certain response functions such as the magnetic susceptibilityin order to answer specific questions relating to the Kondo effect.Finally I would like to comment on the implications of these findings with respect tothe development of artificial light-harvesting systems. Research into photosynthesis ismotivated not only by our general interest in a complete understanding of the biophysicsbut also by the potential applications to our own light-capturing technologies. Naturehas had millions of years to optimise these processes, and reverse engineering the en-ergy transfer mechanism in photosynthesis could help us make our own artificial lightharvesting technologies more efficient. If we are to achieve higher efficiencies in light-capturing technologies, it is possible that this can be achieved by promoting inelasticphonon-assisted tunnelling processes in the exciton transport schemes of these devices.This would mean reducing the inter-molecular spacing as much as possible, as this wouldfacilitate greater exciton wavefunction overlap between molecules, and therefore aug-mented phonon-assisted tunnelling. We have learned from photosynthesis, that even inhot and messy molecular environments, exciton transfer can be very efficient, providedthe inter-molecular distances are kept sufficiently small. In this case, quantum effects areallowed to persist in otherwise non-ideal environments.124Bibliography[1] D. P. DiVincenzo, D. Loss. Quantum computers and quantum coherence. Journalof Magnetism and Magnetic Materials, 200(202-218), 1999.[2] D. A. Ryndyk. Theory of Quantum Transport at Nanoscale. Springer, 2016.[3] Hwang, H and Rossky, P. Electronic decoherence induced by intramolecular vi-brational motions in a betaine dye molecule. J. Phys. Chem. B, 108(6723–6732),2004.[4] Franco, I; Shapiro, M; Brumer, P. Femtosecond dynamics and laser controlof charge transport in trans-polyacetylene. The Journal of Chemical Physics,128(244905), 2008.[5] Gregory S. Engel, Tessa R. Calhoun, Elizabeth L. Read, Tae-Kyu Ahn, TomasMancal, Yuan-Chung Cheng, Robert E. Blankenship, and Graham R. Fleming.Evidence for wavelike energy transfer through quantum coherence in photosyntheticsystems. Nature, 446(7137):782–786, April 2007.[6] H. Lee, Y. C. Cheng, and G. R. Fleming. Coherence dynamics in photosynthesis:Protein protection of excitonic coherence. Science, 316(1462), 2007.[7] Panitchayangkoon, G et al. Long-lived quantum coherence in photosynthetic com-plexes at physiological temperature. PNAS, 107(12766–12770), 2010.[8] Collini, E et al. Coherently wired light-harvesting in photosynthetic marine algaeat ambient temperature. Nature, 463(644), 2010.[9] Fassioli, F et. al. Photosynthetic light harvesting: excitons and coherence. J. R.Soc. Interface, 11(20130901), 2014.125[10] Roseanne J. Sension. Biophysics: Quantum path to photosynthesis. Nature,446(7137):740–741, April 2007.[11] I. Burghardt, V. May, D. A. Micha, and E. R. Bittner, editors. Energy TransferDynamics in Biomaterial Systems. Springer Series in Chemical Physics 93, 2009.[12] R. E. Blankenship. Molecular Mechanisms of Photosynthesis. Blackwell Science,London, 2002.[13] X. Hu and K. Schulten. How nature harvests sunlight. Physics Today, 50(28), 1997.[14] C. Curutchet and B. Mennucci. Quantum chemical studies of light harvesting.Chemical Reviews, 117(294-343), 2016.[15] M. B. Plenio and S. F. Huelga. Dephasing-assisted transport: quantum networksand biomolecules. New J. Phys., 10(11):113019, November 2008.[16] F. Caruso, A. W. Chin, A. Datta, S. F. Huelga, and M. B. Plenio. Highly efficientenergy excitation transfer in light-harvesting complexes: The fundamental role ofnoise-assisted transport. J. Chem. Phys., 131(10):105106, September 2009.[17] A. W. Chin, A. Datta, F. Caruso, S. F. Huelga, and M. B. Plenio. Noise-assistedenergy transfer in quantum networks and light-harvesting complexes. New J. Phys.,12(6):065002, June 2010.[18] A.W. Chin et al. The role of non-equilibrium vibrational structures in electroniccoherence and recoherence in pigment–protein complexes. Nature Physics, 9, 2013.[19] J. Lim et al. Phonon-induced dynamic resonance energy transfer. New Journal ofPhysics, 16, 2014.[20] J. Cao and R.J. Silbey. Optimization of exciton trapping in energy transfer pro-cesses. The Journal of Physical Chemistry A, 113(50), 2009.[21] Ishizaki, A; Fleming, G.R. Unified treatment of quantum coherent and incoher-ent hopping dynamics in electronic energy transfer: Reduced hierarchy equationapproach. THE JOURNAL OF CHEMICAL PHYSICS, 130(234111), 2009.[22] A. Kolli et al. The fundamental role of quantized vibrations in coherent lightharvesting by cryptophyte algae. The Journal of Chemical Physics, 137(174109),2012.126[23] Alexandra Olaya-Castro, Chiu Fan Lee, Francesca Fassioli Olsen, and Neil F. John-son. Efficiency of energy transfer in a light-harvesting system under quantum co-herence. Phys. Rev. B, 78(8):085115, August 2008.[24] Ishizaki, A; Fleming, G.R. On the adequacy of the redfield equation and relatedapproaches to the study of quantum dynamics in electronic energy transfer. J.Chem. Phys, 130(234110), 2009.[25] Sharp, L.Z; Egorova, D; Domcke, W. Efficient and accurate simulations of two-dimensional electronic photonecho signals: Illustration for a simple model of thefenna–matthews–olson complex. J. Chem. Phys, 132(014501), 2009.[26] P. Nalbach, D. Braun and M. Thorwart. Exciton transfer dynamics and quan-tumness of energy transfer in the fenna-matthews-olson complex. Phys. Rev. E,84(041926), 2011.[27] Leonardo A. Pacho´n and Paul Brumer. The physical basis for long-lived electroniccoherence in photosynthetic light harvesting systems. J. Phys. Chem. Lett., 2(21),2011.[28] V.J Emery, A. Luther. Low-temperature properties of the kondo hamiltonian. Phys.Rev. B, 9(1), 1974.[29] Caldeira, A & Leggett, A. Influence of dissipation on quantum tunneling in macro-scopic systems. PRL, 46(4), 1981.[30] Caldeira, A & Leggett, A. Quantum tunnelling in a dissipative system. Annals ofPhysics, 149(374-456), 1983.[31] Leggett, A.J., Chakravarty, S. et al. Dynamics of the dissipative two-state system.Rev. Mod. Phys., 59(1), 1987.[32] Chakravarty, S; Leggett, A. Dynamics of the two-state system with ohmic dissipa-tion. PRL, 52(1), 1984.[33] Dekker, H. Noninteracting-blip approximation for a two-level system coupled to aheat bath. PRA, 35(3), 1987.[34] Weiss, U. Influence of friction and temperature on coherent quantum tunneling.Journal of Low Temperature Physics, 68(3/4), 1987.127[35] Grabert, H; Weiss, U. Quantum tunneling rates for asymmetric double-well systemswith ohmic dissipation. PRL, 54(15), 1985.[36] Fisher, M.P.A. and Dorsey, A.T. Dissipative quantum tunneling in a biased double-well system at finite temperatures. Physical Review Letters, 54(15), 1985.[37] Feynman, R.P.; Vernon, F.L. The theory of a general quantum system interactingwith a linear dissipative system. Annals of Physics, 24(118-173), 1963.[38] Fisher, M.P.A and Zwerger, W. Quantum brownian motion in a periodic potential.Phys. Rev. B, 32(10), 1985.[39] Wilhelm, F.K.; Kleff, S; von Delft, J. The spin-boson model with a structuredenvironment: a comparison of approaches. Chemical Physics, 296(345-353), 2004.[40] Grabert, H. Dissipative quantum tunneling of two-state systems in metals. PRB,46(19), 1992.[41] Joel Gilmore, Ross H. McKenzie. Spin boson models for quantum decoherence ofelectronic excitations of biomolecules and quantum dots in a solvent. J. Phys.:Condens. Matter, 17(1735), 2005.[42] Joel Gilmore, Ross H. McKenzie. Criteria for quantum coherent transfer of excitonsbetween chromophores in a polar solvent. arXiv:quant-ph/0412170, 2006.[43] Gilmore, J and McKenzie, R. Quantum dynamics of electronic excitations inbiomolecular chromophores: Role of the protein environment and solvent. J. Phys.Chem. A, 112(2162-2176), 2008.[44] Jeremy Moix et al. Efficient energy transfer in light-harvesting systems, iii: Theinfluence of the eighth bacteriochlorophyll on the dynamics and efficiency in fmo.J. Phys. Chem. Lett., 2(3045–3052), 2011.[45] Hui Dong, Da-Zhi Xu, Jin-Feng Huang, and Chang-Pu Sun. Coherent excitationtransfer via the dark-state channel in a bionic system. Light Sci. Appl., 1(3):e2,March 2012.[46] Julia Adolphs and Thomas Renger. How proteins trigger excitation energy transferin the FMO complex of green sulfur bacteria. Biophys. J., 91(8):2778–2797, 2006.128[47] T. Holstein. Studies of polaron motion: Part 2. the ’small’ polaron. Ann. Phys.,8(343), 1959.[48] R. Peierls. Quantum Theory of Solids. Clarendon Press, Oxford, 1955.[49] W. P. Su, J. R. Schrieffer, and A. J. Heeger. Solitons in polyacetylene. Phys. Rev.Lett., 42(1698), 1979.[50] D.J. J. Marchand, P.C. E. Stamp, and M.Berciu. Dual coupling effective bandmodel for polarons. Phys. Rev. B, 95(035117), 2017.[51] G. H. Richards, K. E. Wilk, P. M. G. Curmi, H. M. Quiney, and J. A.Davis. Excitedstate coherent dynamics in light-harvesting complexes from photosynthetic marinealgae. J. Phys. B: At. Mol. Opt. Phys., 45(154015), 2012.[52] G.D. Mahan. Many particle physics, 2nd ed. Plenum Press, 1990.[53] N. Wu, K. Sun, Z. Chang, and Y. Zhao. Resonant energy transfer assisted byoff-diagonal coupling. J. Chem. Phys., 136(124513), 2012.[54] A. Nazir. Correlation-dependent coherent to incoherent transitions in resonantenergy transfer dynamics. Phys. Rev. Lett., 102(146404), 2009.[55] Y. Yao, N, Zhou, J, Prior, Y, Zhao. Competition between diagonal and off-diagonalcoupling gives rise to charge-transfer states in polymeric solar cells. Nature Scien-tific Reports, 5(14555), 2015.[56] V. Coropceanu et. al. Charge transport in organic semiconductors. Chem. Rev.,107(926-952), 2007.[57] S.M. Falke, et al. Coherent ultrafast charge transfer in an organic photovoltaicblend. Science, 344(1001), 2014.[58] Wang, T. and Chan, W.-L. Dynamical localization limiting the coherent transportrange of excitons in organic crystals. J. Phys. Chem. Lett., 5(1812C1818), 2014.[59] G. Scholes et. al. Using coherence to enhance function in chemical and biophysicalsystems. Nature, 543(647), 2017.[60] A. Wu¨rger. Strong-coupling theory for the spin-phonon model. Phys. Rev. B.,57(1), 1998.129[61] U Weiss. Quantum Dissipative Sytems. World Scientific, 2008.[62] Tobias Brixner, Jens Stenger, Harsha M. Vaswani, Minhaeng Cho, Robert E.Blankenship, and Graham R. Fleming. Two-dimensional spectroscopy of electroniccouplings in photosynthesis. Nature, 434(7033):625–628, March 2005.[63] D. Karcz et al. Lessons from chlorophylls: Modifications of porphyrinoids towardsoptimized solar energy conversion. Molecules, 19(15938-15954), 2014.[64] V. May and O. Kuhn, editors. Charge and Energy Transfer Dynamics in Molecu-lar Systems. Wiley, 2004.[65] A.S. Davydov. The theory of molecular excitons. Soviet Physics Uspekhi, 7(2),1964.[66] G.D. Scholes. Long-range resonance energy transfer in molecular systems. Annu.Rev. Phys. Chem., 54(57), 2012.[67] V. Helms, editor. Principles of Computational Cell Biology. Wiley, 2008.[68] Peterman, E.J, Dukker, F.M, van Grondelle, R, van Amerongen, H. Chlorophyll aand carotenoid triplet states in light-harvesting complex ii of higher plants. Biophys.J., 69(2670-2678), 1995.[69] B.W. Matthews, R.E. Fenna. Structure of a green bacteriochlorophyll protein. Acc.Chem. Res., 13(9), 1980.[70] Tronrud,D.E.;Wen,J.;Gay,L.;Blankenship,R.E. Thestructural basis for the differ-ence in absorbance spectra for the fmo antenna protein from various green sulfurbacteria. Photosynth. Res, 100(79-87), 2009.[71] A. Lohner et al. Fluorescence-excitation and emission spectroscopy on single fmocomplexes. Nature Scientific Reports, 6(31875), 2016.[72] D. Singh and S. Dasgupta. Coherence and its role in excitation energy transfer infenna- matthews-olson complex. J. Phys. Chem. B, 121(1290-1294), 2017.[73] M. Aghtar et al. The fmo complex in a glycerol-water mixture. J. Phys. Chem. B,117(7157), 2013.130[74] Hayes, D. & Engel, G. S. Extracting the excitonic hamiltonian of the fenna-matthews-olson complex using three-dimensional third-order electronic spec-troscopy. Biophys. J., 100(2043–2052), 2011.[75] C. Olbrich, T. L. C. Jansen, J. Liebers, M. Aghtar, J. Stru¨mpfer, K. Schulten, J.Knoester, and U. Kleinekatho¨fer. From atomistic modeling to excitation transferand two-dimensional spectra of the fmo light-harvesting complex. J. Phys. Chem.B, 115(8609), 2011.[76] S. Shim, P. Rebentrost, S. Valleau, and A. Aspuru-Guzik. Atomistic study of thelong-lived quantum coherences in the fenna-matthews-olson complex. Biophys. J.,102(649), 2012.[77] B. S. Prall, D. Y. Parkinson, G. R. Fleming, M. Yang, and N. Ishikawa. Two-dimensional optical spectroscopy: Two-color photon echoes of electronically cou-pled phthalocyanine dimers. J. Chem. Phys., 120(2537), 2004.[78] S. M. Vlaming, and R. J. Silbey. Correlated intermolecular coupling fluctuationsin photosynthetic complexes. J. Chem. Phys., 136(055102), 2012.[79] D. J. J. Marchand, G. De Filippis, V. Cataudella, M. Berciu, N. Nagaosa, N. V.Prokofev, A. S. Mishchenko and P. C. E. Stamp. Sharp transition for single polaronsin the one-dimensional su-schrieffer-heeger model. Phys. Rev. Lett., 105(266605),2017.[80] B.W. Matthews, R.E. Fenna. Calculation of couplings and energy-transfer path-ways between the pigments of lh2 by the ab initio transition density cube method.J. Phys. Chem. B, 102(5378-5386), 1998.[81] Sharp, L.Z.; Egorova, D.; Domcke, W. . Efficient and accurate simulations of two-dimensional electronic photon-echo signals: Illustration for a simple model of thefenna–matthews–olson complex. The Journal of Chemical Physics, 132(014501),2010.[82] T. Forster. Delocalized excitation and excitation transfer. Mordern QuantumChemistry, Istanbul Lectures, 3(93), 1965.[83] H. Fukagawa, T. Shimizu, Y. Iwasaki . Operational lifetimes of organic light-emitting diodes dominated by fo¨rster resonance energy transfer. Nature: ScientificReports, 7(1735), 2017.131[84] A. Ishizaki, et al. Quantum coherence and its interplay with protein environ-ments in photosynthetic electronic energy transfer. Phys. Chem. Chem. Phys.,12(7319–7337), 2010.[85] A.G. Redfield. On the theory of relaxation processes. IBM J. Res. Dev., 1(19),1957.[86] P. Lambropoulus and D. Petrosyan. Fundamentals of Quantum Optics and Quan-tum Information. Springer, 2007.[87] H.P. Breuer and F. Petruccione. THE THEORY OF OPEN QUANTUM SYS-TEMS. Oxford University Press, 2002.[88] F. Petruccione H. Breuer. The Theory of Open Quantum Systems. Oxford Univer-sity Press, 2007.[89] H.B. Chen et al. Using non-markovian measures to evaluate quantum master equa-tions for photosynthesis. Nature Scientific Reports, 5(12753), 2015.[90] M. Schmidt et al. The eighth bacteriochlorophyll completes the excitation energyfunnel in the fmo protein. J. Phys. Chem. Lett., 2(93-98), 2011.[91] Thyrhaug, E; et al. Exciton structure and energy transfer in the fenna-matthews-olson complex. J. Phys. Chem. Lett., 7(1653-1660), 2016.[92] K. Saito, T. Suzuki, H. Ishikita. Absorption-energy calculations of chlorophyll aand b with an explicit solvent model. Journal of Photochemistry and PhotobiologyA: Chemistry, 358(422-431), 2017.[93] Leggett, A.J. Percolation, localization and supercon- ductivity. NATO ASI SeriesB: Physics, 109(p.1), 1984.[94] Y. Cheng and G. R. Fleming. Dynamics of light harvesting in photosynthesis.Annu. Rev. Phys. Chem., 60(241), 2009.[95] G. D. Scholes and G. R. Fleming. Energy transfer and photosynthetic light harvest-ing. Adventures in Chemical Physics: A Special Volume in Advances in ChemicalPhysics, 132(57), 2006.[96] Guinea, F. Friction and particle-hole pairs. PRL., 53(13), 1984.132[97] Ha¨nngi, P. and Ingold, G.L. Fundamental aspects of quantum brownian motion.Chaos, 15(026105), 2005.[98] Jang, S.J. and Mennucci, B. Delocalized excitons in natural light harvesting com-plexes. Rev. Mod. Phys., 90(035003), 2018.[99] P. Esquinazi, editor. Tunneling Systems in Amorphous and Crystalline Solids.Springer, 1998.[100] S. Hunklinger. Phonons in amorphous solids. Journal de Physique, 12(c9-461),1982.[101] M. F. Thorpe. Phonons in Amorphous Solids. Springer, 1976.[102] M. Aghtar et al. Different types of vibrations interacting with electronic excitationsin phycoerythrin 545 and fenna-matthews-olson antenna systems. J. Phys. Chem.Lett., 5(3131), 2014.[103] Cho,M et al. Exciton analysis in 2d electronic spectroscopy. J. Phys. Chem. B,109(10542-10556), 2005.[104] B. Li et al. The brownian oscillator model for solvation effects in spontaneouslight emission and their relationship to electron trans fer. J. Am. Chem. SOC,116(11039-11047), 1994.[105] E.T.J. Nibbering, D.A. Wiersma, K.Duppen . Ultrafast electronic fluctuation andsolvation in liquids. Chemical Physics, 183(167-185), 1994.[106] L. Onsager. Electric moments of molecules in liquids. J. Am. Chem. Soc., 58(1486),1936.[107] J. Tomasi and M. Persico. Molecular interactions in solution: An overview ofmethods based on continuous distributions of the solvent. Chem. Rev., 94(2027),1994.[108] J. Newman and K.E. Thomas-Alyea. Electrochemical systems. Wiley, 2004.[109] V.I. Novoderezhkin et al. Excitation dynamics in the lhcii complex of higher plants:Modeling based on the 2.72 a˚crystal structure. J. Phys. Chem. B, 109(10493-10504),2005.133[110] E. Read et al. Visualization of excitonic structure in the fenna-matthews-olsonphotosynthetic complex by polarization-dependent two-dimensional electronic spec-troscopy. Biophysical Journal, 95(847–856), 2008.[111] D. Zigmantas et al. Two-dimensional electronic spectroscopy of the b800–b820light-harvesting complex. . Proc. Natl. Acad. Sci. USA, 103(12672–12677), 2006.[112] T. Renger and R. A. Marcus. On the relation of protein dynamics and excitonrelaxation in pigment–protein complexes: An estimation of the spectral densityand a theory for the calculation of optical spectra. J. Chem. Phys., 116(22), 2002.[113] S.Jang, Y.C. Cheng, D.R. Reichman, and J.D. Eaves. Theory of coherent resonanceenergy transfer. J. Chem. Phys., 129(101104), 2008.[114] Olbrich, C et al. Theory and simulation of the environmental effects on fmo elec-tronic transitions. J. Phys. Chem. Lett., 2(, 1771–1776), 2011.[115] Damjanovic, A; Kosztin, I; Schulten, K. Excitons in a photosynthetic light-harvesting system: A combined molecular dynamics/quantum chemistry and po-laron model study. Physical Review E, 65(031919), 2002.[116] Feynman, R.P. Space-time approach to non-relativistic quantum me- chanics. Rev.Mod. Phys., 20(367-387), 1948.[117] A. Atland and B. Simons. Condensed Matter Field Theory. Cambridge UniversityPress, 2010.[118] A. Devaquet. Avoided crossings in photochemistry. Pure and applied chemistry,41(455), 1975.[119] P. M. Radmore and P. L. Knight. Population trapping and dispersion in a three-level system. J. Phys. B: At. Mol. Phys., 15(4):561, February 1982.[120] Radmore, P.M; Knight, P.L. Population trapping and dispersion in a three-levelsystem. J. Phys. B, 15(561-573), 1982.[121] E.N. Economou. Green’s Functions in Quantum Physics. Springer, 2006.[122] Wurger, A. Dissipative tunneling in insulators: Noninteracting blip approximationand beyond. Physical Review Letters, 9(0031-9007), 1997.134[123] Riccardi, M. Solution of cubic and quartic equations. Formalized Mathematics,17(2), 2009.[124] Holmes, G.C. The use of hyperbolic cosines in solving cubic polynomials. TheMathematical Gazette, 86(507), 2002.[125] Nickalls, R.W.D. A new approach to solving the cubic: Cardan’s solution revealed.The Mathematical Gazette, 77(480), 1993.[126] I.M. Gradshteyn, I.S. & Ryzhik. Table of Integrals, Series, and Products. Elsevier,7th ed, 2007.[127] Davis, P.J. Leonard euler’s integral: A historical profile of the gamma function.The American Mathematical Monthly, 66(849-869), 1959.[128] L. Kleinman. Theory of phonon-assisted tunneling in semiconductors. Phys. Rev.,140(A637), 1965.[129] W. P. Su, J. R. Schrieffer, and A. J. Heeger. Solitons in polyacetylene. Phys. Rev.Lett., 42(1698), 1979.[130] R. Silbey, R.W. Munn. General theory of electronic transport in molecular crys-tals. i. local linear electron–phonon coupling. The Journal of Chemical Physics,72(2763), 1980.[131] R.W. Munn, R. Silbey. Theory of electronic transport in molecular crystals. ii.zeroth order states incorporating nonlocal linear electron–phonon coupling. TheJournal of Chemical Physics, 83(1843), 1985.[132] R.W. Munn, R. Silbey. Theory of electronic transport in molecular crystals. iii.diffusion coefficient incorporating nonlocal linear electron–phonon coupling. TheJournal of Chemical Physics, 83(1854), 1985.[133] D. Chen, J. Ye, H. Zhang, Y. Zhao. On the munn-silbey approach to polarontransport with off-diagonal coupling and temperature-dependent canonical trans-formations. J. Phys. Chem. B, 115(5312), 2011.[134] Z. Huang, et al. Polaron dynamics with off-diagonal coupling: beyond the ehrenfestapproximation. Phys. Chem. Chem. Phys, 19(1655), 2017.135[135] M. Abramowitz and I.A. Stegun. Handbook of Mathematical Functions. DoverPublications, 1972.[136] J. Matthews and R.L Walker, editors. Mathematical Methods of Physics (2ndEdition). Addison-Wesley, 1971.136Appendix ACalculation of the electron-phononcorrelation function for DCSB modelHere I present the calculation of the electron-phonon correlation function used to de-scribe the effects of the phonon bath on a central system. In section chapter 5 I requireda calculation of the electron-phonon correlation function 〈Bx+(t)Bz−(t)〉, which involvesa correlation between phonon displacement operators pertaining to both diagonal andnon-diagonal couplings. While following a derivation in Mahan [52] in evaluating thesecorrelation functions, I present the specifics of my deviations rather than simply quotingthe results, as I must be careful to distinguish between the diagonal and non-diagonalcouplings λx and λz respectively. In the traditional literature, the electron correlationfunction calculation usually assumes just diagonal couplings present in the system there-fore I must present the more general derivation case here. Starting from the calculationof the trace over the phonon distributionsFxz(t) =〈Bx+(t)Bz−(0)〉B= eβΩphTr[e−β∑q ωqnqBx+(t)Bz−(0)](A.1)Averaging each phonon state independently137Fxz(t) =∏qfxzq (t)fxzq (t) = eβΩph∞∑nq=0e−βnqωq 〈nq| e−λxq (b†qeiωqt−bqe−iωqt)eλzq(b†q−bq) |nq〉 (A.2)whereeβΩq = ∞∑nq=0e−βnqωq−1 = 1− e−βωq (A.3)is a normalisation prefactor. For brevity we drop the phonon wavevector subscript q forthe time being.fxz(t) = (1− e−βω)∞∑n=0e−βnω 〈n|−λx(b†eiωt−be−iωt) eλz(b†−b) |n〉 (A.4)The state of n bosonic excitations is given by|n〉 = (b†)n√n!|0〉 (A.5)Performing the Feynman-disentangling of operators [52] and applying the BCH formulaeA+B = eAeBe−(1/2)[A,B] we findBx+(t)Bz−(0) = e−λ2x/2−λ2z/2exp(−λxb†eiωt + λxbe−iωt + λzb† − λzb) (A.6)In the interest of getting all the destruction operators on the right hand side and thecreation operators on the left, the center two operators need to be exchangedeλxb(t)eλzb†= eλzb†[e−λzb†eλxb†(t)eλzb†](A.7)Applying the BCH formula to evaluate the expressione−λzb†beλzb†= b+ λz (A.8)allows us to expresseλxb(t)eλzb†= eλzb†exp[λxe−iωt(b+ λz)]= eλzb†eλxλze−iωteλxb(t) (A.9)138So this leaves us for the electron-phonon correlation function for each wavevector qfxz(t) =(1− e−βω)exp(−λ2x2− λ2z2+ λxλze−iωt)×∞∑n=0e−βωnexp 〈n| [(λz − λxeiωt)b† − (λz − λxe−iωt)b] |n〉 (A.10)defining u = λz − λxeiωt, we can expand in a power seriese−ua |n〉 =∞∑l=0(−u)ll!al |n〉 (A.11)Using the properties of the boson annihilation operators acting on the harmonic oscillatorstates gives use−ua |n〉 =n∑l=0(−u)ll![n!(n− l)!]1/2|n− l〉 (A.12)Using the orthogonality of harmonic oscillator states and identifying the Laguerre poly-nomial of order n〈n| e−u∗a†e−ua |n〉 = Ln(|u|2) (A.13)Performing the summation over boson states by identifying the generating function ofLaguerre polynomials [52](1− z)∞∑n=0Ln(|u|2)zn = e|u|2z/(z−1) (A.14)where z = e−βω and −z/(z − 1) = N = 1/(eβω − 1). This leaves us withf qxz(t) = exp(−λ2x(1 +Nq)/2− λ2z(1 +Nq)/2) exp (λxλz(1 +Nq)e−iωt + λxλzeiωtNq)(A.15)Therefore the full electron-phonon correlation function is the product over phonon states139Fxz(t) =∏qf qxz(t)=∏qexp(−u2x(1 +Nq)/2− u2z(1 +Nq)/2) exp (uxuz [eiωqt +Nq (eiωqt + e−iωqt)])(A.16)which in its continuous form (see the next section B) isFxz(t) = exp(−ˆ ∞0dω2ω2(Jx(ω) + Jz(ω))coth(~βω/2))× exp(ˆ ∞0dω√Jx(ω)Jz(ω)ω2(i sinωt+ cos(ωt) coth(~βω/2)))(A.17)which is the result we use in chapter 5.140Appendix BConnection between continuous anddiscrete electron-phonon correlationfunctionIn the text we use two different forms for the electron-phonon propagator and here makethe link between the two. In the literature, the connection between the two is nevermade and usually the reader is expected to assume the jump between the discrete andcontinuous forms of the electron-phonon propagator. Here we start with the continuousformϕ(t) =ˆ ∞0dωJ(ω)ω2(i sinωt− (1− cosωt) coth (~βω/2)) (B.1)Inserting the general form of the spectral function1412ˆ ∞0dω∑qλ2ω2δ(ω − ωq)(i sinωt− (1− cosωt) coth ~βω/2)=∑qu2q[eiωqt − e−iωqt +(eiωqt + e−iωqt − 1)coth ~βωq/2]=∑qu2q[eiωqt − e−iωqt +(eiωqt + e−iωqt − 1)(1 + 2nq)]=∑qu2q[2eiωqt + 2nq(eiωqt + e−iωqt)− (1 + 2nq)]=∑qu2q[nq(eiωqt − 1)+ (1 + nq)(e−iωqt − 1)](B.2)wherenq =1e~βωq − 1 (B.3)and we recognise the finally form as the discrete form of the electron-phonon propagator(see Mahan).142Appendix CFull eigenvalues ofdetuned-3-site-V-systemHere we present the eigenvalue solutions to the HamiltonianH3 = 1 |1〉〈1|+ 2 |2〉〈2|+ 0 |0〉〈0|+ ∆10(|1〉〈0|+ h.c.) + ∆20(|2〉〈0|+ h.c.) (C.1)where 0 = 0 is the ground state. This Hamiltonian represents the most general caseof the 3-site-V system, with the potential for detuned upper levels 1 6= 2. We includethe solutions here for the bare-3-site-V system primarily to demonstrate the complicatednature of the eigenvalues even without a bath. In the main text we mostly deal with thetuned case 1 = 2. The eigenvalues are calculated to beλ1 =1 + 23−3√23N[(1 + 2)2 + 3Ω˜2 − 12]− N3 3√2λ2 =1 + 23− (1 + i√3)3N[(1 + 2)2 + 3Ω˜2 − 12]+(1− i√3)N6 3√2λ3 =1 + 23− (1− i√3)3N[(1 + 2)2 + 3Ω˜2 − 12]+(1 + i√3)N6 3√2(C.2)where143N =[− 9∆2101 + 18∆2201 + 18∆2102 − 9∆2202+[4(−3 (∆210 + ∆220 − 12)− (1 + 2) 2) 3+(−9∆2101 + 18∆2201 + 18∆2102 − 9∆2202 − 231 + 3221 + 3221 − 232) 2]1/2− 231 + 3221 + 3221 − 232] 13(C.3)144Appendix DFluctuation-dissipation theoremThe fluctuation-dissipation theory is a central feature of linear response theory and isapplied in the perturbative model in the text. The theorem relates the relaxation of aweakly perturbed system to the thermal fluctuations in the environment i.e. the responseof the system to a small applied force is equivalent to the response to a spontaneousfluctuation. The main result of the theory relates the power spectrum S(ω) of thefluctuations, to the Fourier transform of the susceptibility χ(ω) (the linear responsefunction)S(ω) = ~Imχ(ω) coth(~βω/2) (D.1)with the random force correlation function given in the time-domain by〈ξ(t)ξ(0)〉 =ˆ ∞0dωJ(ω)(coth(β~ω/2)− i sin(ωt))(D.2)and the power spectrum is found to be [61]S(ω) = J(ω) coth(β~ω/2) (D.3)145Appendix ECubic polynomialsThe central problem in this thesis involves investigating a 3-site system interacting withits environment. Intrinsic to this problem in the various areas of the parameter spaceof interest involves solutions to cubic polynomials which appear as the pole structure ofthe propagators. Here we detail the general form that these solutions take as a referencefor the main section of the thesis. The formalism presented here can be easily checkedagainst those presented in [123, 124, 125]. Consider the equationx3 + ax2 + bx+ c = 0 (E.1)One solves this equation with the substitutionx = t− a3(E.2)which produces the solvable cubict3 + pt+ q = 0 (E.3)wherep =3b− a23, q =2a3 − 9ab+ 27c27(E.4)Limiting cases include;If q = 0 thent3 + pt = 0, t = ±√−p (E.5)146and the original 3-roots arex1 = −a3, x2,3 = ±√−p− a3(E.6)If p = 0 then t = 3√−q and the 3-roots of our original equation arex1 =3√−q − a3, x2,3 =3√−q(−12±√3i/2− a3(E.7)Beyond the simple cases outlined above we now present the general solutions. First wemust define the discriminant that determines the nature of the roots. The discriminantisD =(q2)2+(p3)3(E.8)and the 3-regimes are;If D > 0 then we have one real root and 2-complex conjugate roots and the solutions arex1 = u− v − a3, x2,3 = −12(u− v)± (u+ v)√32i− a3(E.9)whereu = 3√√D − q2, v = 3√√D +q2(E.10)If D < 0 then all 3 roots are real and distinct and the solutions are obtained with thecosine substitutionx1 = 23√r cos(φ3)− a3x2 = 23√r cos((φ+ 2pi)3)− a3x3 = 23√r cos((φ+ 4pi)3)− a3(E.11)wherer =√(−p3)3, φ = arccos(−q2r)(E.12)147If D = 0 then all the roots are real and 2 are equalx1 = 23√−q2− a3x2,3 = − 3√−q2− a3(E.13)148Appendix FSaddle point integration of opticalphonon correlation functionsHere I present the calculation of the correlation function for an optical phonon as requiredin the text. The influence phase is given byϕop = iQ′op +Q′′op =ˆ ∞0dωω2e− (ω−ω0)2ξ2 cos(ωτ)csch(~βω/2) (F.1)This integral is divergent as it stands however the saddle-point approximation can beused to calculate it [136]. This method is valid in the limit of very large exponentialarguments, which in this case corresponds to 1/ξ2 1. As ξ represents the full-width-half-maximum (FWHM) of the Guassian peak, this condition applies to a very narrowpeak, centred around ω0. The saddle point method evaluates real integrals of the formI =ˆ badωe−xf(ω)g(ω) (F.2)where x is very large. Another necessary condition for the use of the saddle-point methodis for the peak to be symmetrical about its center such that f(ω0) = 0, which is valid forthe Guassian lineshape used here. Expanding about ω0f(ω) = f(ω0) +12f ′′(ω0)(ω − ω0)2 + ..., g(ω) = g(ω0) + ... (F.3)which inserted in to the integral equation givesI = e−xf(ω0)ˆ badωe−x2f ′′(ω0)(ω−ω0)2g(ω0) (F.4)149Making the change of variables ω = ω0 + y√2/f ′′(ω0) such thatI = e−xf(ω0)√2f ′′(ω0)ˆ b′a′dye−xy2g(ω0) (F.5)where a′ = (a−ω0)√f ′′(ω0)/2 and correspondingly for b′. The parameter x is very large,so the reparameterised Guassian is also very narrow with center at y = 0 this time. Thismeans we can extend the limits to infinity. Thus,I = e−xf(ω0)√2f ′′(ω0)ˆ ∞−∞dye−xy2g(ω0) (F.6)Making a further change of variable z =√xy yieldsI = e−xf(ω0)√2xf ′′(ω0)ˆ ∞−∞dze−z2g(ω0) (F.7)and the Guassian integral can be performed to giveI = e−xf(ω0)√2pixf ′′(ω0)g(ω0) (F.8)which is the result we use in the text.150Appendix GPure dephasing dynamics in the3-site-V systemI mention briefly in Section 3.4 of the main text how, in the limit of zero tunnelling∆10 = ∆20 = ∆ = 0 (for the 3-site-V system), the system exhibits ‘pure dephasing’. Thiscorresponds to no change in the diagonal elements of the density matrix, but a decay inthe off-diagonal elements despite the absence of tunnelling between states. This is entirelydue to the fluctuating on-site energies induced by the bath, and leads to decoherence inthe system without energy dissipation. As this is an exactly solvable model, with justthree copies of the well known Independent-boson model, we merely comment on it in thetext and do not quote the corresponding calculations for the dynamics. Instead I includethem here for reference. The system is easily solved by transforming the Hamiltonianwith the Lang-Firsov unitary operator U = eS, where S = −∑αq uαq(bq − b†q) |α〉〈α| anduαq = λαq/ωq such that H → UHU † = eSHe−S. Physically, this represents the shiftingof the boson cloud to its new equilibrium position. So S can be thought of as a shiftoperator. Here we find (using the Baker-Campbell-Hausdorff formula)UHU−1 =∑qωqb†qbq −∑quαq (G.1)The dynamics of the coherences are thus〈Pαβ(t)〉 = TrB[Uρs ⊗ ρBU−1UPαβU−1](G.2)where we’ve assumed an initial spin-bath factorised state and inserted the identity 1 =UU−1. Defining Φα = −2i∑q uαq(bq − b†q) we can apply the BCH formula again to state151Pαβ = |α〉〈β| to evaluateUPαβU−1 = |α〉〈β|+[−i∑α|α〉〈α| Φα2, |α〉〈β|]+ ...= |α〉〈β|(1 + i(Φβ/2− Φα/2) + ...)(G.3)eventually one finds (including the biases now)〈Pαβ(t)〉 = Re e−i(α−β)e−Q′′(t)Q′′αβ(t) =ˆdω|Jα(ω)− Jβ(ω)|ω2(1− cos(ωt)) coth(βω2)(G.4)where we’re interested in the probability of occupancy of the off-diagonal states andso inspect the real part. We see that the bath correlation term above corresponds tothe real part of the reparameterised Feynamn-Vernon influence function as discussedin the introduction. It is indeed this term that is responsible for dephasing owing tothe self-interaction of the off-diagonal paths in the density matrix. The coherences areexponentially damped in this regime with the rate depending on the difference of thespectral densities corresponding to the states constituting the coherent superposition.The bath correlation function Q′′ determines the rate Γαβϕ , which is the pure dephasingrate: the rate at which the off-diagonal elements are suppressed without relaxation effectsdue to tunnelling processes. When tunnelling processes are included, the decay of theoff-diagonal density matrix elements will include a combination of the pure dephasingrate and relaxation rate. For now however, without tunnelling included, the decoherencerate is equivalent to the dephasing rate.The spectral densities are unspecified in the above but generally there will be somedependence of the coupling on the position of the wells. Since the dephasing rate isevidently dependent on the difference between spectral densities, it will also depend onthe well separation and as this increases the dephasing rate does accordingly. This makessense physically as for larger well separation we expect the overlap between wavefunctionsto decrease and coherences to become less pronounced.For the Ohmic case (s=1) we have for the spectral density Jα(ω) = ηαωeω/ωc such thatthe frequency integral can be formed in the bath correlation function to get152Q′′αβ(t) =12γαβ ln(1 + ω2c t2) + γαβ ln[βpitsinh(pitβ)](G.5)where γαβ = (ηα − ηβ)2/pi~. In the various limits we find for the dynamics|〈Pαβ(t)〉| = 1[1 + (ωct)2]γαβ/2, T = 0K (G.6)where for short times below the characteristic time-scale of the bath|〈Pαβ(t)〉| = 1− γαβ (ωct)22, ∀ ωct 1, T = 0K (G.7)For intermediate times we observe the power-law governed dephasing|〈Pαβ(t)〉| = (ωct)−γαβ ∀ 1/ωc t, T = 0K (G.8)So we see that at zero temperature we still have dephasing but at a slower than exponen-tial rate; instead the off-diagonal elements decay algebraically. Recall that when γα = γβthe bath decouples from the central system and we expect the decoherence rate in thiscase to go to zero. When we include non-zero temperatures|〈Pαβ(t)〉| = [(~β/pit) sinh(pit/~β)]−γαβ[1 + (ωct)2]γαβ/2(G.9)We can once again inspect for short times. First we look at the low-T case where weexpand in the dimensionless quantity pit/~β << 1, which essentially defines the limit ofsmall bath fluctuations relative to tunnelling energy i.e. for t ∼ 1/∆, kBT/~∆ 1|〈Pαβ(t)〉| =[1[1 + (ωct)2]γαβ/2][(pit~β)−γαβ−1− γαβ6(pit~β)−γαβ+1 ], ∀ pit/~β 1|〈Pαβ(t)〉| =[1[1 + (ωct)2]γαβ/2](pit2~β)exp(−γαβpit/~β), ∀ pit/~β 1 (G.10)For the super-Ohmic case (s=3) we have for the spectral density Jα(ω) = ρα(ω3/ω2c )eω/ωc .Calculating the bath-correlatorQ˜′′ph(t) =ˆ ∞0dω ωe−ω/ωc [1− cos(ωt)] coth(~βω/2) (G.11)153we find for the zero-temperature caseQ˜′′ph(t) =ω2c (ωct)2[3 + (ωct)2][1 + (ωct)2]2 for T = 0K (G.12)For short and long times respectively (to 2nd-order in the respective small parameters)Q˜′′ph(t) = 3ω4c t2 for ∀ T = 0K,ωct 1Q˜′′ph(t) = ω2c(1 + (ωct)2(ωct)2)∀ T = 0K,ωct 1 (G.13)For non-zero temperatures we get a temperature dependent correction and the integralis computed to beQ˜′′ph(t) = −ω2c (ωct)2[3 + (ωct)2][1 + (ωct)2]2 − 14(~β)2[2ψ′(1/ωc~β)− ψ′(1− iωct~ωcβ)− ψ′(1 + iωct~ωcβ)](G.14)which still contains the small parameter ωc~β 1 intrinsic to the continuous-discretesystem truncation procedure and would be inconsistent to relax it now. We findQ˜′′ph(t) =ω2c (ωct)2[3 + (ωct)2][1 + (ωct)2]2 + pi4t2[1 + (ωct)2]2240(~β)4[1 + (ωct)2]2 (G.15)We can explore the same limiting time-domains as beforeQ˜′′ph(t) = 3ω4c t2 +pi4t2240(~β)4for ∀ ~∆/kBT 1, ωct 1Q˜′′ph(t) = ω2c(1 + (ωct)2(ωct)2)+pi4t2240(~β)4∀ ~∆/kBT 1, ωct 1 (G.16)where in both cases we obtain a temperature dependent correction that goes as T 4. Sofor the case of short times and low temperatures, the loss of coherence is algebraic intime154|〈Pαβ(t)〉| ∼ 1− 3γαβω4c t2 −pi4γαβt2(~β)4(G.17)and for long times and/or high temperatures, the loss of coherence is of course exponentialin time.155
You may notice some images loading slow across the Open Collections website. Thank you for your patience as we rebuild the cache to make images load faster.
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Modelling exciton dynamics in light-harvesting molecules
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Modelling exciton dynamics in light-harvesting molecules Ruocco, Leonard 2019
pdf
Notice for Google Chrome users:
If you are having trouble viewing or searching the PDF with Google Chrome, please download it here instead.
If you are having trouble viewing or searching the PDF with Google Chrome, please download it here instead.
Page Metadata
Item Metadata
Title | Modelling exciton dynamics in light-harvesting molecules |
Creator |
Ruocco, Leonard |
Publisher | University of British Columbia |
Date Issued | 2019 |
Description | I investigate the dynamics of multi-state central systems coupled bilinearly to an external oscillator bath within the noninteracting-blip approximation. I focus on both a 3-site configuration, as well as a 2-site model for the central systems of interest. The 2-site model, dubbed the dual-coupling spin-boson (DCSB) model, includes both diagonal and non-diagonal system-bath couplings, whereas the 3-site model considers only diagonal couplings. The bath spectral densities considered in this work include both Ohmic and super-Ohmic forms, as well as single optical phonon peaks. This work is motivated by the recent observance of long-lived quantum coherence effects in the photosynthetic organism known as the Fenna-Matthews-Olson (FMO) complex. The models investigated in this thesis are applied to this system in an attempt to explain its remarkably efficient exciton transfer mechanism, as well as to shed light on the functionality of coherence. The DCSB model is shown to reproduce the rapid exciton transfer times as well as the relatively long coherence times observed in the FMO complex. The non-diagonal system-bath coupling is shown to play a crucial role in this process. This can be attributed to the inelastic phonon-assisted tunnelling (IPAT) mechanism arising from the presence of significant non-diagonal system-bath interactions. Conversely, the 3-site model predicts rapid but incoherent exciton transfer. This can be attributed to the presence of a resonant state in the 3-site architecture, resulting in a relatively slow exciton transfer mode in the system. Therefore efficient exciton transfer requires a careful configuration of the chromophore energy landscape to avoid a resonant 3-site-V configuration. Furthermore, I conclude that coherence effects arising from excitons delocalised across multiple chromophores, promotes IPAT processes arising from non-diagonal system-bath couplings, producing rapid exciton transfer between chromophores. This offers a potential explanation as to the functional role that coherence plays in the energy transfer mechanism of photosynthesis. |
Genre |
Thesis/Dissertation |
Type |
Text |
Language | eng |
Date Available | 2020-02-28 |
Provider | Vancouver : University of British Columbia Library |
Rights | Attribution-NonCommercial-NoDerivatives 4.0 International |
DOI | 10.14288/1.0376448 |
URI | http://hdl.handle.net/2429/68368 |
Degree |
Doctor of Philosophy - PhD |
Program |
Physics |
Affiliation |
Science, Faculty of Physics and Astronomy, Department of |
Degree Grantor | University of British Columbia |
GraduationDate | 2019-05 |
Campus |
UBCV |
Scholarly Level | Graduate |
Rights URI | http://creativecommons.org/licenses/by-nc-nd/4.0/ |
AggregatedSourceRepository | DSpace |
Download
- Media
- 24-ubc_2019_may_ruocco_leonard.pdf [ 3.11MB ]
- Metadata
- JSON: 24-1.0376448.json
- JSON-LD: 24-1.0376448-ld.json
- RDF/XML (Pretty): 24-1.0376448-rdf.xml
- RDF/JSON: 24-1.0376448-rdf.json
- Turtle: 24-1.0376448-turtle.txt
- N-Triples: 24-1.0376448-rdf-ntriples.txt
- Original Record: 24-1.0376448-source.json
- Full Text
- 24-1.0376448-fulltext.txt
- Citation
- 24-1.0376448.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
data-media="{[{embed.selectedMedia}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.24.1-0376448/manifest