We numerically isolate the limits of validity of the Landauer approximation to describe charge transport along molecular junctions in condensed phase environments. To do so, we contrast Landauer with exact time-dependent non-equilibrium Green’s function quantum transport computations in a two-site molecular junction subject to exponentially correlated noise. Under resonant transport conditions, we find Landauer accuracy to critically depend on intramolecular interactions. By contrast, under nonresonant conditions, the emergence of incoherent transport routes that go beyond Landauer depends on charging and discharging processes at the electrode–molecule interface. In both cases, decreasing the rate of charge exchange between the electrodes and molecule and increasing the interaction strength with the thermal environment cause Landauer to become less accurate. The results are interpreted from a time-dependent perspective where the noise prevents the junction from achieving steady-state and from a fully quantum perspective where the environment introduces dephasing in the dynamics. Using these results, we analyze why the Landauer approach is so useful to understand experiments, isolate regimes where it fails, and propose schemes to chemically manipulate the degree of transport coherence.

## I. INTRODUCTION

The Landauer formula for electron transport plays a central role in nanoscience and molecular electronics.^{1,2} It is routinely used to atomistically understand the ability of molecules in junctions to transport charge as it provides a connection between the observable current and the molecular structure.^{3–12} This formula enables the use of conductance measurements made on single molecules bridging the gap between two electrodes to investigate fundamental chemistry and physics at the nanoscale, including switching,^{9,13} rectification,^{14} quantum interference,^{15–17} and chemical reactivity.^{8,18}

### A. The emergence of incoherent transport routes

In the Landauer approach, the current that passes through a nanojunction as a voltage is applied is due to quantum tunneling of electrons and is solely based on considering the steady-state contributions to the current. The approach is based on a single-particle scattering theory and can be seen as the limiting case of the general Jauho–Meir–Wingreen formula^{19} for electron transport where electron–electron interactions are captured in a mean field way. Landauer currents are regarded as coherent in nature as they emerge from quantum tunneling events. Contributions to the net current that go beyond Landauer and that arise due to interactions of the molecular junction with solvent, vibrations, or other thermal environments are often referred to as incoherent transport routes, as they involve events beyond purely electronic quantum tunneling in which the environment assists or suppresses transport routes.^{15,20–24} A well-known example of this mode of transport is thermally activated hopping,^{25–27} in which the charge transiently localizes on molecular sites and its kinetics is determined by thermal rate equations. Other possible contributions include environment-mediated electron back-scattering, tunneling, and excitation of the molecular bridge. Interestingly, because the observable current in molecular junctions is a purely electronic quantity, the setup is now being used to directly interrogate electronic coherences.^{28}

### B. The challenge in modeling and interpreting experiments

Recent theoretical progress has yielded several techniques that are able to capture incoherent transport contributions. This includes the Time-Dependent Non-Equilibrium Green Function (TD-NEGF) method that solves the quantum master equation for molecules in junctions,^{29–31} Landauer–Büttiker probes^{26,32–38} that introduce decoherence via fictitious thermal electrodes that accept charge from the molecule and inject it back at the same energy but with a random quantum mechanical phase, and the introduction of vibrational self-energies to steady-state formulations of quantum transport.^{39–41} These advances complement efforts to capture electronic many-body effects to the current.^{42–44}

Despite this progress, state-of-the-art experiments are often interpreted solely based on Landauer considerations. This is partly because the Landauer approach often provides a qualitatively useful explanation of experimental observations and because it has been efficiently implemented in a variety of quantum mechanical simulation software packages that enable first-principle computations on atomistically detailed models.^{45–48}

Nevertheless, most experiments in molecular junctions are performed at room temperature and in solvent. As such, the experimentally relevant situation is, in principle, beyond Landauer as the thermal environment can introduce or suppress transport routes.

Experimentally assessing the emergence of incoherent transport directly remains challenging. For example, a common strategy is to analyze the dependence of the charge transport with increasing temperature and/or molecular length.^{49–54} Coherent transport routes are often thought to have no temperature dependence and to exponentially decay with molecular length as they emerge from quantum tunneling. By contrast, incoherent transport is thought to be thermally activated and to have a weaker (polynomial) dependence with molecular length. However, it has been noted that both coherent and incoherent routes have a strong temperature and length dependence^{32,55} and an assignment of the transport mechanism requires detailed computations. For example, in molecular junctions, the conductance can exhibit an Arrhenius-type temperature activated behavior even in the coherent limit due to thermal population of energy levels in the metal that are at resonance with molecular levels.

For these reasons, presently, it is unclear when the Landauer approximation provides a quantitative description of quantum transport for molecular junctions in solvent and when do incoherent transport routes dominate. Understanding when the Landauer approximation is accurate and when it fails is essential to design useful modeling strategies for molecular junctions. Furthermore, it is key to establish how molecular structure determines the emergence of coherent, incoherent, and intermediate charge transport mechanisms. This is needed to understand the operation of photovoltaic cells, redox reactions, and charge transfer processes in biological systems^{56–60} and to design useful chemical strategies to control incoherent and coherent processes in molecules.

### C. Ongoing efforts

In the last few years, our groups have made attempts at understanding the utility of Landauer considerations for molecular junctions in thermal environments. This has included analyses based on driven junctions,^{61,62} noisy junctions,^{63} and the incorporation of vibrational dephasing in NEGF schemes.^{41} These efforts complement programs in both our groups understanding and modeling decoherence processes in molecules and developing methods to capture the dynamics of open quantum systems.^{64–71}

Here, we join forces and provide a systematic analysis of the relevance of coherent and incoherent transport routes. We are motivated by the following:

The lack of clarity of why the Landauer approach is so successful in providing a qualitatively useful description of state-of-the-art experiments, even when these experiments are performed in thermal environments where the decoherence is unavoidable. A notable example is the clear observation of quantum interference effects in molecular transport

^{15,16}even at room temperature and in complex self-assembled monolayers^{72,73}where the decoherence is expected to play an important role.The incisive analysis by Kilgour and Segal

^{32}that demonstrated that the signatures that we usually associate with incoherent transport (strong temperature dependence and weak dependence on molecular length) can also emerge in fully coherent mechanisms. Thus, new qualitative rules to interpret experiments need to be developed.Recent experiments in DNA that suggest that it is possible to chemically engineer molecular transport coherence.

^{74}Interpreting this emerging class of experiments requires a detailed understanding of the importance of incoherent transport routes.A lack of theoretical studies systematically investigating the validity of the Landauer approach for molecular junctions in thermal environments even when it is the main equation in the field of molecular electronics.

