Wie simuliert man einen Phasenübergang in der Molekulardynamik?

Ich versuche, den Phasenübergang für Gold und Silizium getrennt mit LAMMPS zu simulieren. Ich habe den Schmelzpunkt für Gold richtig mit dem folgenden Code.

units   metal
atom_style      atomic
boundary        p p p
variable        a equal 4.0782

lattice fcc 4.0782
region  box block 0 10 0 10 0 10
create_box      1 box
create_atoms    1 box
mass    1 196.97

pair_style      eam
pair_coeff      * * Au_u3.eam

minimize        1.0e-8 1.0e-8 1000 100000
min_style       cg

timestep 0.001
velocity all create 300.0 454883 mom yes rot yes dist gaussian

thermo  50000
thermo_style    custom step pe ke etotal temp vol press density atoms

fix 1 all press/berendsen iso 0.0 0.0 100.0
fix 2 all nvt temp 300.00 2400.00 1.0


run     10000000

Wenn ich diesen Code jedoch auf Silizium anwende, erhalte ich nicht die richtigen Ergebnisse. Ich habe das Gefühl, dass ich den Code und die Physik hinter dem Erhalten eines Phasenübergangsdiagramms mit Molekulardynamik nicht verstehe. Ich schätze, meine Bitte hier ist, ob Sie Empfehlungen für einige Artikel oder Bücher haben, die ich lesen sollte, um die Simulation durchführen zu können (Phasenübergang in der Molekulardynamik). Ich habe das Buch "Molecular Dynamics verstehen von Dan Frenkel und Berend Smit" durchgeblättert, aber ich habe immer noch das Gefühl, dass mir etwas fehlt.

bearbeitet: Die Ausgabe ist unten für Gold und Silizium dargestellt.

Gold: Temperatur vs. Energie

siliocn: temparture vd energy

Verschiedene Randbedingungen ausprobiert?
Zeigen Sie uns, was Sie mit Gold und was mit Silizium erreichen. Geben Sie uns weitere Einzelheiten darüber, warum Ihrer Meinung nach Ihre Simulation fehlgeschlagen ist.
(Nebenbei: Ihre nvt-Zeitkonstante beträgt ca 10 × zu groß)
@ valerio92 Die Ausgabe ist oben gezeigt, wo der Phasenübergang in Gold bei etwa 1300 k liegt, was korrekt ist, und ich habe die latente Wärme berechnet und war auch korrekt. Für Silizium war der Schmelzpunkt jedoch nicht korrekt, da er bei etwa 2000 k lag und der korrekte Wert bei etwa 1600 k liegt.

Antworten (1)

Metalle wie Gold unterliegen Phasenübergängen mit einer sehr begrenzten Keimbildungsbarriere, während dies bei Silizium nicht der Fall ist.

Wenn Sie beispielsweise periodische Randbedingungen (keine Oberflächen) und keine anderen Defekte haben, kann Silizium sehr leicht überhitzt werden (zumindest auf der typischen Zeitskala einer MD-Simulation).

Wenn Sie den Schmelzpunkt messen möchten, ist es richtig, zunächst zwei Phasen - kristallin und amorph - in Kontakt miteinander zu erzeugen und dann herauszufinden, bei welcher Temperatur die beiden im Gleichgewicht sind. (Bei niedrigeren Temperaturen wächst die Kristallphase, bei höheren Temperaturen wächst die amorphe Phase).

Wenn es jedoch Ihr Ziel ist, das Nukleationsereignis zu simulieren, dann ist dies ein sehr schwieriges Problem und ein sehr aktives Forschungsgebiet. Sie könnten einfach zu einer höheren Temperatur gehen oder Ihre Zellgröße korrigieren (entfernen Sie die Presse / Berendsen-Korrektur) - das führt zu einem Phasenübergang, aber es ist nicht unbedingt ein realistischer Übergang. Dafür gibt es weitaus kompliziertere Verfahren, bei denen Rare-Event-Sampling-Methoden wie Metadynamik [ 1 ] oder Seeding-Methoden [ 2 ] zum Einsatz kommen.

Ich versuche den Schmelzpunkt zu messen. Ich denke, ich sollte zwei Phasen erstellen, aber können Sie näher darauf eingehen, wie das geht?
@maryam Es sollte periodisch sein, wenn die Hälfte des Kristalls geschmolzen ist; NPT-Ensemble; Führen Sie die Simulation bei mehreren Temperaturen durch (z. B. 1000 K, 1100 K, ...). Notieren Sie für jeden, ob der Kristall gewachsen ist oder die amorphe Phase gewachsen ist. Sie werden eine Schwellentemperatur finden, die diese beiden Klassen trennt – das ist der Schmelzpunkt. Wenn Sie eine genaue Messung der Schmelztemperatur benötigen, können Sie einen Ordnungsparameter (entweder potentielle Energie oder einen geometrischen wie den Zentrosymmetrieparameter) verwenden, um das Wachstum jeder Phase zu messen. ... Wenn Sie eine spezifische Frage stellen, kann ich eine spezifische Antwort geben.
Danke! Ich habe getan, was Sie gesagt haben, und die Ergebnisse, die ich erhalten habe, sind: Das Erhöhen der Temperatur bei jedem Lauf hat keinen Einfluss auf das System. Die amorphe Phase nimmt nicht zu, obwohl ich versucht habe, sie bei 2500 K laufen zu lassen. Ich denke also, dass etwas mit meinem Skript nicht stimmt. Ich würde es wirklich schätzen, wenn ich das Skript mit Ihnen teilen und mir helfen könnte, zu sehen, wo ich falsch gelaufen bin.