Оптимизация — Оптимизация 3D-цикла (C ++)

Я работаю над многосеточным решателем в C ++ и сейчас пытаюсь улучшить производительность последовательного порта. Самая трудоемкая часть в этом — это сглаживание, которое в моем случае является последовательным решением чрезмерной релаксации. Это выглядит следующим образом (я надеюсь, что это само за себя):

int idx;
int strideY = stride_[level][0];
int strideZ = stride_[level][1];
for(int i = 0; i < steps; ++i) {

for(int z = 1; z <= innerGridpoints_[level][2]; ++z) {
for(int y = 1; y <= innerGridpoints_[level][1]; ++y) {
idx = getIndexInner(level, 1,y,z);
for(int x = 1; x <= innerGridpoints_[level][0]; ++x, ++idx) {
grid[idx] = (1. - omega)  * grid[idx] + omega * 1./6. * (grid[idx+1] + grid[idx-1] +
grid[idx + strideY]  + grid[idx - strideY] +
grid[idx + strideZ]  + grid[idx - strideZ] -
spacing_[level] * spacing_[level] * rhs[idx]);
}
}
}
}

Я уже сделал несколько оптимизаций: циклы расположены так, что внутренний цикл дает большинство локальных записей (т. Е. Соседние элементы расположены вдоль измерения x), и предварительный расчет idx (хотя это встроенная функция, она сохранила довольно много время так).
Я также пытался блокировать, то есть не итерировать по всей сетке, а только по маленьким частям, чтобы увеличить местность, но это не имело никакого влияния.
Последняя идея, которую я имею, это попытаться развернуть цикл, но на самом деле я не ожидаю больших улучшений от этого. Я думал, что, возможно, есть некоторые возможные улучшения в доступе к памяти. Любые советы приветствуются 🙂

Только к вашему сведению: размер сетки будет варьироваться от очень маленького до 255х255х255. Кроме того, сетка имеет некоторые границы в каждом измерении, состоящие из небольшого количества строк, то есть итерация не по всей сетке.