^{1,2}The well-documented quantitative disparity between the experimental and computed quantum conductance in molecular junctions^{75}has prevented the use of experiments to assess the validity of Landauer.

At a technical level, this initiative has faced the following challenges:

The computationally convenient Landauer–Büttiker approach to introduce decoherence is phenomenological.

^{26,32–38}Therefore, it is challenging to establish a physical map between this approach and solvent-induced decoherence.Efforts to assess the validity of Landauer through noise-based models of decoherence are limited by the challenge of converging the noisy statistics in a relative large parameter space.

^{63}This is particularly true for nonresonant transport where one has to capture small conductances*G*= 10^{−3}–10^{−6}*G*_{0}(where*G*_{0}= 2*e*^{2}/*h*is the quantum of conductance) in a setting where the noise can generate conductance fluctuations over several orders of magnitude.Rules isolated by us for the validity of Landauer are based on classical noise models of the environmental effects.

^{63}It has been unclear how to formally adapt these insights to a fully quantum setting where the environment is described by a thermal density matrix.Valuable recent advances in including decoherence effects through vibrational self-energies in the NEGF formalism

^{39–41}have not yet yielded intuitive rules for the emergence of incoherent transport that can help guide numerical and experimental progress.

### D. This work

Here, we investigate coherent and incoherent transport routes in a two-level molecular junction coupled to a thermal environment. We favor model computations, instead of a specific example, because they enable the systematic exploration of the parameter space and, thus, the characterization of coherent and incoherent transport routes with varying temperature, molecular interactions, electrode–molecule interactions, and bath correlation times. The computations are based on the Time-Dependent Non-Equilibrium Green’s function method (TD-NEGF) and on the representation of the quantum thermal environment through colored noise. Through them, we isolate experimentally relevant conditions where incoherent transport routes become important and elucidate key differences that emerge between resonant and nonresonant transport. Using these results, we analyze why the Landauer approximation is so useful to understand experiments and propose schemes to chemically manipulate the degree of transport coherence.

Landauer transport computations routinely used to understand experiments often suppose that the junction can be represented by its minimum energy or a few selected conformations. Recent advances in the description of the molecular dynamics in junctions have now enabled us to investigate the transport contributions of all molecular conformations in the thermal ensemble by performing Landauer transport computations on each conformation and computing averages.^{3,4,10,76–78} This procedure yields a wide dispersion in conductances within the thermal ensemble and averages that can significantly differ from those obtained from minimum-energy geometries.^{79}

We seek to isolate transport behavior that cannot be captured by the Landauer approach even when considering ensemble averages.

Noise-based models naturally lead to a dynamical analysis of the emergence of incoherent transport routes. Recent advances in understanding when classical noise models can accurately mimic quantum decoherence phenomena^{80} now enable us to connect this quantum-classical perspective to fully quantum treatments.

## II. MODEL AND METHODS

We consider transport through a molecular junction in the presence of solvent, see Fig. 1(a). The junction is represented as a two-level system with energy difference 2|*γ*| [Fig. 1(b)] that can exchange charge and energy with the electrodes and experience additional decoherence due to the interaction with the solvent. In the “site” or localized molecular orbital representation, this system corresponds to two sites with energy *ϵ*_{1} = *ϵ*_{2} = 0 in tight-binding coupling *γ* among them. Each terminal site is connected to the electrode adjacent to it and exchanges charge at a rate Γ/*ℏ* (*ℏ* = *h*/2*π*, where *h* is the Planck constant). The molecule–solvent interaction modulates the site energies *ϵ*_{i} with a coupling strength determined by the reorganization energy *λ*. This model for the molecule–solvent interaction is well established in decoherence studies and corresponds to the displaced harmonic oscillator model (or spin-boson model) with solvent-induced pure dephasing contributions to the dynamics (where the solvent leads to coherence loss without net exchange of energy) and molecular electronic transitions.

We focus on high temperatures (for which *k*_{B}*T* ≫ *ℏω*, where *ω* are the solvent frequencies) where the solvent can either be represented as a set of quantum harmonic oscillators or, equivalently,^{41,80} as a source of exponentially correlated classical noise on the site energies *ϵ*_{i}(*t*) of strength *σ* and correlation time *τ*_{corr}. We opt to represent it as a source of classical noise, but the two models can be essentially equivalent in certain limits.

To determine incoherent and coherent contributions to transport, we adopt a dynamical perspective and capture the dynamics through the time-dependent (TD) non-equilibrium Green’s function formalism. The TD formalism solves the exact quantum dynamics (the Liouville–von Neumann equation) for molecules in junctions and, therefore, can capture the steady-state Landauer limit (LD) and any dynamical effects that go beyond LD. The decoherence is captured by following the dynamics of an ensemble of noisy trajectories in which the environment leads to stochastic diffusion of the site energies *ϵ*_{i}(*t*) and then taking ensemble averages of the observables. These trajectories are representative of the time-series of molecular snapshots that are encountered in classical molecular dynamics simulations of junction evolution in the presence of solvent.

To isolate incoherent transport contributions, we contrast the exact current *I*_{TD}(*t*) with that obtained by assuming that transport is at steady-state at each instant of time *t* and, thus, well described by the Landauer approach for each molecular snapshot encountered in the noisy dynamics *I*_{LD}(*t*). Specifically, we contrast the average currents obtained by both time and ensemble averaging,

where ⟨⋯⟩ denotes this type of average. The ensemble average is obtained by computing the current for an ensemble of noisy trajectories with different initial conditions and stochastic dynamics. When the agreement is equal to 1, Landauer is exact and incoherent routes do not play a role in the electron transport. When the agreement is larger or smaller than one, incoherent routes suppress or enhance the electron transport.

Details of the model Hamiltonian, the TD and LD method, the statistics of the noisy trajectories, and its equivalence with fully quantum models of decoherence are discussed below.

### A. Hamiltonian

The Hamiltonian for the two-level molecule interacting with the electrodes and solvent is given by

where *H*_{m}(*t*) describes the molecule and any time-dependence introduced by the solvent degrees of freedom in the molecular dynamics, *H*_{e} represents the electrodes, and *H*_{me} represents the molecule–electrode interaction. In the site representation, the molecular Hamiltonian is given as

with instantaneous eigenenergies

In turn, the electrodes are described by

