如何加快绝对损失矩阵的计算速度

• 启用编译器优化@Jesper Juhl。
• 如果可以，使用float类型和float函数。有时候快4倍。对我来说，快了8%。
• 使用restrict让编译器知道引用的数据没有重叠。否则，编译器必须假设L[i] = ...;可能会更改D[]，这将阻止某些优化。
• 对于引用的数据，在可能的情况下使用const。
• 使用一致的索引类型。
• 更改索引增量。公司简介
• 索引类型：对我来说size_t和unsigned差不多。unsigned short速度快5%。

void absolute_loss(const float * restrict D, float * restrict L,
    mwSize cols, mwSize rows) {
  mwSize rows_cols = rows*cols;
  for (mwSize i = 0; i < rows; i++) {
    for (mwSize j = 0; j < rows; j++){
      float aux = 0.0;
      for (mwSize k = 0; k < rows_cols; k += rows) {
        aux = aux + fabsf(D[j + k] - D[i + k]); // Note: fabsf
      }
      L[i] = L[i] + aux;
    }
  }
  for (mwSize i = 0; i < rows; i++) {
    L[i] /= rows;
  }
}

注意事项：
我希望在函数的开始出现以下内容。

for (mwSize i = 0; i < rows; i++) {
  L[i] = 0.0;
}

提示，而不是rows, cols, i, j，使用M, N, m, n来匹配公式。我不确定你是否正确。
利用variable length arrays和示例用法的候选重写：

#include <math.h>

typedef unsigned short mwSize;

// Note re-ordered parameters.
void absolute_loss(mwSize m_rows, mwSize n_cols, //
    float D[restrict m_rows][n_cols], float L[restrict m_rows]) {

  for (mwSize ell = 0; ell < m_rows; ell++) {
    float ell_sum = 0.0;
    for (mwSize n = 0; n < n_cols; n++) {
      float d_ell_n = D[ell][n];
      for (mwSize m = 0; m < m_rows; m++) {
        ell_sum += fabsf(D[m][n] - d_ell_n);
      }
    }
    L[ell] = ell_sum / (float) m_rows;
  }
}

#include <assert.h>
#include <stdio.h>
#include <stdlib.h>
#include <time.h>

int main(void) {
  // Usually M is around 20,000 and N takes values around {10, 30, 144}.
  mwSize m_rows = (mwSize) (rand() % 1000 + (20000 - 1000));
  mwSize n_cols = (mwSize[3]) {10, 30, 144}[rand() % 3];

  float (*D)[m_rows][n_cols] = malloc(sizeof *D);
  assert(D);
  float (*L)[m_rows] = malloc(sizeof *L);
  assert(L);
  for (mwSize m = 0; m < m_rows; m++) {
    for (mwSize n = 0; n < n_cols; n++) {
      (*D)[m][n] = (float) (rand() % 1000 + 1);
    }
  }

  clock_t t0 = clock();
  absolute_loss(m_rows, n_cols, *D, *L);
  clock_t t1 = clock();
  // Print some of L
  for (mwSize ell = 0; ell < m_rows; ell++) {
    printf(" %-7g", (*L)[ell]);
    if (ell > 10) {
      printf("\n");
      break;
    }
  }
  printf("\n%g seconds.\n", (double) (t1 - t0) / CLOCKS_PER_SEC);
  free(L);
  free(D);
}

我的时间：4.906秒。

你的矩阵D的行和列似乎有一个不寻常的排列。您的函数使用索引访问数据，这些索引会跳来跳去，无法充分利用内存缓存。重新排列循环可以使它处理大部分连续的元素，从而显著提高性能。在我的电脑上，这个运行速度几乎是你发布的M=10000，N=30的函数的3倍。

void absolute_loss2(const double * D, double * L, , mwSize cols, mwSize rows) {

  double Dtemp;
  int i, j, k, rowstimesk;
  for (i = 0; i < rows; i++) {
    L[i] = 0.0;
  }
  for (i = 0; i < rows; i++) {
    for  (k = 0; k < cols; k++) {
      rowstimesk = rows * k;
      Dtemp = D[i + rowstimesk];
      for (j = 0; j < rows; j++){
        L[j] += fabs(D[j + rowstimesk] - Dtemp);
      }
    }
  }

  for (i = 0; i < rows; i++) {
    L[i] /= rows;
  }
}

您可能会以牺牲简单性为代价，用SIMD构建一个更快的函数，特别是如果您想保持可移植性的话。不过，这不是一个微不足道的重构。