2.4.3.9.1 Mathematik der Diskretisierung des Wachstums von Keimen
Die Tücken der Diskretisierung als Voraussetzung für die numerische Integration müssen bekannt sein, um sich über ihre Folgen nicht wundern zu müssen.
Da wir es beim Wachstum eines Objekts in einer Umgebung, die über einen Transportprozess das Material liefern soll, immer mit einem Erhaltungssatz ebendiesen Materials (als Masse) zu tun haben, soll der auch am Anfang stehen:
Transport-Typ: Diffusion in der flüssigen Lösung hin zum festen Keim
ms + ml = mg (feste (solid) + flüssige (liquid) = gesamte)
Die Konzentration der Lösung sei mit c bezzeichnet und stellt ein Verhältnis der gelösten Masse zur Gesamtmasse in einem gegebenen diskreten Volumen der Nummer i dar:
ci = mi/m
Ist erstens diese Konzentration klein gegen 1 und zweitens der Dichteverhältnis von Lösungsmittel und gelöstem Stoff nicht klein oder groß gegen 1, so braucht man bei der Anlagerung an den Keim keine Volumenänderung des Gesamtsystems anzunehmen und kann eine konstante Volumen-Diskretisierung im geometrischen Raum vornehmen.
Der analytische Ausdruck für den Umstand, dass die Diffusionsstromdichte j proportional dem Gradienten der Konzentration c sein soll, wird im diskretisierten Fall (im isotropen Fall radial) zu zwei skalaren Gleichungen (Strom senkrecht zur Oberfläche des Volumenelements i, positiv nach innen ins Volumenelement hinein, negativ nach außen hinaus, also unabhängig von der Richtung des Radiusvektors)
ji+ = D*(c(i+1)-c(i))
ji- = D*(c(i-1)-c(i))
Für radiale Volumenelemente gleicher „Dicke“ dr (im 2D-Fall Kreisringe, im 3D-Fall Kugelschalen [im 1D-Fall wären es ebene Schichten konstanter Dicke]) gilt dann, dass ihre Grenzflächen zu den beiden Nachbarn unterschiedlich groß sind, dass also selbst im Falle sich aufhebender Stromdichten j+ und j- ein Masse- und somit Konzentrations -Änderung auftreten würde, welche sich außerdem noch von den Konzentrationsänderungen der beiden Nachbar-Elemente unterscheiden würde, da deren Volumina nicht gleich sind!
Bezeichnet man mit ri den mittleren Radius des i-ten Volumenelements, so ergibt sich, wenn man innen mit i=1 beginnt,
ri = (i-0,5)*dr
und für die Volumina Vi sowie die Grenzflächen Ai+ und Ai- und die sich daraus ergebenden neuen werte für Massen und Konzentrationen
2D
Vi = pi*((ri+0,5dr)²-(ri-0,5dr)²) = pi*(ri²+ri*dr+0,25dr²-ri²+ri*dr-0,25dr²) = 2pi*ri*dr
Ai+ = 2pi*(ri+0,5dr)
Ai- = 2pi*(ri-0,5dr)
dmi = ((ji+ *Ai+) + (ji- *Ai-))*dt (auch die Zeit ist in Takte dt diskretisiert!!)
mi(neu) = mi(alt)+dmi
ci(neu) = mi(neu)/Vi
3D
Vi = 4/3*pi*((ri+0,5dr)³-(ri-0,5dr)³) = 4/3*pi*(3ri²*dr+0,25dr³) = 4pi*ri²*dr
(dr³ ist gegen ri²*dr vernachlässigt worden; der Keim sollte also wenigstens bis i=2 beim Programmstart schon fest sein)
Ai+ = 4*pi*(ri+0,5dr)²
Ai- = 4*pi*(ri-0,5dr)²
dmi = ((ji+ *Ai+) + (ji- *Ai-))*dt
mi(neu) = mi(alt)+dmi
ci(neu) = mi(neu)/Vi
Aus den hier aufgezeigten Zusammenhängen kann man unmittelbar schließen, dass für einen „stationären“ Fall zeitlich konstanten Massezuwachses am zentralen Festkörper durch die Radiusabhängigkeit der Volumenelement-Oberfläche die Stromdichte folgende Radiusabhängigkeit haben muss:
1D: konstant
2D: prop. 1/r
3D: prop. 1/r²
Da die dazugehörige Konzentration das Weg-Integral über ihren Gradienten ist, der seinerseits der Stromdicht porortional ist, ergibt sich folgender Konzentrationsverlauf:
1D: linear
2D: logarithmisch prop. ln(r)
3D: hyperbolisch prop. 1/r
Daraus folgt unmittelbar, dass nur im 3D-Fall eine räumlich quasi-unendliche (unerschöpfliche) Lösung einen stationären Fall erzeugen kann.
Es ist also stets zu erwarten, dass je nach tatsächlicher Größe von Volumen, Konzentration (Übersättigung) , Diffusionskoeffizient und der Kinetik des Phasenübergangs zeitabhängige Konzentrationsfelder und damit zeitabhängige Wachstumsgeschwindigkeiten auftreten müssen. Der „stationäre Fall“ ist nur für kurze Zeit überhaupt näherungsweise möglich!
(Man kann natürlich am äußeren Rand des berechneten Diffusionsvolumens eine Quelle ansetzen, deren Stärke etwa der Senkenstärke entspricht oder die aus vielen kleinen Keimen gebildet wird, die sich konzentrationsabhängig auflösen… Das wäre eine Idee für spätere, weiterführende Modellierungen!)
Für die Diskretisierung des Prozesses gilt es noch die Berücksichtigung der beiden Ränder zu beachten:
Das Teilvolumen Vi, dass die Phasengrenze enthält, ist der „innere“ Rand. Vi gibt keine Masse weiter nach innen ab, eine der obigen Gleichungen enfällt also. Die Aufnahmefähigkeit dieser Schicht unterscheidet sich von den anderen weiterhn dadurch, dass sie durch die radiusabhängige „Übersättigung“ dar äußeren Nachbarschicht gesteuert wird. Ist durch den Masseuwachs die Konzentration der fetsen Phase erreicht, wird diese Eigenschaft an die nächstäußere Schicht weitergegeben. Hier nuss im Algorithmus also eine Abfrage stehen.
Das Teilvolumen Vi, das die äußerste Schicht bildet, hat ebenfalls nur einen Strom, und zwar den nach innen. Diese Schicht wird erst nach soviel Zyklen überhaupt beeinflusst, wie die Schichtanzahl in der Lösung beträgt. Dann reicht der Diffusionstrichter der Konzentration bis an den Rand des Gesamtvolumens. Man erkennt, dass der gesamte Prozess der Masse-Stroms ab jetzt verlangsamt wird, bis er zum Stillstand kommt, weil an der Oberfläche des wachsenden Festkörpers eine Gleichheit von tatsächlich anliegender Konzentration und radiusabhängiger Gleichgewichts-Konzentration entstanden ist.
Da die Gleichgewichtskonzentration mit zunehmendem Radius sinkt, wird die Verringerung der Wachstumsgeschwindigkeit durch zeitlich sinkende Lösungs-Konzentration teilweise aufgehoben.
Genau dieser Effekt führt beim Vorhandensein mehrerer Festkörper in der Lösung zur Bevorzugung der Großen gegenüber den Kleinen.
Ein passabler Ansatz für diese Abhängigkeit sollte sein, dass man eine Gleichgewichtskonzentration coo für die Ebene als Konstante nimmt und darauf eine Hyperbel setzt, die für einen charakteristischen Radius r=r2 zur Verdoppelung der Ebenen-GG-Konzentration führt:
co(r) = coo * (1+r2/r)
Man erkennt an der roten Kurve, ab welchem Radius ein Keim für zwei verschiedene Umgebungskonzentrationen (als konstant angesehen, also mit Konvektion und Massenachschub gedacht) überhaupt wachsen kann (Wert der Übersättigung muss über 0 sein). Da diese Konzentration im Diffusionsfeld aber aus der Rückkopplung über die Massenaufnahme des Keims beeinflusst wird, wird die Modellierung spannend.
Man darf dann in erster Näherung von einer Proportionalität zwischen Übersättigung und Anlagerungs-Stromdichte ausgehen (also schon in räumlichen und zeitlichen Dimensionen gedacht, ganz entsprechend der Diskretisierung).
(Hinweis: Für polyedrische Kristalle wird es wesentlich komplexer, weil es unterschiedliche Bindungsenergien je nach Anzahl und Orientierung der bindenden Gitternachbarn gibt. Auch der Raum der Lösungsdiffusion muss in Würfelelemente zerlegt werden mit entsprechender Summierung der Austausch-Prozesse… Es entsteht außerdem das Problem der morphologischen Stabilität…)
Fazit: Die doppelte Diskretisierung des Raumes und der Zeit muss so fein gestaltet werden, dass bei der numerischen Integration die Fehler möglichst klein bleiben. Die Existenz der Volumen-Ränder innen und außen muss durch eine entsprechende Programm-Struktur als Spezialfälle in den Schleifen berücksichtigt werden. Dabei ist zu beachten: Je feiner der Raum diskretisiert wird, desto mehr Zeit-Takte sind erforderlich, um lokale Änderungen vermittels der Diffusion zum anderen Ende „durchzureichen“. Es gibt also ein von der Rechenleistung abhängiges Optimum für die Feinheit der Diskretisierung und die Dauer der Rechenzeit pro Parametersatz.
Kommentar abgeben