where {*q*} label electrode states with {*ϵ*_{βq}} energies, for the left (*β* = *l*) and right (*β* = *r*) electrodes, and $c\beta q\u2020$ and *c*_{βq} are, respectively, the electron creation and annihilation operators for these states. Site *i* = 1 (*i* = 2) is assumed to interact only with the left (right) electrode such that the molecule–electrode interaction is

where $Viq\beta $ are the couplings between the site *i* and the electrode *β* and H.c. denotes the Hermitian conjugate. We focus on the wide-band limit^{81} where the effective site-electrode coupling is defined by the spectral density

where *ζ*^{β} is the density of states of the *β* electrode. Here, the rate of charge exchange Γ/*ℏ* between the molecule and the electrodes at the two boundaries is taken to be identical, that is, Γ_{1,l} = Γ_{2,r} = Γ.

### B. Stochastic trajectories

The effect of solvent is captured by subjecting the site energies to exponentially correlated noise of strength *σ* and correlation time *τ*_{corr}. That is, by modeling *ϵ*_{i}(*t*) as noisy trajectories. The initial conditions are sampled from an initial Gaussian distribution of the form^{82}

with mean $\u03f5i\u0304$. The fluctuations between sites are assumed to be uncorrelated such that

We investigate the influence of noise of varying *σ* and *τ*_{corr}. Unless stated otherwise, *τ*_{corr} = 40 fs and $\u03f51\u0304=\u03f52\u0304$. Figure 1(d) shows an example of the resulting noisy trajectory of the site energies (*σ* = 0.05 eV).

### C. Landauer vs TD-NEGF

For every set of trajectories {*ϵ*_{i}(*t*)}, the transport was computed using the Landauer (LD) formula and the TD-NEGF formalism (TD) as developed and implemented by Zhang *et al.*^{29} In the Landauer computations, the current *I*_{LD}(*t*) was obtained via

where *e* is the charge of the electron, *f*_{β} is the Fermi function of electrode *β* with chemical potential *μ*_{β} at 300 K, and *T*(*E*;*t*) is the transmission function at energy *E* for the molecular snapshot at time *t*. The transmission function is computed using the NEGF formalism, where *T*(*E*) = Tr[Γ_{l}*G*^{R}(*E*)Γ_{r}*G*^{A}(*E*)] and *G*^{(R/A)} corresponds to the retarded/advanced Green’s function. The time-dependence in *I*_{LD}(*t*) just reflects the fact that *T*(*E*;*t*) depends parametrically on time as the molecular site energies fluctuate over time due to solvent. However, Eq. (10) remains a steady-state approximation to the time-dependent transport.

In turn, the time-dependent current *I*_{TD}(*t*) was obtained by solving the Liouville–von Neumann equation for the reduced single-electron density matrix *ρ*(*t*) in the presence of electrodes,^{29,30,83}

where *H*_{m}(*t*) is the Hamiltonian of the noisy molecule and *φ*_{β} are auxiliary density matrices that describe charge injection and withdrawal processes from the molecular region. The latter can be expressed in terms of the lesser *G*^{<} and greater *G*^{>} Green’s functions of the molecule and the lesser $\Sigma \beta <$ and greater $\Sigma \beta >$ self-energies due to the coupling between the molecule and the electrode *β*. Specifically,^{83}

The first term in the right-hand side of Eq. (12) represents the outcoming rate of electrons from the molecule to the electrode *β*, and the second term represents the incoming rate of electrons from the same electrode to the molecule. Therefore, *φ*_{β} can be interpreted as the net rate of electron transfer between the molecule and electrode *β*. This allows us to compute the transient time-dependent current at the interface as

The net current passing through the junction is computed as the average current flowing through the two electrodes,

The computations, implemented as in Ref. 29, are based on the wide-band limit of the electrode–molecule interaction and employ a Padé expansion of the Fermi distribution functions to analytically solve the energy integrals in the definition of the self-energies. The number of Padé functions (20) and the integration time-step (*τ*_{corr}/2000) used in the Runge–Kutta four propagation of Eq. (11) were checked for convergence.

### D. Ensemble averages

For TD resonant transport, we computed averages over 150 trajectories each 1150 *τ*_{corr} (46 000 fs) in duration. Nonresonant transport requires significantly more effort to achieve convergence. In that case, we employ 500–2000 trajectories of 1150 *τ*_{corr} in duration. The initial segment of the dynamics where the molecule is equilibrating with the electrodes is not taken into account when computing averages. The Landauer transport computations use the same trajectories but require less sampling to achieve convergence. For the LD averages, we employ only the first 100 *τ*_{corr} (4000 fs) of each noisy trajectory.

## III. EQUIVALENCE BETWEEN CLASSICAL NOISE AND QUANTUM DECOHERENCE

Noise is a useful phenomenological model of decoherence since it can lead to a damping of coherences in the density matrix upon ensemble averaging.^{41} Nevertheless, there are known differences between quantum decoherence and these classical noise-induced processes, and a quantitative connection between the two was only recently established.^{80} For the noise in Eq. (9), Rahman *et al.*^{41} numerically demonstrated that the TD gives identical results to a steady-state NEGF approach with an additional vibrational self-energy acting on the sites. Thus, the TD results can be interpreted fully quantumly.

The noise considered here can introduce both dephasing and population exchange between the molecular eigenstates. In Sec. II B, the noise is expressed in the site basis where it introduces fluctuations of the site energies. In turn, in the basis of eigenstates of the pristine molecular system [with energies *E*_{1/2} = *ϵ*_{0} ± *γ* and eigenstates $E1/2=121\xb11$, where $\u03f50=\u03f51\u0304=\u03f52\u0304$],

where $\sigma x=0110$ and $\sigma z=100\u22121$ are Pauli matrices and *η*_{±}(*t*) = (*ϵ*_{1}(*t*) ±*ϵ*_{2}(*t*))/2 is the noise. From Eq. (9), it follows that $\u27e8\eta i(t)\eta j(t\u2032)\u27e9=\delta ij\sigma 2e\u2212(t\u2212t\u2032)/\tau corr$, where *i*, *j* = ±.

In the Markovian limit, Gu and Franco^{80} showed that upon ensemble averaging, this type of Hamiltonian leads to an identical dynamics to that predicted by the Lindblad equation (the general equation for the dynamics of open quantum systems in the presence of Markovian environments) except that it does not capture spontaneous emission processes. The Lindblad transition rates between the two states are $1\u210f2\u222b0tdt\u2032\eta \u2212(t\u2032)\eta \u2212(t)=\sigma 2\tau corr/\u210f2$ for *t* ≫ *τ*_{corr}.

Beyond Markovian dynamics, in the pure dephasing limit, Gu and Franco^{80} showed that exponentially correlated noise of the form in Eq. (9) in the fluctuations of the transition energy *E* with respect to the average Δ*E* = *E*(*t*) − ⟨*E*(*t*)⟩, $\u27e8\Delta E(t)\Delta E(0)\u27e9=\sigma \u03032\u2061exp\u2212t/\tau \u0303$, leads to a dynamics equivalent to that introduced by quantum harmonic environments. The equivalence holds in the high-temperature limit and for Ohmic environments with Lorentz–Drude spectral density

where $\omega 0=1/\tau \u0303$ is the cutoff frequency and $\lambda eff=\sigma \u03032/2kBT$. Below, we show that in the model, $\u27e8\Delta E(t)\Delta E(0)\u27e9=\sigma \u03032\u2061exp\u2212t/\tau \u0303$, where the effective noise strength $\sigma \u0303$ and correlation time will depend on model parameters.

Beyond these two limits where the correspondence between noise-driven models and fully quantum models has been formally established, Aghtar *et al.* have numerically shown^{84} that noise-driven models can provide a semiquantitative description of fully quantum models in the high temperature limit. In this limit, the failure to capture the thermalization of the system, which has been encountered when classical baths are considered,^{85–87} is avoided.

## IV. TRANSPORT BEYOND LANDAUER

Here, we identify physical conditions that lead to the emergence of incoherent transport routes that go beyond Landauer considerations. We focus on both resonant and nonresonant transport and isolate effects that emerge by changing molecular properties, chemical anchor groups, and solvent–molecule interaction strength. We do so in the context of a two-level molecular junction in the presence of decoherence due to solvent and by computing the exact current *I*_{TD}, its Landauer approximation *I*_{LD}, and contrasting the two. Figure 1 includes a scheme of the model. It consists of two sites with tight-binding coupling *γ*, coupling to electrodes Γ, whose site energies *ϵ*_{i}(*t*) fluctuate in time due to interaction with a solvent that introduces noise of strength *σ* and correlation time *τ*_{corr} into the dynamics. The relevant parameters of the model and their interpretation in both time and energy are summarized in Table I.

Parameter . | Energetical interpretation . | Dynamical interpretation . | Physical origin . |
---|---|---|---|

Γ | Molecular level broadening | Rate of electrode–molecule | Electrode–molecule coupling V and |

due to electrodes | charge transfer (Γ/ℏ) | electrodes’ density of states η (Γ = 2πη|V|^{2}) | |

γ | Bare molecular transition | Timescale molecular | Molecular electronic coupling γ |

energy 2|γ| | electron motion ℏ/2γ | ||

σ | Site level broadening | Spectral diffusion of site | Strength of site-environment coupling |

due to solvent | energy levels | λ/temperature ($\sigma =2\lambda kBT$) | |

τ_{corr} | Cutoff frequency of | Bath correlation time | Highest-frequency environmental mode |

environment 1/τ_{corr} | mode coupled to molecular sites |

Parameter . | Energetical interpretation . | Dynamical interpretation . | Physical origin . |
---|---|---|---|

Γ | Molecular level broadening | Rate of electrode–molecule | Electrode–molecule coupling V and |

due to electrodes | charge transfer (Γ/ℏ) | electrodes’ density of states η (Γ = 2πη|V|^{2}) | |

γ | Bare molecular transition | Timescale molecular | Molecular electronic coupling γ |

energy 2|γ| | electron motion ℏ/2γ | ||

σ | Site level broadening | Spectral diffusion of site | Strength of site-environment coupling |

due to solvent | energy levels | λ/temperature ($\sigma =2\lambda kBT$) | |

τ_{corr} | Cutoff frequency of | Bath correlation time | Highest-frequency environmental mode |

environment 1/τ_{corr} | mode coupled to molecular sites |

The details of the model and the method employed to capture the exact current dynamics *I*_{TD}, the Landauer current *I*_{LD}, and their agreement [Eq. (1)] are included in Sec. II. For resonant transport [Fig. 1(b)], the average molecular site energies ⟨*ϵ*_{i}⟩ = 0 eV and the molecular orbitals are inside the *μ*_{l} − *μ*_{r} = 1 eV bias window with *μ*_{l} = 0.5 eV for all *γ*. For nonresonant transport [Fig. 1(c)], ⟨*ϵ*_{i}⟩ = 2 eV and the molecular orbital energies are above the transport window *μ*_{l} − *μ*_{r} = 0.1 eV with *μ*_{l} = 0.05 eV.

Here, we analyze the time-dependence for the LD and TD currents under the influence of external driving. We first focus on simple periodic driving that exemplifies how external forces can drive transport out of steady-state and beyond the regime of validity of Landauer (Sec. IV A). The insights developed from this deterministic driving are then employed in Sec. IV B to analyze the complex thermal dynamics of the nanojunction and the dramatic differences that emerge for resonant and nonresonant transport.

### A. Model periodic driving

Consider a driven molecular junction where the molecular site energies oscillate periodically with a frequency of 40 fs and an amplitude of 0.05 eV, as shown in Figs. 2(a) and 2(b). The periodic switching between the limiting values of the site energies is obtained by means of a periodic sequence of tanh(*st*) functions. This form has the advantage that by varying *s*, one can obtain oscillations in a form close to sinusoidal [*s* = 0.2 fs^{−1}, (a), (c), and (e)] to almost step-like transitions [*s* = 10 fs^{−1}, (b), (d), and (f)]. Figure 2 shows the LD and TD currents for single trajectories, *ϵ*_{i}(*t*), for both resonant [(c) and (d)] and nonresonant [(e) and (f)] transport. To make the differences between LD and TD dramatic, we choose weakly coupled sites (*γ* = 0.001 eV). A case with a larger *γ* = 0.14 eV is included in Fig. S1 of the supplementary material. Overall, the driving takes transport out of steady-state, leading to TD and LD currents that can be very different for both resonant and nonresonant transport. While *I*_{LD} is strictly positive, *I*_{TD} can be positive or negative.

For *resonant* transport, the *I*_{LD} current consists of a series of peaks that coincide with the crossing times among the site energies. At these crossing times, the molecular orbitals are delocalized across the junction, leading to a large *I*_{LD}. For |*ϵ*_{1} − *ϵ*_{2}| ≫ *γ*, each molecular orbital localizes in a single site, leading to a small Landauer current. In turn, during the quantum dynamics, the system tries to recover the Landauer steady-state determined by the driving. However, the changes in the electronic structure due to driving prevent the system to arrive to the steady-state limit, thus leading to a TD current that is qualitatively and quantitatively different from LD. For example, for *s* = 0.2 fs^{−1}, $I\u0304LD/I\u0304TD=0.15$, where $I\u0304$ denotes the current’s temporal average. In addition, the faster the change (cf. *s* = 10 vs 0.2 fs^{−1}), the further away that the current is driven away from the steady-state, further contributing to the difference between LD and TD transport. In fact, for $s=10fs\u22121,I\u0304LD/I\u0304TD=0.02$.

For *nonresonant* transport, the magnitude of the LD current is small with respect to the TD current. The TD current oscillates around zero as the levels switch. These oscillations arise because of charging/discharging effects at the interface that lead to large differences between *I*_{TD} and *I*_{LD} at each instant of time. To demonstrate that these effects arise because of charging/discharging at the interface, we performed computations with *γ* = 0 for which $I\u0304TD=0$. The results shown in Fig. S2 show the characteristic oscillations in *I*_{TD}(*t*) as the molecular site energies oscillate. For *γ* ≠ 0, these charging/discharging contributions to the current do not completely cancel one another, leading to small but important contribution to the overall nonresonant current that is beyond LD. In this case, $I\u0304LD/I\u0304TD\u223c0.001$ for both *s* values studied.

These charging/discharging effects are also present under resonant conditions. To exemplify this, we performed resonant transport computations for *γ* = 0 that are included in Fig. S3. However, these charging/discharging effects are negligible with respect to the overall ballistic resonant current and can be disregarded in this case. This contrasts with nonresonant transport for which the effect is comparable in magnitude to the LD current.

Increasing *γ* to 0.14 eV (Fig. S1) leads to a dynamics where |*ϵ*_{1}(*t*) − *ϵ*_{2}(*t*)| ≪ 2*γ* for all times. In this case, for resonant transport, *I*_{LD}(*t*) and *I*_{TD}(*t*) currents approximately coincide at all times since the driving does not appreciably change the degree of localization of the molecular orbitals, leading to a dynamics that is well approximated by the steady-state behavior. In fact, the average agreement [Eq. (1)] in this case is 0.85–0.9. For nonresonant transport, while *I*_{TD}(*t*) and *I*_{LD}(*t*) differ at each time because of charging/discharging effects, they approximately coincide on average ($I\u0304LD/I\u0304TD=0.99$ for *s* = 0.2 fs^{−1} and 0.6 for *s* = 10 fs^{−1}).

These results, and the observed stark differences between the behavior of resonant and nonresonant transport, provide key elements to interpret coherent and incoherent transport routes in thermal junctions.

### B. Noisy environment

Consider now the case of thermal junctions. Figure 3 shows example trajectories for the LD and TD current during the thermal noise-driven dynamics under resonant (a) and nonresonant (b) conditions for varying *γ* (see Table I). Overall, as in Sec. IV A, the TD and LD current do not generally coincide. This is because the noise prevents the junction from achieving steady-state, can excite the molecule and generate dephasing, and introduces charging/discharging effects at the interface.

For *resonant transport*, to understand the differences and similarities between *I*_{TD}(*t*) and *I*_{LD}(*t*), it is useful to recall the expression for the resonant LD current across a two-level junction [Eqs. (13) and (15) in Ref. 63],

where *c*_{i} are the coefficients in the expansion of the molecular orbitals $\phi $ in terms of the sites $i$, $\phi =\u2211icii$. The steady-state current is determined by the product |*c*_{1}|^{2}|*c*_{2}|^{2}, which is identical for both orbitals. As the noise changes the site energies, it can change the character of the molecular orbitals. When $(\u03f52\u2212\u03f51)2/4\gamma 2\u226a1$, the molecular orbitals are delocalized across the junction, leading to a large LD current as both sites enter into the expansion. By contrast, when $(\u03f52\u2212\u03f51)2/4\gamma 2\u226b1$, the orbitals are mostly localized in one of the sites, leading to a small |*c*_{1}|^{2}|*c*_{2}|^{2} and, thus, LD current. This basic feature can be used to understand the LD and TD currents in Fig. 3(a). Specifically, for *γ* = 0.001 and 0.012 eV, the noise leads to frequent localization and delocalization of the molecular orbitals, leading to dramatic changes in the LD current. The peaks in the LD current coincide with events when $(\u03f52(t)\u2212\u03f51(t))2/4\gamma 2\u226a1$ and the valleys to events when $(\u03f52(t)\u2212\u03f51(t))2/4\gamma 2\u226b1$. In this case, the TD and LD currents differ dramatically because the noise prevents the junction from achieving steady-state as it is changing the degree of delocalization of the molecular orbitals and, thus, the magnitude of the steady-state currents by several orders of magnitude. In the TD dynamics, the junction tries to reestablish the steady-state, but the steady-state changes faster than it can do so, thus leading to large differences between the dynamics of TD and LD. By contrast, for *γ* = 0.288 eV, the noise does not lead to a dramatic change in the degree of delocalization of the molecular orbitals, the current stays close to steady-state, and TD and LD approximately coincide at all times.

The disparities between LD and TD arise when there is a strong change in the degree of delocalization and, thus, of the magnitude of the steady-state current. If the molecular orbitals remain localized or delocalized, the LD current will be a useful approximation to the TD current. These changes in the degree of localization due to the noise require the noise strength *σ* to be comparable in magnitude to the tight-binding couplings *γ*.

For *nonresonant* transport, *I*_{TD}(*t*) and *I*_{LD}(*t*) differ dramatically at each instant of time as the LD current is always positive and small, while the TD current can be positive or negative and exhibits large fluctuations relative to the average. As revealed by the dynamics in Sec. IV A, these large fluctuations arise because of charging/discharging at the interface that is not captured by the Landauer approach.

## V. EXAMINING THERMAL AVERAGES

Experiments in molecular junctions record currents that are time-averaged over microseconds and cannot resolve the time-dependent features discussed in Sec. IV. Thus, the experimentally relevant quantities are the average currents. To capture them, we average the current traces for the thermal junction over time and over an ensemble of initial conditions and contrast LD and TD currents.

### A. Resonant transport: Incoherent contributions are a molecular affair

#### 1. Agreement between LD and TD

##### a. Agreement for varying *γ**.*

Figure 4(a) shows ⟨*I*_{TD}⟩ and ⟨*I*_{LD}⟩ for increasing *γ*. Both LD and TD currents increase monotonically with *γ* and reach a maximum value of *I* = *e*Γ/*ℏ* [this value follows from Eq. (17) in the *γ* → *∞* limit]. The LD current is always smaller than the full current, as the TD can capture additional modes of transport that go beyond the steady-state. That is, Landauer provides a lower limit to the full current in a wide parameter range. In the parameter range consider, incoherent transport can be even the dominant contribution to the overall current. Landauer is an accurate description of the exact current for either large or small *γ*. Incoherent routes become the most important for intermediate values of *γ*, where *γ* ∼ *σ* or *γ* ∼ Γ. We will refer to the intermediate *γ* value that leads to the worst agreement between LD and TD for a given Γ and *σ* as *γ*^{⋆}.

##### b. Agreement for varying *σ**.*

As the noise strength increases [Fig. 4(b)], incoherent routes become increasingly important. The quantity *γ*^{⋆} is also shifted toward larger *γ* as *σ* increases. Therefore, with respect to this intermediate value, increasingly larger or smaller *γ* are required for LD to provide a more accurate representation of the overall current.

##### c. Agreement for varying Γ*.*

Decreasing the rate of charge transfer Γ/*ℏ* between the molecule and the electrode [Fig. 4(c)] makes the incoherent contributions to the overall current increasingly important for all *γ*. The quantity *γ*^{⋆} shifts toward smaller values as Γ decreases.

##### d. Agreement for varying *τ*_{corr}*.*

Decreasing the correlation time of the bath leads to worse agreement between LD and TD for all *γ*; see Fig. 5. The quantity *γ*^{⋆} remains approximately constant as *τ*_{corr} is changed.

#### 2. Dynamical interpretation

From a time-dependent perspective, the key to understand the importance of incoherent transport to the overall current is to recognize that Landauer is a steady-state approximation. Therefore, parameter regimes that decrease the ability of the system to recover such a steady-state lead to the emergence of incoherent transport contributions.^{63,88} For example, decreasing Γ increases the time that it takes for the junction to reestablish steady-state as it decreases the molecule–electrode charge transfer rate, leading to a worse agreement between LD transport and TD transport. Increasing the noise strength (*σ*) leads to larger current fluctuations during the thermal dynamics and, therefore, to the time needed to recover the steady-state, leading to a worse agreement between LD transport and TD transport. Increasing *τ*_{corr} leads to better agreement because the dynamical events at a molecular level are slower, and therefore, the junction has more time to reestablish the steady-state.

The need for sufficiently large Γ for LD to be an accurate representation of the dynamics has long been recognized.^{63} However, Γ is not the only relevant parameter, and the response of the molecule to the thermal noise is also key. To quantify that response, we follow the physical fluctuations of the current across the junction due to the thermal dynamics ⟨Δ*I*_{TD}⟩. Figure 6(a) shows ⟨Δ*I*_{TD}⟩ for varying *γ* and *σ*. We observe that for fixed Γ, the larger ⟨Δ*I*_{TD}⟩, the worse the agreement between TD and LD. This is because ⟨Δ*I*_{TD}⟩ signals how susceptible is the current to the thermal noise.

For *γ* ≳ *γ*^{⋆}, the observed behavior can be understood in terms of the ability of the noise to change the degree of delocalization of the molecular orbitals, as discussed in Sec. IV. Specifically, when *γ* ∼ *σ*, the noisy excursions of the site energies can lead to strong changes in the degree of delocalization, large changes in the steady-state current and in ⟨Δ*I*_{TD}⟩, and a poor agreement between LD and TD. By contrast, for *γ* ≫ *σ*, the noise does not change the degree of delocalization of the molecular orbitals, leading to only modest changes in the steady-state current, the current fluctuations ⟨Δ*I*_{TD}⟩, and a good agreement between LD and TD.

Surprisingly, for *γ* < *γ*^{⋆}, we observe a recovery of the agreement between LD and TD even when the dynamics of the two currents do not coincide [Fig. 3(a)]. That is, the LD and TD only agree on average but not at each instant of time. The reason for this unexpected behavior is that the current fluctuations decrease with *γ* < *γ*^{⋆} because the overall current is small. That is, in this regime, the noise can take the junction out of the steady-state, but *γ* is so small that the overall changes to the current are only modest, making it easier for the junction to reestablish steady-state behavior.

To better understand these three transport regimes (*γ* ≪ *γ*^{⋆}, *γ* ∼ *γ*^{⋆}, and *γ* ≫ *γ*^{⋆}), we quantified the populations in the two sites of the molecular junction during the dynamics; see Figs. 7(a)–7(c). For *γ* ≪ *γ*^{⋆} and *γ* ≫ *γ*^{⋆}, the site populations remain approximately constant during the dynamics. For *γ* ≫ *γ*^{⋆}, the populations of the two sites are identical as the molecular orbitals are delocalized. For *γ* ≪ *γ*^{⋆}, only site 1 is populated as it is connected to the electrode of higher chemical potential and the molecular orbitals are spatially localized. By contrast, for *γ* ∼ *γ*^{⋆}, there is significant dynamical charge exchange between the two levels due to the noise. That is, the emergence of incoherent transport routes coincides with a physical regime where there is sequential transport across the molecule in agreement with hopping transport mechanisms. Isolating observable signatures of this charge dynamics can be the key to develop experimentally accessible routes to monitor incoherent transport routes.

We further analyze this dynamics in terms of Landau–Zener theory where the probability of charge exchange from site 1 to 2 during a crossing event is

where *ν* = |*d*(*ϵ*_{1} − *ϵ*_{2})/*dt*| is the velocity of crossing of the site energies. For the noisy environment, there will be a distribution of velocities at the crossing times. Using as an estimate *ν* = *σ*/*τ*, Fig. 7(a) is in the diabatic limit where there is no charge exchange between sites (*p*_{1→2} = 0.008) and Fig. 7(c) is in the adiabatic limit where all crossings lead to charge exchange *p*_{1→2} = 1.0. By contrast, Fig. 7(b) is in an intermediate regime where some crossings lead to charge exchange *p*_{1→2} = 0.672, leading to rich populations dynamics.

