Ab Initio Car-Parrinello Molecular Dynamics Research Articles

Charge carrier dynamics in conducting polymer PEDOT using abinitio molecular dynamics simulations.

As conducting polymers become increasingly important in electronic devices, understanding their charge transport is essential for material and device development. Various semi-empirical approaches have been used to describe temporal charge carrier dynamics in these materials, but there have yet to be any theoretical approaches utilizing abinitio molecular dynamics. In this work, we develop a computational technique based on abinitio Car-Parrinello molecular dynamics to trace charge carrier temporal motion in archetypical conducting polymer poly(3,4-ethylenedioxythiophene) (PEDOT). Particularly, we analyze charge dynamics in a single PEDOT chain and in two coupled chains with different degrees of coupling and study the effect of temperature. In our model we first initiate a positively charged polaron (compensated by a negative counterion) at one end of the chain, and subsequently displace the counterion to the other end of the chain and trace polaron dynamics in the system by monitoring bond length alternation in the PEDOT backbone and charge density distribution. We find that at low temperature (T = 1K) the polaron distortion gradually disappears from its initial location and reappears near the new position of the counterion. At the room temperature (T = 300K), we find that the distortions induced by polaron, and atomic vibrations are of the same magnitude, which makes tracking the polaron distortion challenging because it is hidden behind the temperature-induced vibrations. The novel approach developed in this work can be used to study polaron mobility along and between the chains, investigate charge transport in highly doped polymers, and explore other flexible polymers, including n-doped ones.

Competitive adsorption between sodium citrate and naphthenic acids on alumina surfaces: Experimental and computational study

The adsorption of naphthenic acids (NAs) on clay surfaces induces an increase in surface hydrophobicity, which is detrimental to crude oil recovery and froth quality. Sodium citrate (Na3Cit), with more functional groups than NAs, is expected to compete with NAs for adsorption on mineral surfaces. The adsorption of NAs on alumina surfaces was quantified in the presence and absence of Na3Cit at pH 8 ± 0.1. An increase in hydrophobicity of the alumina surface was observed after the adsorption of NAs without Na3Cit. However, adding Na3Cit significantly reduced NAs adsorption and protected the surface wettability. The viscoelasticity and electrokinetic properties of adsorbed materials indicate Na3Cit was preferentially covered on alumina surfaces and repelled NAs by its negative charges. The competition between Na3Cit and NAs for adsorption sites was explained using ab initio Car-Parrinello molecular dynamics (CPMD) simulation and periodic plane-wave density functional theory (DFT) to calculate inner-sphere adsorption and outer-sphere adsorption, respectively. Favourable adsorption of sodium citrate than NAs representatives on aluminum (hydr)oxide cluster and the basal and edge planes of gibbsite were observed from the perspective of free energy barrier and adsorption energy, respectively.

Molecular dynamics study of structure and reactions at the hydroxylated Mg(0001)/bulk water interface.

A molecular level understanding of the aqueous Mg corrosion mechanism will be essential in developing improved alloys for battery electrodes, automobile parts, and biomedical implants. The structure and reactivity of the hydroxylated surface is expected to be key to the overall mechanism because (i) it is predicted to be the metastable surface state (rather than the bare surface) under a range of conditions and (ii) it provides a reasonable model for the outer corrosion film/water interface. We investigate the structure, interactions, and reactivity at the hydroxylated Mg(0001)/water interface using a combination of static Density Functional Theory calculations and second-generation Car-Parrinello ab initio molecular dynamics. We carry out detailed structural analyses into, among other properties, near-surface water orientations, favored adsorption sites, and near-surface hydrogen bonding behavior. Despite the short timescale (tens of ps) of our molecular dynamics run, we observe a cathodic water splitting event; the rapid timescale for this reaction is explained in terms of near-surface water structuring lowering the reaction barrier. Furthermore, we observe oxidation of an Mg surface atom to effectively generate a univalent Mg species (Mg+). Results are discussed in the context of understanding the Mg corrosion mechanism: For example, our results provide an explanation for the catalytic nature of the Mg corrosion film toward water splitting and a feasible mechanism for the generation of the univalent Mg species often proposed as a key intermediate.

Car-Parrinello molecular dynamics study of CuF, AgF, CuPF6 and AgPF6 in acetonitrile solvent and Cluster-Continuum calculation of the solvation free energy of Cu(I), Ag(I) and Li(I)

