Numerische Simulation in Python-Invasive Percoration
Einführung
Python kann wirklich alles, daher gibt es viele Möglichkeiten, es zu verwenden, aber ich möchte vorstellen, wie ich es verwende. Das Paket, das ich oft benutze, ist numpy. Es wird von einer numerischen Berechnung begleitet. Ohne das bin ich wirklich in Schwierigkeiten. Der Rest ist matplotlib und scipy. Dann benutze ich Tkinter, um den Dialog anzuzeigen. Ich denke, dass es ziemlich viele Leute gibt, die diesen Bereich nutzen, daher denke ich, dass es einfach ist, Informationen zu sammeln.
Abschließend möchte ich einen invasiven Perkorationscluster erstellen und ein ganzes Programm zur Berechnung seiner fraktalen Dimension erstellen. Standardmäßig sind numpy, scipy und matplotlib nicht enthalten, daher müssen Sie sie zuerst einschließen.
Als ich es eilig hatte, im Labor in einen Laptop zu steigen
windows7 - Python, SciPy, Installieren von matplotlib (Windows) --Qiita:
Es war wahnsinnig praktisch. Danke für deine Hilfe.
Wenn Sie es in Ubuntu setzen, sollten Sie apt gehorsam verwenden.
Für das Modell des folgenden Programms verzeihen Sie mir bitte, indem Sie meinen Bericht kopieren und einfügen ...
Über invasive Perkoration
Ein dynamischer Prozess, der als invasive Perkoration bekannt ist, wird verwendet, um die Form der Öl-Wasser-Grenzfläche zu modellieren, die auftritt, wenn Wasser in ein poröses Medium gedrückt wird, das Öl enthält. Die Idee ist, mit Wasser so viel Öl wie möglich zu gewinnen. Bei diesem Prozess wachsen Wassercluster über den am wenigsten resistenten Weg in das Öl hinein. Stellen Sie sich ein L × 2L-Gitter vor und nehmen Sie an, dass der linke Rand anfänglich von Wasser (Eindringlingen) besetzt ist. Die Stärke des Widerstands gegen den Eindringling wird jedem Punkt im Gitter mit einer einheitlichen Zufallszahl zwischen 0 und 1 zugewiesen und bleibt während des Invasionsprozesses fest. Die den Gitterpunkten aller Angreifer am nächsten gelegenen Punkte sind die umgebenden Punkte. Zu jeder Teilungszeit wird der umgebende Punkt mit der kleinsten Zufallszahl vom Eindringling besetzt, und das Öl (Verteidiger) an diesem Punkt wird ersetzt. Der Eindringlingscluster wächst, bis ein Pfad gebildet wird, der das linke und das rechte Ende des Gitters verbindet. Periodische Randbedingungen werden für die Ober- und Unterkante verwendet, um den Effekt der Grenze zu minimieren, und alle Größen werden nur im L × L-Bereich der Mitte des Gitters gemessen. Die wichtigsten interessierenden Größen sind der Prozentsatz der von den Invasoren belegten Gitterpunkte und die Wahrscheinlichkeit P (r), dass Gitterpunkte mit zufälligen Werten zwischen r + dr besetzt sind.
Nachdem ich nun das Bild des Modells erklärt habe, möchte ich auf den Inhalt des Codes eingehen.
Codebeschreibung
Zunächst von den Grundlagen
Die Gitterpunkte werden wie eine Matrix unter Verwendung des Numpy-Arrays dargestellt (wobei es als self.lattice geschrieben wird).
numpy.random.random([x, y])
Erstellt ein Array von Zufallszahlen mit den Elementen 0 bis 1 in einer Matrix aus x Zeilen und y Spalten. In diesem Modell wird angenommen, dass der Wert der Zufallszahl, die jedem gegeben wird, den Widerstand darstellt. Als nächstes werden wir, wie wir es in der while-Klausel tun, wie im Kommentar geschrieben, die Gitterpunkte mit dem minimalen Widerstandswert in der Reihenfolge von den Punkten um die besetzten Gitterpunkte besetzen.
Wo ich es mir ausgedacht habe
#Stellen Sie die Koordinaten des Gitters mit dem kleinsten Wert an den umgebenden Punkten auf mm ein
mm = min([(v,k) for k,v in nnsite.items()])[1]
Der Schlüssel und der Wert des Wörterbuchs nnsite werden von nnsite.items () herausgenommen, und der Schlüssel mit dem Minimalwert wird durch Vergleichen der Größe des ersten Werts des übergebenen Taples erhalten, der die Spezifikation von min ist.
class TopWindow:
def quit(self):
self.root.destroy()
sys.exit()
def show_window(self, title="title", *args):
self.root = Tk()
self.root.title(title)
frames = []
for i, arg in enumerate(args):
frames.append(Frame(self.root, padx=5, pady=5))
for k, v in arg:
Button(frames[i],text=k,command=v).pack(expand=YES, fill='x')
frames[i].pack(fill='x')
f = Frame(self.root, padx=5, pady=5)
Button(f,text='quit',command=self.quit).pack(expand=YES, fill='x')
f.pack(fill='x')
self.root.mainloop()
Geben Sie die Funktion an, die der Schaltflächenbeschriftung in dem Satz von Taples zugeordnet ist, die dem Argument zugewiesen wurden. Beim Anruf
run = (('run', pushed), ('save canvas to sample.eps', pr))
run2 = (('calculate D', b4_pushed),)
top.show_window("Invasion Percolation", run, run2)
Machen.
Ich kann nicht alle anderen Details erklären, daher werde ich vorstellen, welche Art von Bildschirm als Ergebnis der Simulation erhalten wird und welche Art von Ausführungsergebnis erhalten wird.
Auf diese Weise wird ein Dialogfeld angezeigt. Klicken Sie auf Ausführen, um einen Percoration-Cluster zu generieren, und berechne_D gibt Ihnen das obige Diagramm. Die Bedeutung wird weggelassen.
invasion_percolation.py
#! /usr/bin/env python
# -*- coding:utf-8 -*-
#
# written by ssh0, August 2014.
from Tkinter import *
import sys
import numpy as np
import scipy.optimize as optimize
import matplotlib.pyplot as plt
class Percolation:
def __init__(self, Lx=40, Ly=20):
self.sub = None
self.Lx = Lx
self.Ly = Ly
def perc_cluster(self):
self.lattice = np.random.random([self.Lx+2, self.Ly])
Lx = self.Lx
Ly = self.Ly
#Das linke Ende sind alle besetzten Standorte
self.lattice[:1, :] = 1.5
self.lattice[Lx+1:, :] = 0
if self.sub is None or not self.sub.winfo_exists():
lattice = self.lattice
ne = [(0, -1), (0, 1), (-1, 0), (1, 0)]
#Halten Sie die Koordinaten der umgebenden Punkte im Wörterbuchformat
nnsite = {(1, y):lattice[1, y] for y in range(Ly)}
percolate = False
while len(nnsite) != 0 and percolate == False:
#Stellen Sie die Koordinaten des Gitters mit dem kleinsten Wert an den umgebenden Punkten auf mm ein
mm = min([(v,k) for k,v in nnsite.items()])[1]
lattice[mm] = 1
del nnsite[mm]
#Listen Sie die Koordinaten von Punkten um mm in nn auf(Wenden Sie eine periodische Randbedingung in y-Richtung an)
nn = [(mm[0] + nx, (mm[1] + ny)%Ly) for nx, ny in ne
if lattice[mm[0] + nx, (mm[1] + ny)%Ly] < 1]
#Schließt nn aus, die bereits in nnsite enthalten sind--> nnlist
nnlist = list(set(nn) - set(nnsite.keys()))
#Fügen Sie nnsite neue periphere Punkte hinzu
for n in nnlist:
nnsite[n] = lattice[n]
#Wenn sich der belegte Gitterpunkt am rechten Ende befindet, wird er als versickert betrachtet.
if mm[0] == Lx:
percolate = True
self.lattice = lattice[1:-1, :]
return self.lattice
def draw_canvas(self, rect, Lx, Ly):
default_size = 640 # default size of canvas
r = int(default_size/(2*Lx))
fig_size_x = 2*r*Lx
fig_size_y = 2*r*Ly
margin = 10
sub = Toplevel()
sub.title('invasion percolation')
self.canvas = Canvas(sub, width=fig_size_x+2*margin,
height=fig_size_y+2*margin)
self.canvas.create_rectangle(margin, margin,
fig_size_x+margin, fig_size_y+margin,
outline='black', fill='white')
self.canvas.pack()
c = self.canvas.create_rectangle
site = np.where(rect==1)
for m, n in zip(site[0], site[1]):
c(2*m*r+margin, 2*n*r+margin,
2*(m+1)*r+margin, 2*(n+1)*r+margin,
outline='black', fill='black')
def get_fractal_dim(self, trial=20, Lmin=20, Lmax=40, Lsample=10):
#Ein ganzzahliger Wert zwischen Lmin und Lmax, damit der Abstand gleichmäßig ist, wenn so viel wie möglich protokolliert wird.
L = np.array([int(i) for i
in np.logspace(np.log10(Lmin), np.log10(Lmax), Lsample)])
M_L = []
for l in L:
self.Lx = l*2
self.Ly = l
m_L = 0
for i in range(trial):
lattice = self.perc_cluster()
#Gesamtzahl der belegten Standorte im zentralen L × L-Raster
m_L += np.sum(lattice[int(l/2)+1:l+int(l/2),:]==1)
M_L.append(m_L/float(trial))
print "L = %d, M_L = %f" % (l, M_L[-1])
M_L = np.array(M_L)
def fit_func(parameter0, L, M_L):
log = np.log
c1 = parameter0[0]
c2 = parameter0[1]
residual = log(M_L) - c1 - c2*log(L)
return residual
parameter0 = [0.1, 2.0]
result = optimize.leastsq(fit_func, parameter0, args=(L, M_L))
c1 = result[0][0]
D = result[0][1]
print "D = %f" % D
def fitted(L, c1, D):
return np.exp(c1)*(L**D)
fig = plt.figure("Fractal Dimension")
ax = fig.add_subplot(111)
ax.plot(L, M_L, '-o', label=r"$M(L)$")
ax.plot(L, fitted(L, c1, D), label="fit func: D = %f" % D)
ax.set_xlabel(r'$\ln L$', fontsize=16)
ax.set_ylabel(r'$\ln M(L)$', fontsize=16)
ax.set_xscale('log')
ax.set_yscale('log')
ax.set_ymargin(0.05)
fig.tight_layout()
plt.legend(loc='best')
plt.show()
class TopWindow:
def quit(self):
self.root.destroy()
sys.exit()
def show_window(self, title="title", *args):
self.root = Tk()
self.root.title(title)
frames = []
for i, arg in enumerate(args):
frames.append(Frame(self.root, padx=5, pady=5))
for k, v in arg:
Button(frames[i],text=k,command=v).pack(expand=YES, fill='x')
frames[i].pack(fill='x')
f = Frame(self.root, padx=5, pady=5)
Button(f,text='quit',command=self.quit).pack(expand=YES, fill='x')
f.pack(fill='x')
self.root.mainloop()
if __name__ == '__main__':
Lx = 40
Ly = 20
top = TopWindow()
per = Percolation(Lx, Ly)
count = 1
def pr():
global count
d = per.canvas.postscript(file="figure_%d.eps" % count)
print "saved the figure to a eps file"
count += 1
def pushed():
per.perc_cluster()
per.draw_canvas(per.lattice, Lx, Ly)
def b4_pushed():
trial = 100; Lmin = 20; Lmax = 100; Lsample = 10
per.get_fractal_dim(trial, Lmin, Lmax, Lsample)
run = (('run', pushed), ('save canvas to sample.eps', pr))
run2 = (('calculate D', b4_pushed),)
top.show_window("Invasion Percolation", run, run2)