#### 3. Energetical interpretation

It is instructive to consider the comparison between LD and TD from a purely quantum perspective. From Eq. (15), the role of the bath is to introduce transitions and dephasing between the two molecular eigenstates. Since both levels are equally occupied in resonant transport, noise-induced transitions are suppressed, suggesting that the dynamics can be interpreted through a pure dephasing limiting model. In this limit, exponentially correlated noise of the eigenenergies is equivalent to the decoherence induced by a quantum harmonic environment with a Lorentz–Drude spectral density with effective reorganization energy $\lambda eff=\sigma \u03032/2kBT$ and cutoff frequency $\omega 0=1/\tau \u0303$ [see Eq. (16)]. To demonstrate that the eigenenergies (as opposed to the sites) are subject to exponentially correlated noise, Fig. 8(a) shows the time autocorrelation function ⟨Δ*E*(*t*)Δ*E*(0)⟩ for the energy gap fluctuations Δ*E*. As can be seen, in fact, $\u27e8\Delta E(t)\Delta E(0)\u27e9=\sigma \u03032\u2061exp\u2212t/\tau \u0303$ with effective noise strength $\sigma \u0303$ and correlation time $\tau \u0303$ that depend on model parameters. Here, we will use this equivalence to understand the results from an energetical perspective.