The CuF species in acetonitrile solution is a potential intermediate in copper catalyzed fluorination. In the process, different species such as Cu+, Ag+, CuF and AgF could be involved. The objective of this work is to provide relevant data on the structure of these species in solution phase, as well as the corresponding solvation thermodynamics to better understand the equilibrium for formation of CuF. To attain this goal, ab initio Car-Parrinello molecular dynamics method and cluster-continuum quasichemical approach were used. It was found that Cu+ and Ag+ ions are tetracoordinated by acetonitrile molecules, being possible to observe the release of one coordinated acetonitrile to the solution to form a tricoordinate structure, suggesting a dissociative mechanism for ligand exchange. The CuF species was found to form a predominant tricoordinate complex, with copper coordinated to two acetonitrile molecules. The AgF species have presented more labile ligand feature, with silver ion interconverting between bi- and tricoordinate structure, involving coordination to one and two acetonitrile molecules, respectively. The solvation free energy calculated for Cu+ and Ag+ are in good agreement with available experimental data. The present results indicate that CuF can be easily formed from dissolved AgF and Cu+ species in acetonitrile solution. Consequently, CuF could play an important role in copper catalyzed fluorination.

Using Car-Parrinello simulations and microscopic order descriptors to reveal two locally favored structures with distinct molecular dipole moments and dynamics in ambient liquid water

Water is essential for life and technological applications, mainly for its unique thermodynamic and dynamic properties, often anomalous or counterintuitive. These anomalies result from the hydrogen-bonds fluctuations, as evidenced by studies for supercooled water. However, it is difficult to characterize these fluctuations under ambient conditions. Here, we fill this knowledge gap thanks to the Car-Parrinello ab initio molecular dynamics (MD) simulation technique. We calculate the local structural order parameter ζ, quantifying the coordination shells separation, and find two locally-favored structures or states: High-ζ and Low-ζ. On average, High-ζ molecules have a tetrahedral arrangement, with four hydrogen bonds, and the first and the second coordination shell well separated. The Low-ζ molecules are less connected, partially merging the first and the second shells. The appearance of isosbestic points in the radial distribution functions and the collective density fluctuations at different length scales and timescales reveal that the two-state model, consistent with available experimental data for supercooled water, also holds under ambient conditions, as we confirm by analyzing the vibrational spectrum of both types of water molecules. Significant consequences of the structural differences between the two states are that High-ζ molecules have a dipole moment 6 % higher than Low-ζ. At the same time, Low-ζ structures are more disordered and with more significant angular fluctuations. These differences are also reflected in the dynamics under ambient conditions. The Low-ζ molecules decorrelate their reorientation faster than High-ζ and merge their coordination shells within 0.2 ps, while the High-ζ preserve the shell separation for longer times. Our analysis shows that first-principle calculations make predictions under ambient conditions calling for new and definitive experiments confirming the two-state model.

Behaviors of sodium and calcium ions at the borosilicate glass-water interface: Gaining new insights through an ab initio molecular dynamics study.

We study reactivity and leaching at the calcium sodium borosilicate (CNBS)-water interface by means of a Car-Parrinello ab initio molecular dynamics simulation over a simulation time of 100ps. With an emphasis on the comparison between the behaviors of Ca2+ and Na+ cations at the CNBS glass-water interface, different mechanism events during the trajectory are revealed, discussed, and correlated with other density functional theory calculations. We show that Na+ ions can be released in solution, while Ca2+ cannot leave the surface of CNBS glass. This release is correlated with the vacancy energy of Ca2+ and Na+ cations. Here, we found that the CNBS structure with the Na+ cation vacancy is energetically more favorable than the structure with the Ca2+ cation vacancy. The calcium adsorption site has been shown to have a greater affinity for water than can be found in the case of the sodium site, demonstrating that affinity may not be considered a major factor controlling the release of cations from the glass to the solution.

QEHeat: An open-source energy flux calculator for the computation of heat-transport coefficients from first principles

