निष्पादन समय कम करने के लिए शुद्ध NumPy में लूप के लिए एक रिवाइटिंग

Question

I recently asked about trying to optimise a Python loop for a scientific application, और मेरे लिए an excellent, smart way of recoding it within NumPy which reduced execution time by a factor of around 100 प्राप्त हुआ!

हालांकि, B मूल्य की गणना वास्तव में, कुछ अन्य छोरों के भीतर नेस्ट क्योंकि यह पदों की एक नियमित ग्रिड पर मूल्यांकन किया जाता है है। क्या इस प्रक्रिया को बंद करने के लिए एक समान स्मार्ट NumPy पुनर्लेख है?

मुझे संदेह है कि इस भाग के प्रदर्शन लाभ को कम चिह्नित किया जाएगा, और नुकसान संभवतः यह होगा कि गणना की प्रगति पर उपयोगकर्ता को वापस रिपोर्ट करना संभव नहीं होगा, परिणाम न लिखे जा सकते हैं गणना फ़ाइल के अंत तक आउटपुट फ़ाइल, और संभवतः एक बड़े कदम में ऐसा करने से स्मृति प्रभाव पड़ता है? क्या इनमें से किसी को बाधित करना संभव है?

import numpy as np 
import time 

def reshape_vector(v): 
    b = np.empty((3,1)) 
    for i in range(3): 
     b[i][0] = v[i] 
    return b 

def unit_vectors(r): 
    return r/np.sqrt((r*r).sum(0)) 

def calculate_dipole(mu, r_i, mom_i): 
    relative = mu - r_i 
    r_unit = unit_vectors(relative) 
    A = 1e-7 

    num = A*(3*np.sum(mom_i*r_unit, 0)*r_unit - mom_i) 
    den = np.sqrt(np.sum(relative*relative, 0))**3 
    B = np.sum(num/den, 1) 
    return B 

N = 20000 # number of dipoles 
r_i = np.random.random((3,N)) # positions of dipoles 
mom_i = np.random.random((3,N)) # moments of dipoles 
a = np.random.random((3,3)) # three basis vectors for this crystal 
n = [10,10,10] # points at which to evaluate sum 
gamma_mu = 135.5 # a constant 

t_start = time.clock() 
for i in range(n[0]): 
    r_frac_x = np.float(i)/np.float(n[0]) 
    r_test_x = r_frac_x * a[0] 
    for j in range(n[1]): 
     r_frac_y = np.float(j)/np.float(n[1]) 
     r_test_y = r_frac_y * a[1] 
     for k in range(n[2]): 
      r_frac_z = np.float(k)/np.float(n[2]) 
      r_test = r_test_x +r_test_y + r_frac_z * a[2] 
      r_test_fast = reshape_vector(r_test) 
      B = calculate_dipole(r_test_fast, r_i, mom_i) 
      omega = gamma_mu*np.sqrt(np.dot(B,B)) 
      # write r_test, B and omega to a file 
    frac_done = np.float(i+1)/(n[0]+1) 
    t_elapsed = (time.clock()-t_start) 
    t_remain = (1-frac_done)*t_elapsed/frac_done 
    print frac_done*100,'% done in',t_elapsed/60.,'minutes...approximately',t_remain/60.,'minutes remaining' 

Answer 1

एक स्पष्ट बात आप कर सकते हैं के साथ लाइन

r_test_fast = reshape_vector(r_test) 

की जगह है

r_test_fast = r_test.reshape((3,1)) 

शायद प्रदर्शन में कोई बड़ा अंतर नहीं होगा, लेकिन किसी भी मामले में यह समझ में आता है पहिया को पुनर्निर्मित करने के बजाय numpy buildins का उपयोग करने के लिए।

आम तौर पर, जैसा कि आपने शायद अब तक देखा है, numpy को अनुकूलित करने वाली चाल को अंजीर पूर्ण-सरणी संचालन की सहायता से या कम से कम स्लाइस के साथ पाइथन कोड में प्रत्येक तत्व पर पुनरावृत्ति के बजाय एल्गोरिदम व्यक्त करना है। इस तरह के "वेक्टरेशन" को रोकने के लिए क्या होता है, तथाकथित लूप-आधारित निर्भरता है, यानी लूप जहां प्रत्येक पुनरावृत्ति पिछले पुनरावृत्ति के परिणाम पर निर्भर होती है। अपने कोड पर संक्षेप में देख रहे हैं, आपके पास ऐसी कोई चीज़ नहीं है, और यह आपके कोड को ठीक से व्यवस्थित करना संभव होना चाहिए।

संपादित करें: एक समाधान

मैं सत्यापित नहीं किया है कि यह सही है, लेकिन आप इसे कैसे दृष्टिकोण की एक विचार देना चाहिए।

सबसे पहले, cartesian() function, which we'll use लें। तब

 

def calculate_dipole_vect(mus, r_i, mom_i): 
    # Treat each mu sequentially 
    Bs = [] 
    omega = [] 
    for mu in mus: 
     rel = mu - r_i 
     r_norm = np.sqrt((rel * rel).sum(1)) 
     r_unit = rel/r_norm[:, np.newaxis] 
     A = 1e-7 

     num = A*(3*np.sum(mom_i * r_unit, 0)*r_unit - mom_i) 
     den = r_norm ** 3 
     B = np.sum(num/den[:, np.newaxis], 0) 
     Bs.append(B) 
     omega.append(gamma_mu * np.sqrt(np.dot(B, B))) 
    return Bs, omega 


# Transpose to get more "natural" ordering with row-major numpy 
r_i = r_i.T 
mom_i = mom_i.T 

t_start = time.clock() 
r_frac = cartesian((np.arange(n[0])/float(n[0]), 
        np.arange(n[1])/float(n[1]), 
        np.arange(n[2])/float(n[2]))) 
r_test = np.dot(r_frac, a) 
B, omega = calculate_dipole_vect(r_test, r_i, mom_i) 

print 'Total time for vectorized: %f s' % (time.clock() - t_start) 
 

खैर, मेरी परीक्षण में, इस तथ्य पाश आधारित दृष्टिकोण मैं से शुरू की तुलना में थोड़ा धीमी है। बात यह है कि, प्रश्न में मूल संस्करण में, यह आकार के सरणी (20000, 3) पर पूरे सरणी संचालन के साथ पहले से ही सदिश हो गया था, इसलिए कोई भी आगे वेक्टरेशन वास्तव में और अधिक लाभ नहीं लाता है। वास्तव में, यह ऊपर के रूप में प्रदर्शन को खराब कर सकता है, शायद बड़े अस्थायी सरणी के कारण।

Answer 2

आप profile अपने कोड, तो आप उस चलने का समय के 99% तो यह वास्तव में पाशन निष्पादन समय में एक उल्लेखनीय कमी नहीं देंगे के लिए समय को कम करने calculate_dipole में है देखेंगे। यदि आप इसे तेज़ी से बनाना चाहते हैं तो आपको अभी भी calculate_dipole पर ध्यान केंद्रित करने की आवश्यकता है। मैंने इस पर calculate_dipole के लिए अपने साइथन कोड की कोशिश की और समग्र समय में 2 के कारक में कमी आई। साइथन कोड को बेहतर बनाने के अन्य तरीके भी हो सकते हैं।