While the effective correlation time $\tau \u0303$ of the Δ*E*(*t*) dynamics is independent of *γ* [see Fig. 8(b)], the effective strength of the noise $\sigma \u0303=\u27e8\Delta E(0)\Delta E(0)\u27e9$ is modulated by the coupling between sites. As Eq. (16) indicates, increasing $\sigma \u0303$ increases the coupling strength between environmental modes and the molecule and thus leads to stronger effects due to the thermal environment. As shown in Fig. 8(c), increasing *γ* results in lower effective $\sigma \u0303$, which favors the agreement between the LD and TD currents. While the effective noise strength is maximum for *γ* → 0, in this regime, the magnitude of the current fluctuations is small [see Fig. 6(a)], leading to a recovery of the ⟨*I*_{LD}⟩/⟨*I*_{TD}⟩ agreement.

In turn, when Γ decreases, the disagreement between LD and TD increases as the frequency bandwidth of the environment contains increasingly important contributions at frequencies that cannot be resolved by the electrodes.

### B. Nonresonant transport: Incoherent transport routes are all about charging/discharging at the interface

#### 1. Agreement between LD and TD

##### a. Agreement for varying *γ**.*

Figure 4(d) shows ⟨*I*_{TD}⟩ and ⟨*I*_{LD}⟩ for varying *γ*. Both LD and TD currents increase ∝ *γ*^{2}. Contrary to the resonant case, the values of ⟨*I*_{TD}⟩ and ⟨*I*_{LD}⟩ do not reach a plateau with increasing *γ* in the range considered. The average TD is always larger than the LD counterpart. However, the difference between them remains approximately constant. As a result of this, for small values of *γ*, in which the currents are small, the importance of the incoherence contributions to the transport become dominant and the agreement ⟨*I*_{LD}⟩/⟨*I*_{TD}⟩ goes to zero as *γ* → 0 limit, as shown in Fig. 4(e).

