Using finite element codes as a numerical platform to run molecular dynamics simulations

F. Niederhöfer, J. Wackerfuß

A mathematically rigorous methodology for embedding the governing equations of molecular dynamics in the formalism of the finite element method is presented. Only one generalized finite element type is needed to cover all different types of existing interatomic potentials. The finite element type is simply specified by two parameters characterizing the type of the interatomic potential to be considered. Built on this formulation a partitioned Runge–Kutta method—summarizing a wide range of explicit and implicit, single- and multi-stage, lower and higher order time integration schemes—is embedded in a unified manner. The required finite element residual vector and the related Jacobian matrix are stated explicitly. The related FE-mesh coincides with the neighborhood lists used in standard molecular dynamics enabling the use of common tools. The range, versatility and performance of the proposed finite element formulation have been demonstrated by means of several numerical examples.

Heutzutage stellen Molekulardynamik-Simulationen in vielen Bereichen der Materialwissenschaften sowie der Bio- und Nanotechnologie ein wichtiges Werkzeug zur Charakterisierung bekannter und zur Entwicklung neuer Molekülstrukturen dar. Die Effizienz von derartigen numerischen Simulationen hängt u.a. von der Wahl des zugrundeliegenden Zeitintegrationsverfahrens ab. Bei dieser Wahl sind unterschiedliche Faktoren zu berücksichtigen (Gesamtanzahl der Atome, Anzahl der Zeitschritte, geforderte Genauigkeit, Speicherressourcen u.v.m.). Auf der anderen Seite gilt es Einflüsse zu berücksichtigen, die von der Wahl des Zeitintegrationsverfahrens unabhängig sind (z.B. Rechenungenauigkeiten).

Ein Aspekt dabei ist die Empfindlichkeit der numerischen Ergebnisse bezüglich minimaler Änderungen der Anfangswerte. Die folgende Animation soll diesen Sachverhalt exemplarisch illustrieren.

Dargestellt sind die Ergebnisse von 3 unterschiedlichen Molekulardynamik-Simulationen für einen um 10 % gegenüber der Gleichgewichtslage „aufgepumpten“ Buckyball. Zur Zeitintegration wurde das Verfahren von Störmer-Verlet mit einer Zeitschrittweite von Δt = 0.025 Femtosekunden (bewusst klein!) und insgesamt 1000 Zeitschritten gewählt. Die 3 Simulationen unterscheiden sich dabei lediglich durch die Genauigkeit der Anfangsdaten: In der linken Simulation sind diese auf 4, in der mittleren auf 8 und in der rechten auf 12 Nachkommastellen exakt angegeben. Der optische Vergleich zeigt, dass der Zeitpunkt des Übergangs von einer rein radialen in eine „chaotische“ Bewegung von der Genauigkeit der Anfangsdaten abhängt. Dieses Verhalten ist unabhängig von dem gewählten Zeitschrittverfahren.

In einem weiteren Simulationsbeispiel wird ein (8,0) carbon nanotube (CNT) der Länge ca. 4 Nanometer in axiale Schwingung versetzt. Das CNT, das am rechten Rand unterschiedlich gelagert ist, wird am linken Rand durch eine vorgegebene (axiale) Auslenkung angeregt. Mit einer Zeitschrittweite von Δt = 0.25 Femtosekunden wurden insgesamt 10000 Zeitschritte gerechnet. Die Videosequenz zeigt deutlich eine von links nach rechts durch das CNT propagierende Longitudinalwelle, die am rechten Rand reflektiert wird. Diese Reflektion leitet einen sukzessiven Übergang zu einer unbestimmten Gitterschwingung ein.

Publikation

Niederhöfer, F & Wackerfuß, J 2012, 'High-order time integration methods in molecular dynamics', PAMM Proceedings in Applied Mathematics and Mechanics, Volume 12, S. 47-48.

Wackerfuß, J., Niederhöfer, F., 2018. Using finite element codes as a numerical platform to run molecular dynamics simulations. Computational Mechanics 18, 1–30, https://doi.org/10.1007/s00466-018-1594-5.