Analyse af smelter og væsker fra Ab Initio Molecular Dynamics Simulations med UMD-pakken

Summary

Smelter og væsker er allestedsnærværende vektorer for massetransport i naturlige systemer. Vi har udviklet en open source-pakke til at analysere ab initio molekylær-dynamik simuleringer af sådanne systemer. Vi beregner strukturelle (limning, klyngedannelse, kemisk speciation), transport (diffusion, viskositet) og termodynamiske egenskaber (vibrationelle spektrum).

Abstract

Vi har udviklet en Python-baseret open source-pakke til at analysere resultaterne stammer fra ab initio molekylær-dynamik simuleringer af væsker. Pakken er bedst egnet til applikationer på naturlige systemer, som silikat og oxid smelter, vandbaserede væsker, og forskellige superkritiske væsker. Pakken er en samling af Python scripts, der omfatter to store biblioteker, der beskæftiger sig med filformater og med krystallografi. Alle scripts køres på kommandolinjen. Vi foreslår et forenklet format til lagring af atombaner og relevante termodynamiske oplysninger om simuleringerne, som gemmes i UMD-filer, der står for Universal Molecular Dynamics. UMD-pakken gør det muligt at beregning af en række strukturelle, transport- og termodynamiske egenskaber. Begyndende med parfordelingsfunktionen definerer den bindingslængder, opbygger en interatomar konnektivitetsmatrix og bestemmer til sidst den kemiske speciation. Bestemmelse af den kemiske arts levetid gør det muligt at køre en fuld statistisk analyse. Derefter beregner dedikerede scripts de gennemsnitlige firkantede forskydninger for atomer såvel som for den kemiske art. Den implementerede selvkorrelationsanalyse af atomhastighederne giver diffusionskoefficienterne og vibrationsspektret. Den samme analyse, der anvendes på stress giver viskositet. Pakken er tilgængelig via GitHub hjemmeside og via sin egen dedikerede side af ERC IMPACT projektet som open-access pakke.

Introduction

Væsker og smelter er aktive kemiske og fysiske transportvektorer i naturlige miljøer. De forhøjede satser for atom diffusion favorisere kemiske udvekslinger og reaktioner, den lave viskositet kombineret med varierende opdrift favoriserer stor masseoverførsel, og krystal-smelte tæthed relationer favorisere lagdeling inde planetariske organer. Fraværet af et periodisk gitter, typiske høje temperaturer, der kræves for at nå den smeltede tilstand, og vanskeligheden ved at slukke gør den eksperimentelle bestemmelse af en række indlysende egenskaber, som tæthed, diffusion og viskositet, ekstremt udfordrende. Disse vanskeligheder gør alternative beregningsmetoder stærke og nyttige værktøjer til at undersøge denne klasse af materialer.

Med fremkomsten af computerkraft og tilgængeligheden af supercomputere anvendes der i øjeblikket to store numeriske atomistiske simuleringsteknikker til at studere den dynamiske tilstand af et ikke-krystallinsk atomistisk system, Monte ^Carlo1 og molekylær dynamik (MD)^1,2. I Monte Carlo simuleringer er konfigurationsrummet tilfældigt udtaget; Monte Carlo metoder viser lineær skalering i parallelisering, hvis alle prøveudtagning observationer er uafhængige af hinanden. Kvaliteten af resultaterne afhænger af kvaliteten af generatoren af tilfældige tal og prøveudtagningens repræsentativitet. Monte Carlo metoder viser lineær skalering i parallelisering, hvis prøveudtagningen er uafhængig af hinanden. I molekylær dynamik (MD) prøves konfigurationsrummet efter tidsafhængige atombaner. Startende fra en given konfiguration, er atomare baner beregnes ved at integrere newtonske ligninger af bevægelse. De interatomare kræfter kan beregnes ved hjælp af model interatomare potentialer (i klassisk MD) eller ved hjælp af første principper metoder (i ab initio, eller første principper, MD). Kvaliteten af resultaterne afhænger af længden af banen og dens evne til ikke at blive tiltrukket af lokale minima.

Molekylære dynamik simuleringer indeholder en overflod af information, alle relateret til den dynamiske adfærd i systemet. Termodynamiske gennemsnitlige egenskaber, som intern energi, temperatur og tryk, er ret standard at beregne. De kan udvindes fra outputfilen eller -filerne af simuleringerne og beregnes i gennemsnit, mens mængder, der er direkte forbundet med atomernes bevægelse samt deres gensidige relation, skal beregnes efter udvinding af atompositioner og hastigheder.

Derfor har en stor indsats været dedikeret til at visualisere resultaterne, og forskellige pakker er tilgængelige i dag, på forskellige platforme, open source eller ej [^Ovito3, ^VMD4, ^Vesta5, ^Travis6 osv.]. Alle disse visualiseringsværktøjer beskæftiger sig effektivt med interatomare afstande, og som sådan tillader de effektiv beregning af parfordelingsfunktioner og diffusionskoefficienter. Forskellige grupper, der udfører store molekylære dynamik simuleringer har proprietære software til at analysere forskellige andre egenskaber som følge af simuleringer, undertiden i shareware eller andre former for begrænset adgang til samfundet, og undertiden begrænset i omfang og brug til nogle specifikke pakker. Sofistikerede algoritmer til at udtrække oplysninger om interatomar limning, geometriske mønstre og termodynamik er udviklet og implementeret i nogle af disse pakker3,4,5,6,7 osv.

