 |
 |
 |
Als jonge wetenschapper hield ik me bezig met computationele statistische thermodynamica. Dat is nogal een mondvol, dus laat me het even uitleggen. Enkele jaren geleden schreef ik een inleiding tot de klassieke thermodynamica waarin ik uitlegde dat deze theorie over bulkmaterie zich niet bezighoudt met de structuur van materie. Het werkt prima voor stoommachines of onderzeeërs, maar volstaat niet voor het beschrijven van processen op zeer kleine schaal. Daarvoor hebben we de statistische thermodynamica. Briljante geesten als Maxwell, Boltzmann en Liouville hebben het raamwerk voor deze theorie ontwikkeld. Kort samengevat verklaart de statistische thermodynamica de waarschijnlijkheid dat minuscule deeltjes zich in een bepaalde configuratie bevinden, en dat de waarde van een bulkeigenschap het gemiddelde is van die eigenschap in alle haalbare configuraties, gewogen naar hun waarschijnlijkheid.
Het enige wat we zouden hoeven te doen, is exacte vergelijkingen afleiden om alles te berekenen wat er te weten valt over gasvormige, vloeibare of vaste materie. Helaas is het in de meeste gevallen wiskundig onmogelijk om deze vergelijkingen op papier exact op te lossen. In computationele statistische thermodynamica berekenen we numeriek de resultaten van deze onoplosbare vergelijkingen. Het idee is om een (groot) aantal deeltjes te observeren, hun beweging een tijdje te volgen en gaandeweg de eigenschappen te middelen. Twee veelgebruikte methoden hiervoor zijn moleculaire dynamica of Monte Carlo simulaties. Soms worden deze methoden ook wel “computationele microscopen” genoemd. In plaats van een wiskundige formule krijg je een tabel met getallen. Theoretici gebruiken de tabellen om hun wiskundige modellen te ontwikkelen, en ingenieurs gebruiken ze in hun werk.
Monte Carlo was de eerste techniek die ik gebruikte voor mijn afstudeeronderzoek. Simpel gezegd is het idee achter Monte Carlo-simulaties in de statistische thermodynamica als volgt. Je plaatst de deeltjes in een plausibele beginconfiguratie. Vervolgens pas je deze aan, bijvoorbeeld door een deeltje een stukje te verplaatsen, twee deeltjes van plaats te verwisselen of het volume van het hele systeem te wijzigen. Dat levert een nieuwe configuratie op met een eigen waarschijnlijkheid van bestaan. Nu werp je een dobbelsteen om te beslissen of je de nieuwe configuratie als de nieuwe accepteert of terugkeert naar de oude, afhankelijk van de configuratiewaarschijnlijkheden. Herhaal deze aanpassingen vele malen, en je genereert uiteindelijk de meest waarschijnlijke configuraties voor je set deeltjes. Bereken het gemiddelde van wat je wilt weten in deze configuraties, en je komt uit op een goede benadering van de theoretische resultaten.
Mijn afstudeerbegeleider aan de Universiteit van Nijmegen, een scherpzinnige theoretisch natuurkundige, bood me een afstudeeronderwerp aan waarbij ik een geheel nieuw Monte Carlo-programma zou schrijven en daarmee de thermodynamica van vast-vloeistofgrensvlakken zou bestuderen. De resultaten van dergelijke onderzoeken hebben bijvoorbeeld betrekking op industriële kristalgroei, wat het onderwerp technisch relevant maakt. De opdracht klonk als muziek in mijn oren - ik zou intensief gebruikmaken van het universitair rekencentrum en als chemicus zou ik samenwerken met theoretisch natuurkundigen. Het was echter een opmerkelijk dat op de afdeling scheikunde in Nijmegen alleen werk op het gebied van moleculaire quantummechanica als theoretische chemie werd beschouwd.
Ik heb al die tijd een kopie bewaard van de broncode van het bovengenoemde programma, dat een uitstekende testcase bleek te zijn voor mijn eigen Fortran-compiler. Het bleek dat het Monte Carlo-programma wederom uitstekend werkte. Zoals verwacht zijn de prestaties aanzienlijk verbeterd in vergelijking met het mainframe waarop het in de jaren tachtig draaide. In een uur kan ik nu doen waar ik indertijd weken over deed, niet alleen vanwege de verbeterde verwerkingssnelheid, maar ook omdat ik geen capaciteit meer hoef te delen. Destijds werden zware rekentaken buiten kantooruren in batch uitgevoerd, en moest ik geduldig in de wachtrij staan met mijn collega’s uit de theoretische chemie, kristallografie, hoge-energiefysica, astrofysica, enzovoort.
Natuurlijk kon ik het niet laten om, puur voor de lol, nog eens terug te kijken naar een deel van het onderzoek dat ik destijds heb gedaan, nu met behulp van moderne ray-tracing-software om prachtige afbeeldingen van de deeltjesconfiguraties te maken. Hier wil ik een paar observaties laten zien van een mengsel van xenon en neon bij een zeer lage temperatuur. Edelgassen zijn populaire modelsystemen in simulaties, omdat de interacties tussen atomen vrij goed begrepen worden en relatief eenvoudig zijn. Heel kort gezegd is dit hoe twee edelgasatomen op elkaar inwerken - ze stoten elkaar af op korte afstand, trekken elkaar aan op middellange afstand, maar negeren elkaar op grote afstand, terwijl een derde atoom in de buurt geen invloed heeft op de interactie tussen de eerste twee.
Ik heb de simulaties uitgevoerd alsof we nog in de jaren tachtig leefden en ben uitgegaan van een aantal atomen dat destijds net haalbaar zou zijn geweest, als ik toegang had gehad tot de supercomputers uit die tijd. Daarom gaan we uit van iets meer dan 10.000 xenonatomen en hetzelfde aantal neonatomen, in een rechthoekige doos. Op het mainframe van onze universiteit lag het haalbare aantal atomen ongeveer een orde van grootte lager. Om een bulksysteem te modelleren, wordt de rechthoekige doos met de atomen periodiek in alle drie de dimensies voortgezet. Dat is een benadering van een oneindig groot volume, hoewel periodieke randvoorwaarden leiden tot eindig-systeemeffecten; elke goede simulatiestudie houdt daar rekening mee.
In de praktijk worden de meeste laboratoriumexperimenten bij constante druk uitgevoerd. We zullen in onze simulaties iets vergelijkbaars doen en simulaties bij verschillende temperaturen, maar bij dezelfde druk, met elkaar vergelijken. Je zult begrijpen dat de interacties tussen de deeltjes bijdragen aan de druk – het duwen en trekken tussen atomen veroorzaakt fysieke spanning. Een mogelijke aanpak is om de deeltjes "spanningsvrij" te maken, dat wil zeggen het volume zodanig aan te passen dat de bijdragen aan de druk door deeltjesinteracties elkaar gemiddeld opheffen. Dit betekent dat we de druk in de simulaties afstemmen op de ideale-gas druk.
Onze eerste simulatie vindt plaats bij een temperatuur van 25K, het smeltpunt van neon. Hieronder zie je een momentopname van de simulatie: de witte bolletjes zijn xenonatomen en de oranje bolletjes zijn neonatomen. Gesmolten neon is voor de duidelijkheid weergegeven als grote druppels. Er zijn nog geen xenonatomen opgelost in de neonvloeistof. Xenon smelt bij 161K, en links zien we dat xenon inderdaad bevroren is met een vlak kristaloppervlak. Rechts zien we dezelfde configuratie, waarbij de neonatomen het kristaloppervlak omringen. Merk op dat deze atomen aan het xenonkristal lijken vast te zitten – een vloeistof dicht bij een grensvlak blijkt structuur aan te nemen. Om zich als een "vloeistof" te gedragen, moet een deeltje zich vrij in elke richting kunnen bewegen. In de buurt van een vaste wand is die vrijheid beperkt.
Hieronder zie je een momentopname van de configuratie bij 120K. Nu zien we inderdaad dat enkele xenonatomen zijn opgelost in de neonvloeistof. We zien ook dat het kristaloppervlak niet langer vlak is. Een meer diepgaande analyse van deze simulaties toont aan dat de buitenste lagen van een kristal smelten bij een lagere temperatuur dan het smeltpunt van het bulk-kristal. Destijds was dit een interessante ontdekking, aangezien toevallig tegelijkertijd "oppervlaktesmelten" experimenteel werd aangetoond op echte loodkristallen bij AMOLF. Ons model is minder complex dan een loodoppervlak, maar toch hebben we hetzelfde fenomeen in onze simulaties waargenomen.
Na mijn afstuderen kreeg ik een functie aangeboden bij de faculteit Technische Natuurkunde van de Universiteit Twente, waar ik opnieuw samenwerkte met theoretische natuurkundigen die ook altijd op zoek waren naar meer rekenkracht. Een slimme persoon bouwde voort op mijn simulatiewerk in Nijmegen, en uiteindelijk werd er een artikel over gepubliceerd. Na mijn academische loopbaan raakte ik steeds verder verwijderd van intensieve berekeningen, om daar de afgelopen jaren weer naar terug te keren dankzij de opkomst van industriële toepassingen van kunstmatige intelligentie. Het opnieuw leven inblazen van veertig jaar oude broncode bracht echter veel goede herinneringen terug aan mijn eerste substantiële programmeer- en onderzoekswerk.
¹ J.P. van der Eerden, A. Roos, J.M. van der Veer. Surface roughening versus surface melting on Lennard-Jones crystal surfaces. Journal of Crystal Growth 99(1–4) 77-82 [1990].
