Atom-centered electric multipole moments can be extremely useful in chemistry, as they enable the systematic mapping of a complex electrostatic problem to a simpler model. However, since they do not correspond to physical observables, there is no unique way to define them. In this study, we present an extension of the dynamically generated RESP charges (D-RESP) method, referred to as xDRESP, where atom-centered multipoles are computed from mixed quantum mechanics/molecular mechanics molecular dynamics simulations. We compare the ability of xDRESP charges to reproduce the electrostatic potential, as well as molecular multipoles, against the performance of fixed point-charge models commonly used in force fields. Moreover, we highlight cases where xDRESP atomic multipoles can provide valuable information about chemical systems, such as indicating when polarization plays a significant role, and chemical reactions in which xDRESP atomic multipoles can be used as an on-the-fly analysis tool to track changes in electron density.

This study presents the polarizable embedding cluster perturbation framework (PE-CP) and derives the working equations needed to perform calculations in a singles and perturbative doubles excitation space, denoted PE-CPS(D). The method is implemented in the Penguin program through interfaces to the external libraries PyFraMe and CPPE. To obtain the cluster amplitudes and multipliers, which constitute two coupled sets of equations in PE-CCSD, a standard scheme of outer and inner iterations is usually applied. In the PE-CPS(D) model, only the singles excitation space requires this iterative procedure to determine the wave function parameters; the perturbative doubles contributions can be evaluated directly. This feature highlights one of the principal strengths of the PE-CPS(D) approach. Finally, we present numerical tests comparing PE-CPS(D) with PE-CCSD, to which it formally converges. These results demonstrate that PE-CPS(D) provides an efficient route to high-level coupled-cluster accuracy for large quantum-mechanical systems embedded in a classical polarizable environment.

Channelrhodopsin-2 (ChR2) is a light-gated ion channel widely used in optogenetics, a technique that enables precise control of neuronal activity by genetically engineering light-sensitive proteins into cell membranes. This protein exists in dimeric form, with each monomer containing a retinal Schiff base (RSB) moiety covalently bonded that undergoes trans-cis isomerization upon light absorption. However, the limited penetration depth of visible light in biological tissues motivates the use of multiphoton-absorption techniques, which enhance tissue penetration, improve focality, and reduce phototoxicity, thereby offering a promising alternative for optogenetic applications. In this paper, we present a fully atomistic multiscale methodology for computing the one-, two-, and three-photon absorption spectra of ChR2, where the protein, lipid bilayer, and solvent are explicitly considered throughout the workflow. This methodology integrates classical molecular mechanics (MM) molecular dynamics (MD), quantum mechanics/molecular mechanics (QM/MM)-MD, and fragment-based polarizable embedding (PE) to derive environment-specific PE potentials from the explicit protein-lipid-solvent environment. The final step in the methodology is to use these potentials to compute accurate spectra via PE-time-dependent density functional theory (PE-TD-DFT). Validation against experimental one-photon absorption spectra demonstrates excellent agreement. For the first time, we report the theoretical two- and three-photon absorption in ChR2, albeit without direct experimental comparison. We compare the multiphoton absorption (MPA) spectra where the two RSB moieties are sampled using classical MD and QM/MM-MD, respectively. The resulting spectral differences are attributed to variations in key structural parameters that we analyze and document.

We present a workflow, benchmarks, and applications to provide a roadmap for simulating harmonic IR and Raman spectra for solute-solvent systems by employing a polarizable-embedding quantum-mechanics (PE-QM) approach. This multiscale modeling scheme divides the system into a central core region described by quantum-mechanical methods and an environment region described through the fragment-based polarizable embedding (PE) model. The workflow involves generating representative structures, calculating properties, and postprocessing data. Benchmark calculations quantify errors introduced by some of the key approximations used in our approach and discuss its strengths and weaknesses. Finally, we apply the workflow to acetone in three different solvents, comparing simulated spectra to experimental results to further evaluate our approach and identify potential weaknesses. Accurate simulations of solute-solvent systems are an important step toward modeling more complex molecular systems with a fragment-based PE approach.

