Real-time feedback from iterative electronic structure calculations
RReal-time feedback from iterative electronic structure calculations
Real-time feedback from iterative electronic structure calculations
Alain C. Vaucher, Moritz P. Haag, and Markus Reiher a) ETH Z¨urich, Laboratorium f¨ur Physikalische Chemie, Vladimir-Prelog-Weg 2, CH-8093 Z¨urich,Switzerland (Dated: 4 October 2018)
Real-time feedback from iterative electronic structure calculations requires to mediate between the inherentlyunpredictable execution times of the iterative algorithm employed and the necessity to provide data in fixedand short time intervals for real-time rendering. We introduce the concept of a mediator as a componentable to deal with infrequent and unpredictable reference data to generate reliable feedback. In the context ofreal-time quantum chemistry, the mediator takes the form of a surrogate potential that has the same localshape as the first-principles potential and can be evaluated efficiently to deliver atomic forces as real-timefeedback. The surrogate potential is updated continuously by electronic structure calculations and guaranteesto provide a reliable response to the operator for any molecular structure. To demonstrate the application ofiterative electronic structure methods in real-time reactivity exploration, we implement self-consistent semi-empirical methods as the data source and apply the surrogate-potential mediator to deliver reliable real-timefeedback.Keywords: interactive quantum chemistry, real-time quantum chemistry, real-time feedback, electronic struc-ture calculations
I. INTRODUCTION
Knowledge of the atomic rearrangements occurring dur-ing chemical reactions is crucial for understanding reac-tivity. It allows designing new reactions and optimizingexisting ones. In recent years, approaches toward inter-active reactivity studies have been introduced.We developed Haptic Quantum Chemistry and Real-Time Quantum Chemistry to actively interact withsome molecular system under consideration and perceiveits immediate (quantum-mechanical) response. This im-mersion into the molecular world is enabled by immedi-ate feedback of quantum-chemical calculations, which isof paramount importance to guide operators during the(interactive) exploration of molecular systems. Interac-tivity enables chemists to understand reactivity more in-tuitively and more efficiently than with traditional tools.In other approaches toward the interaction with molec-ular systems, the emphasis is on the possibility to steersimulations in real time. For example, in the field ofinteractive molecular dynamics, Stone et al. allowed op-erators to drive classical molecular dynamics simulationstoward events that would happen too rarely or not atall otherwise . This field was recently extended to abinitio molecular dynamics by Luehr et al. Another ap-proach was implemented by Bosson et al., who applieda non-iterative semi-empirical quantum-chemical methodthat instantaneously optimizes a molecular structure setup with their molecular editor Samson . Such struc-ture relaxation approaches to quickly generate reason-able molecular structures (often based on classical forcefields) are becoming a standard tool for structure gener- a) Electronic mail: [email protected] ation attempts in such graphical user interfaces (cf. foranother example the
Avogadro program ).Real-time quantum chemistry allows exploring theBorn–Oppenheimer potential energy surface interactivelyin real time to gain intuitive insight into reactivity. Suchreal-time reactivity explorations rely on two types offeedback. Firstly, operators can actively manipulate themolecular structure, for example, with a haptic devicethat allows them to move atoms while experiencing thequantum-chemical force acting on them. The haptic feed-back provides an immediate and intuitive understandingof which parts of the potential energy surface are acces-sible (at a given temperature or energy) and indicateshow a molecular system is prepared to react. Secondly,the operator can directly observe the effect of structuremanipulations on the molecular system as a whole. Thisvisual feedback relies on a continuously running real-timegeometry optimization that drives the molecular systemto nearby local minima. Both types of feedback rely onquantum-mechanical forces calculated in real time by fastelectronic structure methods.For an optimal immersion into reactivity explorations,such real-time feedback needs to be reliable in differentaspects. From the operator’s perspective, feedback mustalways be provided independently of the current state ofthe system, even if its exact calculation is not possibleat all. With respect to accuracy, feedback must give aqualitatively correct description of molecular behavior or,if such a description is unavailable, it must not lead theoperator to inaccessible areas of configuration space.Recently, we demonstrated the application of real-timereactivity studies described by the non-self-consistentdensity-functional tight-binding (DFTB) method. DFTB is a non-iterative method featuring constantexecution times and therefore delivers forces at a con-stant frequency. a r X i v : . [ phy s i c s . c o m p - ph ] O c t eal-time feedback from iterative electronic structure calculations 2Non-iterative methods such as DFTB are, however, tooapproximate to describe complex molecular systems reli-ably. Hence, for most systems the application of self-consistent field (SCF) methods is required. However,the application of SCF methods introduces an additionallayer of complexity due to their iterative nature. With it-erative methods, it cannot be guaranteed that the forceswill be delivered at a constant frequency or that they willbe delivered at all (if convergence cannot be achieved).If, for a given structure, the execution time of the elec-tronic structure optimization is incompatible with thereal-time requirement, it will be necessary to freeze themolecular structure until the optimization converges asdone by Luehr et al. in the field of steered ab initio molec-ular dynamics. However, in the context of real-time re-activity explorations, immersion and hence interactivitywould be jeopardized by freezing haptic and visual feed-back.To maintain interactivity in real-time reactivity explo-rations based on iterative methods, the feedback can-not (and need not) reflect the exact quantum-mechanicalforces directly because they might be provided only in-termittently. In this work, we present a strategy thataddresses this issue and allows performing interactivequantum-chemical reactivity studies with SCF methods.
II. TIMESCALES FOR DATA GENERATION ANDREAL-TIME FEEDBACK
In interactive chemical reactivity studies, informationabout the behavior of a molecular system duringoperator-induced structure manipulations is conveyed tothe operator by real-time rendering of data generatedby electronic structure calculations. Rendering in thiscontext entails the complete process of transforming rawdata into data perceptible by human senses.To maintain interactivity, the rendering process must beperformed in real time. Depending on the sense ad-dressed by the rendering, different real-time requirementsmust be fulfilled. Addressing the visual and haptic sen-sory systems impose different requirements for the up-date rate that is to be maintained. Achieving certainminimum update rates translates into maximum timeintervals that are allowed for the data generation andrendering process.The time interval available for the calculation of a re-sult is therefore defined by the senses addressed by thefeedback mechanism. For a setting where feedback con-sists of rendering molecular forces and visualizing struc-tural changes upon manipulations by the operator, Ta-ble I gives an overview of the imposed time limits andthe types of calculation to be performed in the availabletime. visual hapticFrequency 60 Hz 1 kHzTime interval 16 ms 1 msCalculation structure update force calculationTABLE I. Senses addressed in the haptic quantum chemistrysetting and their implications for quantum-chemical calcula-tions.
Basically all reliable electronic structure methods requirean iterative optimization of the electronic wave functionor density (e.g., an SCF optimization of the orbitals).The iterative nature of these methods, however, leadsto unpredictable execution times. Assuming that all pa-rameters of a calculation are fixed and only the molec-ular structure varies from calculation to calculation, asis the case in explorations of chemical reactivity, it can-not be known a priori whether and when a calculationwill converge. This also holds for fast methods that usu-ally converge in few iterations but might not converge atall for far-from-equilibrium structures generated by theoperator.The finite execution time of electronic structure calcu-lations introduces, at a given instant of time, shifts inconsistency between the molecular structure displayed tothe operator, the structure used for the data generation,and the structure to which the feedback corresponds. Inpractice, the electronic structure can be optimized onlyfor a fraction of the molecular structures visited as soonas the electronic structure calculation cannot be per-formed at 60 Hz. The feedback is always given witha certain delay because of the finite execution time ofthe electronic structure optimization and because of thetime elapsed between the completion of the last electronicstructure optimization and the instant of time at whichthe feedback is needed. This is represented schematicallyin Fig. 1.
Structure sequence ttt
CalculationFeedback
FIG. 1. Changing molecular structure (depicted by coloredshapes) of a real-time exploration. Top: Real-time evolutionof the molecular structure as displayed to the operator. Themolecular structure is updated with 60 Hz. Middle: Elec-tronic structure calculations are performed only for some ofthe visited structures. Bottom: The feedback relies on molec-ular structures for which the calculations finished. The verti-cal dashed lines represent instants of time when an electronicstructure calculation finished and another one starts.
The delay will not be perceptible if the electronic struc-ture calculations are very fast, but it can lead to arti-eal-time feedback from iterative electronic structure calculations 3facts in the operator’s perception as soon as their exe-cution time increases. In Fig. 1, a new electronic struc-ture calculation starts only when the previous one is fin-ished. Starting parallel electronic structure calculationsfor consecutive molecular structures can reduce the delayby updating the data on which the feedback relies morefrequently, but the delay will never vanish.The inherently unpredictable execution times and thenon-vanishing delays even in the most well-behaved sit-uations prohibit the direct application of iterative elec-tronic structure calculations as a data source for real-timefeedback as required in interactive applications.
III. MEDIATOR STRATEGY
In a setting that allows an operator to interact with vir-tual objects, the response of the system to the operator’sactions is presented by a feedback component. This feed-back component is usually fed with data by some datageneration component—usually some sort of calculation.Providing real-time feedback imposes strict deadlines forthe provision of the data. Hence, the data flow from gen-erator component to feedback component can only be di-rect if the data generation is performed in a predictableand constant time interval. If, however, the update ratesof data generation and feedback component differ, filter-ing techniques must be utilized to render the incomingdata stream. data generation userfeedback mediator
FIG. 2. The three main components necessary for interac-tive systems with real-time feedback and unpredictable dataavailability.
If the data generator cannot guarantee the provision ofnew data in fixed time intervals, an additional compo-nent is necessary to mediate between the unpredictabledata output stream and the strict real-time requirementsimposed by the operator feedback (see Fig. 2). The me-diator consumes the data stream coming from the gen-erating calculation and, based on it, provides real-timedata for the operator feedback. The mediator must ful-fill three important requirements: • It must guarantee to deliver data for the feedbackin fixed and constant time intervals (whose value isset by the human sense addressed). • It must provide reliable feedback in any situation. • The data output must be based on the most recentdata provided by the data generator.The mediator in turn must rely on the ability of thedata generator to provide a time ordering of the datastream. This means that the raw data must be assignedto a certain stage in the exploration so that the medi-ator is able to provide an adequate response. The timeordering will automatically be achieved if the data gener-ator performs the calculations consecutively. In the caseof parallel calculations, the data generator must provideadditional information on what results belong to whichinstant of time.A trivial implementation of the mediator would be toprovide the last data received to the operator feedbackcomponent until it receives new data from the generator.However, as the most recent data will no longer corre-spond to the situation the operator is experiencing (mis-match of property data for a user-modified new struc-ture), interactivity will be severely hampered or will belost entirely. Therefore, this trivial implementation doesnot fulfill the second requirement for a mediator.
A. The Surrogate Function
When the feedback is a function of the current configura-tion of the system, a mediator will be required if the timeneeded for the evaluation of the function is large or un-predictable. In such cases, a mediator fulfilling the aboverequirements can be implemented in terms of a surrogatefunction that can be evaluated rapidly. We will refer tothe unknown but correct function as ’reference’ or ’exact’function.The reference function is still sampled as often as possibleand delivers data allowing for updating the shape of thesurrogate function. The surrogate function can be eval-uated at the high frequencies required by the feedbackcomponent.If the reference function is a potential or energy and thereal-time feedback relies on forces, we will call the sur-rogate function a ’surrogate potential’. The surrogatepotential is chosen to approximate the local shape of thereference potential. The force feedback is evaluated fromthe surrogate potential, as illustrated in Fig. 3.For the surrogate potential as the central component ofthe mediator, it is not important how new reference cal-culations are started as long as the results are strictlyordered in time. The mediator only needs to know whichdata is most recent. The calculations can therefore beperformed consecutively (a new calculation starts onlywhen the previous one finished) or in parallel (a new cal-culation can start before the previous one finished).The requirement on the mediator to always yield reliabledata that prevents the operator to explore unreasonableeal-time feedback from iterative electronic structure calculations 4 raw data generation(unpredictable) mediator force-feedback (strict every millisecond)
FIG. 3. The reference potential (red) is sampled (blue arrows)to generate the surrogate potentials (black), which are in turnsampled to provide the real-time feedback data. Unused partsof the surrogate potentials due to the arrival of new referencedata are indicated in gray. areas of configuration space translates into the require-ment for the surrogate potentials to be always boundedfrom below. Such surrogate potentials ensure conserva-tive forces that will drive the operator toward a guess fora local minimum of the reference potential, which usu-ally also corresponds to the part of the reference potentialwell reproduced by the surrogate potential.It must be emphasized that the surrogate potential isnot an interpolation scheme as it uses only local infor-mation from the most recent calculation to predict thereference potential until new data arrives. An interpo-lation of the accumulated exploration history would inmost cases not improve the prediction as most of thetime the operator’s manipulations will explore new areasof configuration space that cannot be foreseen.
B. Obtaining the Surrogate Potential for Feedback fromElectronic Structure Calculations
For real-time quantum chemistry, the reference potentialcan, for instance, be approximated by surrogate poten-tials of quadratic form, which allow for an efficient cal-culation of the forces needed for the structure relaxationand for the haptic force rendering.Quadratic potentials are characterized by the general ex-pression V sur ( x ) = V + a T ( x − x )+ 12 ( x − x ) T B ( x − x ) (1)where x = { x , . . . , x m } is a coordinate vector for theatomic (i.e., nuclear) positions, V a constant potentialshift, a an m -dimensional vector and B an m × m matrix. x is a reference position and, in the context of real-time quantum chemistry, will correspond to the configurationat which the reference potential was sampled last. Fig. 4illustrates the quadratic approximation of a given poten-tial in two-dimensional configuration space ( m = 2). FIG. 4. Left: Arbitrary potential energy surface in twodimensions and the surrogate potential approximating it.Right: Close-up representation.
Quadratic potentials allow for a straightforward evalua-tion of the forces, F sur ( x ) = − ∇ V sur ( x ) = − a − B ( x − x ) , (2)provided that B can be efficiently calculated. F sur is avector composed of atomic forces and represents the basisfor real-time force feedback.To reproduce the reference potential around the atomicpositions as accurately as possible, the surrogate poten-tial is chosen to reproduce the local curvature of thepotential around the atomic nuclei. Given a molecularstructure for which the electronic wave function was op-timized, the coefficients of the surrogate potential around x are derived from the electronic energy E and itsderivatives as: V = E ( x ) , (3) a i = ∂E ( x ) ∂x i (cid:12)(cid:12)(cid:12)(cid:12) x = x , (4) B ij = ∂ E ( x ) ∂x i ∂x j (cid:12)(cid:12)(cid:12)(cid:12) x = x . (5)Applying Eqs. (3)–(5) will occasionally generatequadratic potentials with maxima or saddle points, whichwill not be bounded from below. Such potentials arecharacterized by B not being positive definite, i.e., B then possesses one or several negative eigenvalues. Tochange them into potentials with a global minimum, i.e.,into bounded potentials, can be achieved by replacing allnegative eigenvalues of B by their absolute values:1. Diagonalize B to yield B diag = U − BU = { B diag ,ij } , B diag ,ij = b ij = (cid:40) i (cid:54) = jb i i = j eal-time feedback from iterative electronic structure calculations 52. In B diag , replace all negative eigenvalues by theirabsolute value, b (cid:48) i = abs( b i ), to generate the matrix B (cid:48) diag .3. Back-transform to a new matrix B (cid:48) = U B (cid:48) diag U − with positive eigenvalues b (cid:48) i .The gradient at the reference position x is preserved andthe original stationary point is changed into a minimum.The effect of converting an unbounded quadratic poten-tial into a bounded quadratic potential is illustrated inFig. 5 for two dimensions. FIG. 5. The reference quadratic potential, possessing amaximum (left) or a saddle point (right), is changed into thepotential printed in light gray, which strictly posseses a min-imum. The reference point around which the conversion isperformed is represented by a red dot.
For a molecular system consisting of N atoms, it is possi-ble to apply quadratic surrogate potentials in two ways.First, one can approximate the reference potential bya single surrogate potential in m = 3 N dimensions forall degrees of freedom of the system. In this case, B is the Hessian matrix. Second, the reference potentialcan be approximated by a combination of N quadraticpotentials in m = 3 dimensions, each depending on oneatomic position. The second option corresponds to ap-proximating the Hessian matrix by its block-diagonal andwill therefore be less accurate. For real-time quantumchemistry, this approximation is tolerable because struc-ture modifications are mainly local. The second option iscomputationally more efficient as it involves the calcula-tion of 6 N second derivatives instead of (3 N )(3 N + 1) / N × N matrix, but of N × m = 3 allows for a fast evaluation dependingon the coordinates of the manipulated atom only. There-fore, a force feedback of 1 kHz can be guaranteed evenif the electronic structure calculations finish at a lowerrate. This results in an excellent responsiveness uponmanipulation.The surrogate potentials also provide all necessary forcesfor the structure response (e.g., steepest-descent struc-ture relaxation). If only the first derivatives of the po-tential energy were used for the relaxation, atoms could easily move so far during the time required for the elec-tronic structure optimization that explicitly calculatedforces and molecular structures hardly correspond. Uponstructure changes, the surrogate potentials limit the mo-tion of atoms until the next electronic structure calcu-lation converges because the quadratic potentials can beguaranteed to feature one single minimum. A schematicrepresentation of the effect of the application of surrogatepotentials in the context of structure relaxation is shownin Fig. 6. xV FIG. 6. Trajectory of an atomic nucleus following a steepest-descent relaxation in real time. The black line represents thereference potential for the atom, which is in parts unknown inthe course of an interactive exploration. The structural evo-lution according to the steepest-descent algorithm is basedon the surrogate potentials shown by dashed parabolas. Themotion of the atom (green arrows) occurs toward the mini-mum of the current surrogate potential. When the surrogatepotential is updated (blue vertical arrows), the motion canproceed again.
The implementation of surrogate potentials in real-timequantum chemistry is represented schematically in Fig. 7.Three independent loops (Fig. 7, middle) run continu-ously. In the first one, electronic structure calculations(of unpredictable execution time) are performed succes-sively. This loop fetches the current molecular structurebefore each electronic structure calculation starts andgenerates the surrogate potentials from the energy andits derivatives when it finishes. The second loop runs at60 Hz and is responsible for the structure relaxation ofthe system based on the steepest-descent algorithm. Thethird loop runs at 1 kHz. It tracks the structure ma-nipulations by the operator and renders a force feedbackthrough the haptic device.eal-time feedback from iterative electronic structure calculations 6
Energy and derivativesSurrogatepotentials MolecularstructureElectronicstructurecalculationStructuralrelaxation60 HzHapticloop1 kHz UserinteractionFeedback
FIG. 7. Schematic flow of force generation in real-time quan-tum chemistry. The loops are represented by curved red ar-rows and the information flow by orange arrows.
IV. NUMERICAL RESULTS
We extended our real-time quantum chemistryframework within the SAMSON molecular edi-tor to reactivity studies based on iterative electronicstructure methods. As SCF methods, we implementedthe semi-empirical Parametrized Method 6 (PM6) andself-consistent DFTB variants . This implementationallowed us to assess the application of mediators in theform of surrogate potentials. In the following, we focuson PM6 as the data source. The PM6 reference potentialis approximated by a combination of three-dimensionalsurrogate potentials, one for each atom of the molecularsystem under consideration. The surrogate potentialsare generated from the PM6 energy, first derivativesand frozen-density second derivatives. Frozen-densitysecond derivatives are obtained by assuming a constantdensity matrix in the expression for the analytic secondderivatives. The generation of potentials from exactanalytical second derivatives has not been consideredbecause their calculation involves the solution of thecoupled perturbed Hartree–Fock equations and iscomputationally significantly more expensive than thecalculation of the gradients, even for semi-empiricalmethods. To assess the reliability of surrogate potentials for real-time feedback, the real-time exploration of a [1,5] hydrideshift reaction with variable side chain length, shown inFig. 8, was studied with and without surrogate poten-tials.The variable size of the side chain allows us to studydifferent molecular sizes with different average executiontimes of the electronic structure optimization. For thedata production, the hydrogen atom displayed in bold in
FIG. 8. [1,5] hydride shift reaction. The side chain lengthdepends on the value of the integer n . Fig. 8 was moved in a straight line at constant speed.
A. Haptic Force
To evaluate the effect of surrogate potentials on the hap-tic force, the forces acting on the hydrogen atom duringthe reaction were recorded. The structures visited dur-ing the reaction were stored for the calculation of the ex-act PM6 forces after the exploration. Exact PM6 forcesand real-time forces were then decomposed into paral-lel and perpendicular components with respect to thedirection of motion of the hydrogen atom. The resultswith and without surrogate potentials for the model re-action with n = 0 , ,
16 are presented in Fig. 9. For thestructure relaxation, the steepest-descent algorithm with γ = 0 . / Hartree was applied ( γ is defined in Eq.(6) below).For small molecular systems with short calculation times(left column in Fig. 9), the exact PM6 forces are repro-duced almost exactly at all times. The differences result-ing from the introduction of surrogate potentials becomemore pronounced for longer execution times of the elec-tronic structure optimization (middle and right columnsin Fig. 9). Without surrogate potentials, the appliedforces correspond to earlier structures, which results in atemporal shift that becomes larger when the SCF proce-dure needs more iterations. With a quadratic surrogatepotential, the surrogate forces are piecewise linear sincethe motion of the hydrogen atom is constant. Impor-tantly, the surrogate potential allows for an immediateforce feedback when the motion of the hydrogen atom isinitiated, as shown in Fig. 10. This is a consequence ofthe increased responsiveness of force feedback upon sud-den changes of exploration direction, made possible bysurrogate potentials.As mentioned above, the purpose of the surrogate po-tential is to deliver reliable forces in unknown areas ofconfiguration space and not to exactly match the refer-ence potential. Although the forces appear to be, onaverage, less accurate when the surrogate potential isapplied, Fig. 9 shows that this requirement is fulfilledeven when the surrogate potential cannot be updated ata high frequency. Furthermore, the provided forces areconservative as they are stronger than the reference whenthey are opposing the manipulation by the operator andeal-time feedback from iterative electronic structure calculations 7 F o r c e / ( − a . u . ) F o r c e / ( − a . u . ) F o r c e / ( − a . u . ) F o r c e / ( − a . u . ) F o r c e / ( − a . u . ) F o r c e / ( − a . u . ) FIG. 9. Parallel (red) and perpendicular (blue) forces calculated in real time for the hydrogen atom involved in the hydrideshift reaction of Fig. 8. The exact PM6 forces are depicted by black dashed lines. The motion of the hydrogen atom wasinitiated at time t = 0 s. Top row: No surrogate potentials are applied. Bottom row: Surrogate potentials are applied. Leftcolumn: n = 0; Middle column: n = 8; Right column: n = 16, with n being the number of double bonds in the side chain forthe model reaction in Fig. 8. F o r c e / ( − a . u . ) F o r c e / ( − a . u . ) F o r c e / ( − a . u . ) F o r c e / ( − a . u . ) F o r c e / ( − a . u . ) F o r c e / ( − a . u . ) FIG. 10. Enlargement of the graphs in Fig. 9 displaying the rendered forces shortly after the motion of the hydrogen isinitiated. Notation as in Fig. 9. eal-time feedback from iterative electronic structure calculations 8smaller than the reference when they are in the samedirection.
B. Structural relaxation
Our framework allows for structural relaxation by apply-ing the steepest-descent algorithm with a frequency of60 Hz. Note that other algorithms can be applied as welland lead to similar conclusions. In a steepest-descentstep, the position x i ( t n ) of each atom i evolves along theatomic force F i ( t n ): x i ( t n +1 ) = x i ( t n ) + γ F i ( t n ) (6)The application of surrogate potentials has a direct effecton the maximal value of γ that can be chosen and there-fore on the maximal relaxation speed. Table II shows,for the model reaction of Fig. 8 with n = 0 , , ,
24, themaximal value of γ that can be employed without lead-ing to artifacts such as artificial oscillations of atomicpositions or tearing apart the molecular system. n = 0 n = 8 n = 16 n = 24w/o V sur .
35 0 .
45 0 .
15 0 . V sur .
35 1 .
45 1 .
50 1 . γ (in (bohr) / Hartree) applica-ble for the model reaction of Fig. 8 with and without surrogatepotentials. The structural relaxation was performed at 60 Hzand n is the number of double bonds in the side chain for themodel reaction of Fig. 8. Without surrogate potentials, the maximal value of γ needs to be adapted to the system size. The larger themolecular structures, the longer the execution time ofthe electronic structure optimization is and, at constantstructure-relaxation frequency, the more steps are to beperformed based on the same forces. Applying surrogatepotentials decouples γ and the size of the system, andtherefore allows for larger values of γ . In this case, themaximal values for γ rather depend on the local shape ofthe potential, i.e., of the second derivative of the surro-gate potential, and on the frequency at which the struc-ture relaxation is performed. By default, our frameworkwill employ γ = 0 . / Hartree independent of themolecular system size if surrogate potentials are applied.
V. CONCLUSION
A solution is presented to mediate between the inher-ently unpredictable execution times of iterative electronicstructure calculations and the time requirements of real-time feedback. In the context of real-time quantumchemistry, such feedback relies on quantum-mechanical forces and must be generated at high frequencies thatcannot be guaranteed when iterative electronic structuremethods are applied.The strategy followed in this work is to calculate theatomic forces from surrogate potentials that can easilybe evaluated for molecular structures for which the elec-tronic structure is still unknown. For consecutive molec-ular structures visited during real-time reactivity explo-rations, one or several surrogate potentials are generatedfrom the electronic energy and its first and second deriva-tives to approximate the reference potential energy sur-face. Surrogate potentials such as bounded quadraticpotentials allow for a fast evaluation of forces compatiblewith real-time requirements.To calculate forces at such high frequencies offers imme-diate benefits for real-time feedback in the form of on-the-fly structure relaxation and haptic force feedback.In the context of structure relaxation, surrogate poten-tials avoid the need for special measures such as freez-ing the system when electronic structure calculations donot converge or are very time-consuming — the molec-ular structure will evolve toward the minimum of thesurrogate potential and remain there until the surrogatepotential is updated. Consequently, the structural evolu-tion can be performed faster (translated into larger stepsin the steepest-descent algorithm) since the surrogate po-tentials guarantee to lead to a minimum. In particular,the step size of the structure relaxation algorithm doesnot need to be adjusted for different system sizes or cal-culation times.For haptic feedback, the ability to update forces at 1 kHzimproves substantially on the responsiveness of forcefeedback upon structure manipulations by the operator.Furthermore, in the case of non-converging calculations,the surrogate potential will drive the operator back to itsminimum and hinder large manipulations from his partwhen insufficient information about the reference forcesis available.The presented procedure is independent of the sys-tem size and is also advantageous when the executiontime of the electronic structure calculation is constantand predictable. Despite the additional computationaltime needed for the calculation of (approximate) secondderivatives, the application of surrogate potentials provesnecessary and beneficial for interactive reactivity studiesbased on iterative methods.
ACKNOWLEDGMENTS
This work was generously supported by ETH ResearchGrant ETH-20 15-1 and ETH Pioneer Fellowship GrantPIO-11-14-2.eal-time feedback from iterative electronic structure calculations 9 Haag, M. P. and Reiher, M.,
Faraday Discuss. , 2014, , 89–118. Marti, K. H. and Reiher, M.,
J. Comput. Chem. , 2009, , 2010–2020. Haag, M. P.; Marti, K. H. and Reiher, M.,
ChemPhysChem ,2011, , 3204–3213. Haag, M. P. and Reiher, M.,
Int. J. Quantum Chem. , 2013, ,8–20. Stone, J. E.; Gullingsrud, J. and Schulten, K. In
Proceedings ofthe 2001 Symposium on Interactive 3D Graphics , I3D ’01, pages191–194, New York, NY, USA, 2001. ACM. Luehr, N.; Jin, A. G. B. and Mart´ınez, T. J.,
J. Chem. TheoryComput. , 2015, , 4536–4544. Bosson, M.; Richard, C.; Plet, A.; Grudinin, S. and Redon, S.,
J. Comput. Chem. , 2012, , 779–790. S. Redon et al. Hanwell, M. D.; Curtis, D. E.; Lonie, D. C.; Vandermeersch, T.;Zurek, E. and Hutchison, G. R.,
J. Cheminform. , 2012, , 1–17. Haag, M. P.; Vaucher, A. C.; Bosson, M.; Redon, S. and Reiher,M.,
ChemPhysChem , 2014, , 3301–3319. Porezag, D.; Frauenheim, T.; K¨ohler, T.; Seifert, G. andKaschner, R.,
Phys. Rev. B , 1995, , 12947–12957. Seifert, G.; Porezag, D. and Frauenheim, T.,
Int. J. QuantumChem. , 1996, , 185–192. Bosson, M.; Grudinin, S.; Bouju, X. and Redon, S.,
J. Comput.Phys. , 2012, , 2581–2598. Stewart, J. J. P.,
J. Mol. Model. , 2007, , 1173–1213. Elstner, M.; Porezag, D.; Jungnickel, G.; Elsner, J.; Haugk, M.;Frauenheim, T.; Suhai, S. and Seifert, G.,
Phys. Rev. B , 1998, , 7260–7268. Gaus, M.; Cui, Q. and Elstner, M.,
J. Chem. Theory Comput. ,2011, , 931–948. Frisch, M.; Scalmani, G.; Vreven, T. and Zheng, G.,
Mol. Phys. ,2009,107