2016-09-29 2 views
2

Für meine Computational Physics Klasse müssen wir die Madelung Constant für NaCl berechnen. Mein Code hierfür verwendet drei verschachtelte for-Schleifen und läuft daher sehr langsam. Ich habe mich gefragt, ob es eine Möglichkeit gibt, Arrays oder andere Methoden zu verwenden, um die Geschwindigkeit der Berechnung zu erhöhen. DankWie beschleunigt man die dreidimensionale Summe für die Madelung Constant?

from math import sqrt 
l= int(input("The lattice size is :")) 
M = 0.0 
for i in range(-L,L+1): 
    for j in range(-L,L+1): 
     for k in range(-L,L+1): 
      if not (i==j==k==0): 
       M += ((-1)**(i+j+k+1))/sqrt(i*i +j*j +k*k) 

print("The Madelung constant for NaCl with lattice size",l,"is",M) 
+0

nicht Ihren Computational Physics Klasse lehren Sie numpy? –

+0

Ja, aber ich verstehe nicht, wie man es für dieses Problem benutzt. Der Professor hat kaum diskutiert, Arrays zu benutzen, außer zu sagen, dass sie hier verwendet werden können. – Tom

Antwort

1

Da Sie in einem Kommentar darauf hingewiesen haben, dass Sie numpy verwenden kann, schlage ich vor, dass zu tun. Sie können ein 3D-Gitter für Ihre Ganzzahlen konstruieren und jeden Ausdruck gleichzeitig berechnen, um so Ihre Berechnung zu vektorisieren. Sie müssen nur für den Einzelfall achten, wo jede ganze Zahl 0, zum Beispiel ist numpy.where mit:

import numpy as np 

ran = np.arange(-L,L+1) 
i,j,k = np.meshgrid(ran,ran,ran) 
M = np.where((i!=0) | (j!=0) | (k!=0), 
      (-1)**(i+j+k+1)/np.sqrt(i**2+j**2+k**2), 
      0).sum() 

ran sind ein numpy Array mit den gleichen Elementen wie in range() (wenn gegossen in eine Liste in Python 3) . meshgrid konstruiert dann drei 3D-Arrays, die zusammen den 3D-Raum überspannen, in dem Sie Ihre Summe ausführen müssen.

Beachten Sie, dass für große Domänen dieser Ansatz viel höheren Speicherbedarf hat. Dies ist normalerweise bei der Vektorisierung der Fall: Sie können CPU-Zeit auf Kosten eines erhöhten Speicherbedarfs sparen.

+0

Vielen Dank. – Tom

+0

@Tom Ich schlage vor, dass Sie sich Divakars neue Antwort ansehen. Die Verwendung von 'Ogrid' sollte viel effizienter sein, da diese riesigen 3D-Arrays nicht im Speicher gespeichert sind. Seine Lösung, um den einzelnen singulären Punkt zu fixieren, ist auch viel einfacher und sicherlich effizienter. –

2

Hier ist ein Ansatz mit offenen Maschen mit np.ogrid statt tatsächlichen Maschen zu schaffen, die auf this other post beruhen muss ziemlich effizient sein -

# Create open meshes for efficiency purposes 
I,J,K = np.ogrid[-L:L+1,-L:L+1,-L:L+1] 

# Perform all computations using those open meshes in a vectorized manner 
all_vals = ((-1)**(I+J+K+1))/np.sqrt(I**2 +J**2+ K**2) 

# Corresponding to "not (i==j==k==0)", which would be satisfied by 
# just one combination and that will occur at location (L,L,L), 
# let's set it as zero as a means to sum later on without adding for that elem 
all_vals[L,L,L] = 0 
M_out = all_vals.sum() 
+0

Ich vergesse immer wieder, dass 'Ogrid' existiert ... sehr nette Antwort :) –

+1

@AndrasDeak Ja sogar ich habe es nicht viel benutzt. Aber diese offenen Netze schreien nur, weil sie 'broadcasting' nutzen wollen und deshalb sollte ich mehr darauf schauen;) – Divakar

Verwandte Themen