Pobierz Modele wzrostu kryształów stałych i więcej Egzaminy w PDF z Transport tylko na Docsity! Materiały do wykładu Modele wzrostu kryształów stałych Marek Izdebski Instytut Fizyki PŁ 2016 i Spis treści Temat 1. Termodynamiczne podstawy równowagi fazowej i krystalizacji. ........................................1 1.1. Równowaga i quasi-równowaga ...............................................................................................1 1.2. Funkcje stanu.............................................................................................................................1 1.3. Równanie stanu .........................................................................................................................2 1.4. Pierwsza zasada termodynamiki ...............................................................................................2 1.5. Druga zasada termodynamiki....................................................................................................3 1.6. Odwracalność procesu i entropia ..............................................................................................3 1.7. Potencjał chemiczny..................................................................................................................4 1.7.1. Definicja potencjału chemicznego .....................................................................................4 1.7.2. I zasada termodynamiki uwzględniająca wymianę masy ..................................................4 1.7.3. Równanie Gibbsa-Duhema ................................................................................................6 1.8. Potencjały termodynamiczne ....................................................................................................6 1.9. Równowaga...............................................................................................................................7 1.10. Równowaga chemiczna...........................................................................................................8 Literatura do tematu 1 ......................................................................................................................8 Temat 2. Równowaga fazowa i diagramy fazowe ...............................................................................9 2.1. Reguła faz Gibbsa .....................................................................................................................9 2.2. Przyczyny krystalizacji ...........................................................................................................10 2.3. Równanie Clausiusa-Clapeyrona ............................................................................................10 2.4. Klasyfikacja przemian fazowych ............................................................................................11 2.5. Układy jednoskładnikowe .......................................................................................................12 2.6. Układy dwuskładnikowe .........................................................................................................14 2.6.1. Roztwór doskonały...........................................................................................................15 2.6.2. Roztwór niedoskonały......................................................................................................18 2.6.3. Roztwór regularny, w którym ∆hXS < 0 ...........................................................................19 2.6.4. Roztwór regularny, w którym ∆hXS > 0 ...........................................................................20 2.6.5. Punkt inwariantny ............................................................................................................22 2.6.6. Faza przejściowa ..............................................................................................................26 2.6.7. Topnienie kongruentne i niekongruentne.........................................................................27 Literatura do tematu 2 ....................................................................................................................28 Temat 3. Termodynamiczne modele blokowe wzrostu kryształów...................................................29 3.1. Termodynamiczny model Jacksona [1]...................................................................................29 3.2. Termodynamiczny model Temkina [2]...................................................................................32 3.3. Literatura podstawowa do tematu 3 ........................................................................................34 Temat 4. Zarodkowanie kryształów. ..................................................................................................35 4.1. Klasyczna teoria zarodkowania homogenicznego ..................................................................35 4.2. Zarodkowanie heterogeniczne.................................................................................................37 4.3. Literatura podstawowa do tematu 4 ........................................................................................38 Temat 5. Komputerowe symulacje wzrostu kryształów metodą Monte Carlo ..................................39 5.1. Idea Metody Monte Carlo .......................................................................................................39 5.2. Generatory liczb losowych......................................................................................................39 5.3. Symulacje wzrostu kryształów - skala zadania i typowe założenia ........................................40 5.4. Częstotliwość zachodzenia kreacji i anihilacji........................................................................41 5.5. Dyfuzja powierzchniowa.........................................................................................................44 5.6. Algorytmy symulacji...............................................................................................................44 5.7. Szybkość i czas wzrostu ściany kryształu ...............................................................................46 5.8. Wyniki symulacji MC .............................................................................................................47 5.9. Czy podejścia kinetyczne i termodynamiczne i są równoważne?...........................................52 5.10. Literatura podstawowa do tematu 5 ......................................................................................53 3 1.5. Druga zasada termodynamiki Istnieje wiele sformułowań II zasady termodynamiki. Sformułowanie 1 Spontaniczna tendencja układu do dążenia do stanu równowagi nie może być odwrócona bez zamiany w tym samym czasie jakiejś formy zorganizowanej energii (pracy) w energię niezorganizowaną (ciepło). Przykład Rozważmy przepływ ciepła. Doświadczenie pokazuje, że ciepło przepływa tylko od ciał cieplejszych do zimniejszych. Jeżeli użyjemy np. lodówki, żeby odwrócić ten proces, to musimy także wykonać pracę, która w otoczeniu lodówki prowadzi do zamiany energii zorganizowanej w niezorganizowaną. W większym układzie, który obejmuje lodówkę, chłodzone ciało i otoczenie lodówki razem ze źródłem energii zorganizowanej i ogrzewanym ciałem, tendencja do dążenia do stanu równowagi pozostaje zachowana. Sformułowanie 2 (wg. Clausiusa) Niemożliwe jest zbudowanie urządzenia, które pracując cyklicznie, nie daje innego efektu tylko przepływ ciepła od ciała zimniejszego do cieplejszego. 1.6. Odwracalność procesu i entropia Odwracalność procesu może być zdefiniowana na kilka sposobów: Sformułowanie 1 Jeżeli układ wymienia ciepło w powtarzalnym cyklu, to cykl ten jest odwracalny gdy 0=∆ ∑ i i i T Q lub 0=∫ T dQ . (1.6a) Przykładem układu pracującego w cyklu odwracalnym jest silnik Carnota, w którym: 0 2 3 1 1 =∆+∆ T Q T Q , (1.6b) gdzie ∆Q1 > 0 jest ciepłem pobranym w procesie izotermicznym ze źródła ciepła o temperaturze T1 (1 faza cyklu Carnota), zaś ∆Q3 < 0 jest ciepłem oddanym do chłodnicy o temperaturze T2 w procesie izotermicznym (3 faza cyklu Carnota). Każdy inny cykl odwracalny można przedstawić jako kombinację cykli Carnota. Ponieważ realizacja takiego procesu wymaga nieskończonego czasu, w praktyce nie istnieją procesy idealnie odwracalne. Definicja entropii w termodynamice klasycznej Równanie (1.6a) wprowadza nową funkcję stanu zwaną entropią S, której różniczka T Q S revd d = , (1.7) gdzie „rev” oznacza proces odwracalny (ang. reversible). Entropia jest ekstensywną zmienną stanu. To oznacza, że entropia może być wyznaczona tylko na podstawie znajomości stanu początkowego i końcowego, także dla dowolnego procesu nieodwracalnego, jednakże tylko proces odwracalny umożliwia wyznaczenie zmiany entropii na podstawie wymiany ciepła ∫=∆ T Q S revd . (1.8) 4 Przykład Podczas przemiany fazowej, takiej jak topnienie, ciepło utajnione ∆Q jest dostarczane do układu w ustalonej temperaturze T. Stąd wynika, że entropia cieczy jest większa niż entropia fazy stałej o czynnik ∆S = ∆Q/T. Ten wzrost entropii jest związany z obniżeniem strukturalnego porządku. Sformułowanie 2 W procesie odwracalnym entropia jest zachowana i przenoszona z jednego zbiornika ciepła do innego. Przyrost entropii układu nieizolowanego jest wówczas spowodowany tylko dopływem ciepła z otoczenia ∫=∆ T Q S revd . (1.9a) W procesie nieodwracalnym przenoszona jest mniejsza entropia, a różnica wytwarzana w procesie przenoszenia zwiększa entropię układu ∫≥∆ T Q S revd . (1.9b) Sformułowanie 3 Entropia układu adiabatycznego nie może nigdy zmaleć ∫≥∆ T dQ S . (1.10) Entropia jest miarą nieporządku w danym układzie. Entropia jest także miarą maksymalnej wymienialności ciepła na pracę. Jeżeli znamy zmianę entropii ∆S związana z przepływem ciepła przez układ, to możemy obliczyć bezpośrednio tą część ciepła, która jest niedostępna do konwersji na pracę, gdyż jest tracona w rezerwuarze ciepła o niższej temperaturze T2 Qniedost. ≥ T2 ∆S. (1.11) Pozostała część energii wewnętrznej, która może być wykorzystana do wykonania pracy w danym procesie, nazywana jest energią swobodną. Definicja entropii w termodynamice statystycznej Entropia statystyczna układu według Boltzmanna S = kB ln(W), (1.12) gdzie: kB = 1.38·10−23 J/K jest stałą Boltzmanna, W jest liczbą sposobów liczba sposobów, na jakie makroskopowy stan termodynamiczny układu może być zrealizowany poprzez stany mikroskopowe. 1.7. Potencjał chemiczny 1.7.1. Definicja potencjału chemicznego Potencjał chemiczny i-tego składnika układu opisuje zmianę energii wewnętrznej układu U gdy liczba moli tego składnika w układzie zmienia się o bardzo małą wartość dni ij nnMVSi i n U ≠ ∂ ∂=µ ,,, . (1.13) 1.7.2. I zasada termodynamiki uwzględniająca wymianę masy Jeżeli do układu wprowadzamy kilka substancji, to pierwsza zasada termodynamiki w postaci uwzględniającej jednoczesną zmianę masy i zmianę objętości układu o objętość dodanej materii, przybiera postać 5 ∑ +++−= i iii dnPvuPdVQU K)(dd . (1.14) gdzie … oznacza inne formy pracy nad układem, np. praca magnetyczna, elektryczna, powierzchniowa, ui - cząstkowa molowa energia wewnętrzna i-tego składnika ij nnPTi i n U u ≠ ∂ ∂= ,, , (1.15) oraz vi - cząstkowa molowa objętość i-tego składnika ij nnPTi i n V v ≠ ∂ ∂= ,, . (1.16) W przypadku procesów odwracalnych bez wymiany masy dQ = TdS, natomiast całkowita zmiana entropii podczas procesu odwracalnego z jednoczesną wymianą ciepła i masy ∑+= i ii ns T Q S d d d . (1.17) Tak więc wymiana ciepła podczas procesu odwracalnego z wymianą masy wynosi ∑−= i ii nsTSTQ ddd , dla procesu odwracalnego, (1.18) gdzie si to cząstkowa molowa entropia ij nnPTi i n S s ≠ ∂ ∂= ,, . (1.19) natomiast w przypadku procesu nieodwracalnego wymiana ciepła jest mniejsza. Warto zauważyć, że cząstkowe molowe wielkości ui, vi, si nie są specyficzne dla danych składników, lecz zależą także od koncentracji innych składników. Przykładowo istnieją takie układy, gdzie vi może być nawet ujemne. Analogiczne wielkości u, v, s dla czystych substancji są specyficzne dla danego materiału. Podstawiając dQ z równania (1.18) do (1.14) otrzymujemy kolejną postać I zasady termodynamiki K+−++−= ∑ i iiii nTsPvuVPSTU d)(ddd , (1.20a) co można także zapisać K+µ+−= ∑ i ii nVPSTU dddd , (1.20b) gdzie K−−+=µ iiii TsPvu (1.21) jest molowym potencjałem chemicznym. Równanie (1.20) możemy zapisać najogólniej oznaczając przez Yk zestaw intensywnych funkcji stanu −P, H, E, σ,… i przez Xk odpowiedni zestaw zmiennych ekstensywnych V, M, P, A,… ∑∑ µ++= i ii k kk nXYSTU dddd . (1.22) Ze względu na ekstensywny charakter U, S, Xk, ni możliwe jest przeprowadzenie całkowania związku (1.22), co prowadzi do wyznaczenia absolutnej wartości energii wewnętrznej: ∑∑ µ++= i ii k kk nYXTSU . (1.23) 8 )2()1( ii µ=µ . (1.39) Spontaniczny transport materii z fazy (2) do (1) zachodzi gdy )2()1( ii µ<µ . (1.40) Uwaga: Błędem jest zakładanie, że procesy dyfuzyjne zachodzą zawsze w kierunku malejącej koncentracji. Przykładowo w niektórych roztworach gradient temperatury powoduje rozdzielenie molekuł: cięższe molekuły gromadzą się w chłodniejszym obszarze. Prawo Ficka wiążące strumień dyfuzji z gradientem koncentracji ma więc zastosowanie tylko wtedy gdy inne parametry wchodzące w skład definicji potencjału chemicznego (1.21) są ujednolicone w całym układzie. Równowaga termiczna Rozważmy układ odizolowany złożony z jednego składnika w dwóch fazach. Z warunku osiągania przez entropię maksymalnej wartości wynika warunek równych temperatur faz )2( )2( )1( )1( )2()1( dd ddd0 T Q T Q SSS +=+== . (1.41) Ponieważ w układzie izolowanym dQ(1) = –dQ(2) otrzymujemy T(1) = T(2). (1.42) Równowaga mechaniczna Analogicznie z warunku osiągania minimum przez energię swobodną wynika warunek równych ciśnień. dF = dF(1) + dF(2) = –S(1)dT(1) – P(1)dV(1) + µ(1)dn(1) – S(2)dT(2) – P(2)dV(2) + µ(2)dn(2) = 0 (1.43) i ponieważ w układzie izolowanym dV(1) = – dV(2), dn(1) = – dn(2) oraz wcześniej wykazaliśmy, że µ(1) = µ(2) oraz T1 = T2 ⇒ T = const. dF = dF(1) + dF(2) = – (P(1) – P(2))dV(1) = 0, P(1) = P(2). (1.44) 1.10. Równowaga chemiczna Do tej pory rozważaliśmy równowagę tylko w układach, w których nie zachodzą reakcje chemiczne. Teraz podamy (bez dowodu), że warunek równowagi pomiędzy substratami i produktami danej reakcji odpowiadający minimum funkcji Gibbsa G jest następujący: 0=µν∑ i ii , (1.45) gdzie i jest numerem substancji, zaś νi są współczynnikami stechiometrycznymi reakcji. Przykładowo dla reakcji H2 + Cl2 ↔ 2HCl współczynniki dla substratów (zanikających reagentów) są ujemne ν1 = ν2 = −1, zaś współczynnik dla produktu jest dodatni ν3 = 2. Problemem przy próbie wykorzystania warunku (1.45) jest zależność µi od P i V. Literatura do tematu 1 [1] F. Rosenberger, "Fundamentals of Crystal Growth I" (Springer 1979). 9 Temat 2. Równowaga fazowa i diagramy fazowe 2.1. Reguła faz Gibbsa Rozważmy układ złożony z C składników (ang. components) rozdzielonych pomiędzy P homogeniczne fazy (ang. phases). Układ jest w równowadze jeżeli spełnia następujące warunki: )()2()1( )( 2 )2( 2 )1( 2 )( 1 )2( 1 )1( 1 )()2()1( )()2()1( P CCC P P P P ppp TTT µ==µ=µ µ==µ=µ µ==µ=µ === === K K K K K K , (2.1) gdzie indeksy górne oznaczają numer fazy, a indeksy dolne numer składnika. Rozważmy liczbę zmiennych niezależnych opisujących dany układ Koncentracja C-tego składnika w jednej fazie wynika jednoznacznie z koncentracji pozostałych C−1 składników. Stąd skład wszystkich P faz jest opisany jednoznacznie przez P(C−1) zmiennych niezależnych. Po uwzględnieniu temperatury i ciśnienia, otrzymujemy P(C−1) + 2 zmiennych niezależnych opisujących stan układu, który może nie być w równowadze. Układ (2.1) narzuca dodatkowo (P−1)C związków pomiędzy potencjałami chemicznymi, pozostawiając P(C−1) + 2 − (P−1)C zmiennych niezależnych. Reguła faz Gibbsa Liczba zmiennych niezależnych opisujących jednoznacznie stan układu w równowadze 2+−= PCf . (2.2a) Liczba f nazywana jest liczbą stopni swobody, tzn. jest to liczba niezależnych intensywnych zmiennych stanu, które można zmieniać bez zmiany liczby faz. Przykład 1 W układzie złożonym z samej wody w stanie ciekłym f = 2, co oznacza, że możemy niezależnie zmieniać np. P i T. Przykład 2 W układzie jednoskładnikowym dwufazowym otrzymujemy f = 1. Wystarczy podać jedną z dwóch wielkości, ciśnienie albo temperaturę, żeby określić jednoznacznie stan układu. Przykład 3 W układzie jednoskładnikowym trzyfazowym f = 0, co oznacza, że możliwy jest tylko jeden stabilny stan układu zwany „punktem potrójnym”. W tym stanie żadna zmienna nie może zostać zmieniona. W przypadku możliwości zajścia reakcji chemicznych liczba składników C powinna zostać zmniejszona o liczbę R różnych reakcji możliwych pomiędzy składnikami obecnymi w układzie: 2)( ++−= RPCf . (2.2b) Liczba R obejmuje tylko następujące reakcje: nie mogą być zapisane jako sekwencja innych możliwych reakcji, 10 zachodzą z mierzalną szybkością w rozważanym czasie. Przykłady będą podane jako zadania na ćwiczeniach. 2.2. Przyczyny krystalizacji Rozważmy układ zawierający jedną substancję w fazie ciekłej i stałej. Krystalizacja, według warunku )2()1( ii µ<µ (1.40), wymaga obniżenia potencjału fazy stałej poniżej potencjału fazy ciekłej. Z równania Gibbsa-Duhema (1.26) K−⋅−⋅−+−=−−=µ EpHm ddddddd PvTsYxTs (2.3) wynika, że obniżenie potencjału chemicznego można osiągnąć na różne sposoby: 1) Obniżenie temperatury - potencjały chemiczne cieczy i ciała stałego zależą od temperatury w inny sposób (Rys. 2.1) Rys. 2.1. Temperaturowa zależność potencjału chemicznego fazy stałej µs (ang. solid) i ciekłej µl (ang. liquid). Tm (ang. melting) oznacza temperaturę równowagi. Źródło: Fig. 3.1. [1]. 2) Zmiana ciśnienia. Znak tej zmiany zależy od znaku zmiany objętości podczas zamrażania: a) vl > vs (większość substancji), b) vl < vs (H20, Bi, Si, Ge). 3) Równanie stanu, np. gazu doskonałego, zawiera 3 zmienne P, T, V z czego do opisu stanu wybrano dwie niezależne P i T. Kontrola trzeciej zależnej zmiennej V może być jednak także wykorzystana do spowodowania krystalizacji. 4) Poddanie układu ciecz-ciało stałe działaniu pola elektrycznego lub magnetycznego – zazwyczaj jednak zmiany potencjałów chemicznych wywołane w ten sposób są za małe. 5) W układach wieloskładnikowych pojawiają się nowe możliwości kontroli, np. potencjał chemiczny jednego składnika może zmienić się na skutek dodania innego składnika, który nie podlega krystalizacji. 2.3. Równanie Clausiusa-Clapeyrona Rozważmy układ jednoskładnikowy dwufazowy. Równowaga w układzie podczas zmiany warunków zostanie zachowana gdy )2()1( dd µ=µ . (2.4) Jeżeli wszystkie inne intensywne zmienne stanu pozostają ustalone, to dla małych zmian ciśnienia P i temperatury T P P T T P P T T TPTP dddd )2()2()1()1( ∂ µ∂+ ∂ µ∂= ∂ µ∂+ ∂ µ∂ . (2.5) 13 Rys. 2.4. Typowy diagram fazowy we współrzędnych T-P dla układu jednoskładnikowego. Pominięto stany metastabilne. W czystym materiale może istnieć tylko jedna faza gazowa i tylko jedna faza ciekła, gdyż możliwe jest tylko jedno losowe uporządkowanie nieodróżnialnych cząsteczek. Możliwe są natomiast różne formy krystaliczne tej samej czystej substancji, zwane odmianami alotropowymi w przypadku pierwiastków lub odmianami polimorficznymi substancji chemicznej. Uwaga: niektóre pierwiastki mogą występować w fazie gazowej i/lub ciekłej w postaci cząsteczek złożonych z różnej liczby atomów, np. tlen o cząstkach O2 i O3. Mieszanin takich cząstek nie możemy traktować jak układu jednoskładnikowego. Przemiany polimorficzne mogą zachodzić jako rezultat zmian temperatury i ciśnienia. Warunkiem stabilnego współistnienia dwóch odmian polimorficznych jest, analogicznie jak w przypadku faz o różnym stanie skupienia, równość ich potencjałów chemicznych (rys. 2.5). Rys. 2.5. Obszary stabilności dwóch stałych odmian polimorficznych we współrzędnych µ-T. Źródło: Fig. 3.6. [1]. Szczególnie skomplikowany diagram fazowy czystej substancji wykazuje bizmut (rys. 2.6). Rys. 2.6. Diagram fazowy bizmutu we współrzędnych T-P. Źródło: Fig. 3.8. [1]. 14 Przykładem polimorfizmu/alotropowości o dużym znaczeniu przemysłowym jest przemiana ferromagnetycznego żelaza α w paramagnetyczne żelazo γ. 2.6. Układy dwuskładnikowe Rozważmy układ dwóch składników, między którymi nie zachodzą reakcje chemiczne. Diagram fazowy układu binarnego można przedstawić w postaci powierzchni równowagi faz we współrzędnych X-T-P, gdzie X jest względną koncentracją jednego ze składników i przyjmuje wartości 0…1. Takie podejście prowadzi jednak do skomplikowanej reprezentacji nawet dla szczególnie prostych przypadków (rys. 2.7). Znacznie prostsze są dwuwymiarowe wykresy fazowe, które powstają jako izobaryczne przecięcie przez wykres trójwymiarowy. Ponadto, na większości takich wykresów pokazywane są tylko obszary fazy stałej i ciekłej, gdzie zmiany ciśnienia mają znacznie mniejszy wpływ na wartości potencjałów chemicznych niż zmiany temperatury. Rys. 2.7. Uproszczony trójwymiarowy diagram fazowy układu binarnego. Źródło: Fig. 3.9. [1]. W rzeczywistości potrzebne są dwie powierzchnie do opisu dwóch współistniejących faz (czego nie widać na rys. 2.7), gdzie obszar pomiędzy powierzchniami reprezentuje mieszaninę. Mieszanie składników może zachodzić w trzech formach: 1) mieszanina mechaniczna (ang. mechanical mixture) – obszary jednolitej substancji są na tyle duże, że część atomów na granicy faz jest zaniedbywalnie mała, 2) roztwór (ang. solution) – mieszanina składników w skali atomowej, 3) mieszanina mechaniczna roztworów – mieszanina obszarów znacznie większych od rozmiarów atomowych, które zbudowane są z różnych roztworów. O tym która z form jest stabilna dla danego składu mieszaniny X oraz wartości P i T decyduje to, która forma ma najmniejszą wartość funkcji Gibbsa. Molowa funkcja Gibbsa mieszaniny mechanicznej jest wielkością addytywną względem czystych składników * BB * AABBAABBAA m )( µ+µ=+−+= XXsXsXThXhXg , (2.16) gdzie hA i hB są molowymi entalpiami składników A i B, sA i sB są molowymi entropiami, µ* A i µ* B reprezentują potencjały chemiczne czystych substancji, zaś XA i XB są ułamkami molowymi BA A A nn n X + = i BA B B nn n X + = . (2.17) W przypadku roztworów funkcja Gibbsa nie jest addytywna, gdyż rzeczywiste entalpie molowe, entropie molowe i potencjały chemiczne składników zależą od składu mieszaniny. Ponadto mieszanie składników A i B w roztworze prowadzi do pojawienia się entropii konfiguracyjnej mieszania. 15 . 2.6.1. Roztwór doskonały Roztwór doskonały (ang. ideal solution) to taki roztwór, przy powstawaniu którego: nie występuje zmiana entalpii (tzn. hi jak w substancji czystej), nie jest pobierane ani oddawane ciepło (tzn. si jak w substancji czystej), nie ma zmiany objętości. Oznacza to, że z roztworem doskonałym mamy do czynienia wtedy, gdy oddziaływania pomiędzy cząstkami rozpuszczalnika są podobne do oddziaływań substancji rozpuszczanej (izotopy). Roztwór idealny stosuje się do prawa Raoulta. Prawo Raoulta Prężność pary nasyconej i-tego składnika nad roztworem Pi jest równa prężności nad czystą substancją ciekłą Pi 0 pomnożoną przez ułamek molowy składnika w roztworze Xi iii XPP 0= . (2.18) Molowa funkcja Gibbsa dla roztworu idealnego m mids sTgg ∆−= , (2.19) gdzie gm jest molową funkcją Gibbsa mieszaniny mechanicznej (2.16), ∆sm jest entropią mieszania (wyprowadzona na ćwiczeniach) )]ln()ln([ BBAAm XXXXRs +−=∆ . (2.20) Jeżeli wyrazimy molową funkcję Gibbsa poprzez potencjały chemiczne analogicznie jak we wzorze (2.16) dotyczącym mieszaniny BBAA ids µ+µ= XXg , (2.21) to, jak wynika ze wzorów (2.16a), (2.19) i (2.20), potencjały chemiczne µi substancji w mieszaninie zależą od ułamków molowych )ln(* iii XRT+µ=µ . (2.22) Wniosek Entropia mieszania ∆sm (2.20) roztworu idealnego jest zawsze dodatnia dla dowolnego składu Xi ∈ [0; 1]. Stąd wynika, że gids < gm, czyli roztwór doskonały jest zawsze stabilniejszy niż mieszanina mechaniczna (rys. 2.8). Rys. 2.8. Molowa funkcja Gibbsa mieszaniny mechanicznej gm i roztworu idealnego gids oraz entropia mieszania jako funkcja składu roztworu. Źródło: Fig. 3.10. [1]. Dla ustalonego ciśnienia gm zgodnie ze wzorem (2.16) maleje ze wzrostem temperatury, ponadto czynnik T∆sm we wzorze (2.19) powoduje, że zależność gids(XB) staje się coraz bardziej wygięta ze wzrostem T (rys. 2.9). 18 Wykres równowagi fazowej w kształcie soczewki został wyprowadzony teoretycznie przy założeniu doskonałości roztworu. Znane są jednak rzeczywiste układy dwuskładnikowe o bardzo podobnych diagramach fazowych, np. german-krzem, miedź-nikiel. Rozważmy równowagowe krzepnięcie roztworu o składzie początkowym X0 przy stałym strumieniu oddawanego ciepła. Po osiągnięciu temperatury Ti zaczyna powstawać kryształ o składzie początkowym Xi (rys. 2.12a), co odpowiada skokowej zmianie nachylenia krzywej chłodzenia (rys. 2.12b). W miarę krzepnięcia coraz większej części układu skład roztworu ciekłego przesuwa się w kierunku Xf, natomiast skład całej fazy stałej w kierunku X0. Po osiągnięciu temperatury Tf cała faza stała ulega zakrzepnięciu, co odpowiada kolejnej zmianie nachylenia krzywej chłodzenia. Opisany proces odbywa się przy założeniu, że dyfuzja składnika B z fazy ciekłej do skrystalizowanej już masy na bieżąco wyrównuje skład całej fazy stałej. Ponieważ dyfuzja w fazie stałej jest procesem powolnym, konieczna jest bardzo powolna krystalizacja żeby uniknąć gradientu koncentracji składników w fazie stałej. Gwałtowne schładzanie wykorzystuje się do rozdziału składników podczas krystalizacji. Model roztworu doskonałego wyklucza topnienie kongruentne, tzn. topienie zachodzące bez zmiany składu, z wyjątkiem przypadku topnienia jednego czystego składnika. Digramy fazowe bardzo zbliżone do soczewkowego (jak na rys. 2.12a) posiadają np. układy Ge-Si oraz Cu-Ni. 2.6.2. Roztwór niedoskonały Molowa funkcja Gibbsa dla roztworu niedoskonałego XS m ms gsTgg ∆+∆−= , (2.24) gdzie gm jest molową funkcją Gibbsa binarnej mieszaniny mechanicznej (2.16), ∆sm jest entropią mieszania molekularnego w roztworze idealnym (2.20), oraz ∆gXS opisuje nadmiar energii swobodnej (ang. excess free energy) na skutek obecności drugiego składnika w roztworze XSXSXS sThg ∆−∆=∆ . (2.25) W roztworach idealnych ∆gXS = 0. Liczne rzeczywiste roztwory można zadowalająco opisać jako roztwory regularne, w których ∆sXS = 0 i wzór (2.24) redukuje się do postaci XS m ms hsTgg ∆+∆−= . (2.26) Dodatkowa zmiana entalpii przy mieszaniu ∆hXS wynika z zamiany części wiązań między atomami w czystych substancjach A oraz B na nierównoważne wiązania typu AB i może być obliczana na różne sposoby. Zakładając ustaloną liczbę wiązań w całym roztworze, na każde utworzenie jednego wiązania typu AB przypada średnio 0,5 zerwanego wiązania AA i 0,5 zerwanego wiązania BB (rys. 2.13). Stąd dla 1 mola roztworu ( )[ ]BBAA2 1 ABAB XS HHHLh +−=∆ , (2.27) gdzie LAB jest liczbą wiązań typu AB w 1 molu roztworu, HXY są energiami 1 wiązania. A B Rys. 2.13. Zamiana atomów pierwotnie wbudowanych w czyste substancje A i B. Powstaje 8 wiązań AB kosztem zerwania 4 wiązań AA i 4 wiązań BB. 19 Zakładając całkowicie losowe rozmieszczenie atomów typu A i B, co znane jest jako przybliżenie Bragga-Williamsa lub przybliżenie zerowego rzędu, możemy łatwo stwierdzić, że prawdopodobieństwo znalezienia atomu B w pewnym wybranym miejscu jest równe ułamkowi molowemu XB, natomiast prawdopodobieństwo znalezienia na jednym wybranym sąsiednim miejscu atomu A wynosi XA = 1 – XB. Stąd liczba wiązań AB na 1 mol roztworu ),1( BBAAB XXzNL −= (2.28) gdzie: z - liczba koordynacyjną, NA – stała Avogradro. Ze wzorów (2.27) i (2.28) otrzymujemy zmianę entalpii mieszania w postaci Ω−=∆ )1( BB XS XXh , (2.29) gdzie Ω jest parametrem oddziaływania stałym dla danych substancji A i B. Podstawiając wzory (2.16), (2.20) i (2.29) do (2.26) i uwzględniając, że XA = 1 − XB otrzymujemy ostatecznie )]ln()1ln()1[()1()1( BBBBBB * BB * AB s XXXXTRXXXXg +−−+Ω−+µ+µ−= . (2.30) Zależność funkcji gs danej wzorem (2.30) od XB przedstawiono na rys. 2.14. 0 0.2 0.4 0.6 0.8 1 ← g XB → µA* µB* Ω > 0 Ω = 0 Ω < 0 T = const. gm Rys. 2.14. Zależność molowej energii Gibbsa gs danej wzorem (2.30) od składu roztworu XB przy różnych wartościach parametru Ω. 2.6.3. Roztwór regularny, w którym ∆hXS < 0 Rozważmy układ złożony z roztworu regularnego z fazie ciekłej oraz stałej, w którym entalpia mieszania ∆hXS jest ujemna dla obu faz. W tym przypadku zarówna entropia mieszania ∆Sm jak i entalpia mieszania ∆hXS powodują wygięcie zależności gs(XB) w dół. Ponieważ jednak możliwe są różne relacje ∆hXS(l) < ∆hXS(s), albo ∆hXS(l) = ∆hXS(s), albo ∆hXS(l) < ∆hXS(s) (2.31) możemy otrzymać trzy typy diagramów fazowych, które przedstawiono na rys. 2.15. W roztworach regularnych, w których ∆hXS(l) ≠ ∆hXS(s), możliwe jest spełnienie równości g(s)(XB) = g(l)(XB) także dla pośredniego składu roztworu. Równość ta odpowiada ekstremum linii solidusu i likwidusu na diagramie fazowym T~XB (rys. 2.15a i c). W punkcie ekstremum faza stała i ciekła o takim samym składzie mogą współistnieć w stanie równowagi. W przypadku gdy ∆hXS(l) ≈ ∆hXS(s) diagram fazowy dla roztworu regularnego ma kształt soczewkowaty (rys. 2.15b), podobny do kształtu diagramu dla roztworu idealnego. Warto zauważyć, że w przypadku ujemnej entalpii mieszania ∆hXS nie występują obszary niemieszalności składników. 20 a) ∆hXS(l) < ∆hXS(s), b) ∆hXS(l) ≈ ∆hXS(s), c) ∆hXS(l) > ∆hXS(s), Rys. 2.15. Schematyczne diagramy fazowe dla ∆hXS < 0. Źródło: Fig. 3.17. [1]. 2.6.4. Roztwór regularny, w którym ∆hXS > 0 W przypadku gdy entalpia mieszania jest dodatnia (także parametr Ω > 0), entalpia ∆hXS we wzorze (2.26) konkuruje z ujemnym wyrażeniem −T∆sm. Wykres zależności gs(XB) danej wzorem (2.30) może więc przyjmować różne kształty w zależności od temperatury: w wysokich temperaturach, np. T1 na rys. 2.16, wyraz −T∆sm jest dominujący i wykres zależność gs(XB) jest wygięty w dół, wraz ze spadkiem temperatury konkurujące wyrazy stają się porównywalne i wykres zależności gs(XB) staje się coraz silniej pofalowany (porównać T2 i T3 na rys. 2.16). 0 0.2 0.4 0.6 0.8 1 ← g XB → µA* µB* T3 Ω = const. > 0 gm T2 > T3 T1 > T2 a b Rys. 2.16. Wpływ temperatury T na krzywe zależności gs(XB) dla roztworu regularnego (jedna faza) przy Ω > 0. Jeżeli wykres zależność gs(XB) leży powyżej stycznej do wykresu w dwóch punktach, np. punkty a i b dla temperatury T3 na rys. 2.16, to stanem najkorzystniejszym energetycznie jest mieszanina mechaniczna dwóch roztworów o składach danych punktami styczności a i b. Rozważając cały przedział temperatur, w którym istnieje dana faza, można skonstruować diagram przedstawiony na rys. 2.17, na którym solwus rozgranicza obszar równowagowego istnienia jednego roztworu od obszaru współistnienia dwóch różnych roztworów w tej samej fazie. W tym drugim przypadku względne ilości roztworów mogą być wyznaczone przy użyciu reguły dźwigni (2.23). Niemieszalność roztworów rzadko występuje w fazie ciekłej, dla której entropia mieszania jest zwykle wystarczająco duża by zapobiec wygięciu do góry wykresu zależności gs(XB). 23 Rys. 2.20. Układ eutektyczny, T1 < T2 < … < T6. a) Krzywe funkcji Gibbsa dla fazy ciekłej (l) i dwóch faz stałych (α) i (β). b) diagram fazowy. Źródło: Fig. 3.22. [1]. Rys. 2.21. Krzepnięcie układu eutektycznego dla trzech przykładowych składów. a) Krzywe chłodzenia. b) diagram fazowy z uproszczoną linią solwusu Źródło: Fig. 3.23. [1]. (LS – liquid solution, SS – solid solution). Linia solwusu na rysunku 2.21b została przybliżona prostokątem. W rzeczywistych układach linia ta może być nachylona, co prowadzi do zmian w składzie fazy stałej omówionych wcześniej na przykładzie rys. 2.19a. Linia reakcji eutektycznej może sięgać na diagramie fazowym aż do składu odpowiadającego jednemu czystemu składnikowi (np. w układzie Al-Si, rys. 2.22) lub rozciągać się w całym zakresie składów 0…100%. czas 1 2 3 te m pe ra tu ra LS SSα+SSβ LS LS SSα +SSβ SSα LS +SSβ LS +SSα a) a) b) b) 24 Rys. 2.22. Diagram fazowy układu Al-Si i przykładowa struktura fazy stałej zawierającej 80% Al i 20% Si (Si - większe czarne ziarna). Źródło: Fig. 3.25. [1]. Jeżeli w układzie eutektycznym jeden z czystych składników ma dwie odmiany alotropowe, np. składnik A na rys 2.23, to oprócz przemiany eutektycznej L ↔ Sβ + Sγ w temperaturze T3, możliwa jest także analogiczna przemiana eutektoidalna całkowicie w fazie: Sβ ↔ Sα + Sγ. Rys. 2.23. Diagram fazowy z przemianą eutektyczną w temp. T3 i przemianą eutektoidalną w temp. T5. Źródło: Fig. 3.26. [1]. W układach, w których występuje niemieszalność w fazie ciekłej, zachodzi przemiana monotektyczna L(1) L(2) + S(1). chłodzenie ogrzewanie 25 W układach, w których występuje przemiana montektyczna, zazwyczaj możliwe są także inne przemiany inwariantne, jak pokazano na rys. 2.24. W zakresie temperatur TM…TC dwie ciecze α i β tworzą osobne warstwy, natomiast w temperaturach T > TC ciecze te tworzą jeden roztwór. Podczas chłodzenia składu monotektycznego XM w temperaturze TM następuje rozdzielenie fazy ciekłej Lα na Lβ i Sα. Dalsze chłodzenie przesuwa skład fazy Lβ w kierunku większego udziału składnika A oraz skład fazy Sα w kierunku większego udziału składnika B. W momencie osiągnięcia temperatury eutektycznej TE skład ostatnich kropli fazy Lβ wynosi XE. Rys. 2.24. Diagram fazowy z przemianą monotektyczną w temp. TM i eutektyczną w temp. TE. Źródło: Fig. 3.27. [1]. Kolejnym typem trójfazowej przemiany inwariantnej jest przemiana perytektyczna Układy perytektyczne powstają przy dużej różnicy temperatur topnienia dwóch czystych składników, np. MnO o temp topnienia TmA = 1785°C i Fe0 TmB = 1365°C (rys. 2.25). W układzie perytektycznym minimum krzywej funkcji Gibbsa g(XB) dla fazy ciekłej L przypada na osi XB poza minimami dla dwóch faz stałych α i β (dla porównania: w układzie eutektycznym minimum dla fazy ciekłej leży pomiędzy minimami dla dwóch faz stałych – rys. 2.20). Rys. 2.25. Diagram fazowy z przemianą perytektyczną w temp. Tp. Źródło: Fig. 3.28. [1]. Rozważmy teraz równowagowe chłodzenie stopów o różnych składach Xi oznaczonych na 2.25: L(1) + S(1) S(2). chłodzenie ogrzewanie 28 Rys. 2.29. Topnienie inkongruentne fazy przejściowej γ w układzie dwuskładnikowym. Źródło: Fig. 3.35. [1]. Literatura do tematu 2 [1] F. Rosenberger, "Fundamentals of Crystal Growth I" (Springer 1979). Topnienie inkongruentne 29 Temat 3. Termodynamiczne modele blokowe wzrostu kryształów W modelach blokowych wzrostu kryształów wszystkie zjawiska zachodzące na powierzchni kryształu zostały sprowadzone do przyłączania i odłączania identycznych bloków fazy stałej. Bloki te mogą znaleźć się tylko w węzłach prostej sieci krystalicznej – zwykle prostej sieci regularnej lub tetragonalnej. W termodynamicznych modelach blokowych rozważa się minimum energii swobodnej granicy faz. Modele znane z literatury można podzielić ze względu na ograniczenia nakładane na grubość granicy faz ciało stałe-faza macierzysta: 1. Model Jacksona – granica o grubości jednej warstwy bloków, 2. Model Mutaftschieva – granica dwuwarstwowa, 3. Model Temkina – bez ograniczeń dla grubości granicy faz. W kinematycznych modelach blokowych rozważa się sekwencje elementarnych zjawisk opisanych poprzez średnie częstotliwości zachodzenia. Takie modele są głównie podstawą komputerowych symulacji Monte Carlo. 3.1. Termodynamiczny model Jacksona [1] W modelu Jacksona rozważana jest tylko jedna ściana (001) kryształu z układu regularnego lub tetragonalnego. Warunki wzrostu kryształu są opisane dwoma parametrami: 1). Współczynnik Temkina zdefiniowany wzorem , 4 BTk Φ=α (3.1) gdzie 4 – liczba wszystkich bocznych sąsiadów stałych i ciekłych danego bloku, kB – stała Boltzmanna, T – temperatura, Φ – średnia zmiana energii na jedno wiązanie gdy wymieniany jest blok z obszaru całkowicie stałego z blokiem z obszaru całkowicie ciekłego ( ) , 2 1 sfffss Φ−Φ+Φ=Φ (3.2) gdzie Φss, Φff, Φsf są negatywnymi energiami wiązania bloków odpowiednio solid-solid, fluid-fluid, solid-fluid. Zrywane są: 4 wiązania solid-solid, 4 wiązania fluid-fluid, powstaje: 8 wiązań solid-fluid. Rys. 3.1. Przemieszczenie bloku stałego pierwotnie wbudowanego w obszar całkowicie stały. 2). Współczynnik chropowacenia kinetycznego - opisuje siłę napędową krystalizacji TkB sf µ−µ=β , (3.3) gdzie: µf – potencjał chemiczny fazy ciekłej/gazowej (ang. fluid), µs – potencjały chemiczne fazy stałej (ang. solid). 30 Można wykazać, że w roztworze doskonałym β = ln(1 + σ) , (3.4) gdzie σ jest względnym przesyceniem fazy macierzystej. Dla niewielkich przesyceń β ≈ σ . (3.5) Definicja Roztwór doskonały (lub mieszanina doskonała) to taki roztwór (mieszanina), który powstając w izotermiczno-izobarycznym procesie dyfuzji nie doznaje zmiany objętości i energii wewnętrznej. W modelu Jacksona jako termodynamiczny stan odniesienia (ang. reference state) rozważa się jednowarstwową granicę faz złożoną z N = Ns + Nf bloków: Ns bloków w obszarze w 100% stałym oraz Nf bloków w obszarze w 100% ciekłym/gazowym. Energia swobodna takiego układu F r = µs Ns + µf Nf . (3.6) Uwaga: w tym rozdziale potencjały chemiczne µs oraz µf odnoszą się do zmiany funkcji Gibbsa przy zmianie liczby cząstek o 1. Używany jest także potencjał chemiczny odnoszący się do zmiany liczby cząstek o 1 mol (np. rozdział 1). Gdy stałe i ciekłe/gazowe jednostki ulegają wymieszaniu energia swobodna jest równa F = F r + ∆Umix − T ∆Smix , (3.7) gdzie zakład się stałość F r podczas mieszania (stałość potencjałów chemicznych), ∆Umix jest energią mieszania, ∆Smix jest entropią mieszania. Energia mieszania ∆Umix wynika z powstawania wiązań typu solid-fluid w miejsce rozrywanych wiązań solid-solid oraz fluid-fluid ,4sfsfmix NkTNU α=Φ=∆ (3.8) gdzie Nsf jest liczbą wiązań pomiędzy blokami solid-fluid w warstwie granicznej. Przybliżenie Bragga-Williamsa Rozkład statystyczny Nsf/N w zależności od parametrów α i β nie został analitycznie znaleziony. Cechą wspólną wszystkich trzech modeli Jacksona, Mutaftschieva oraz Temkina jest wykorzystanie przybliżenia Bragga-Williamsa (zwanego też przybliżeniem zerowego rzędu) w celu oszacowania Nsf/N. Zakłada się całkowicie losowe wymieszanie bloków stałych oraz ciekłych/gazowych. W konsekwencji prawdopodobieństwo znalezienia bloku stałego w pewnym wybranym miejscu nie zależy od wypełnienia miejsc sąsiednich i jest równe ułamkowi bloków stałych wśród wszystkich bloków w warstwie xs = Ns/N. Stąd wynika, że Nsf = 4 N xs (1 − xs) . (3.9) Przybliżenie to prowadzi do zawyżonego oszacowania liczby Nsf w stosunku do wyników symulacji Monte Carlo, wskazujących na preferowanie wiązań między blokami jednego rodzaju. Entropia w termodynamicznym stanie odniesienia wynosi S r = 0 . (3.10) Entropię po wymieszaniu bloków stałych i ciekłych/gazowych wyznaczymy ze wzoru Boltzmanna Smix ≈ kB ln(g), (3.11) gdzie g jest liczbą sposobów na które można rozmieścić N = Ns + Nf bloków w warstwie 33 )]!([ )!( )!( 11 1 ++ + − = nnn n n xxNxN xN g . (3.22) Entropia mieszania dla całej granicy faz ∏ +∞= −∞= ++ − ==∆ n n nnn n xxNxN xN kgkS )]!([ )!( )!( )ln( 11 BBmix . (3.23) Wykorzystując wzór Stirlinga ln N! ≈ N ln N − N oraz uwzględniając warunki brzegowe x−∞ = 1 i x+∞ = 0 po przekształceniach otrzymuje się (pośrednie kroki przekształceń podano w pracy [2]): ( )∑ +∞ −∞= ++ −−−=∆ n nnnn xxxxkTNST 11mix ln)( . (3.24) Po podstawieniu wzorów (3.20), (3.21) i (3.24) do (3.19) otrzymujemy całkowitą zmianę energii swobodnej Gibbsa podczas chropowacenia powierzchni początkowo gładkiej ( )∑∑∑∑ ∞ −∞= ++ ∞ −∞= ∞ =−∞= −−+−α+ −−β=∆ n nnnn n nn n n n n xxxxxxxx TNk G 11 1 0 B ln)()1()1( , (3.25) gdzie α jest współczynnikiem Temkina, zaś β współczynnikiem chropowacenia kinetycznego. Energia ∆G jest funkcją nieskończonej liczby parametrów xn. W celu sprawdzenia, kiedy granica faz będzie stabilna, należy znaleźć minimum ∆G względem wszystkich xn. Warunkiem koniecznym istnienia minimów jest: 0 )/( = ∂ ∆∂ nx NkTG dla wszystkich xn . (3.26) Układ równań (3.26) został rozwiązany tylko na drodze numerycznej [2]. Krzywa β(α) przedstawiona na rys. 26.4 rozdziela dwa charakterystyczne obszary: 1). Pod krzywą β(α) energia ∆G posiada jedno minimum względem każdego xn → przy β > 0 możliwy jest tylko warstwowy mechanizm wzrostu przez dwuwymiarowe zarodkowanie lub wzrost spirali wokół dyslokacji śrubowej. 2). Nad krzywą β(α) brak minimum ∆G względem xn → brak barier termodynamicznych dla swobodnego wzrostu kryształu – ciągły mechanizm wzrostu. 0,00001 0,0001 0,001 0,01 0,1 1 1 2 3 4 α β brak minimum G , powierzchnia niestabilna powierzchnia stabilna Rys. 3.4. Liczba minimów funkcji ∆G(xs) danej wzorem (3.25) w zależności od wartości parametrów α i β. 0 0,2 0,4 0,6 0,8 1 1,2 1 1,5 2 2,5 3 3,5 4 α β brak bariery dla krystalizacji, wzrot ciągły, chropowacenie kinetyczne Jackson Temkin Rys. 3.5. Porównanie przewidywań wg.modeli Jacksona i Temkina. Zniesienie ograniczeń na grubość granicy faz w modelu Temkina prowadzi do zauważalnej zmiany krzywej α(β) rozdzielającej obszary charakterystycznych mechanizmów wzrostu w 34 porównaniu do analogicznej krzywej w modelu Jacksona, jednakże kształt obu krzywych jest zbliżony (rys. 26.5). 3.3. Literatura podstawowa do tematu 3 [1] P. Bennema, KONA Powder and Particle, nr 10 (1992), str. 25-40. [2] P. Bennema, G.H. Gilmer, w: Crystal Growth: an Introduction, ed. P. Hartman, North-Holland Publ., (1973), str. 263-327. Kopie prac [1] i [2] można pobrać ze strony internetowej: www.if.p.lodz.pl/marek.izdebski/MC/ 35 Temat 4. Zarodkowanie kryształów. Zarodki kryształów powstają na skutek fluktuacji wielkości fizycznych, np. ciśnienia, składu, itp. Zarodkowanie nazywamy: homogenicznym, gdy tworzący się zarodek powstaje w kontakcie z macierzystą fazą, heterogenicznym, gdy stabilny zarodek pozostaje na uprzednio istniejącej powierzchni (powierzchnia obcego ciała, rosnącego kryształu lub warstwy). Teorie zarodkowania: teorie klasyczne, teorie atomistyczne, inne podejścia oparte na rozwiązaniach mikroskopowych równań kinetycznych. Zarodkami krytycznymi nazywamy takie zarodki, które mają równe prawdopodobieństwo wzrostu i rozpadu. 4.1. Klasyczna teoria zarodkowania homogenicznego Zarodkowanie homogeniczne jest zarodkowaniem występującym rzadko, wymagającym wielkich wartości przechłodzenia. Powstawaniu zarodków w fazie pary lub fazie ciekłej towarzyszy zmiana funkcji Gibbsa układu ∆G ∆G = −∆GV + ∆GS , (4.1) gdzie ∆GS jest zmianą związana z utworzeniem powierzchni zarodka ∆GS = α S, (4.2) α jest energią powierzchniową właściwą, S – powierzchnia granicy rozdziału faz. ∆GV jest pracą objętościową przy kondensacji gazu ∆GV = ∆gV V, (4.3) gdzie różnica między entalpią swobodną jednostkowej objętości kryształu i fazy macierzystej molmol )( VVV n V G g ckckV V µ∆=µ−µ=µ−µ=∆=∆ , (4.4) gdzie w przypadku wzrostu z roztworu lub pary ∆µ można powiązać ze względnym przesyceniem σ zależnością ∆µ = RT ln(1 + σ) ≈ RT σ (ostatnia równość tylko dla małych wartości przesyceń σ). W przypadku wzrostu z roztopu równowadze fazy stałej i ciekłej w temperaturze TE (ang. equilibrium) odpowiada ∆gV = 0, natomiast w innych temperaturach .VVV sThg ∆−∆=∆ (4.5) Przyjmując, że entalpia krystalizacji ∆hV i entropia krystalizacji ∆sV na jednostkę objętości nie zależą od temperatury (co jest poprawne dla niezbyt dużych przechłodzeń ∆T) możemy podstawić ET h s V V ∆=∆ (4.6) otrzymując z równania (4.5) , EE E T T h T TT hg VVV ∆∆=−∆=∆ (4.7) gdzie wielkość ∆hV oznacza tutaj utajone ciepło topnienia na jednostkę objętości. Jeżeli energia granicy zarodek-faza macierzysta nie zależy od jej orientacji, to ze względów energetycznych najkorzystniejszy jest kulisty kształt zarodka i wówczas απ+∆π−=∆ 23 4 3 4 rgrG V . (4.8) 38 )cos2()cos1( 3 )3( 3 2 32 θ+θ−π=−π= r hr h V , (4.19) pole powierzchni podstawy Sp = π ρ2 = π r2 (1 − cos2θ) (4.20) oraz pole powierzchni czaszy kuli Scz = 2π rh = 2π r2(1 − cos θ). (4.21) Warunkiem istnienia zarodka jest równowaga napięć powierzchniowych αAP = αBP + αAB cos θ, (4.22) gdzie αAP, αBP i αAB są właściwymi energiami swobodnymi [N/m] = [J/m2] pomiędzy powierzchniami rozdziału odpowiednio: faza pary-podłoże, kondensat-podłoże i kondensat-faza pary. Składowa powierzchniowa entalpii swobodnej utworzenia zarodka ∆GS = Scz αAB + Sp (αBP − αAP). (4.23) Stąd, po podstawieniu wzorów (4.20), (4.21), αBP − αAP = −αAB cos θ z warunku równowagi (4.22) otrzymujemy ∆GS = 2π r2(1 − cos θ) αAB − π r2(1 − cos2θ)cos(θ) αAB = = π r2(1 − cos θ)2 (2 + cos θ) αAB. (4.24) Składowa objętościowa entalpii swobodnej utworzenia zarodka )cos2()cos1( 3 2 3 θ+θ−π∆=∆=∆ r gVgG VVV . (4.25) Całkowita entalpia swobodna ∆G = ∆GS + ∆GV = π (1 − cos θ)2 (2 + cos θ) (r2 αAB − r3 ∆gV /3). (4.26) Zakładając stałość kąta kontaktowego θ otrzymujemy ( ) ,02)cos2()cos1( )( 2 AB 2 =∆−αθ+θ−π= ∂ ∆∂ Vgrr r G (4.27) skąd wynika promień krytyczny zarodka Vg r ∆ α= AB k 2 (4.28) oraz, po odstawieniu do (4.26), energia utworzenia zarodka o promieniu krytycznym .)cos2()cos1( 3 4 2 3 2 k V AB g G ∆ αθ+θ−π=∆ (4.29) Porównując tą energię z analogicznym wzorem (4.11) dotyczącym zarodkowania homogenicznego otrzymujemy ,)( homo k het k GG ∆θχ=∆ (4.30) gdzie funkcja )cos2()cos1( 4 1 )( 2 θ+θ−π=θχ (4.31) jest miara katalitycznej siły podłoża, mianowicie: jeżeli αAP ≥ αBP + αAB ⇒ θ = 0°, χ(θ) = 0, ∆Gk = 0 i łatwo powstają płaskie krople, jeżeli αBP ≥ αAP + αAB ⇒ θ = 180°, χ(θ) = 1 i zarodkowanie na podłożu staje się mało prawdopodobne. 4.3. Literatura podstawowa do tematu 4 [1] Józef Żmija, Podstawy teorii wzrostu kryształów, PWN Warszawa 1987. 39 Temat 5. Komputerowe symulacje wzrostu kryształów metodą Monte Carlo 5.1. Idea Metody Monte Carlo W metodzie Monte Carlo (MC) decyzje dotyczące przebiegu symulacji zapadają na podstawie wartości liczbowych otrzymanych z generatora liczb (pseudo)losowych. Metoda MC jest stosowana głównie do modelowania procesów zbyt złożonych, by można było przewidzieć wyniki w sposób analityczny, np.: - symulacje zderzeń wysokoenergetycznych cząstek z jądrem atomowym, - obliczanie całek oznaczonych z funkcji o nieznanych funkcjach pierwotnych, - symulacje ruchu samochodowego w dużym mieście. Przykład - liczba π jako pole S okręgu o promieniu 1 Losujemy współrzędne punktów (x, y) przy użyciu generatora o rozkładzie równomiernym i zakresie 0…1. N N S kole w22= . Rys. 5.1. Przykładowy wynik losowania dla N=9. Dokładność wyznaczenia wartości π rośnie ze wzrostem liczby losowań N. 5.2. Generatory liczb losowych 1). Generatory programowe W bibliotekach standardowych dostarczanych z językami programowania najczęściej dostępne są generatory multiplikatywno-addytywne obliczające kolejną wartość pseudolosową xi na podstawie poprzedniej xi−1 według wzoru xi = (a xi−1 + b) mod m, (5.1) gdzie a, b i m są odpowiednio dobranymi stałymi. Stała a jest wybierana jako duża liczba pierwsza. Przykładowo, generator liczb pseudolosowych całkowitych z przedziału 0…(215−1) dostarczany jako funkcja rand() w bibliotece standardowej kompilatora Borland C++ 5.5 wykonuje następujące działania: ,2mod)2/( ,mod1(22695477 1516 1 ii ii xRND xx = 2 )+= 32 − (5.2) gdzie RNDi są wynikami zwracanymi przez funkcję rand(). N S 1 4 2 4 3 2,667 4 3 5 3,2 6 3,333 7 2,857 8 3 9 3,111 +∞ 3,142 0 1 1 –1 –1 x y 7 w kole x2 + y2 ≤ 1 2 poza kołem x2 + y2 ≥ 1 40 2). Generatory sprzętowe Rys. 5.2. Schemat blokowy typowego sprzętowego generatora liczb losowych. 5.3. Symulacje wzrostu kryształów - skala zadania i typowe założenia Przeprowadzenie symulacji 3D wzrostu kryształu o makroskopowych rozmiarach (np. sześcian 1 mm × 1 mm × 1 mm) obejmującej pełny czas wzrostu przekracza możliwości dostępnych obecnie komputerów. W krysztale NaCl: • Odległość między najbliższymi atomami 2,8·10−10 m • Liczba atomów w krysztale 1 mm × 1 mm × 1 mm (1 mm / 2,8·10−10 m)3 ≈ 4,5·10+19 • Liczba atomów na powierzchni ściany 1 mm × 1 mm (1 mm / 2,8·10−10 m)2 ≈ 1,3·10+13 Typowy komputer osobisty: • Szybkość procesora (na 1 rdzeń) 3·109 cykli/s → ok. 105 zjawisk/s ≈ 9·10+9 zjawisk/dzień • Pamięć operacyjna 8 GB ≈ 8·10+9 B W celu umożliwienia przeprowadzenia pierwszych symulacji na przełomie lat 60-tych i 70-tych XX wieku przyjęto szereg założeń upraszczających przebieg symulacji [1,2]: 1. Identyczne jednostki-molekuły fazy stałej mogą zajmować miejsca tylko w węzłach idealnej prymitywnej sieci regularnej lub tetragonalnej (model blokowy). 2. Każda jednostka przestrzeni należy do fazy stałej (ang. solid) albo fazy ciekłej/gazowej (ang. fluid), tzn. nie ma miejsc pustych ani miejsc pośredniego typu. 3. Rozpatrywana jest powierzchnia tylko jednej ściany kryształu, zwykle o wskaźnikach (001). 4. Pod każdym blokiem stałym znajduje się blok stały. Założenie to znane jest jako SOS (ang. „solid on solid”). W konsekwencji wypełnienie dwóch kolejnych warstw blokami stałymi musi spełniać nierówność: nn CC ≤+1 , (5.3) gdzie Cn = Nns /N; Nns – liczba bloków stałych w n-tej warstwie, N – stała liczba wszystkich bloków w jednej warstwie. Ponadto izolowany blok całkowicie otoczony cieczą nie może należeć do kryształu; nie są też możliwe nawisy bloków stałych. generator szumu białego przetwornik analogowo-cyfrowy liczby losowe … bity najmniej znaczące 43 Według Binsbergena [3] zachodzenie zjawisk w roztworze jest utrudnione z powodu występowania lepkości, z którą związana jest energia swobodna ∆Gη aktywacji przepływu lepkiego lub rotacji bloku (rys. 5.6). Całkowita energia aktywacji dla przejścia od pewnego stanu o energii swobodnej ∆Gj−1 do stanu o energii ∆Gj wynosi ( )12 1 −η ∆−∆+∆= jja GGGE , (5.10) gdzie energie ∆Gj−1 oraz ∆Gj wynikają z rozważań nad wiązaniami zrywanymi oraz tworzonymi (rys. 5.7). Rys. 5.6. Energia swobodna i-tej konfiguracji w przypadku wzrostu kryształu z roztworu [3]. Częstotliwości kreacji i anihilacji w miejscu przy granicy faz otoczonym przez i sąsiadów bocznych można wyrazić wzorami (wzory pochodzące z [3] w zmienionym zapisie wg. [1]): ( ) β+−α−⋅=+ 2 2 4 expt ifki , (5.11) ( ) β−−α⋅=− 2 2 4 expt ifki . (5.12) Wyraz ft można interpretować jako liczbę prób zajścia zjawiska w jednostce czasu związaną z energią swobodną aktywacji przepływu lepkiego ∆Gη ∆ −= η kT G h kT f expt . (5.13) Φsf Φsf Φff Φss 1) Przed kreacją 2) Po kreacji zerwanie 10 wiązań stałe-ciekłe 5 nowych wiązań ciekłe-ciekłe 5 nowych wiązań stałe-stale energia aktywacji przepływu lepkiego ∆Gη Rys. 5.7. Bariera energetyczna dla kreacji przy wzroście kryształu z roztworu. 44 5.5. Dyfuzja powierzchniowa W modelu Gilmera i Bennemy dyfuzja powierzchniowa jest realizowana jako złożenie anihilacji z kreacją na jednej z 4-ech sąsiednich pozycji. Częstotliwość dyfuzji kij z miejsca otoczonego przez i sąsiadów do miejsca mającego j sąsiadów bocznych można więc wyrazić [1] + jiij kkrk −= , (5.14) gdzie r jest stałą o wymiarze czasu, którą można obliczyć wykorzystując równanie dyfuzji Einsteina τ = 2a Ds , (5.15) gdzie: a – stała sieci krystalicznej, Ds – stała dyfuzji powierzchniowej, τ = 1/k00 jest średnim czasem miedzy przeskokami dyfuzji powierzchniowej po gładkiej powierzchni. Średnia droga dyfuzji powierzchniowej Xs wynosi css DX τ= , (5.16) gdzie τc = 1/k0 − jest średnim czasem życia jednostki na gładkiej powierzchni. Po podstawieniu i = j = 0 do (5.14) i wykorzystaniu związków (5.15) i (5.16) otrzymujemy 2 000 00 1 == ++− a X kkk k r s (5.17) Związki (5.14)-(5.17) obowiązują dla dowolnych wzorów opisujących częstotliwości − ik oraz + jk , w tym także dla wzorów (5.11) i (5.12) dotyczących wzrostu z roztworu. Można więc uzupełnić model Binsbergena [3], który w swej oryginalnej postaci nie uwzględniał zjawiska dyfuzji powierzchniowej. 5.6. Algorytmy symulacji Spotykane w literaturze algorytmy komputerowych symulacji wzrostu kryształów metodą Monte Carlo można podzielić na trzy grupy: 1) Algorytmy bez pustych losowań, w którym decyzje o zdarzeniach wybranych do realizacji podejmowane są przez porównania liczb (pseudo)losowych z przedziałami o szerokościach zmieniających się odpowiednio do zmian stanu powierzchni kryształu [4]. 2) Algorytmy z pustymi losowaniami, w którym wybór zdarzenia odbywa się przez porównanie wylosowanych liczb z przedziałami o ustalonych szerokościach, natomiast dostosowanie rozkładu prawdopodobieństwa realizowanych zdarzeń do bieżącego stanu powierzchni kryształu uzyskuje się dzięki dodatkowej decyzji czy wybrane już zjawisko zostanie zrealizowane [1]. 3) Algorytm „waiting list” – zdarzenia wylosowane dla wszystkich możliwych miejsc na powierzchni kryształu umieszcza się na liście uporządkowanej wg. wylosowanych czasów ich zajścia. Po realizacji zdarzenia pierwszego jest ono usuwane z listy i losowane ponownie dla późniejszej chwili czasu. Ponadto wszystkie inne zdarzenia zależne muszą być usunięte z listy i wylosowane ponownie z uwzględnienie zmian stanu powierzchni. Algorytm „waiting list” nie będzie szerzej omawiany w niniejszym skrypcie. W algorytmie bez pustych losowań (rys. 5.8) wykorzystuje się tablicę przechowującą częstotliwości k1 … kn wszelkich możliwych procesów dla wszystkich dozwolonych miejsc. Przykładowo w modelu Gilmera i Bennemy w wybranym miejscu na powierzchni kryształu może zajść sześć różnych procesów: kreacja, anihilacja oraz dyfuzje powierzchniowe w czterech kierunkach. Tak więc dla kryształu o wymiarach 100 × 100 jednostek wzrostu konieczne jest przechowywanie 60000 częstotliwości. Po zrealizowaniu jednego procesu trzeba zaktualizować 45 częstotliwości ki dla kilkudziesięciu procesów związanych z danym miejscem i miejscami sąsiednimi oraz zaktualizować sumę K. Uzyskanie dużej wydajności programu wymaga starannego wyboru takiej struktury danych, która gwarantuje optymalny kompromis między szybkością aktualizacji i szybkością wyszukiwania danych. Wybór zdarzenia 0 K k1 k2 … … kn start RND km … Wykonaj wybrane zjawisko Aktualizuj częstotliwości ki Rys. 5.8. Ogólny algorytm symulacji metodą Monte Carlo bez pustych losowań. W praktyce podejmowanie decyzji (oznaczone tutaj jako jeden krok) może być rozdzielone na kilka etapów. W algorytmie z pustymi losowaniami (rys. 5.9) o wyborze zjawiska dyfuzji, kreacji albo anihilacji decyduje porównanie liczby losowej z przedziału 0…IYmax z przedziałami o stałych szerokościach [1]: max maxmaxmax max D )()()( )( IY kkk k L ijii ij ++ = −+ , (5.18) max maxmaxmax max K )()()( )( IY kkk k L ijii i ++ = −+ + , (5.19) max maxmaxmax max A )()()( )( IY kkk k L ijii i ++ = −+ − . (5.20) Na kolejnym etapie symulacji zapada decyzja czy wylosowane zdarzenie może zajść. W odpowiednim przedziale LD, LK albo LA wydzielany jest podprzedział o długości LD,ij , LA,i albo LK,i albo proporcjonalnej do częstotliwości zachodzenia zjawiska w wylosowanym miejscu D ij ij ijD L k k L max , )( = ; K i i iK L k k L max , )( + + = ; A i i iA L k k L max , )( − − = . (5.21) Jeśli liczba losowa trafi w ten podprzedział, to powierzchnia symulowanego kryształu jest odpowiednio modyfikowana; jeśli nie - od razu rozpoczyna się nowy cykl symulacji. Ze względu na prostą implementację tego algorytmu jego wydajność może być konkurencyjna przy niezbyt dużych wartościach parametru α. Liczba wygenerowanych zjawisk w porównaniu do liczby losowań maleje jednak szybko ze wzrostem wartości α. 48 Rys. 5.10. Wzrost ściany (001) z roztworu dla α = 1,5; β = 0,10; Xs = 1a; tablica 50 × 50. Rys. 5.11. Wzrost ściany (001) z roztworu dla α = 3,2; β = 0,10; Xs = 1a; tablica 50 × 50. Rys. 5.12. Wzrost ściany (001) z roztworu dla α = 5,0; β = 0,35; Xs = 1a; tablica 50 × 50. 49 0,00 0,02 0,04 0,06 0,08 0,10 0 30 60 90 120 150 180 210 240 długość boku tablicy (r '- s zy bk oś ć w zr os tu ) 1 /2 β = 0.42 β = 0.50 β = 0.38 Rys. 5.13. Wpływ długości boku tablicy opisującej ścianę (001) kryształu na otrzymaną bezwymiarową szybkość wzrostu ściany r' = r/(d·ft) dla α = 6,0; Xs = 1a; wzrost z roztworu; powierzchnia bez defektów [4]. Symulacje MC dają możliwość uwzględnienia defektów struktury kryształu poprzez specjalny sposób liczenia bocznych sąsiadów bloków w wybranych lokalizacjach z dodatkowym uskokiem o ±1 w kierunku normalnym do ściany kryształu. Wprowadzenie takiego uskoku tylko na wybranym odcinku na powierzchni ściany kryształu odpowiada zespołowi dwóch dyslokacji spiralnych połączonych prostoliniowym stopniem, który nazywany jest źródłem Franka-Reada. Przy odpowiednio dużych wartościach parametru α wokół dyslokacji narastają spiralne stopnie (rys. 5.14). Możliwe jest znalezienie takich kombinacji parametrów α i β, przy których znaczenie ma więcej niż jeden mechanizm wzrostu, np. wzrost spirali wokół dyslokacji oraz dwuwymiarowe zarodkowanie (rys. 5.15). Mikrostruktura powierzchni przewidywana przez symulacje MC znajduje potwierdzenie w obserwacjach kryształów pod mikroskopem sił atomowych (rys. 5.16). a) α = 7,0; β = 0,5; Xs = 0, b) α = 9,0; β = 0,9; Xs = 0. Rys. 5.14. Wzrost ściany (001) z roztworu, tablica 250 × 250, powierzchnia ze źródłem Franka-Reada. Wzrost ściany odbywa się przez wzrost spirali wokół dyslokacji. wzrost ciągły wzrost warstwowy 50 a) α = 5,5; β = 0,3; Xs = 0, b) α = 7,0; β = 1,0; Xs = 0. Rys. 5.15. Wzrost ściany (001) z roztworu, tablica 250 × 250, powierzchnia ze źródłem Franka-Reada. Do wzrostu ściany przyczynia się zarówno wzrost spirali wokół dyslokacji, jak i dwuwymiarowe zarodkowanie. a) Spirala wokół dyslokacji oraz zarodki 2D. b) Wzrost przez dwuwymiarowe zarodkowanie. Rys. 5.16. Ściana rosnącego kryształu pod mikroskopem AFM (atomic force microscope – mikroskop sił atomowych). Źródło: K. Maiwa, M. Plomp, W.J.P. van Enckevort, P. Bennema, „AFM obserwation of barium nitrate {111} and {100} faces: spiral growth and two-dimensional nucleation growth”, J. Crystal Growth, vol. 186 (1998), str. 214-223. W przypadku ciągłego mechanizmu wzrostu ściany kryształu (niskie wartości α) szybkość wzrostu ściany rośnie niemal liniowo ze wzrostem parametru β (rys. 5.17.a). Związek parametru β z przesyceniem roztworu lub pary podano wcześniej we wzorze (3.4). Jeżeli na powierzchni rosnącego kryształu nie występują żadne defekty sieci krystalicznej, to dla odpowiednio dużych wartości α możliwy jest tylko wzrost przez dwuwymiarowe zarodkowanie, przy czym wzrost ten obserwowany jest dopiero po przekroczeniu pewnej granicznej wartości przesycenia (rys. 5.17.b). Obecność stopni na powierzchni ściany kryształu, np. spiralnych stopni narastających wokół dyslokacji śrubowej, umożliwia wzrost ścian kryształu także przy małych przesyceniach. Wzrost przez przemieszczanie się stopni odpowiada liniowej zależności szybkości wzrostu r' od β na rys. 5.18. Przy dużych wartościach β zależność ta staje się silnie nieliniowa, co oznacza, że dominującym mechanizmem wzrostu staje się dwuwymiarowe zarodkowanie.