Multiscale simulations are essential techniques in computational chemistry, providing insights into complex phenomena across extended temporal and spatial scales. With a particular interest in the dynamics of such processes, we developed MiMiC, a framework for efficient multiscale molecular dynamics simulations suited for high-performance computing. One of its key characteristics is a flexible design where external specialized programs handle individual subsystems. This article reviews the core features and some recent advancements in MiMiC, particularly the integration of OpenMM and CP2K.

MiMiC is a framework for modeling large-scale chemical processes that require treatment at multiple resolutions. It does not aim to implement single-handedly all methods required to treat individual subsystems, but instead, it relegates this task to specialized computational chemistry software while it serves as an intermediary between these external programs and computes the interactions between the subsystems. MiMiC minimizes issues typically associated with molecular dynamics performed with multiple programs by adopting a multiple-program multiple-data paradigm combined with a loose-coupling model. In this work, we present the addition of a new client program, CP2K, to the MiMiC ecosystem. Moreover, to align the implementation of MiMiC with its modular philosophy, we performed a major refactoring of the entire framework. This endeavor unlocks its full flexibility and reduces any future efforts for introducing new programs to a minimum. Furthermore, by thorough timing analysis, we verify that the introduced changes do not affect the performance of MiMiC or CP2K, and neither are they a source of significant computational overheads that would be detrimental to simulation efficiency. Finally, we demonstrate the benefits of the framework's modular design, by performing a QM/MM MD simulation combining CP2K with previously interfaced OpenMM, with no additional implementation effort required.

MiMiC is a flexible and efficient framework for multiscale simulations in which different subsystems are treated by individual client programs. In this work, we present a new interface with OpenMM to be used as an MM client program and we demonstrate its efficiency for QM/MM MD simulations. Apart from its high performance, especially on GPUs, and a wide selection of features, OpenMM is a highly flexible and easily extensible program, ideal for the development of novel multiscale methods. Thanks to the open-ended design of MiMiC, the OpenMM-MiMiC interface will automatically support any new QM client program interfaced with MiMiC for QM/MM and, with minimal changes needed, new multiscale methods implemented, opening up new research directions beyond electrostatic embedding QM/MM.

G protein-coupled receptors are key drug targets due to their role in cellular signaling. Among them, bistable Rhodopsins such as the Jumping Spider Rhodopsin 1 (JSR1), are promising for optogenetic applications, but their transduction mechanisms remain poorly understood. In this study, we used microsecond equilibrium molecular dynamics simulations, network analysis, and machine learning to investigate allosteric communication paths between the retinal chromophore and the intracellular G protein-binding site in JSR1. We analyzed structural differences in three functional states with retinal chromophores in 9-cis, 11-cis, and all-trans configurations. Results revealed that Trp290 is crucial for transmitting the movements of the retinal after isomerization to the G protein-binding site during JSR1 activation as well as residues along TM6 helix. Overall, these findings advance our understanding of bistable Rhodopsins and their potential in light-driven technologies.

We present a workflow, benchmarks, and applications to provide a roadmap for simulating harmonic IR and Raman spectra for large solute-solvent systems by employing a polarizable-embedding quantum-mechanics (PE-QM) approach. This multiscale modeling scheme divides the system into a central core region described by quantum-mechanical methods and an environment region described through the fragment-based polarizable embedding (PE) model. The workflow involves generating representative structures, calculating properties, and post-processing data. Benchmark calculations quantify errors introduced by some of the key approximations used in our approach and discuss its strengths and weaknesses. Finally, we apply the workflow to acetone in three different solvents, comparing simulated spectra to experimental results to further evaluate our approach and identify potential weaknesses. Accurate simulations of solute-solvent systems are an important step toward modeling more complex molecular systems with a fragment-based PE approach.

The complexity of biological systems and processes, spanning molecular to macroscopic scales, necessitates the use of multiscale simulations to get a comprehensive understanding. Quantum mechanics/molecular mechanics (QM/MM) molecular dynamics (MD) simulations are crucial for capturing processes beyond the reach of classical MD simulations. The advent of exascale computing offers unprecedented opportunities for scientific exploration, not least within life sciences, where simulations are essential to unravel intricate molecular mechanisms underlying biological processes. However, leveraging the immense computational power of exascale computing requires innovative algorithms and software designs. In this context, we discuss the current status and future prospects of multiscale biomolecular simulations on exascale supercomputers with a focus on QM/MM MD. We highlight our own efforts in developing a versatile and high-performance multiscale simulation framework with the aim of efficient utilization of state-of-the-art supercomputers. We showcase its application in uncovering complex biological mechanisms and its potential for leveraging exascale computing.

