Я не привык к векторизации, но мне это кажется очень подозрительным:
chi[j]=vchi[0];
chi[j+1]=vchi[1];
chi[j+2]=vchi[2];
chi[j+3]=vchi[3];
И на самом деле, замена его на то, что очень похоже на правильную функцию для работы, а именно _mm256_store_pd()
, похоже, решает проблему.
Ваша функция теперь может выглядеть так (с несколькими стилистическими исправлениями)
void dd_m_dd(double *ahi, double *bhi, double *chi, int m, int n) {
#pragma omp parallel for
for (int j = 0; j < m*n; j+=4) {
__m256d vbhi = _mm256_broadcast_sd(&bhi[j]);
__m256d vahi = _mm256_load_pd(&ahi[j]);
__m256d vchi=vahi*vbhi;
_mm256_store_pd( &chi[j], vchi );
}
}
Другая проблема заключается в том, что вы не применяете правильное выравнивание ваших указателей ... Переписывание распределений, как это, просто исправляет это:
double *xhi=(double *)aligned_alloc(256, sizeof(double) * m*1);
double *yhi=(double *)aligned_alloc(256, sizeof(double) * m*1);
double *z=(double *)aligned_alloc(256, sizeof(double) * m*1);