Signifikans
Syre-base-reaktioner er blandt de vigtigste kemiske processer. Alligevel mangler vi en enkel måde at beskrive denne klasse af reaktioner som en funktion af atomkoordinaterne. Når H+ og dets konjugerede anion OH- er opløst i vand, har H+ og dets konjugerede anion OH- nemlig en meget fluxional struktur, som er vanskelig at fastlægge. Her løser vi dette problem ved at anlægge det synspunkt, at syre-base-reaktioner beskrives som en ligevægt mellem den opløste substans og hele opløsningsmidlet. Dette gør det muligt at identificere generelt anvendelige deskriptorer. Som følge heraf er det nu muligt at udføre kvantitativ forbedret samplingssimulering af syre-base-reaktioner i vand og i andre miljøer som f.eks. zeolithhulrum eller ved overflader.
Abstract
Syre-base-reaktioner er allestedsnærværende i naturen. Forståelse af deres mekanismer er afgørende inden for mange områder, fra biokemi til industriel katalyse. Desværre giver eksperimenter kun begrænset information uden megen indsigt i den molekylære adfærd. Atomistiske simuleringer kan supplere eksperimenterne og kaste et værdifuldt lys over mikroskopiske mekanismer. De store frie energibarrierer, der er forbundet med protonedissociation, gør det imidlertid obligatorisk at anvende forbedrede prøveudtagningsmetoder. Her udfører vi en ab initio-molekylær dynamik (MD)-simulering og forbedrer prøveudtagningen ved hjælp af metadynamik. Dette er blevet muligt ved at indføre deskriptorer eller kollektive variabler (CV’er), som er baseret på et konceptuelt anderledes syn på syre-base-ligevægte. Vi afprøver med succes vores fremgangsmåde på tre forskellige vandige opløsninger af eddikesyre, ammoniak og bikarbonat. Disse er repræsentative for sur, basisk og amfoterisk adfærd.
- syre-base
- metadynamik
- kollektive variabler
- forbedret prøveudtagning
Syre-base-reaktioner spiller en central rolle inden for mange grene af kemien. Uorganiske komplekseringsreaktioner, proteinfoldning, enzymatiske processer, polymerisering, katalytiske reaktioner og mange andre transformationer inden for forskellige områder er følsomme over for ændringer i pH. Forståelse af pH’s rolle i disse reaktioner indebærer, at man har kontrol over deres reaktivitet og kinetik.
Den afgørende betydning af pH har stimuleret indsamling af en stor mængde data om syre-base ligevægte. Disse måles typisk i gas- og kondenserede faser ved hjælp af spektroskopiske og potentiometriske teknikker. Der er imidlertid praktiske begrænsninger for nøjagtigheden af disse metoder, især i kondenserede faser (1). Desuden er det meget vanskeligt at uddrage et mikroskopisk billede af de involverede processer fra eksperimentelle data. Det er derfor ikke overraskende, at syre-base ligevægt har været genstand for intens teoretisk aktivitet (1⇓⇓⇓⇓⇓⇓⇓⇓⇓⇓-12).
Syreniveauet for en kemisk art i vand kan udtrykkes i form af pKa, den negative logaritme af syredissociationskonstanten. Der findes to måder at beregne disse værdier på, den ene statisk og den anden dynamisk.
Den mest standardiserede tilgang er den statiske, hvor opløsningsfasens frie energier og dermed pKa’er fås ved at lukke en Born-Haber-cyklus bestående af gasfase- og solvationsfri energier (1, 3⇓⇓⇓⇓-7). Selv om den statiske tilgang er yderst vellykket i mange tilfælde, har den statiske tilgang nogle begrænsninger. Der skal vælges en solvationsmodel, og kontinuumopløsningsmiddelmodeller har en begrænset nøjagtighed. Dette gælder især i systemer som zeolitter eller proteiner, der er karakteriseret ved uregelmæssige hulrum, hvor en implicit beskrivelse af opløsningsmidlet er en udfordring. Det er klart, at der ikke kan opnås dynamiske oplysninger ved en sådan fremgangsmåde. Desuden kan der være konkurrerende reaktioner, som der ikke kan tages hensyn til, medmindre de eksplicit indgår i modellen.
I princippet kan disse begrænsninger ophæves i en mere dynamisk tilgang baseret på molekylær dynamik (MD)-simuleringer, hvor opløsningsmiddelmolekylerne behandles eksplicit. Hvis man havde ubegrænset computertid, ville sådanne simuleringer kunne udforske alle mulige veje og tildele den relative statistiske vægt til de forskellige tilstande. Desværre forpurrer tilstedeværelsen af kinetiske flaskehalse denne mulighed ved at fange systemet i metastabile tilstande, da de forskellige protoniseringstilstande er adskilt af store barrierer. Desuden brydes og dannes der kemiske bindinger i syre-base-reaktioner. Dette kræver anvendelse af ab initio MD, hvor de interatomare kræfter beregnes i farten ud fra teorier om elektronisk struktur. Dette gør beregningen dyrere og reducerer yderligere den tidsskala, der kan udforskes.
For at overvinde denne vanskelighed bliver brugen af forbedrede samplingmetoder (13), der fremskynder udforskningen af konfigurationsrummet, obligatorisk. En meget populær klasse af forbedrede prøveudtagningsmetoder er baseret på identifikation af de frihedsgrader, der er involveret i den langsomme reaktion af interesse. Disse frihedsgrader kaldes normalt kollektive variabler (CV’er) og udtrykkes som eksplicitte funktioner af atomkoordinaterne R. Samplingen forbedres derefter ved at tilføje en bias, der er en funktion af de valgte CV’er (14⇓-16). Endvidere har udformningen af et korrekt sæt af gode CV’er også en dybere betydning. Vellykkede CV’er indfanger på en kondenseret måde problemets fysik, identificerer dets langsomme frihedsgrader og fører til en nyttig modelistisk beskrivelse af processen.
I kemiske standardreaktioner er dette relativt enkelt, da der kan tildeles reaktanter og produkter veldefinerede strukturer (17⇓-19). Dette er ikke tilfældet for syre-base-reaktioner, hvor en proton tilføjes til eller trækkes fra den opløste substans. Når denne proces har fundet sted, opløses vandioner (H+ eller/og OH-), og deres struktur bliver ubestemmelig. Faktisk kan vandioner hurtigt diffundere i mediet via en Grotthuss-mekanisme (20). De er blevet meget fluktuationelle, og identiteten af de atomer, der indgår i deres struktur, ændres løbende. Det er således vanskeligt at indfange disse arters natur i en eksplicit analytisk funktion af R. I betragtning af relevansen af syre-base-reaktioner er der imidlertid gjort mange forsøg på at definere disse enheder (8⇓⇓⇓⇓-12). Desværre har disse CV’er en ad hoc-natur, og selv om de er vellykkede i dette eller hint tilfælde, kan de ikke anvendes generelt.
For at opbygge generelle og nyttige CV’er foretager vi to konceptuelle trin. Det ene er at se på syre-base-processen som en reaktion, der kun involverer nogle få dele, nemlig hele opløsningsmidlet og de reagerende rester i det solventerede molekyle. Når der f.eks. kun er én type dissocierende rest, tænker vi på syre-base-ligevægten som en reaktion af typenA+H2NON⇌Bq0+H2N+q1ONq1,hvor N er antallet af vandmolekyler, A og B er henholdsvis et generisk syre-base-molekyle i opløsning og dets konjugerede art, q0 og q1 er hele tal, der kan antage værdierne +1 og -1 alt efter artens syre-base-adfærd, og q1+q0=0.
Dette indebærer, at vi ikke ser på opløsningsmidlet som et sæt af molekyler, der konkurrerer om at reagere med syre-base-arten. Vi betragter snarere opløsningsmidlet i sin helhed som en af de to addukter. Dette synspunkt er især relevant i polære opløsningsmidler som vand, der er kendetegnet ved stærkt strukturerede netværk. I dette tilfælde ændrer tilstedeværelsen af et overskud eller et underskud af protoner lokalt netværksstrukturen, og denne forvrængning forplanter sig langs hele netværket.
Siden Wicke og Eigen (21) og Zundel og Metzger (22) har forskerne kæmpet med, hvor mange molekyler der skal medregnes i definitionen af forstyrrelsen (23⇓-25). Da der ikke findes fysiske parametre, der kan give et klart og utvetydigt svar på dette spørgsmål, omgår idéen om at betragte opløsningsmidlet som en helhed dette problem. Opløsningsmidlet er således ikke blot et medium med en passiv rolle, men betragtes som et ensemble af molekyler, der tilsammen bidrager til dannelsen af det konjugerede syre-basepar. Dette synspunkt er meget tættere på det oprindelige synspunkt foreslået af Brønsted og Lowry, hvor reaktionen kan ses som en simpel udveksling af et hydrogenkation mellem et syre-basepar.
For at reaktionen kan finde sted, skal forstyrrelsens centrum bevæge sig væk fra den opløste stofmængde. Det andet vigtige trin er således at overvåge forstyrrelsens centrum. På grund af Grotthuss-lignende mekanismer bevæger forstyrrelsen sig langs netværket. Dette kan føre til forskellige definitioner af forstyrrelsens centrum. Hvis vi imidlertid tesselliserer hele rummet ved hjælp af Voronoi-polyedre med centrum på vandsyreatomer, kan vi utvetydigt henføre hvert hydrogenatom til et og kun et af disse polyedre. Det sted, hvis Voronoi-polyhedron indeholder et unormalt antal protoner, tages som centrum for forstyrrelsen (fig. 1).
To eksempler på opdeling af rummet. (Til venstre) Vi viser en konvektionel tilgang, hvor afstanden fra iltatomet bruges til at definere dets omgivelser. Der kan tydeligt ses kunstige overlejringer. (Til højre) Voronoi-tesselationen lider ikke af disse mangler.
Dette synspunkt giver metoden en meget generel karakter, der gør den anvendelig på ethvert syre-basesystem, uden at det er nødvendigt på forhånd at fastsætte de reagerende par. Det er således muligt at udforske alle de relevante protoneringstilstande, selv i systemer, der består af mere end ét syre-basepar.
Denne generelle tilgang gør det muligt at definere CV’er uden at skulle pålægge specifikke strukturer eller vælge identiteten af de involverede atomer. Vi tester vores metode ved at udføre metadynamiske simuleringer i et svagt syretilfælde (eddikesyre), i en svag base (ammoniak) og i en amfoterisk art (bikarbonat), der er valgt som benchmarks på grund af deres sammenlignelige styrke, men forskellige syre-baseadfærd.
Metoder
Som diskuteret ovenfor introducerer vi to CV’er, en relateret til protoneringstilstanden og en anden, der lokaliserer ladningsdefekter og måler deres relative afstand. Begge disse CV’er har brug for en robust definition for at tildele hydrogenatomerne til det respektive syre-base-sted. For at opnå dette resultat opdeler vi hele rummet i Voronoi-polyedre med centrum på de syre-base-steder i, der er placeret på Ri. Disse steder omfatter alle de atomer, der er i stand til at bryde og danne bindinger med en syreproton. Standardfordelingen af Voronoi-rummet beskrives ved et sæt indeksfunktioner wi(r) centreret på de forskellige Ris, således at wi(r)=1, hvis det i-te atom er det nærmeste til r, og wi(r) = 0 i modsat fald. For at kunne anvendes i forbedrede prøveudtagningsmetoder skal CV’er være differentierbare. Med henblik herpå indfører vi en glat version af indeksfunktionerne, wis(r). Disse er defineret ved hjælp af softmax-funktionerwis(r)=e-λ|Ri-r|∑me-λ|Rm-r|,hvor i og m løber over alle syre-base-stederne, og λ styrer den stejlhed, hvormed kurverne aftager til 0, dvs. funktionens selektivitet. Med et passende valg af λ opnår denne definition det ønskede resultat som vist i fig. 2. På denne måde tildeles et hydrogenatom i en position Rj til polyederet centreret på stedet i med vægten wi(Rj). Derefter er det samlede antal hydrogenatomer, der er tildelt det i’te syre-base-sted,Wi=∑j∈Hwis(Rj),hvor summeringen på j løber over alle hydrogenatomer.
Smooth tessellation af et 2D-rum med celler centreret på de tre iltatomer i vandmolekylet. De flade blå områder repræsenterer den del af rummet, hvor funktionen antager en værdi på 1, og de gule områder repræsenterer grænserne mellem cellerne. Denne overflade er opnået med en værdi af λ=4.
Man kan til hvert syre-base-sted knytte en referenceværdi Wi0, der tæller antallet af bundne hydrogenatomer i den neutrale tilstand. Forskellen mellem den øjeblikkelige værdi af hydrogenatomer og referenceværdien erδi=Wi-Wi0. Når δi er forskellig fra nul, vil det signalere, om det i-te sted har fået eller mistet en proton. I tilfælde af vands oxygenatomer har en hydroniumion en δi=+1, mens en hydroxidion har δi=-1.
Vi grupperer derefter syre-base-stederne i arter. For eksempel vil antallet af arter Ns i tilfælde af den enkleste aminosyre glycin i vandig opløsning være lig med 3. Alle vandsyreatomer tilhører en art; derefter tæller man i en anden art de to carboxylsyreatomer, og endelig betragter man som den tredje art nitrogenatomet i aminogruppen.
I ånden af dette arbejde tæller vi det samlede overskud eller defekt af protoner, der er knyttet til hver art: qk=∑i∈kδi.Dette indebærer, at vi ikke er interesseret i den specifikke identitet af det reagerede sted, men om den k-te art i sin helhed er steget (qk=+1), er faldet (qk=-1) eller ikke har ændret sit antal protoner. Hvis vi betragter en opløst stof med kun én reaktiv del, kan hver mulig tilstand i systemet beskrives ved en af de tre 2D-vektorer (0,0), (-1,1) eller (1,-1).
I det generelle tilfælde kan hver protoneringstilstand beskrives ved en vektor q→=(q0,q1,…qNs-1) med en dimension lig med antallet af uækvivalente reaktive steder, Ns. En mere udtømmende forklaring findes i SI Appendix.
Til brug ved forbedret prøveudtagning skal disse vektorer udtrykkes som en skalarfunktion f=f(q→), således at f for hvert fysisk relevant q→ opnår værdier, der kan skelne mellem de forskellige overordnede protoneringstilstande. Der er uendeligt mange måder at konstruere en skalar på ud fra en vektor. Muligvis er det enkleste valg at skrive f(q→)=X→⋅q→ og, for at skelne mellem forskellige protoneringstilstande, at vælge X→=(20,21,22,…2Ns-1).
Dette fører til følgende definition for CV’et, der bruges til at beskrive systemets protoneringstilstand,sp=∑k=0Ns-12k⋅qk,hvor k er de indekser, der bruges til at mærke de respektive reaktive stedgrupper. I SI Appendix udarbejdes et eksempel i detaljer. CV’et gøres naturligvis kontinuerligt ved brug af Wi i beregningen af δi, der er nødvendig for at evaluere qk i Eq. 5.
Det andet CV er en summering af afstandene mellem alle syre-base-stederne multipliceret med deres partielle ladning δi,sd=∑i,m>i-rim⋅δi⋅δm,hvor indeksene i og m løber over alle syre-base-stederne, der tilhører forskellige k-grupper, og rim er afstanden mellem de to atomer. På denne måde giver kun det syre-basepar, der har udvekslet en proton, et bidrag forskelligt fra nul. Ligning 7 er kun gyldig, når der er et enkelt konjugeret syre-basepar til stede. På grund af forspænding kan det dog lejlighedsvis forekomme, at der dannes flere syre-basepar. For at undgå at udtage prøver af disse meget usandsynlige hændelser anvender vi en begrænsning på antallet af par. Yderligere oplysninger findes i SI-tillægget.
Resultater
Vi har anvendt vores metode på tre vandige opløsninger af eddikesyre, ammoniak og bicarbonat som repræsentationer af henholdsvis en svag syre, en svag base og en amfotere forbindelse. Opsætningerne i alle tre simuleringer er identiske bortset fra identiteten af de solvente molekyler. Dette sikrer, at resultatet afspejler den forskellige kemi i disse tre systemer, og at der ikke er nogen bias på grund af startbetingelsen.
Hver simulation af systemerne blev udført med Born-Oppenheimer MD-simuleringer kombineret med godt tempereret metadynamik (14, 26) ved hjælp af CP2K-pakken (27) patchet med PLUMED 2 (28) og stærkt begrænset og passende normeret funktionel (29) for xc-energien, Exc. Se SI Appendix for nærmere oplysninger.
I Fig. 3 plotter vi de frie energioverflader (FES) som en funktion af sp og sd. Disse FES’er gengiver levende den forventede opførsel. De har alle et minimum ved sp=0, der svarer til den tilstand, hvor der ikke er nogen ladninger til stede i opløsningsmidlet. I eddikesyre-FES (fig. 3A) er der et andet minimum tæt på sp=-1, som afspejler dens syreadfærd. I modsætning hertil viser ammoniak-FES (fig. 3B) et andet minimum tæt ved sp=1. Formen af ammoniak- og eddikesyre-FES’erne er omtrent beslægtet med en spejlsymmetri, der afspejler deres kontrastfyldte adfærd. På samme måde afspejler den bikarbonatsymmetriske FES (fig. 3C) dens amfoteriske karakter.
(A-C) Frienergioverflader langs sp og sd for eddikesyre (A), ammoniak (B) og bikarbonat (C) i vandig opløsning. Farvebjælker angiver den frie energi udtrykt i kJ⋅mol-1-enheder. CV sd er udtrykt i angstromer.
Når det konjugerede par er dannet, begynder sd at antage positive værdier svarende til adskillelse og diffusion af det konjugerede par. Sammenlignet med den ikke-dissocierede tilstand, hvor kun sd=0 er tilladt, viser tilstande, hvor der er et konjugeret par til stede, en langstrakt form af bassinerne langs denne variabel. Dette skyldes hydronium- og hydroxidionernes diffusionsadfærd i opløsningen, som gør et kontinuum af afstande tilgængeligt. Desuden kan vi langs denne CV observere en barriere omkring 1,5, der svarer til bruddet af den kovalente binding mellem hydrogenatomet og syre-base-stedet.
Konklusioner
Den generelle anvendelighed af denne metode til systemer med forskellige naturer er et vigtigt skridt i deres forståelse og beskrivelse. Skemaet kan udvides til at omfatte kvantekerneeffekter ved hjælp af path integral molecular dynamics (30). Dette ville være af kvantitativ betydning, da f.eks. pKa-værdier påvirkes af deuterering. Desuden gør fraværet af antagelser eller pålæg om reaktive kandidater eller reaktionsveje det muligt at udvide denne metode til systemer af stigende kompleksitet, som ikke kan behandles med traditionelle metoder. Eksempler på spørgsmål, der nu kan besvares, er tautomeriske ligevægte i biokemiske processer og syreadfærd i zeolitter og på overfladen af oxider, der udsættes for vand.
Anerkendelser
Beregningerne blev udført på ETH Euler-klyngen og på Mönch-klyngen på Swiss National Supercomputing Center. Denne forskning blev støttet af Den Europæiske Union Grant ERC-2014-AdG-670227/VARMET.
Publiceret under PNAS-licens.