In this work, we present the development of a fully-polarizable KS-DFT/AMOEBA embedding scheme for delocalized basis sets such as plane-waves and real-space grids. The augmented problem of electron spill-out inherent to a polarizable QM/MM implementation with plane-wave basis sets is addressed and the periodicity for the MM subsystem is taken into account, as implemented in the Tinker-HP software. We discuss the software design and how computational efficiency is enabled through the interoperable multiscale simulation framework MiMiC. The implementation is validated on QM/MM energies for dimer systems and a quantitative assessment of molecular dipoles of embedded solutes in order to estimate the magnitude of errors related to the damping parameters used in the model.

MiMiC is a framework for performing multiscale simulations in which loosely coupled external programs describe individual subsystems at different resolutions and levels of theory. To make it highly efficient and flexible, we adopt an interoperable approach based on a multiple-program multiple-data (MPMD) paradigm, serving as an intermediary responsible for fast data exchange and interactions between the subsystems. The main goal of MiMiC is to avoid interfering with the underlying parallelization of the external programs, including the operability on hybrid architectures (e.g., CPU/GPU), and keep their setup and execution as close as possible to the original. At the moment, MiMiC offers an efficient implementation of electrostatic embedding quantum mechanics/molecular mechanics (QM/MM) that has demonstrated unprecedented parallel scaling in simulations of large biomolecules using CPMD and GROMACS as QM and MM engines, respectively. However, as it is designed for high flexibility with general multiscale models in mind, it can be straightforwardly extended beyond QM/MM. In this article, we illustrate the software design and the features of the framework, which make it a compelling choice for multiscale simulations in the upcoming era of exascale high-performance computing.

The partial Hessian approximation is often used in vibrational analysis of quantum mechanics/molecular mechanics (QM/MM) systems because calculating the full Hessian matrix is computationally impractical. This approach aligns with the core concept of QM/MM, which focuses on the QM subsystem. Thus, using the partial Hessian approximation implies that the main interest is in the local vibrational modes of the QM subsystem. Here, we investigate the accuracy and applicability of the partial Hessian vibrational analysis (PHVA) approach as it is typically used within QM/MM, i.e., only the Hessian belonging to the QM subsystem is computed. We focus on solute-solvent systems with small, rigid solutes. To separate two of the major sources of errors, we perform two separate analyses. First, we study the effects of the partial Hessian approximation on local normal modes, harmonic frequencies, and harmonic IR and Raman intensities by comparing them to those obtained using full Hessians, where both partial and full Hessians are calculated at the QM level. Then, we quantify the errors introduced by QM/MM used with the PHVA by comparing normal modes, frequencies, and intensities obtained using partial Hessians calculated using a QM/MM-type embedding approach to those obtained using partial Hessians calculated at the QM level. Another aspect of the PHVA is the appearance of normal modes resembling the translation and rotation of the QM subsystem. These pseudotranslational and pseudorotational modes should be removed as they are collective vibrations of the atoms in the QM subsystem relative to a frozen MM subsystem and, thus, not well-described. We show that projecting out translation and rotation, usually done for systems in isolation, can adversely affect other normal modes. Instead, the pseudotranslational and pseudorotational modes can be identified and removed.

We introduce periodic boundary conditions (PBCs) for the induced electrostatics in the polarizable embedding (PE) model through the minimum-image convention (MIC). It is a simple yet effective approach that includes a more physically accurate description of the polarization throughout the molecular system. Using PE with MIC (PE-MIC), we shed new light on the limitations of commonly employed cutoff models, such as the droplet model, when used in PE calculations. Specifically, we investigate the effects of the unphysical polarization at the outer boundary by comparing induced dipoles and the associated electrostatic potentials, as well as some optical properties of solute–solvent and biomolecular systems. We show that the magnitude of the inaccuracies caused by the unphysical polarization depends on multiple parameters: the nature of the quantum subsystem and of the environment, the cutoff model and distance, and the calculated property.

