7 Στοχαστικά Συστήματα – Μέθοδοι Monte Carlo

7.5 Παράδειγμα – Μοριακό Μοντέλο Σκληρών Σφαιρών

Εδώ παρουσιάζουμε ένα πρώτο παράδειγμα εφαρμογής αλγορίθμων Monte Caro βασισμένων σε Μαρκοβιανές Αλυσίδες, και πιο συγκεκριμένα του αλγορίθμου Metropolis Monte Carlo.

Έστω N σφαιρικά σωματίδια (π.χ. άτομα ή μόρια), με συντεταγμένες

𝒓=(r1,r2,,rN)

σε ένα κουτί. Θεωρούμε ότι τα σωματίδια αλληλεπιδρούν ανά δύο, με ένα δυναμικό ζεύγους (pair potential) αλληλεπίδρασης της μορφής:

U(𝒓)=i=1,j>iNU(rij),U(rij)={,ανrijd0,ανrij>d (7.21)

όπου rij είναι το μέτρο της απόστασης μεταξύ δύο σωματιδίων i,j και d η διάμετρος κάθε σωματιδίου.

Αν εξετάσουμε το πρόβλημα σε 3 διαστάσεις (η πιο γενική μορφή, αλλά επίσης το πρόβλημα μπορεί να μελετηθεί σε 2 ή ακόμη και 1 διάσταση), με xi,yi, και zi, αντίστοιχα τις x,y,z συντεταγμένες, στο Καρτεσιανό επίπεδο, του μορίου i τότε η απόσταση μεταξύ δύο σωματιδίων i και j είναι:

tij=(xi-xj)2+(yi-yj)2+(zi-zj)2.

Πρακτικά, το ότι η ενέργεια είναι άπειρη για αποστάσεις μικρότερες ή ίσες του d σημαίνει ότι τα σωματίδια δεν μπορούν να βρεθούν σε αποστάσεις μικρότερες ή ίσες με d, ενώ μηδενική ενέργεια για αποστάσεις μεγαλυτερες του d, σημαίνει ότι σε αυτές τις αποστάσεις τα σωματίδια δεν αλληλεπιδρούν. Για το λόγο αυτό το παραπάνω μοντέλο ονομάζεται μοντέλο «σκληρών σφαιρών» (hard spheres model).

Παρατήρηση 7.9.

Αν θεωρήσουμε ότι τα παραπάνω σωματίδια αντιπροσωπεύουν ένα μοριακό σύστημα στο κανονικό στατιστικό σύνολο NVT (σταθερός αριθμός σωματιδίων N, όγκου του συστήματος V και θερμοκρασίας Τ) τότε η κατανομή την οποία υπακούουν είναι η κατανομή Boltzmann (δες Εξ. (7.11)):

π(𝒓)=exp[-βU(𝒓)]exp[-βU(𝒓)]d𝒓=1Zexp[-βU(𝒓)]

όπου β είναι η αντίστροφη θερμική ενέργεια, β=1/(kBT), όπου kB είναι η σταθερά Boltzmann και T η θερμοκρασία του συστήματος.

Σχήμα 7.7: Δυναμικό αλληλεπίδρασης μοντέλου σκληρών σφαιρών.

Στο Σχήμα 7.7 φαίνεται η μορφή του δυναμικού σκληρών σφαιρών. Προσέξτε ότι το παραπάνω δυναμικό καταρχήν δεν είναι ένα ρεαλιστικό δυναμικό για συγκεκριμένα μοριακά συστήματα. Σε μακροσκοπικό επίπεδο θα λέγαμε ότι τα σωματίδια συμπεριφέρονται σαν «μπάλες μπιλιάρδου».

Αντίθετα υπάρχουν πιο ρεαλιστικά δυναμικά (π.χ. δυναμικό Lennard – Jones) τα οποία εμφανίζουν πιο ρεαλιστική περιοχή άπωσης και έλξης. Ωστόσο το δυναμικό σκληρών σφαιρών είναι μια πολύ καλή πρώτη προσέγγιση για μεγάλο εύρος μοριακών συστημάτων, όπως αερίων, απλών μορίων, κολλοειδών, κλπ.

Ο αλγόριθμος Metropolis Monte Carlo για το παραπάνω μοντέλο σωματιδίων τύπου σκληρών σφαιρών αποτελείται από τα παρακάτω βήματα:

  • Βήμα 1: Έστω ότι τα σωματίδια βρίσκονται σε θέσεις r. Τυχαία επιλογή ενός σωματιδίου, π.χ. το σωματίδιο i με συντεταγμένες xi, yi, και zi.

  • Βήμα 2: (Κανόνας Μετάβασης) Προτεινόμενη νέα θέση μέσω τυχαίας μετατόπισης της θέσης του σωματιδίου σε κάθε διεύθυνση, δηλαδή:

    ritemp=ri+Δr,ή  (xitemp,yitemp,zitemp)=(xi+Δx,yi+Δy,zi+Δz)

    όπου οι ποσότητες Δx,Δy και Δz είναι τυχαίες μεταβλητές από μια συγκεκριμένη κατανομή, π.χ. την τυπική κανονική, δηλαδή ΔxN(0,1),ΔyN(0,1),ΔzN(0,1).

  • Βήμα 3: Κατόπιν υπολογίζουμε την ποσότητα:

    exp[-βΔU]=exp[-β(U(rtemp)-U(r))]
  • Βήμα 4: Αποδοχή της νέας θέσης για το σωματίδιο i, με πιθανότητα ΔU, δηλαδή δημιουργία τυχαίου αριθμού u, από ομοιόμορφη κατανομή στο διάστημα (0,1),uU(0,1) και:

    rinew={𝒓itemp,ανuexp[-ΔU]𝒓i,ανu>exp[-ΔU]
Παρατήρηση 7.10.

Προσέξτε ότι η παραπάνω ποσότητα ΔU, στο μοντέλο των σκληρών σφαιρών, είναι ή 0. Συνεπώς πρακτικά το βήμα 4 δεν χρειάζεται καθώς η προτεινόμενη μετάβαση είναι αποδεκτή κάθε φορά που ΔU=0 και μη αποδεκτή αν ΔU=.

Με τον παραπάνω αλγόριθμο μπορούμε να δημιουργήσουμε μια Μαρκοβιανή αλυσίδα Τ.Μ. των θέσεων των σωματιδίων 𝒓, η οποία, σε μεγάλους χρόνους, υπακούει την επιθυμητή κατανομή Boltzmann. Κατόπιν οι τιμές της Μαρκοβιανής αλυσίδας μπορούν να χρησιμοποιηθούν για να υπολογίσουν αναμενόμενες (μέσες) τιμές για διάφορες ποσότητες – συναρτήσεις των 𝒓, όπως για παράδειγμα την ενέργεια, την πίεση ως συνάρτηση τη θερμοκρασίας T.