Moleküldynamiksimulation: fluktuierende Dipolmodellimplementierung

Ich führe eine Molekulardynamiksimulation für Siliciumdioxid durch. Vor einiger Zeit habe ich mich dem fluktuierenden Dipolmodell zugewandt und nach vielen Bemühungen habe ich immer noch Probleme damit, es umzusetzen.

Kurz gesagt sind alle Sauerstoffatome im System polarisierbar und ihre Dipolmomente hängen von ihrer Position in Bezug auf alle anderen Atome im System ab. Insbesondere verwende ich TS-Potential ( Ссылка ), wo Dipole sind zu jedem Zeitschritt iterativ gefunden.

Dies bedeutet, dass ich bei der Bewertung von Kräften, die auf Atome einwirken, diese potentielle Energieabhängigkeit von Koordinaten berücksichtigen muss. Vorher verwendete ich einfache paarweise potentielle Modelle, und so würde ich mein Programm auf die Berechnung von Kräften unter Verwendung von analytischen Formeln einstellen, die durch Differenzieren des potentiellen Energieausdrucks erhalten wurden.

Jetzt bin ich ratlos: Wie setzt man dieses neue Potenzial um? In allen Artikeln, die ich gefunden habe, geben sie nur Formeln, aber nicht den Algorithmus. Wie ich es sehe, muss ich bei der Berechnung von Kräften, die auf ein bestimmtes Atom wirken, die Änderung des Dipols dieses Atoms, die Änderung der Dipole aller benachbarten Atome und dann die Änderung der Dipole von noch mehr Atomen berücksichtigen und so weiter, so wie sie voneinander abhängen. Schließlich liegt es an dieser Interdependenz, dass die Dipole bei jedem Zeitschritt iterativ gefunden werden. Natürlich kann ich die Kräfte nicht iterativ für jedes Atom berechnen, weil die Rechenkomplexität des Algorithmus viel zu hoch wäre. Sollte ich einige einfache Funktionen verwenden, um den Wechsel der Dipole zu erklären? Das sieht auch nicht nach einer guten Idee aus, denn Dipole werden iterativ mit hoher Genauigkeit berechnet, und dann, wo es wirklich darauf ankommt (Rechenkräfte), würden wir grobe Funktionen verwenden?

Wie implementiere ich dieses Modell? Ist es auch möglich, Kräfte analytisch zu berechnen, wie ich es vorher getan habe, oder ist es notwendig, sie unter Verwendung der endlichen Differenzenformel für die Ableitung zu berechnen?

Ich habe die Antwort auf meine Frage in der Literatur nicht gefunden, aber wenn Sie einen Artikel, eine Website oder ein Buch kennen, wo dieses Material hervorgehoben ist, leiten Sie mich bitte zu dieser Quelle.

Danke für Ihre Zeit!

================================================== ====================================

UPDATE:

Danke für Ihre Antwort. Leider war das nicht meine Frage. Ich habe nicht gefragt, wie man Dipole berechnet, sondern wie man Kräfte berechnet, vorausgesetzt, dass diese Dipole mit der Bewegung beträchtlich variieren.

Ich habe versucht, die Kräfte auf einfache Art und Weise zu berechnen (ohne Dipole Interdependenz über ihre Abstände zu berücksichtigen, berechne einfach die Dipole auf jeder Stufe, dann berechne ich die Kräfte, als wären diese Dipole statisch), aber die Ergebnisse waren physikalisch nicht richtig.

Um die Situation zu analysieren, habe ich eine Simulation eines Systems bestehend aus nur zwei Atomen: Si und O. Sie haben entgegengesetzte Ladungen, und so oszillieren sie. Und die Energiezeitabhängigkeitsgrafik sieht so aus:

Die obere Kurve stellt die kinetische Energie dar, die mittlere stellt die potentielle Energie dar, ohne die Dipol-Wechselwirkung zu berücksichtigen, und die untere stellt die potentielle Energie des Systems dar, wobei die Dipol-Wechselwirkung berücksichtigt wurde. p>

Sie können aus dieser Grafik deutlich erkennen, dass das System tut, was es nicht tun sollte: die mögliche Steigung hinaufzusteigen. Also entschied ich, dass dies auf die Tatsache zurückzuführen ist, dass ich die Koordinationsabhängigkeit des Dipolmoments nicht berücksichtigt habe. Zum Beispiel berechnen wir zu einem bestimmten Zeitpunkt die Kräfte, und sie sind so gerichtet, dass sie beide Atome aufeinander zu bewegen. Aber wenn wir sie aufeinander zu bewegen (auch nur geringfügig), ändert sich das Dipolmoment, und wir finden heraus, dass wir tatsächlich eine höhere potentielle Energie als zuvor hatten! Im nächsten Zeitschritt ist die Situation gleich.

Es stellt sich also die Frage, wie dieser Effekt berücksichtigt werden kann, denn die wenigen Möglichkeiten, die ich mir vorstellen kann, sind entweder zu rechenintensiv oder viel zu grob.