SIMD के साथ कॉलम-वार अधिकतम का अनुकूलन
मेरे पास यह फ़ंक्शन है जहां मैंने अपने कोड में महत्वपूर्ण समय बिताया है, और यदि संभव हो तो मैं वेक्टराइजेशन-सिमड-कंपाइलर इंट्रिंसिक्स द्वारा इसे अनुकूलित करना चाहूंगा।
यह अनिवार्य रूप से स्तंभों पर एक मैट्रिक्स पर मूल्य और अधिकतम का स्थान पाता है, और उन्हें संग्रहीत करता है:
- val_ptr: इनपुट मैट्रिक्स: स्तंभ-प्रमुख (फोरट्रान-शैली) n_rows-by-n_cols (जहां आमतौर पर n_rows >> n_cols)
- Opt_pos_ptr: लंबाई का int वेक्टर n_rows जहां अधिकतम की स्थिति संग्रहीत करने के लिए। शून्य से भरे प्रवेश पर।
- max_ptr: लंबाई n_rows का फ्लोट वेक्टर जहां अधिकतम स्टोर करने के लिए। Val_ptr के पहले कॉलम की प्रतियों से भरी प्रविष्टि पर
- फ़ंक्शन को समानांतर लूप में कहा जाएगा
- स्मृति क्षेत्र को अतिव्यापी नहीं होने की गारंटी है
- मुझे वास्तव में अधिकतम भरे जाने की आवश्यकता नहीं है, वर्तमान में इसका उपयोग केवल पुस्तक-रखने और स्मृति आवंटन से बचने के लिए किया जाता है
- मैं विंडोज 10. पर MSVC, C ++ 17 का उपयोग करता हूं। आधुनिक इंटेल सीपीयू चलाने के लिए
कोड, जहां टेम्पलेट का प्रकार फ्लोट या डबल होना है:
template <typename eT>
find_max(const int n_cols,
const int n_rows,
const eT* val_ptr,
int* opt_pos_ptr,
eT* max_ptr){
for (int col = 1; col < n_cols; ++col)
{
//Getting the pointer to the beginning of the column
const auto* value_col = val_ptr + col * n_rows;
//Looping over the rows
for (int row = 0; row < n_rows; ++row)
{
//If the value is larger than the current maximum, we replace and we store its positions
if (value_col[row] > max_ptr[row])
{
max_ptr[row] = value_col[row];
opt_pos_ptr[row] = col;
}
}
}
}
मैंने अब तक क्या कोशिश की:
- मैंने इनर लूप के लिए ओपनएमपी समानांतर का उपयोग करने की कोशिश की, लेकिन केवल बहुत बड़ी पंक्तियों पर कुछ लाता है, जो मेरे वर्तमान उपयोग से थोड़ा बड़ा है।
- अगर आंतरिक लूप में अगर #pragma omp simd को काम करने से रोकता है, और मैं इसके बिना इसे फिर से लिखने में सक्षम नहीं था।
जवाब
आपके द्वारा पोस्ट किए गए कोड नमूने के आधार पर, ऐसा लगता है कि आप एक ऊर्ध्वाधर अधिकतम मूल्य की गणना करना चाहते हैं, जिसका अर्थ है कि आपके मामले में "कॉलम" क्षैतिज हैं। C / C ++ में तत्वों के क्षैतिज अनुक्रम (अर्थात जहां दो आसन्न तत्वों में स्मृति में एक तत्व की दूरी होती है) को सामान्य रूप से पंक्तियों और ऊर्ध्वाधर कहा जाता है (जहां दो आसन्न तत्वों में स्मृति में पंक्ति आकार की दूरी होती है) - कॉलम। नीचे दिए गए मेरे जवाब में मैं पारंपरिक शब्दावली का उपयोग करूंगा, जहां पंक्तियाँ क्षैतिज हैं और स्तंभ लंबवत हैं।
इसके अलावा, संक्षिप्तता के लिए मैं एक संभावित प्रकार के मैट्रिक्स तत्व पर ध्यान केंद्रित करूंगा - float। मूल विचार के लिए समान है double, मुख्य अंतर प्रति वेक्टर तत्वों की संख्या और _ps/ _pdआंतरिक चयन है। मैं doubleअंत में एक संस्करण प्रदान करूंगा ।
विचार यह है कि आप समानांतर स्तंभों के लिए लंबवत अधिकतम का उपयोग कर सकते हैं समानांतर _mm_max_ps/ _mm_max_pd। पाए गए अधिकतम की स्थिति को रिकॉर्ड करने के लिए, आप वर्तमान तत्वों के साथ पिछले अधिकतम की तुलना कर सकते हैं। तुलना का परिणाम एक मुखौटा है, जहां तत्व सभी हैं जहां अधिकतम अपडेट किया जाता है। उस मास्क का उपयोग यह चुनने के लिए किया जा सकता है कि किस स्थिति को अद्यतन करने की आवश्यकता है।
मुझे ध्यान देना चाहिए कि नीचे एल्गोरिथ्म मानता है कि यह महत्वपूर्ण नहीं है कि किसी कॉलम में एकाधिक समान अधिकतम तत्व हैं, तो अधिकतम तत्व की स्थिति दर्ज की गई है। इसके अलावा, मुझे लगता है कि मैट्रिक्स में NaN मान नहीं है, जो तुलना को प्रभावित करेगा। इस पर और बाद में।
void find_max(const int n_cols,
const int n_rows,
const float* val_ptr,
int* opt_pos_ptr,
float* max_ptr){
const __m128i mm_one = _mm_set1_epi32(1);
// Pre-compute the number of rows that can be processed in full vector width.
// In a 128-bit vector there are 4 floats or 2 doubles
int tail_size = n_rows & 3;
int n_rows_aligned = n_rows - tail_size;
int row = 0;
for (; row < n_rows_aligned; row += 4)
{
const auto* col_ptr = val_ptr + row;
__m128 mm_max = _mm_loadu_ps(col_ptr);
__m128i mm_max_pos = _mm_setzero_si128();
__m128i mm_pos = mm_one;
col_ptr += n_rows;
for (int col = 1; col < n_cols; ++col)
{
__m128 mm_value = _mm_loadu_ps(col_ptr);
// See if this value is greater than the old maximum
__m128 mm_mask = _mm_cmplt_ps(mm_max, mm_value);
// If it is, save its position
mm_max_pos = _mm_blendv_epi8(mm_max_pos, mm_pos, _mm_castps_si128(mm_mask));
// Compute the maximum
mm_max = _mm_max_ps(mm_value, mm_max);
mm_pos = _mm_add_epi32(mm_pos, mm_one);
col_ptr += n_rows;
}
// Store the results
_mm_storeu_ps(max_ptr + row, mm_max);
_mm_storeu_si128(reinterpret_cast< __m128i* >(opt_pos_ptr + row), mm_max_pos);
}
// Process tail serially
for (; row < n_rows; ++row)
{
const auto* col_ptr = val_ptr + row;
auto max = *col_ptr;
int max_pos = 0;
col_ptr += n_rows;
for (int col = 1; col < n_cols; ++col)
{
auto value = *col_ptr;
if (value > max)
{
max = value;
max_pos = col;
}
col_ptr += n_rows;
}
max_ptr[row] = max;
opt_pos_ptr[row] = max_pos;
}
}
सम्मिश्रण के कारण उपरोक्त कोड को SSE4.1 की आवश्यकता है। आप का एक संयोजन के साथ उन की जगह ले सकता _mm_and_si128/ _ps, _mm_andnot_si128/ _psऔर _mm_or_si128/ _ps, जिस स्थिति में आवश्यकताओं SSE2 को कम कर दी जाएगी। विशेष आंतरिक पर अधिक विवरण के लिए इंटेल इंट्रिंसिक्स गाइड देखें , जिसमें निर्देश सेट एक्सटेंशन शामिल हैं जिनकी उन्हें आवश्यकता है।
NaN मूल्यों के बारे में एक नोट। यदि आपके मैट्रिक्स में NaN हो सकते हैं, तो _mm_cmplt_psपरीक्षण हमेशा गलत होगा। के रूप में _mm_max_ps, यह आम तौर पर ज्ञात नहीं है कि यह क्या लौटेगा। maxpsअनुदेश कि रिटर्न के लिए आंतरिक तब्दील हो अपनी दूसरी (स्रोत) संकार्य अगर ऑपरेंड में से किसी एक NaN है तो शिक्षा का ऑपरेंड की व्यवस्था के द्वारा आप या तो व्यवहार प्राप्त कर सकते हैं। हालांकि, यह प्रलेखित नहीं है कि _mm_max_psआंतरिक का कौन सा तर्क निर्देश के किस संचालक का प्रतिनिधित्व करता है, और यह भी संभव है कि संकलक विभिन्न मामलों में अलग-अलग संघ का उपयोग कर सकता है। देखें इस अधिक जानकारी के लिए इस सवाल का जवाब।
सही व्यवहार को सुनिश्चित करने के लिए। NaN आप इनलाइन असेंबलर का उपयोग maxpsऑपरेंड के सही क्रम को लागू करने के लिए कर सकते हैं । दुर्भाग्य से, यह x86-64 लक्ष्य के लिए MSVC के साथ एक विकल्प नहीं है, जिसे आपने कहा था कि आप उपयोग कर रहे हैं, इसलिए इसके बजाय आप _mm_cmplt_psइस तरह से एक दूसरे मिश्रण के लिए परिणाम का पुन: उपयोग कर सकते हैं :
// Compute the maximum
mm_max = _mm_blendv_ps(mm_max, mm_value, mm_mask);
यह परिणामस्वरूप अधिकतम मानों में NaN को दबा देगा। यदि आप इसके बजाय NaN रखना चाहते हैं, तो आप NaN का पता लगाने के लिए दूसरी तुलना का उपयोग कर सकते हैं:
// Detect NaNs
__m128 mm_nan_mask = _mm_cmpunord_ps(mm_value, mm_value);
// Compute the maximum
mm_max = _mm_blendv_ps(mm_max, mm_value, _mm_or_ps(mm_mask, mm_nan_mask));
यदि आप व्यापक वैक्टर ( __m256या __m512) का उपयोग करते हैं और एक छोटे से कारक द्वारा बाहरी लूप को अनियंत्रित करते हैं, तो आप संभवतः ऊपर एल्गोरिथ्म के प्रदर्शन को और बेहतर बना सकते हैं , ताकि आंतरिक लूप के प्रत्येक पुनरावृत्ति पर कम से कम कैश डेटा की पंक्ति रेखा लोड हो जाए।
यहाँ कार्यान्वयन के लिए एक उदाहरण है double। यहां ध्यान देने वाली महत्वपूर्ण बात यह है कि क्योंकि doubleवेक्टर के प्रति केवल दो तत्व होते हैं और वेक्टर के प्रति अभी भी चार स्थान हैं, हमें doubleएक बार में दो वैक्टर को संसाधित करने के लिए बाहरी लूप को अनियंत्रित करना होगा और फिर दो मास्क को तुलनाओं से संपीड़ित करना होगा। 32-बिट पोजीशन को मिश्रित करने के लिए पिछले मैक्सिमम।
void find_max(const int n_cols,
const int n_rows,
const double* val_ptr,
int* opt_pos_ptr,
double* max_ptr){
const __m128i mm_one = _mm_set1_epi32(1);
// Pre-compute the number of rows that can be processed in full vector width.
// In a 128-bit vector there are 2 doubles, but we want to process
// two vectors at a time.
int tail_size = n_rows & 3;
int n_rows_aligned = n_rows - tail_size;
int row = 0;
for (; row < n_rows_aligned; row += 4)
{
const auto* col_ptr = val_ptr + row;
__m128d mm_max1 = _mm_loadu_pd(col_ptr);
__m128d mm_max2 = _mm_loadu_pd(col_ptr + 2);
__m128i mm_max_pos = _mm_setzero_si128();
__m128i mm_pos = mm_one;
col_ptr += n_rows;
for (int col = 1; col < n_cols; ++col)
{
__m128d mm_value1 = _mm_loadu_pd(col_ptr);
__m128d mm_value2 = _mm_loadu_pd(col_ptr + 2);
// See if this value is greater than the old maximum
__m128d mm_mask1 = _mm_cmplt_pd(mm_max1, mm_value1);
__m128d mm_mask2 = _mm_cmplt_pd(mm_max2, mm_value2);
// Compress the 2 masks into one
__m128i mm_mask = _mm_packs_epi32(
_mm_castpd_si128(mm_mask1), _mm_castpd_si128(mm_mask2));
// If it is, save its position
mm_max_pos = _mm_blendv_epi8(mm_max_pos, mm_pos, mm_mask);
// Compute the maximum
mm_max1 = _mm_max_pd(mm_value1, mm_max1);
mm_max2 = _mm_max_pd(mm_value2, mm_max2);
mm_pos = _mm_add_epi32(mm_pos, mm_one);
col_ptr += n_rows;
}
// Store the results
_mm_storeu_pd(max_ptr + row, mm_max1);
_mm_storeu_pd(max_ptr + row + 2, mm_max2);
_mm_storeu_si128(reinterpret_cast< __m128i* >(opt_pos_ptr + row), mm_max_pos);
}
// Process 2 doubles at once
if (tail_size >= 2)
{
const auto* col_ptr = val_ptr + row;
__m128d mm_max1 = _mm_loadu_pd(col_ptr);
__m128i mm_max_pos = _mm_setzero_si128();
__m128i mm_pos = mm_one;
col_ptr += n_rows;
for (int col = 1; col < n_cols; ++col)
{
__m128d mm_value1 = _mm_loadu_pd(col_ptr);
// See if this value is greater than the old maximum
__m128d mm_mask1 = _mm_cmplt_pd(mm_max1, mm_value1);
// Compress the mask. The upper half doesn't matter.
__m128i mm_mask = _mm_packs_epi32(
_mm_castpd_si128(mm_mask1), _mm_castpd_si128(mm_mask1));
// If it is, save its position
mm_max_pos = _mm_blendv_epi8(mm_max_pos, mm_pos, mm_mask);
// Compute the maximum
mm_max1 = _mm_max_pd(mm_value1, mm_max1);
mm_pos = _mm_add_epi32(mm_pos, mm_one);
col_ptr += n_rows;
}
// Store the results
_mm_storeu_pd(max_ptr + row, mm_max1);
// Only store the lower two positions
_mm_storel_epi64(reinterpret_cast< __m128i* >(opt_pos_ptr + row), mm_max_pos);
row += 2;
}
// Process tail serially
for (; row < n_rows; ++row)
{
const auto* col_ptr = val_ptr + row;
auto max = *col_ptr;
int max_pos = 0;
col_ptr += n_rows;
for (int col = 1; col < n_cols; ++col)
{
auto value = *col_ptr;
if (value > max)
{
max = value;
max_pos = col;
}
col_ptr += n_rows;
}
max_ptr[row] = max;
opt_pos_ptr[row] = max_pos;
}
}