##### b. Agreement for varying *σ**.*

As shown in Fig. 4(e), the transport incoherent routes become more important when increasing the strength of the noise *σ*. In all cases, the agreement between the LD and TD currents goes to zero for small values of the coupling between sites *γ*. For *γ* ≫ *σ* [data points to the right of the red marks in Fig. 4(e)], the agreement asymptotically approaches 1.

##### c. Agreement for varying Γ*.*

Increasing the rate of electrode–molecule charge transform favors coherent transport routes for all *γ*; see Fig. 4(f). Landauer becomes an increasingly poorer approximation as the molecule–electrode coupling decreases.

#### 2. Dynamical interpretation

Under nonresonant conditions, the charging and discharging of the sites due to their interaction with the electrodes as they thermally fluctuate take the system out of steady-state, resulting in incoherent contributions to the transport that can be dominant. To interpret these effects, we analyze how changing the junction model and bath parameters affects the ability of the system to recover steady-state transport.

The spectral density Γ determines the charge transfer rate at the electrode–molecule interface. Increasing Γ enhances the ability of the electrodes to time resolve these charging/discharging effects and reestablish the steady-state, leading to a better agreement between the LD and TD currents [see Fig. 4(f)].

The amplitude of the sites energy fluctuations increases with *σ*. This transiently brings the sites closer to the electrodes’ Fermi energy, leading to a larger charge transfer without changing Γ. As a result, the charging and discharging processes are more effective in taking the system out of the steady-state for large values of *σ* and the agreement between the LD and TD currents gets worse [see Fig. 4(e)].