We give a detailed presentation of the theory and numerical implementation of an expression for the adiabatic energy flux in extended systems, derived from density-functional theory. This expression can be used to estimate the heat conductivity from equilibrium ab initio molecular dynamics, using the Green-Kubo linear response theory of transport coefficients. Our expression is implemented in an open-source component of the Quantum ESPRESSO suite of computer codes for quantum mechanical materials modeling, which is being made publicly available. Program summaryProgram title:QEHeatCPC Library link to program files:https://doi.org/10.17632/c6wxkvy4z3.1Licensing provisions: GPLv3Programming language: FortranNature of problem: The computation of thermal transport coefficients via equilibrium molecular dynamics and the Green-Kubo theory of linear response requires the definition of a heat-flux describing the instantaneous flow of energy. When considering predictive first-principles methods, a definition of the heat-flux compatible with density-functional theory is required [1]. The evaluation of such a heat flux requires an extension of state-of-the-art atomic simulation codes.Solution method: This work describes in detail the numerical implementation of the adiabatic energy current derived in Refs. [1,2] and makes it available to the users of the Quantum ESPRESSO suite of computer codes [3]. Used in conjunction with the cp.x code, to perform Car-Parrinello ab initio molecular dynamics, and the SPORTRAN post-processing tool for data analysis [4–6], the program allows to estimate heat transport coefficients in extended systems entirely from first principles. The new code provides as well to developers a modular and easily extendable framework to evaluate time derivatives of electronic properties (e.g. electronic densities or potentials) via a finite difference approach.

Oxygen depletion in FLASH ultra-high-dose-rate radiotherapy: A molecular dynamics simulation.

We present a first-principles molecular dynamics (MD) simulation and expound upon a mechanism of oxygen depletion hypothesis to explain the mitigation of normal tissue injury observed in ultra-high-dose-rate (UHDR) FLASH radiotherapy. We simulated damage to a segment of DNA (also representing other biomolecules such as RNA and proteins) in a simulation box filled with O and molecules. Attoseconds physical interactions (ionizations, electronic, and vibrational excitations) were simulated by using the Monte Carlo track structure code Geant4-DNA. Immediately after ionization, ab initio Car-Parrinello molecular dynamics (CPMD) simulation was used to identify which O and molecules surrounding the DNA molecule were converted into reactive oxygen species (ROS). Subsequently, the femto- to nanosecond reactions of ROS were simulated by using MD with reactive force field (ReaxFF), to illustrate ROS merging into new types of non-reactive oxygen species (NROS) due to strong coupling among ROS. A coarse-grained model was constructed to describe the relevant collective phenomenon at the macroscopic level on ROS aggregation and formation of NROS agglomerates consistent with the underlying microscopic pathways obtained from MD simulations. Time-dependent molecular simulations revealed the formation of metastable and transient spaghetti-like complexes among ROS generated at UHDR. At the higher ROS densities produced under UHDR, stranded chains (i.e., NROS) are produced, mediated through attractive electric polarity forces, hydrogen bonds, and magnetic dipole-dipole interactions among hydroxyl radicals. NROS tend to be less mobile than cellular biomolecules as opposed to the isolated and sparsely dense ROS generated at conventional dose rates (CDR). We attribute this effect to the suppression of biomolecular damage induced per particle track. At a given oxygen level, as the dose rate increases, the size and number of NROS chains increase, and correspondingly the population of toxic ROS components decreases. Similarly, at a given high dose rate, as the oxygen level increases, so do the size and number of NROS chains until an optimum level of oxygen is reached. Beyond that level, the amount of oxygen present may be sufficient to saturate the production of NROS chains, thereby reversing the sparing effects of UHDRs. We showed that oxygen depletion, hypothesized to lead to lower normal-tissue toxicity at FLASH dose rates, takes place within femto- to nanoseconds after irradiation. The mechanism is governed by the slow dynamics of chains of ROS complexes (NROS). Under physoxic (≈ 4-5% oxygen) conditions (i.e., in normal tissues), NROS are more abundant than in hypoxic conditions (e.g., <0.3% in parts of tumors), suggesting that biomolecular damage would be reduced in an environment with physoxic oxygen levels. Hence irradiation at UHDRs would be more effective for sparing physoxic normal tissues but not tumors containing regions of hypoxia. At much higher levels of oxygen (e.g., >10-15%), oxygen depletion by UHDRs may not be sufficient for tissue sparing.

Prebiotic Soup Components Trapped in Montmorillonite Nanoclay Form New Molecules: Car-Parrinello Ab Initio Simulations.