Polarizable embedding (PE) refers to classical embedding approaches, such as those used in quantum mechanics/molecular mechanics (QM/MM), that allow mutual polarization between the quantum and classical regions. The quality of the embedding potential is critical to provide accurate results, e.g., for spectroscopic properties and dynamical processes. High-quality embedding-potential parameters can be obtained by dividing the classical region into smaller fragments and deriving the parameters from ab initio calculations on the fragments. For solvents and other systems composed of small molecules, the fragments can be individual molecules, while a more complicated fragmentation procedure is needed for larger molecules, such as proteins and nucleic acids. One such fragmentation strategy is the molecular fractionation with conjugate caps (MFCC) approach. As is widely known, hydrogen bonds play a key role in many biomolecular systems, e.g., in proteins, where they are responsible for the secondary structure. In this work, we assess the effects of including hydrogen-bond fragmentation in the MFCC procedure [MFCC(HB)] for deriving the embedding-potential parameters. The MFCC(HB) extension is evaluated on several molecular systems, ranging from small model systems to proteins, directly in terms of molecular electrostatic potentials and embedding potentials and indirectly in terms of selected properties of chromophores embedded in water and complex protein environments.

Multiscale simulations have been established as a powerful tool to calculate and predict excitation energies in complex systems such as photoreceptor proteins. In these simulations the chromophore is typically treated using quantum mechanical (QM) methods while the protein and surrounding environment are described by a classical molecular mechanics (MM) force field. The electrostatic interactions between these regions are often treated using electrostatic embedding where the point charges in the MM region polarize the QM region. A more sophisticated treatment accounts also for the polarization of the MM region. In this work, the effect of such a polarizable embedding on excitation energies was benchmarked and compared to electrostatic embedding. This was done for two different proteins, the lipid membrane-embedded jumping spider rhodopsin and the soluble cyanobacteriochrome Slr1393g3. It was found that the polarizable embedding scheme produces absorption maxima closer to experimental values. The polarizable embedding scheme was also benchmarked against expanded QM regions and found to be in qualitative agreement. Treating individual residues as polarizable recovered between 50% and 71% of the QM improvement in the excitation energies, depending on the system. A detailed analysis of each amino acid residue in the chromophore binding pocket revealed that aromatic residues result in the largest change in excitation energy compared to the electrostatic embedding. Furthermore, the computational efficiency of polarizable embedding allowed it to go beyond the binding pocket and describe a larger portion of the environment, further improving the results.

MiMiC is a highly flexible, extremely scalable multiscale modeling framework. It couples the CPMD (quantum mechanics, QM) and GROMACS (molecular mechanics, MM) codes. The code requires preparing separate input files for the two programs with a selection of the QM region. This can be a tedious procedure prone to human error, especially when dealing with large QM regions. Here, we present MiMiCPy, a user-friendly tool that automatizes the preparation of MiMiC input files. It is written in Python 3 with an object-oriented approach. The main subcommand PrepQM can be used to generate MiMiC inputs directly from the command line or through a PyMOL/VMD plugin for visually selecting the QM region. Many other subcommands are also provided for debugging and fixing MiMiC input files. MiMiCPy is designed with a modular structure that allows seamless extensions to new program formats depending on the requirements of MiMiC.

We present a fully self-consistent polarizable embedding(PE) modelthat does not suffer from unphysical boundary polarization. This isachieved through the use of the minimum-image convention (MIC) inthe induced electrostatics. It is a simple yet effective approachthat includes a more physically accurate description of the polarizationthroughout the molecular system. Using PE with MIC (PE-MIC), we shednew light on the limitations of commonly employed cutoff models, suchas the droplet model, when used in PE calculations. Specifically,we investigate the effects of the unphysical polarization at the outerboundary by comparing induced dipoles and the associated electrostaticpotentials, as well as some optical properties of solute-solventand biomolecular systems. We show that the magnitude of the inaccuraciescaused by the unphysical polarization depends on multiple parameters:the nature of the quantum subsystem and of the environment, the cutoffmodel and distance, and the calculated property.

The introduction of halogen atoms in the donor molecules in organic solar cells leads to a decrease in the reorganization energy, which in turn results in reduced non-radiative voltage losses and an improved open-circuit voltage in the devices.