Increasing *γ* leads to a slight increase in the incoherent contributions to the current as it makes the charge transfer between sites more effective [see Eq. (18)]. However, this increment is small [2 × 10^{−4} nA in Fig. 4(d)] with respect to the magnitude of the net current [5 × 10^{−3} nA in Fig. 4(d)]. Therefore, the incoherent contributions to the transport remain approximately constant when varying *γ*. Their importance decreases as the net Landauer current increases with *γ* [Figs. 4(e) and 4(f)]. By contrast, as *γ* goes to zero, the net current drops and the incoherent mechanisms that are activated by the charging and discharging of the sites are the only remaining contributions to the current, making the agreement ⟨*I*_{LD}⟩/⟨*I*_{TD}⟩ go to zero. As a consequence of the insensitivity of the incoherent transport routes to changes in the coupling between sites, the average fluctuations of the current are approximately constant when changing *γ* (see Fig. 6). In addition, the dynamics of the sites’ population under nonresonant conditions [Fig. 7(b)] is insensitive to varying the coupling between sites as both sites are unpopulated.

#### 3. Energetical interpretation

Since the molecular sites are unpopulated in this nonresonant model [see Figs. 7(d)–7(f)], dephasing and molecular transitions due to the thermal environment are muted. This is consistent with the observation that the incoherent contributions to the transport are small with respect to the net current and remain essentially unaltered by varying *γ*. Nonetheless, the environmental effects at the site-electrode interface remain and increase with *σ*. This results in worse agreements between the LD and TD current when increasing *σ* [see Fig. 4(e)].

In turn, as Γ decreases, the range of bath frequencies contributing to the charging/discharging that the electrodes can resolve decreases, leading to the deterioration of the ⟨*I*_{LD}⟩/⟨*I*_{TD}⟩ agreement [see Fig. 4(f)].

## VI. OUTLOOK

The Landauer formula is central to the computation of charge transport across molecular junctions using atomistically detailed models.^{1,2,47,48} This approach assumes that only coherent tunneling mechanisms solely determined by the molecular electronic levels contribute to the charge transport. However, realistic molecular junctions are always subject to the effects of thermal environments that can trigger incoherent transport routes.