The catalytic effects of complex minerals or meteorites are often mentioned as important factors for the origins of life. To assess the possible role of nanoconfinement within a catalyst consisting of montmorillonite (MMT) and the impact of local electric field on the formation efficiency of the simple hypothetical precursors of nucleic acid bases or amino acids, we performed ab initio Car–Parrinello molecular dynamics simulations. We prepared four condensed-phase systems corresponding to previously suggested prototypes of a primordial soup. We monitored possible chemical reactions occurring within gas-like bulk and MMT-confined four simulation boxes on a 20-ps time scale at 1 atm and 300 K, 400 K, and 600 K. Elevated temperatures did not affect the reactivity of the elementary components of the gas-like boxes considerably; however, the presence of the MMT nanoclay substantially increased the formation probability of new molecules. Approximately 20 different new compounds were found in boxes containing carbon monoxide or formaldehyde molecules. This observation and an analysis of the atom–atom radial distribution functions indicated that the presence of Ca2+ ions at the surface of the internal MMT cavities may be an important factor in the initial steps of the formation of complex molecules at the early stages of the Earth’s history.

Hydrogen-Bonded Complexes of p-Hydrobenzoic Acid and Its Derivatives with a Polar Cosolvent in Supercritical Carbon Dioxide

The results of numerical modeling of p-hydroxybenzoic acid, methylparaben, and propylparaben solutions in methanol-modified supercritical carbon dioxide were presented. The composition, structure, and formation energy of the solute–cosolvent hydrogen-bonded complexes were determined. The duration of the existence of hydrogen bonds via individual paraben atoms was evaluated based on the data of Car–Parrinello ab initio molecular dynamics simulation, and the potential of mean force of interaction with methanol was calculated taking into account the medium. It was shown that, despite the presence of several potential centers of specific interaction, p-hydroxybenzoic acid and its derivatives form hydrogen bonds mainly via carboxyl and hydroxyl protons, respectively.

Impact of finite-temperature and condensed-phase effects on theoretical X-ray absorption spectra of transition metal complexes.

The impact of condensed-phase and finite-temperature effects on the theoretical X-ray absorption spectra of transition metal complexes is assessed. The former are included in terms of the all-electron Gaussian and augmented plane-wave approach, whereas the latter are taken into account by extensive ensemble averaging along second-generation Car-Parrinello ab initio molecular dynamics trajectories. We find that employing the periodic boundary conditions and including finite-temperature effects systematically improves the agreement between our simulated X-ray absorption spectra and experimental measurements. © 2018 Wiley Periodicals, Inc.

Explicit Aqueous Solvation Treatment of Epinephrine from Car-Parrinello Molecular Dynamics: Effect of Hydrogen Bonding on the Electronic Absorption Spectrum.

The electronic absorption spectrum of the neurotransmitter epinephrine (EPN) in water solution is studied, combining ab initio Car-Parrinello molecular dynamics (CPMD) with a quantum mechanical approach within the framework of the time-dependent density functional theory (TDDFT) scheme. By selecting 52 uncorrelated snapshots, the excitation modes were calculated at the LC-ωPBE/6-31+G(d) level of theory, using an optimal range-separation parameter ω, determined by means of the gap-tuning scheme in the presence of the solvent molecules. By comparing with static approaches (vacuum and implicit solvation), we show here that explicit solvation treatment dramatically enhances the photophysical properties of the EPN, especially because of the more realistic dynamic description of the molecular geometry. The agreement between the simulated and experimental spectra is demonstrated when TDDFT calculations are performed with the optimally tuned version of the DFT hybrid, not only improving the relative intensities of the absorption bands but also the λmax position. These results highlight that accounting for the nuclear motions, that is, thermal effects (of both chromophore and solvent molecules), is imperative to predict experimental absorption spectra. In this paper, we have addressed the critical importance of explicit solvation effects on the photophysics of the EPN, raking in performance when the simulation is performed based on first-principles molecular dynamics such as CPMD.

Solid–Vapor Reaction Growth of Transition‐Metal Dichalcogenide Monolayers

AbstractTwo‐dimensional (2D) layered semiconducting transition‐metal dichalcogenides (TMDCs) are promising candidates for next‐generation ultrathin, flexible, and transparent electronics. Chemical vapor deposition (CVD) is a promising method for their controllable, scalable synthesis but the growth mechanism is poorly understood. Herein, we present systematic studies to understand the CVD growth mechanism of monolayer MoSe2, showing reaction pathways for growth from solid and vapor precursors. Examination of metastable nanoparticles deposited on the substrate during growth shows intermediate growth stages and conversion of non‐stoichiometric nanoparticles into stoichiometric 2D MoSe2 monolayers. The growth steps involve the evaporation and reduction of MoO3 solid precursors to sub‐oxides and stepwise reactions with Se vapor to finally form MoSe2. The experimental results and proposed model were corroborated by ab initio Car–Parrinello molecular dynamics studies.