Her foreslår vi UMD-pakken - en open source-pakke skrevet i Python for at analysere produktionen af molekylære dynamics simuleringer. UMD-pakken giver mulighed for beregning af en lang række strukturelle, dynamiske og termodynamiske egenskaber (figur 1). Pakken er tilgængelig via GitHubs hjemmeside (https://github.com/rcaracas/UMD_package) og via en dedikeret side (http://moonimpact.eu/umd-package/) af ERC IMPACT-projektet som en open-access-pakke.

For at gøre det universelt og lettere at håndtere, er vores tilgang først at udtrække alle de oplysninger, der er relateret til den termodynamiske tilstand og atomare baner fra outputfilen i den faktiske molekylære dynamikkørsel. Disse oplysninger gemmes i en dedikeret fil, hvis format er uafhængigt af den oprindelige MD-pakke, hvor simuleringen blev kørt. Vi navngiver disse filer "umd" filer, som står for Universal Molecular Dynamics. På denne måde kan vores UMD-pakke nemt bruges af enhver ab initio gruppe med enhver software, alt sammen med en minimal indsats for tilpasning. Det eneste krav for at bruge den nuværende pakke er at skrive den relevante parser fra output af den pågældende MD-software i umd-filformatet, hvis dette endnu ikke er eksisterende. For øjeblikket leverer vi sådanne parsere til ^VASP8 og QBox9-pakkerne.

Figur 1: Rutediagram over UMD-biblioteket.
Fysiske egenskaber er i blå, og store Python scripts og deres muligheder er i rødt. Klik her for at se en større version af dette tal.

Umd-filerne er ASCII-filer. typisk udvidelse er "umd.dat", men ikke obligatorisk. Alle analysekomponenter kan læse ASCII-filer i UMD-formatet, uanset det faktiske navneudvidelse. Men nogle af de automatiske scripts designet til at udføre hurtige store statistikker over flere simuleringer specifikt kigge efter filer med umd.dat udvidelse. Hver fysisk egenskab udtrykkes på én linje. Hver linje starter med et nøgleord. På denne måde formatet er meget tilpasningsdygtige og giver mulighed for nye egenskaber, der skal føjes til umd-filen, alt imens bevare dens læsbarhed i hele versioner. De første 30 linjer i umd-filen af simuleringen af pyrolit ved 4,6 GPa og 3000 K, der anvendes nedenfor i diskussionen, er vist i figur 2.

Figur 2: Begyndelsen af umd-filen, der beskriver simuleringen af flydende pyrolit ved 4,6 GPa og 3000 K.
Overskriften efterfølges af beskrivelsen af hvert øjebliksbillede. Hver egenskab er skrevet på én linje, der indeholder navnet på den fysiske egenskab, værdien eller værdien og enhederne, alle adskilt af mellemrum. Klik her for at se en større version af dette tal.

Alle umd-filer indeholder en header, der beskriver indholdet af simuleringscellen: antallet af atomer, elektroner og atomtyper samt detaljer for hvert atom, såsom dets type, kemisk symbol, antal valenselektroner og dets masse. En tom linje markerer slutningen af sidehovedet og adskiller den fra hoveddelen af umd-filen.

Derefter er hvert trin i simuleringen detaljeret. For det første angives de øjeblikkelige termodynamiske parametre, hver på en anden linje, med angivelse af (i) navnet på parameteren, som energi, belastninger, tilsvarende hydrostatisk tryk, tæthed, volumen, gitterparametre osv., (ii) dens værdi (r) og (iii) dens enheder. En tabel, der beskriver atomer kommer næste. En overskriftslinje giver de forskellige mål, f.eks. kartesiske positioner, hastigheder, gebyrer osv. Derefter er hvert atom beskrevet på en linje. Ved grupper på tre, svarende til de tre x, y, z akser, posterne er: de reducerede positioner, de kartesiske positioner foldet ind i simuleringscellen, de kartesiske positioner (der korrekt tager hensyn til det faktum, at atomer kan krydse flere enhedsceller under en simulering), atomare hastigheder og atomstyrkerne. De sidste to indgange er skalarer: opladning og magnetisk øjeblik.

To store biblioteker sikrer, at hele pakken fungerer korrekt. Det umd_process.py bibliotek beskæftiger sig med umd filer, som læsning og udskrivning. Det crystallography.py bibliotek beskæftiger sig med alle de oplysninger, der er relateret til den faktiske atomare struktur. Den underliggende filosofi i crystallography.py bibliotek er at behandle gitteret som et vektorielt rum. Enhedscelleparametrene sammen med deres orientering repræsenterer basisvektorerne. "Rummet" har en række skalaregenskaber (specifik volumen, densitet, temperatur og specifikt antal atomer), termodynamiske egenskaber (intern energi, tryk, varmekapacitet osv.) og en række spændte egenskaber (stress og elasticitet). Atomer udfylder dette rum. "Lattice" klassen definerer dette ensemble, sammen med forskellige få korte beregninger, som specifikke volumen, tæthed, opnå den gensidige gitter fra den direkte, osv. Klassen "Atomer" definerer atomer. De er karakteriseret ved en række skalaregenskaber (navn, symbol, masse, antal elektroner osv.) og en række vektoregenskaber (positionen i rummet, enten i forhold til det vektorbaserede grundlag, der er beskrevet i gitterklassen, eller i forhold til universelle kartesiske koordinater, hastigheder, kræfter osv.). Bortset fra disse to klasser indeholder crystallography.py biblioteket en række funktioner til at udføre en række tests og beregninger, såsom atomafstande eller celle multiplikation. Elementernes periodiske system medtages også som en ordbog.

De forskellige komponenter i umd-pakken skriver flere outputfiler. Som hovedregel er de alle ASCII-filer, alle deres poster er adskilt af faner, og de er lavet så selvforklarende som muligt. For eksempel angiver de altid tydeligt den fysiske ejendom og dens enheder. Umd.dat filer fuldt ud overholder denne regel.

Protocol

1. Analyse af molekylær-dynamik kører

BEMÆRK: Pakken er tilgængelig via GitHubs hjemmeside (https://github.com/rcaracas/UMD_package) og via en dedikeret side (http://moonimpact.eu/umd-package/) af ERC IMPACT-projektet som en open-access-pakke.

1. Uddrag hvert enkelt sæt fysiske egenskaber ved hjælp af et eller flere dedikerede Python-scripts fra pakken. Kør alle scripts på kommandolinjen. de anvender alle en række flag, som er så konsekvente som muligt fra et script til et andet. Flagene, deres betydning og standardværdierne er alle opsummeret i tabel 1.

Flag	Betydning	Script ved hjælp af det	Standardværdi
-h	Kort hjælp	al
-f	UMD-filnavn	al
-I	Termiske trin, der skal kasseres	al	0
-I	Inputfil, der indeholder de interatomare bindinger	artsdannelse	bonds.input
-s	Prøveudtagning af hyppigheden	msd, speciation	1 (hvert skridt overvejes)
-a	Liste over atomer eller anioner	artsdannelse
-c	Liste over kationer	artsdannelse
- L	Bond længde	artsdannelse	2
-t	Temperatur	vibrationer, rheologi
-v	Diskretisering af bredden af prøveudtagningsvinduet i banen til den gennemsnitlige firkantede forskydningsanalyse	Msd	20
-z	Diskretisering af starten af prøveudtagningsvinduet i banen til den gennemsnitlige firkantede forskydningsanalyse	Msd	20

Tabel 1: De fleste almindelige flag, der bruges i UMD-pakken, og deres mest almindelige betydning.

2. Start med at omdanne outputtet af MD-simuleringen, der udføres i en første-principkode, som ^VASP8 eller ^QBox9, til en UMD-fil.
   1. Hvis MD-simuleringerne blev udført i VASP, skal du skrive følgende på kommandolinjetypen:
      VaspParser.py -f -i
      hvor –f flag definerer navnet på VASP OUTCAR-filen, og -i termiskisering længde.
      BEMÆRK: Det første trin, defineret af -jeg giver mulighed for at kassere de første trin i simuleringer, som repræsenterer thermalization. I et typisk molekylær-dynamikløb repræsenterer den første del af beregningen termikiseringen, dvs. den tid, det tager systemet for alle atomer at beskrive en gaussisk-lignende fordeling af temperaturen, og for hele systemet at udvise udsving i temperatur, tryk, energi osv. omkring ligevægtsværdier. Denne termiskgørelsesdel af simuleringen bør ikke tages i betragtning ved analyse af væskens statistiske egenskaber.
3. Transformer . umd-filer i . xyz filer til at lette visualisering på forskellige andre pakker, som ^VMD4 eller ^Vesta5. På kommandolinjetypen:
   umd2xyz.py -f -i -s
   hvor –f definerer navnet på . umd fil, -i definerer thermalization periode, der skal kasseres, og -s hyppigheden af prøveudtagningen af bane gemt i . umd-fil . Standardværdierne er –i 0 -s 1, dvs.
4. Tilbagefør umd-filen til VASP-type POSCAR-filer ved hjælp af umd2poscar.py-scriptet. snapshots af simuleringerne kan vælges med en foruddefineret frekvens. På kommandolinjetypen:
   umd2poscar.py -f -i -l -s
   hvor – l repræsenterer det sidste trin, der skal omdannes til POSCAR-fil. Standardværdierne er -i 0 -l 10000000 -s 1. Denne værdi af -l er stor nok til at dække en typisk hele bane.

2. Udfør strukturanalysen

1. Kør gofrs_umd.py scriptet for at beregne PDF-funktionen (Pair Distribution Function) _gᴀʙ(r) for alle par af atomtyperne A og B (figur 3). Outputtet er skrevet i én ASCII-fil, tabulatorsepareret, med filtypenavnet gofrs.dat. På kommandolinjetypen:
   gofrs_umd.py -f -s < Sampling_Frequency > -d -i
   BEMÆRK: Standarderne er Sampling_Frequency (hyppigheden for prøveudtagning af bane) = 1 trin; DiscretizationInterval (til plotning af g(r)) = 0,01 Å; InitialStep (antal trin i begyndelsen af den bane, der kasseres) = 0. Den radiale PDF, _gᴀʙ(r) er det gennemsnitlige antal atomer af type B i en afstand d_ᴀʙ inden for en sfærisk skal af radius r og tykkelse dr centreret om atomer af type A (figur 3):
   
   
   med ρ atomtætheden, _NA og _NB antallet af atomer af type A og B, og δ(r−_rᴀʙ) deltafunktionen, som er lig med 1, hvis atomer A og B ligger i en afstand mellem r og r +dr. Abscissa af det første maksimum af _gᴀʙ(r) giver den højeste sandsynlighedsbindingslængde mellem atomer af type A og B, som er tættest på en gennemsnitlig obligationsafstand, som vi kan bestemme. Det første minimum afgrænser omfanget af det første koordineringsområde. Derfor giver den integrerede del af PDF-filen op til det første minimum det gennemsnitlige koordineringsnummer. Summen af Fourier omdanner af _gᴀʙ(r) for alle par af atomare typer A og B giver diffraktion mønster af væsken, som opnås eksperimentelt med et diffraktometer. Men i virkeligheden, som ofte den høje orden koordinering sfærer mangler fra _gᴀʙ(r), diffraktion mønster kan ikke opnås i sin helhed.

Figur 3: Bestemmelse af parfordelingsfunktion.
a) For hvert atom af en art (f.eks. rødt) tælles alle atomer af de koordinerende arter (f.eks. grå og/eller rød) som en funktion af afstanden. (b) Den resulterende afstandsfordelingsgraf for hvert øjebliksbillede, som på dette stadium kun er en samling deltafunktioner, beregnes derefter i gennemsnit over alle atomer og alle snapshots og vægtes af den ideelle gasfordeling for at generere (c) den parfordelingsfunktion, der er kontinuerlig. Det første minimum af g(r) er radius for den første koordinationssfære, der anvendes senere i speciationsanalysen. Klik her for at se en større version af dette tal.

2. Uddrag de gennemsnitlige interatomare bindingsafstande som radier i de første koordinationssfærer. Til dette skal du identificere placeringen af det første maksimum af _gᴀʙ(r) funktioner: plot gofrs.dat-filen i et regnearksprogram og søg efter maxima og minima for hvert par atomer.
3. Identificer radius for den første koordinationssfære som det første minimum af PDF-filen , _gᴀʙ(r), ved hjælp af regnearkssoftware. Dette er grundlaget for hele den strukturelle analyse af væsken; PDF giver den gennemsnitlige limning status for atomer i væsken.
4. Uddrag afstandene til den første minima, dvs abscissa, og skriv dem i en separat fil, kaldet for eksempel bonds.input. Du kan også køre et af umd-pakkens analyze_gofr scripts for at identificere maxima- og minima-funktionerne i _gᴀʙ(r). På kommandolinjetypen:
   analyze_gofr_semi_automatic.py
5. Klik på placeringen af det maksimale og det mindste af _gᴀʙ(r), der vises i den graf, der åbnes af programmet. Scriptet scanner automatisk den aktuelle mappe, identificerer alle gofr.dat filer og udfører analysen for hver enkelt af dem. Klik igen på det maksimale og det mindste i vinduet, hver gang scriptet har brug for en uddannet første gæt.
6. Åbn og se på den automatisk genererede fil kaldet bonds.input , der indeholder de interatomare obligationsafstande.

3. Udfør speciationsanalysen

1. Beregne topologien af limning mellem atomer, ved hjælp af begrebet tilslutningsmuligheder inden for grafteori: atomer er knudepunkter og interatomare obligationer er stierne. Det speciation_umd.py script skal bruge de interatomare bindingsafstande, der er defineret i filen bonds.input .
   BEMÆRK: Konnektivitetsmatrixen er konstrueret på hvert trin: To atomer, der ligger i en afstand, der er mindre end radiusen af deres tilsvarende første koordinationssfære, anses for at være bundet, dvs. Forskellige atomare netværk er bygget ved at behandle atomer som knudepunkter i en graf, hvis forbindelser er defineret af dette geometriske kriterium. Disse netværk er atomarterne, og deres ensemble definerer atomspektationen i den pågældende væske (figur 4).

Figur 4: Identifikation af atomklyngerne.
Koordinationen polyhedra er defineret ved hjælp af interatomare afstande. Alle atomer i en afstand, der er mindre end en bestemt radius, anses for at være bundet. Her svarer tærsklen til den første koordinationssfære (de lyse røde cirkler), defineret i figur 1. Polymerisering og dermed kemiske arter er fremstillet af netværk af de bundne atomer. Bemærk den centrale _Red1Grey2 klynge, som er isoleret fra de andre atomer, som danner en uendelig polymer. Klik her for at se en større version af dette tal.

2. Kør speciation script for at opnå tilslutning matrix og opnå koordinering polyhedra eller polymerisering. På kommandolinjetypen:
   speciation_umd.py -f -s -i -l -c -a -m -r
   hvor -i-flaget giver filen med de interatomare bindingsafstande, som f.eks. blev produceret i det foregående trin. Du kan også køre scriptet med en enkelt længde for alle obligationer, der er defineret af -l-flaget.
   BEMÆRK: -c flag angiver de centrale atomer, og -en flag ligands. Både centrale atomer og ligands kan være af forskellige typer; i dette tilfælde skal de adskilles af kommaer. -m flaget giver den mindste tid, en art skal leve for at komme i betragtning i analysen. Som standard er denne minimumstid nul, og alle forekomster tælles med i sidste ende.
   1. Kør speciation_umd.py scriptet med flaget -r 0, som prøver tilslutningsgrafen på første niveau for at identificere koordinationspolyhedraen. For eksempel kan et centralt atom, der er angivet som en kation , være omgivet af en eller flere anioner (figur 4). Speciation script identificerer hver eneste af koordineringen polyhedra. Det vejede gennemsnit af alle koordinationspolyhedraen giver koordinationsnummeret, identisk med det, der opnås ved integrationen af PDF-filen. På kommandolinjetypen:
      speciation_umd.py -f -i -c -a -r 0
      BEMÆRK: Gennemsnitlige koordinationstal i væsker er fraktioneret antal. Denne fraktionerethed kommer fra den gennemsnitlige egenskab ved koordineringen. Definitionen baseret på speciation giver en mere intuitiv og informativ repræsentation af væskens struktur, hvor de relative proportioner af de forskellige arter, dvs.
   2. Kør speciation_umd.py script med flaget -r 1, som prøver tilslutningsgrafen på alle dybdeniveauer for at opnå polymeriseringen. Netværket gennem atomgrafen har en vis dybde, da atomer er bundet længere væk til andre bindinger (f.eks. i sekvenser af skiftende kationer og anioner) (figur 4).
3. Åbn de to filer . popul.dat og . stat.dat fortløbende; disse udgør output af speciation script. Hver klynge er skrevet på en linje og angiver dens kemiske formel, det tidspunkt, hvor den dannede, det tidspunkt, hvor den døde, dens levetid, en matrix med listen over de atomer, der danner denne klynge. Plot levetiden for hver atomklynge af alle de kemiske arter, der findes i simuleringen som findes i .popul.dat fil (Figur 5).
4. Plot populationen analyse med overflod af hver art, som findes i . stat.dat fil. Denne analyse, både absolut og relativ, svarer til de faktiske statistikker for koordinationspolyhedraen for sagen -r 0; for polymerisering, med -r 1 dette skal behandles omhyggeligt som nogle normalisering over det relative antal atomer kan være nødvendigt at anvende. Overfloden svarer til den integrerede i løbet af levetiden. Den. stat.dat fil viser også størrelsen af hver klynge, dvs. hvor mange atomer danner det.

4. Beregningsdiffusionskoefficienter

1. Uddrag de gennemsnitlige kvadratiske forskydninger (MSD) af atomer som en funktion af tiden til at opnå selv-diffusivitet. Standardformlen for MSD er:
   
   hvor præfaktorerne er renormalisationer. Med MSD-værktøjet er der forskellige måder at analysere de dynamiske aspekter af væskerne på.
   BEMÆRK: T er den samlede tid for simuleringen og N _α er antallet af atomer af typen α. Den første gang _t0 er vilkårlig og spænder over den første halvdel af simuleringen. _Ninit er antallet af første gange. τ er bredden af det tidsinterval, som MSD beregnes over. dens maksimale værdi er halvdelen af simuleringens tidstid. I typiske MSD-implementeringer starter hvert vindue i slutningen af det foregående. Men en sparsom prøveudtagning kan fremskynde beregningen af MSD, uden at ændre hældningen af MSD. Til dette starter det i-th-vindue på tidspunktet _t0(i), men (i+1)-th-vinduet starter på tidspunktet _t0(i) + τ + v, hvor værdien af v er brugerdefineret. På samme måde øges vinduets bredde i diskrete trin, der er defineret af brugeren, som sådan: τ(i) = τ(i-1) + z. Værdierne z ("vandret trin") og v ("lodret trin") er positive eller nul; Standarden for begge er 20.
2. Beregn MSD ved hjælp af serien af msd_umd scripts. Deres output udskrives i en . msd.dat fil, hvor MSD af hver atomtype, atom eller klynge udskrives på en kolonne som en funktion af tiden.
   1. Beregne den gennemsnitlige MSD for hver atomtype. MSD beregnes for hvert atom og derefter gennemsnit for hver atomtype. Outputfilen indeholder én kolonne for hver atomtype. På kommandolinjetypen:
      msd_umd.py -f -z -v -b
   2. Beregne MSD af hvert atom. MSD beregnes for hvert atom og derefter gennemsnit for hver atomtype. Outputfilen indeholder en kolonne for hvert atom i simuleringen og derefter en kolonne for hver atomtype. Denne funktion giver mulighed for at identificere atomer, der diffus i to forskellige miljøer, som væske og gas, eller to væsker. På kommandolinjetypen:
      msd_all_umd.py -f -z -v -b
   3. Beregne den kemiske arts msd. Brug populationen af klynger, der er identificeret med speciationsscriptet, og som udskrives i . popul.dat fil. Msd beregnes for hver enkelt klynge. Outputfilen indeholder én kolonne for hver klynge. For at undgå at overveje store polymerer skal du sætte en grænse for klyngens størrelse; dens standard er 20 atomer. På kommandolinjetypen:
      msd_cluster_umd.py -f -p -s -b -c
      BEMÆRK: Standardværdierne er: –b 100 –s 1 –c 20.
3. Plot msd ved hjælp af et regneark-baseret software (Figur 6). I en loglogrepræsentation af MSD i forhold til tid skal du identificere hældningsændringen. Adskil den første del, normalt kort, som repræsenterer det ballistiske regime, dvs. bevarelsen af atomernes hastighed efter kollisioner. Den anden længere del repræsenterer diffusivt regime, dvs. spredning af atomers hastighed efter kollisioner.
4. Beregne diffusionskoefficienterne fra msd-hældningen som:
   
   hvor Z er antallet af frihedsgrader (Z = 2 for diffusion i plan, Z = 3 for diffusion i rummet), og t er tidstrinnet.

5. Tidskorrelationsfunktioner

1. Beregning af tidskorrelationsfunktionerne som et mål for systemets inerti ved hjælp af den generelle formel:
   
   A kan være en række tidsafhængige variabler, såsom atomare positioner, atomare hastigheder, belastninger, polarisering osv., der hver giver - via Green-Kubo-forbindelserne12,13 - forskellige fysiske egenskaber, nogle gange efter en yderligere transformation.
2. Analysere de atomare hastigheder for at opnå det vibrationelle spektrum af væsken og alternative udtryk for de atomare selvdiffusionskoefficienter.
   1. Kør vibr_spectrum_umd.py script til at beregne atomare hastighed-hastighed auto-korrelation (VAC) funktion for hver atomare type og udføre sin hurtige Fourier transformering. På kommandolinjetypen:
      vibr_spectrum_umd.py -f -t
      hvor –t er den temperatur, der skal defineres af brugeren. Scriptet udskriver to filer: . vels.scf.dat fil med VAC-funktionen for hver atomtype og . vibr.dat fil med vibrationelle spektrum nedbrudt på hver atomart og den samlede værdi.
   2. Åbn og læs vels.scf.dat. Plot VAC-funktionen ud fra vels.scf.dat-filen ved hjælp af regnearkslignende software.
   3. Behold den rigtige del af Fourier VAC. Dette er, hvad der giver vibrationelle spektrum, som en funktion af frekvens:
      
      hvor m er atommasserne.
   4. Plot vibrationelle spektrum fra vibr.dat fil ved hjælp af regneark-lignende software (Figur 7). Identificer den finite værdi ved ω=0, der svarer til væskens diffusive karakter og de forskellige toppe af spektret ved begrænset frekvens. Identificer hver atomtypes deltagelse i vibrationsspektret.
      BEMÆRK: Nedbrydningen på atomtyper viser, at forskellige atomer har forskellige ω= 0 bidrag, svarende til deres diffusionskoefficienter. Den generelle form af spektret er meget glattere med færre funktioner end for en tilsvarende solid.
   5. Ved skallen skal du læse den integrerede over det vibrationelle spektrum, som giver diffusionskoefficienterne for hver atomart.
      BEMÆRK: Termodynamiske egenskaber kan opnås ved integration fra vibrationsspektret, men resultaterne skal anvendes med forsigtighed på grund af to tilnærmelser: integrationen er gyldig inden for den kvasi-harmoniske tilnærmelse, som ikke nødvendigvis holder ved høje temperaturer; og den gaslignende del af spektret, der svarer til diffusionen, skal kasseres. Integrationen bør så kun ske over den gitterlignende del af spektret. Men denne adskillelse kræver normalt flere yderligere efterbehandlingstrin og ^{beregninger14}, som ikke er omfattet af den nuværende UMD-pakke.
3. Kør viscosity_umd.py scriptet for at analysere selvkorrelationen af komponenterne stress tensor at estimere viskositeten af smelten. På kommandolinjetypen:
   viscosity_umd.py -f -i -s -o -l
   BEMÆRK: Denne funktion er sonderende, og eventuelle resultater skal tages med forsigtighed. For det første skal viskositetens konvergens grundigt kontrolleres med hensyn til simuleringens længde.
   1. Udlede viskositeten af væsken fra selvkorrelationen af stress ^tensor15 som:
      
      hvor V og T er henholdsvis volumen og temperatur, _κB er Boltzmann konstant og σ _ij ij off-diagonal komponent af stress-tensor, udtrykt i kartesiske koordinater.
   2. Brug en mere passende pasform til at opnå et mere robust skøn over viskositeten15,16 og undgå støj fra den stress-tensor autokorrelationsfunktion, der kan opstå som følge af simuleringernes begrænsede størrelse og begrænsede varighed. For den automatiske korrelation funktion af stress tensor, skal du bruge følgende funktionelle ^form15,16, der giver gode resultater:
      
      hvor A, B, _τ1, τ2 og ω er fit parametre. Efter integration bliver udtrykket for viskositeten:
      

6. Termodynamiske parametre, der stammer fra simuleringerne.

1. Kør averages.py for at udtrække gennemsnitsværdierne og spredningen (som standardafvigelse) for tryk, temperatur, tæthed og intern energi fra umd-filerne . På kommandolinjetypen:
   averages.py -f -s
   med –s 0 som standard.
2. Beregn den statistiske fejl i gennemsnittet ved hjælp af blokeringsmetoder.
   BEMÆRK: Der er forskellige varianter af denne metode. Efter arbejdet i Allen og Tildesley2, er det almindeligt at gennemsnittet over sekvenser af tidsblokke, af stadig længere længde, og anslå standardafvigelsen med hensyn til det aritmetiske ^gennemsnit17. Konvergens kan opnås inden for grænsen af mange- og længe nok blokstørrelser, når prøveudtagningen ikke er korreleret. Selv om den faktiske tærskelværdi for konvergensen normalt skal vælges manuelt.
   1. Brug halveringsmetoden18: Start med den første dataprøve, ved hvert trin κ, halvere antallet af prøverne ved i gennemsnit over hver to tilsvarende på hinanden følgende prøver fra det foregående trin κ−1:
      
   2. Kør fullaverages.py scriptet for at udføre den komplette statistiske analyse, herunder middelværdiens fejl. På kommandolinjetypen:
      fullaverages.py
      BEMÆRK: Scriptet automatiseres til det punkt, hvor der søges efter alle .umd.dat filer i den aktuelle mappe og udfører analysen for dem alle. Standardindstillingerne er –s 0 –u 0. For -u 0 output er minimal, og for -u 1 output er i fuld, med flere alternative enheder udskrevet. Dette script kræver grafisk understøttelse, da det skaber et grafisk billede til kontrol af konvergensen til estimering af fejlen på middelværdien.

Representative Results

Pyrolite er en model multi-komponent silikat smelte (0,5Na2O 2CaO 1,5Al2O3 4FeO 30MgO 24SiO2), der bedst tilnærmer sammensætningen af bulk silikat Jorden-den geokemiske gennemsnit eller hele vores planet, bortset fra sin jern-baserede ^kerne19. Den tidlige Jord var domineret af en række store ^{smeltehændelser20}, den sidste kunne have opslugt hele planeten, efter dens kondens til protolunar ^disk21. Pyrolite repræsenterer den bedste caimant til sammensætningen af sådanne planetariske skala magma oceaner. Derfor har vi grundigt undersøgt pyrolitens fysiske egenskaber i 3.000\u20125.000 K temperaturområdet og 0\u2012150 GPa-trykområdet fra ab initio molekylær-dynamics simuleringer i VASP-implementeringen. Disse termodynamiske forhold karakteriserer helt Jordens mest ekstreme magmahavsforhold. Vores undersøgelse er et glimrende eksempel på en vellykket brug af UMD-pakken til hele den dybdegående analyse af ^smelter22. Vi beregnede fordelingen og de gennemsnitlige obligationslængder, vi sporede ændringerne i kation-ilt koordinering, og sammenlignede vores resultater med tidligere eksperimentelle og beregningsmæssige undersøgelser af amorfe silikater af forskellige kompositioner. Vores dybdegående analyse hjalp nedbryde standard koordinering numre i deres grundlæggende bestanddele, skitsere tilstedeværelsen af eksotiske koordinering polyhedra i smelten, og udtrække levetid for alle koordinering polyhedra. Det skitserede også betydningen af prøveudtagning i simuleringer med hensyn til både længden af banen og også antallet af atomer, der er til stede i det modellerede system. Hvad angår efterbehandling, er UMD-analysen uafhængig af disse faktorer, men de bør tages i betragtning ved fortolkningen af resultaterne fra UMD-pakken. Her viser vi et par eksempler på, hvordan UMD-pakken kan bruges til at udtrække flere karakteristiske træk ved smelterne, med en applikation til smeltet pyrolit.

Si-O-parfordelingsfunktionen fra gofrs_umd.py scriptet viser, at radiusen af den første koordinationssfære, som er det første minimum af g(r)-funktionen, ligger omkring 2,5 angstromer ved T = 3000 K og P = 4,6 GPa. Det maksimale af g (r) ligger ved 1.635 Å-dette er den bedste tilnærmelse til bøjning længde. Den lange hale skyldes temperaturen. Ved hjælp af denne grænse som Si-O-bindingsafstand viser speciationsanalysen, at _SiO4-enheder, som kan vare i op til et par picosekunder, dominerer smelten (figur 5). Der er en vigtig del af smelten, der viser delvis polymerisering, som afspejles af tilstedeværelsen af dimers som _Si2O7, og trimere som _Si3Ox enheder. Deres tilsvarende levetid er i rækkefølgen af picosecond. Polymerer med højere orden har alle betydeligt kortere levetid.

Figur 5: Si-O-kemiske arters levetid.
Speciationen er identificeret i en multikomponentsmeltning ved 4,6 GPa og 3000 K. Etiketterne markerer _SiO3, _SiO4 og _SiO5 monomerer og de forskellige _SixOy polymerer. Klik her for at se en større version af dette tal.

De forskellige værdier for de lodrette og vandrette trin, defineret af –z og -v-flagene ovenfor, giver forskellige udsnit af MSD (figur 6). Selv store værdier af z og v er nok til at definere skråningerne og dermed diffusionskoefficienterne for de forskellige atomer. Gevinsten i tide til efterbehandlingen er bemærkelsesværdig, når man går til store værdier af z og v. Msd tilbyder et meget stærkt valideringskriterium for kvaliteten af simuleringerne. Hvis diffusionsdelen af MSD ikke er tilstrækkelig lang, er det et tegn på, at simuleringen er for kort og ikke når væsketilstanden i statistisk forstand. Minimumskravet for den diffusive del af muskel- og skeletbesvær afhænger i høj grad af systemet. Man kan kræve, at alle atomer ændre deres websted mindst én gang i strukturen af smelten for at det skal betragtes som en ^fluid10. Et glimrende eksempel med anvendelser i planetarisk videnskab er komplekse silikat smelter ved højt tryk tæt på eller endda under deres liquidus ^linje11. Si-atomer, de store netværksdannende kationer, skifter websteder efter mere end to dusin picosekunder. Simuleringer, der er kortere end denne tærskel, vil være betydeligt underudtagning af det mulige konfigurationsareal. Men da de koordinerende anioner, nemlig O-atomer, bevæger sig hurtigere end de centrale Si-atomer, kan de kompensere for en del af Si's langsomme mobilitet. Som sådan kunne hele systemet faktisk dække en bedre prøveudtagning af konfigurationsrummet end antaget kun fra Si-forskydningerne.

Figur 6: Middelkvalert forskydninger (MSD).
MSD er illustreret for et par atomare typer af en multi-komponent silikat smelte. Prøveudtagningen med forskellige vandrette og lodrette trin, z og v, giver ensartede resultater. Faste cirkler: -z 50 –v 50. Åbne cirkler: -z 250 –v 500. Klik her for at se en større version af dette tal.

Endelig giver de atomare VAC-funktioner smeltenes vibrationelle spektrum. Figur 7 viser spektret ved samme tryk- og temperaturforhold som ovenfor. Vi repræsenterer bidragene fra Mg, Si og O atomer, samt den samlede værdi. Ved nul frekvens er der en endelig værdi af spektret, hvilket svarer til smeltens diffusionskarakter. Udvinding af de termodynamiske egenskaber fra vibrationsspektret skal fjerne denne gaslignende diffusive karakter fra nul, men også for korrekt at tage hensyn til dets forfald ved højere frekvenser.

Figur 7: Pyrolitens vibrationelle spektrum smelter.
Den virkelige del af Fourier omdanne af atomare hastighed-hastighed selvkorrelation funktion giver vibrationelle spektrum. Her beregnes spektret for en multikomponents silikatsmeltning. Væsker har en ikke-nul gas-lignende diffusiv karakter ved nul frekvens. Klik her for at se en større version af dette tal.

Discussion

UMD-pakken er designet til at fungere bedre med ab initio simuleringer, hvor antallet af snapshots typisk er begrænset til titusinder til hundredtusindvis af snapshots, med et par hundrede atomer pr. enhedscelle. Større simuleringer kan også kan leveres, forudsat at den maskine, hvor efterbehandlingskørslerne har nok aktive hukommelsesressourcer. Koden adskiller sig ved de mange forskellige egenskaber, den kan beregne, og ved sin open source-licens.

Umd.dat filer er passende til ensembler, der bevarer antallet af partikler uændret i hele simuleringen. UMD-pakken kan læse filer, der stammer fra beregninger, hvor simuleringsboksens form og volumen varierer. Disse dækker de mest almindelige beregninger, som NVT og NPT, hvor antallet af partikler, N, temperatur T, volumen, V og / eller tryk, P, holdes konstant.

For den tid begynder parret distribution funktion samt alle de scripts, der har brug for at estimere de interatomare afstande, ligesom speciation scripts, arbejde kun for ortogonale enhed celler, hvilket betyder for kubiske, tetragonale og ortopædiske celler, hvor vinklerne mellem akserne er 90 °.

De vigtigste udviklingslinjer for version 2.0 er fjernelse af ortogonalitetsbegrænsningen for afstande og tilføjelse af flere funktioner til speciationscripts: at analysere individuelle kemiske bindinger, at analysere de interatomare vinkler og at implementere den anden koordineringssfære. Med hjælp fra eksternt samarbejde arbejder vi på at portere koden til en GPU for hurtigere analyse i større systemer.

Materials

Name	Company	Catalog Number	Comments
getopt library			open-source
glob library			open-source
matplotlib library			open-source
numpy library			open-source
os library			open-source
Python software	The Python Software Foundation	Version 2 and 3	open-source
random library			open-source
re library			open-source
scipy library			open-source
subprocess library			open-source
sys library			open-source