Here, we have established well-defined limits in which the electronic currents across a molecular junction subject to a thermal environment can be quantitatively captured via the Landauer steady-state approximation. For this, we calculated the exact time-dependent non-equilibrium Green’s function current along a model two-site molecular junction, in which the site energies are subject to correlated noise, and contrasted it with that obtained from the Landauer approach. The analysis pertains to cases in which there is not back action of the current in the thermal environment. Therefore, it does not capture incoherent routes that imply molecular reorganization due to charging as those that have been identified in electrochemical events.^{27} The analysis also does not take into account additional effects that can arise when electron-correlation plays an important role in transport.^{42–44,47,48}

Under resonant conditions, the Landauer accuracy critically depends on how the intramolecular interactions compare with the noise strength. By contrast, under nonresonant conditions, the emergence of incoherent transport routes depends on charging and discharging processes at the electrode–molecule interface largely independent of intramolecular properties. In both cases, decreasing the rate of charge exchange between the electrodes and molecule (Γ/*ℏ*) and increasing the interaction strength with the thermal environment (*σ*) cause Landauer to become less accurate.

We interpreted the results from a dynamical and an energetical perspective. From a dynamical perspective, Landauer fails when the noise takes transport out of steady-state. This occurs when the electrode cannot resolve the timescale of molecular events (i.e., when *ℏ*/Γ ≫ *τ*_{corr}) or, for resonant transport, when the noise appreciably changes the degree of delocalization of the transport-determining molecular orbitals (i.e., when *σ* ∼ *γ*).

We used the equivalence between the exponentially correlated noise on the eigenenergies and the dephasing induced by a quantum harmonic environment with a Lorentz–Drude spectral density to interpret the results from an energetical perspective. We showed that as *γ* increases, the effective reorganization energy of the bath decreases, leading to improved agreement between LD and TD. Γ/*ℏ* defines the range of bath frequencies that can be resolve by the electrodes, and the agreement improves with increasing Γ.

An interesting challenge for future theoretical analyses would be to isolate intuitive conditions for the validity of Landauer in thermal environments starting from rigorous NEGF formulations^{89,90} of this many-body problem.

### A. Experimental relevance

It is instructive to consider realistic parameters to understand the relevance of these computational observations in experiments. The values of Γ can be determined experimentally by fitting the Landauer formula to temperature- or voltage-dependent currents.^{91,92} The reported values are surprisingly low in the 10^{−3}–1 meV range. For example, Γ = 0.003 meV for tercyclohexylidenes with sulfur anchor groups and Au electrodes^{91} and 0.17–0.55 meV for functionalized difurylethenes in Au junctions.^{92} For such low values of Γ, incoherent transport routes are expected for environments with *τ*_{corr} ≲ 1.2 ps (for Γ = 0.55 meV) and *τ*_{corr} ≲ 0.22 ns (for Γ = 0.003 meV). In turn, using a tight-binding model of a molecular junction, Zelovich *et al.*^{93} found that values of Γ ≈ 7–9 meV yield appropriate currents within their model. For these values of Γ, environments with *τ*_{corr} < 73–94 fs are expected to exhibit important incoherent transport.

Simulations of charge transport in DNA^{62,74,94} have assigned values of *σ* ∼ 0.1 eV, *γ* ∼ 0.03 eV, Γ = 0.001 eV, and *τ*_{corr} ∼ 200 fs, thus making *σ*/*γ* ≫ 1 and Γ/*ℏτ*_{corr} ≪ 1, which clearly suggests that the incoherent transport routes are important and, therefore, Landauer would fail as an approximation to compute the charge transport.

### B. Why is Landauer useful?

Based on this analysis, incoherent transport routes should be important in a wide variety of chemical systems. Nevertheless, the Landauer formula has proven to be extremely useful to interpret and predict experimental trends in molecular electronics.^{95,96} Why do incoherent contributions to the transport appear to be of low relevance in experiments?

The key to understand why Landauer remains useful is to realize that most experiments are done under low bias voltages to prevent the damage of the molecular junction due to heating. This constrains the transport to be under nonresonant conditions as most of the transporting molecular levels are out of the bias windows. As shown here, in this regime, the incoherent transport routes are mostly modulated by the charging and discharging in the electrode–molecule interface and their contribution to the net current remains approximately constant with varying molecular parameters. Thus, while this effect can have a quantitative influence on the current, it will keep trends such as those observed with varying molecular length qualitatively correct.

In addition, measured transport typically takes place through chemical bonds (as opposed to through space). This guarantees that the intramolecular couplings are strong (*γ* ≫ *σ*), bringing the transport to a regime in which the coherent routes are dominant. Exceptions to this would be situations in which the electron donor–acceptor pairs are spatially separated and through-space transport contributions are dominant.

### C. Toward the control of transport coherence

The parameters that determine the degree of transport coherence can be controlled by chemical or mechanical means. For instance, the spectral function Γ depends on the coupling between molecular sites and electrode levels. The strength of these couplings can be chemically tuned by changing the nature of the anchor groups in the molecular junction. Anchor groups that bond covalently to the electrodes are expected to favor coherent transport routes. By contrast, anchor groups that interact weekly with the electrodes (like in physisorption) are expected to trigger contributions from incoherent transport routes. Furthermore, for a given anchor group, it is possible to modify Γ by changing the applied bias voltage^{91} or by weakening the electrode–molecule interactions through mechanical pulling of the junction.

The analysis opens up the possibility to mechanically engineer transport coherence. Specifically, the coupling between sites (*γ*) can be mechanically manipulated. When the transport takes place through molecular electron donor–acceptor complexes that are not covalently bound, it is possible to mechanically change their coupling through pulling^{13,97} and thus tune and stabilize molecular complexes that can sustain coherent, incoherent, and intermediate transport routes. The strength of the noise can be controlled simply by changing the temperature or the chemical nature of the thermal bath.

Overall, by establishing well-defined limits in which coherent, incoherent, and intermediate transport mechanisms operate in molecular junctions subject to thermally fluctuating environments, this analysis offers the technical means to chemically and mechanically control transport coherence.

## SUPPLEMENTARY MATERIAL

See the supplementary material for charge transport computations under periodic driving for strongly coupled and uncoupled two-level model systems.

## ACKNOWLEDGMENTS

This work was supported by the National Science Foundation under Grant Nos. CHE-1553939 and CHE-2102386.

## AUTHOR DECLARATIONS

### Conflict of Interest

The authors have no conflicts to disclose.

## DATA AVAILABILITY

The data that support the findings of this study are available from the corresponding author upon reasonable request.