Solid-Vapor Reaction Growth of Transition-Metal Dichalcogenide Monolayers.

Two-dimensional (2D) layered semiconducting transition-metal dichalcogenides (TMDCs) are promising candidates for next-generation ultrathin, flexible, and transparent electronics. Chemical vapor deposition (CVD) is a promising method for their controllable, scalable synthesis but the growth mechanism is poorly understood. Herein, we present systematic studies to understand the CVD growth mechanism of monolayer MoSe2 , showing reaction pathways for growth from solid and vapor precursors. Examination of metastable nanoparticles deposited on the substrate during growth shows intermediate growth stages and conversion of non-stoichiometric nanoparticles into stoichiometric 2D MoSe2 monolayers. The growth steps involve the evaporation and reduction of MoO3 solid precursors to sub-oxides and stepwise reactions with Se vapor to finally form MoSe2 . The experimental results and proposed model were corroborated by ab initio Car-Parrinello molecular dynamics studies.

Efficient "on-the-Fly" Calculation of Raman Spectra from Ab-Initio Molecular Dynamics: Application to Lithium-Sulfur Batteries

We present a novel computational method to accurately calculate Raman spectra from first principles. Together with an extension of the second-generation Car-Parrinello ab-initio molecular dynamics (AIMD) method to propagate maximally localized Wannier functions along with the nuclei, a speed-up of one order of magnitude can be obtained. Due to expression of the total polarizability as a sum of Wannier polarizabilities, the method also enable us to find partial Raman activities corresponding to local vibrations of different fragments in the system. We have used this method to investigate the Raman spectra of (Li2S4)\U0001d45b, \U0001d45b=1,4,8 clusters which are believed to be the last intermediates in the S8→Li2S transition in lithium-sulfur (Li-S) batteries, during the discharge cycle. However, it is not fully established whether or not any other products coexist with the Li2S crystal in the discharged state. In this study, we have observed a clear evidence of Li2S4→Li2S2 transition by investigating the systematic changes in the simulated Raman spectra as the cluster size increases. In line with recent experiments, we found that the Raman-active sulfur-sulfur stretching mode at 440 cm-1 can be considered as a signature of covalent bonding between two and more sulfur atoms per formula. We have also demonstrated that the transition is mainly due to the strong electrostatic interactions between polar Li2S4 monomers which result in energy lowering by arranging the local dipole moments in anti-parallel fashion. This arrangement result in the formation of crystal-like Li2S2 structures at the interfaces of Li2S4 monomers. Furthermore, through minimum-energy conformation search combined with AIMD simulations, we have shown that upon carbonization of polyacrylonitrile in early discharge stages of Li-S batteries and further exposure of the product to octasulfur, S2 groups form a bridge between two consecutive nitrogen atom on the edge of carbon fiber. Moreover, our Raman simulations reveal a characteristic Raman activity for these -N-S-S-N- groups in 300 - 400 cm-1 frequency range, in form of two distinct peaks which researchers could not experimentally assign before. Finally, we will discuss our recent studies on structure of complex cathode materials using our approach together with experimental measurements.

Alloy counter electrodes in dye-sensitized solar cells

Dye-sensitized solar cells (DSSCs) have gained widespread attention of researchers in recent years due to their unique advantages, such as low cost of production, ease of manufacturing and relatively high power conversion efficiency. Traditionally, a DSSC has several components including photoanode material, dyes, electrolytes and the counter electrode (CE), which must cooperate intimately to transform incident light into current. As a crucial part of DSSCs, designing CE with excellent properties is a key to improve the power conversion efficiency of DSSCs. As the good CE material, noble metal Pt material exhibiting high catalytic activity and electrical conductivity has obtained high power conversion efficiency in DSSCs. However, some problems of Pt, especially high scarcity in the nature and the corrosivity of I–/I3– electrolyte, limit its further application in DSSCs. Meanwhile, lots of Pt-free materials involving carbon materials, transition metal compounds, polymers, and hybrids have been developed as CE materials to solve the challenges of the Pt-based DSSCs, and achieved favorable power conversion efficiency. Recently, another novel type of alloy materials, which include Pt-based and Pt-free alloys, show superior catalytic activity and have achieved significant power conversion efficiency when applied as CEs in DSSCs. The alloy materials have the merits of simple preparation process and low cost and are the outstanding candidate to replace Pt in DSSCs. The present review gives a summary of the alloy materials as CE in DSSCs. We stress the advantages and disadvantages of alloy CE materials, and put forward ideas for preparation methods, DSSCs’ designing, ab-initio Car-Parrinello molecular dynamics simulations, first-principle density functional theory calculations and stability evaluations in future research for DSSCs featuring alloy CEs. Also, some crucial issues in alloy-based DSSCs and the further applications of the alloy materials in other field of catalysis are discussed.

Ab Initio MD Simulations of the Brønsted Acidity of Glutathione in Aqueous Solutions: Predicting pKa Shifts of the Cysteine Residue.

The tripeptide glutathione (GSH) is one of the most abundant peptides and the major repository for nonprotein sulfur in both animal and plant cells. It plays a critical role in intracellular oxidative stress management by the reversible formation of glutathione disulfide with the thiol-disulfide pair acting as a redox buffer. The state of charge of the ionizable groups of GSH can influence the redox couple, and hence the pKa value of the cysteine residue of GSH is critical to its functioning. Here we report ab initio Car-Parrinello molecular dynamics simulations of glutathione solvated by 200 water molecules, all of which are considered in the simulation. We show that the free-energy landscape for the protonation-deprotonation reaction of the cysteine residue of GSH computed using metadynamics sampling provides accurate estimates of the pKa and correctly predicts the shift in the dissociation constant values as compared with the isolated cysteine amino acid.

Ab Initio Molecular Dynamics Simulations of Amino Acids in Aqueous Solutions: Estimating pKa Values from Metadynamics Sampling.

Changes in the protonation and deprotonation of amino acid residues in proteins play a key role in many biological processes and pathways. Here, we report calculations of the free-energy profile for the protonation-deprotonation reaction of the 20 canonical α amino acids in aqueous solutions using ab initio Car-Parrinello molecular dynamics simulations coupled with metadynamics sampling. We show here that the calculated change in free energy of the dissociation reaction provides estimates of the multiple pKa values of the amino acids that are in good agreement with experiment. We use the bond-length-dependent number of the protons coordinated to the hydroxyl oxygen of the carboxylic and the amine groups as the collective variables to explore the free-energy profiles of the Bronsted acid-base chemistry of amino acids in aqueous solutions. We ensure that the amino acid undergoing dissociation is solvated by at least three hydrations shells with all water molecules included in the simulations. The method works equally well for amino acids with neutral, acidic and basic side chains and provides estimates of the multiple pKa values with a mean relative error, with respect to experimental results, of 0.2 pKa units.

Structure and Dynamics of the Instantaneous Water/Vapor Interface Revisited by Path-Integral and Ab Initio Molecular Dynamics Simulations.

The structure and dynamics of the water/vapor interface is revisited by means of path-integral and second-generation Car-Parrinello ab initio molecular dynamics simulations in conjunction with an instantaneous surface definition [Willard, A. P.; Chandler, D. J. Phys. Chem. B 2010, 114, 1954]. In agreement with previous studies, we find that one of the OH bonds of the water molecules in the topmost layer is pointing out of the water into the vapor phase, while the orientation of the underlying layer is reversed. Therebetween, an additional water layer is detected, where the molecules are aligned parallel to the instantaneous water surface.

Accurate and inexpensive prediction of the color optical properties of anthocyanins in solution.

The simulation of the color optical properties of molecular dyes in liquid solution requires the calculation of time evolution of the solute absorption spectra fluctuating in the solvent at finite temperature. Time-averaged spectra can be directly evaluated by combining ab initio Car-Parrinello molecular dynamics and time-dependent density functional theory calculations. The inclusion of hybrid exchange-correlation functionals, necessary for the prediction of the correct transition frequencies, prevents one from using these techniques for the simulation of the optical properties of large realistic systems. Here we present an alternative approach for the prediction of the color of natural dyes in solution with a low computational cost. We applied this approach to representative anthocyanin dyes: the excellent agreement between the simulated and the experimental colors makes this method a straightforward and inexpensive tool for the high-throughput prediction of colors of molecules in liquid solvents.