Типы параллельных данных C++26 и алгоритмы

std::vector<int> v = /* большой массив данных */;
std::sort(std::execution::par, v.begin(), v.end());

// Пример использования SSE-интринсиков (до C++26)
__m128 a = _mm_load_ps(data);
__m128 b = _mm_load_ps(data + 4);
__m128 result = _mm_add_ps(a, b);
_mm_store_ps(output, result);

#pragma omp parallel for
for (int i = 0; i < N; i++) {
    result[i] = heavyComputation(data[i]);
}

tbb::parallel_for(0, N, [&](int i) {
    result[i] = heavyComputation(data[i]);
});

// OpenMP
#pragma omp parallel for
for (int i = 0; i < v.size(); i++) {
    v[i] = process(v[i]);
}
 
// C++17
std::transform(std::execution::par, v.begin(), v.end(), v.begin(), process);

std::vector<double> prices(10000);
std::for_each(std::execution::par, prices.begin(), prices.end(), [](double& price) {
    // ОПАСНО: параллельная запись в соседние элементы может вызвать проблемы из-за
    // возможной переаллокации вектора одним из потоков
    price = calculateNewPrice();
});

// Как разбить данные? Как балансировать нагрузку?
// Стандарт умалчивает, это решает реализация
std::sort(std::execution::par, v.begin(), v.end());

// Обычный код для всех платформ
double dotProduct(const double* a, const double* b, size_t n) {
    double sum = 0.0;
    for (size_t i = 0; i < n; ++i) {
        sum += a[i] * b[i];
    }
    return sum;
}
 
// Специализированный код для процессоров с поддержкой AVX
double dotProductAVX(const double* a, const double* b, size_t n) {
    // Куча низкоуровневого кода с интринсиками
    __m256d sum = _mm256_setzero_pd();
    for (size_t i = 0; i < n; i += 4) {
        __m256d va = _mm256_loadu_pd(a + i);
        __m256d vb = _mm256_loadu_pd(b + i);
        sum = _mm256_add_pd(sum, _mm256_mul_pd(va, vb));
    }
    // И еще больше низкоуровневого кода для получения результата
    // ...
}

// Пример базового использования simd-типа
std::simd<float, 8> a{1.0f, 2.0f, 3.0f, 4.0f, 5.0f, 6.0f, 7.0f, 8.0f};
std::simd<float, 8> b{8.0f, 7.0f, 6.0f, 5.0f, 4.0f, 3.0f, 2.0f, 1.0f};
std::simd<float, 8> c = a + b; // Векторизованное сложение!

std::simd<int, 4> a{10, 20, 30, 40};
std::simd<int, 4> b{15, 15, 15, 15};
std::simd_mask<int, 4> mask = (a > b); // Результат: {false, true, true, true}
std::simd<int, 4> result = std::where(mask, a, b); // {15, 20, 30, 40}

#define IDX(x, y, z, width, height) ((z) * (width) * (height) + (y) * (width) + (x))
// Использование: pixel = data[IDX(x, y, channel, width, height)];

// Одномерный массив в памяти
std::vector<float> data(width * height * depth);
 
// Интерпретация как трехмерного массива
std::mdspan<float, std::extents<size_t, std::dynamic_extent, std::dynamic_extent, std::dynamic_extent>> 
    volume(data.data(), depth, height, width);
 
// Доступ как к трехмерному массиву
float value = volume[z, y, x];

// Row-major макет (по умолчанию)
std::mdspan<float, std::extents<size_t, height, width>, std::layout_right> img_cpp(data.data());
 
// Column-major макет
std::mdspan<float, std::extents<size_t, height, width>, std::layout_left> img_fortran(data.data());

// Выделение блока 10x10 из большого изображения
auto block = submdspan(image, std::pair{start_y, start_y + 10}, std::pair{start_x, start_x + 10});

// Адаптер для нормализации данных
template<class T>
struct NormalizingAccessor {
  T scale;
  T operator()(T value) const { return value / scale; }
  void set(T& ref, T value) const { ref = value * scale; }
};
 
// Использование с mdspan
auto view = mdspan<float, extents<size_t, h, w>, layout_right, 
                   NormalizingAccessor<float>>(data.data(), {255.0f});

// Паттерн для параллельной обработки
std::array<size_t, 2> strides = {2*n, 2};
auto view = mdspan<int, extents<size_t, n/2, n/2>, layout_stride>
            (data.data(), mapping(extents(n/2, n/2), strides));

// Пример использования simd_vector
std::experimental::simd_vector<float, 16> vec(1000); // 1000 элементов, оптимизированных для SIMD-доступа
std::experimental::native_simd<float> simd_val{1.0f};
vec.fill(simd_val); // Заполнение всех элементов значением 1.0f через SIMD

std::experimental::simd_vector<float, 16> a(1000, 1.0f);
std::experimental::simd_vector<float, 16> b(1000, 2.0f);
std::experimental::simd_vector<float, 16> c(1000);
 
std::transform(a.begin(), a.end(), b.begin(), c.begin(), 
               [](auto x, auto y) { return x + y; });
// Трансформация будет автоматически векторизована!

// Пример операции gather
std::vector<int> indices = {5, 10, 15, 20};
std::experimental::simd<float, 4> values = simd_gather(data_array, indices.data());

// Старый подход
for (size_t i = 0; i < data.size(); ++i) {
    result[i] = std::min(data1[i], data2[i]);
}
 
// Новый подход с C++26
auto simd_result = stdx::min(simd_data1, simd_data2);

auto [minimums, maximums] = stdx::minmax(simd_data1, simd_data2);

// Ограничение значений в диапазоне [0, 255]
auto clamped = stdx::clamp(pixels, stdx::simd<int, 8>(0), stdx::simd<int, 8>(255));

auto red_channel = submdspan(image, std::full_extent, std::full_extent, 0);
auto min_values = stdx::reduce(red_channel, stdx::simd<float>(std::numeric_limits<float>::max()), 
                              [](auto a, auto b) { return stdx::min(a, b); });

stdx::fixed_size_simd<float, 8> a = /* ... */;
stdx::fixed_size_simd<float, 8> b = /* ... */;
stdx::fixed_size_simd<float, 8> result = a + b;

stdx::fixed_size_simd<int, 8> values = /* ... */;
stdx::fixed_size_simd<int, 8> threshold{42};
stdx::fixed_size_simd_mask<int, 8> mask = values > threshold;
values = where(mask, values, stdx::fixed_size_simd<int, 8>{0});

std::vector<int> node_ids = {23, 57, 129, 42, 78, 15, 91, 111};
stdx::fixed_size_simd<float, 8> node_values = 
    stdx::simd_gather<float, 8>(node_data.data(), node_ids.data());

// Было: обычный цикл
float rms = 0.0f;
for (size_t i = 0; i < samples.size(); ++i) {
    rms += samples[i] * samples[i];
}
rms = std::sqrt(rms / samples.size());
 
// Стало: SIMD-оптимизированная версия
using simd_type = stdx::native_simd<float>;
simd_type sum{0.0f};
for (size_t i = 0; i < samples.size(); i += simd_type::size()) {
    simd_type batch = stdx::simd_load<simd_type>(samples.data() + i);
    sum += batch * batch;
}
float rms = std::sqrt(stdx::reduce(sum) / samples.size());

// Выровненная загрузка (быстрее)
auto aligned_data = stdx::simd_aligned_allocator<float>().allocate(size);
auto simd_val = stdx::simd_load<simd_type>(aligned_data);
 
// Невыровненная загрузка (медленнее, но более гибкая)
auto simd_val = stdx::simd_load<simd_type>(regular_data);

std::sort(std::execution::par_unseq, data.begin(), data.end());

auto custom_policy = std::execution::simd.with(
  std::execution::prefer_gpu,
  std::execution::chunk_size(1024)
);
 
std::transform(custom_policy, input.begin(), input.end(), output.begin(), 
              [](auto x) { return std::sqrt(x); });

// Параллельная редукция с использованием SIMD
auto result = std::reduce(std::execution::par_unseq, 
                        data.begin(), data.end(),
                        stdx::native_simd<float>{0.0f},
                        std::plus<>());

std::execution::prefer_parallel_size(8) // Предпочитать разделение на 8 частей

// Расчет экспоненциального скользящего среднего (EMA)
std::vector<double> calculateEMA(const std::vector<double>& prices, int period) {
    std::vector<double> ema(prices.size());
    double multiplier = 2.0 / (period + 1.0);
    
    ema[0] = prices[0]; // Первое значение - просто цена
    for (size_t i = 1; i < prices.size(); ++i) {
        ema[i] = (prices[i] - ema[i-1]) * multiplier + ema[i-1];
    }
    
    return ema;
}

std::vector<double> calculateEMASIMD(const std::vector<double>& prices, int period) {
    std::vector<double> ema(prices.size());
    double multiplier = 2.0 / (period + 1.0);
    
    // Рассчитываем множители для каждой позиции
    std::vector<double> factors(prices.size());
    factors[0] = 1.0;
    for (size_t i = 1; i < prices.size(); ++i) {
        factors[i] = factors[i-1] * (1.0 - multiplier);
    }
    
    // Теперь используем SIMD для параллельных вычислений
    using simd_t = stdx::native_simd<double>;
    for (size_t i = 0; i < prices.size(); i += simd_t::size()) {
        simd_t sum{0.0};
        for (size_t j = 0; j <= i; j += simd_t::size()) {
            size_t k = i - j;
            if (j <= i) {
                simd_t price_batch = stdx::simd_load<simd_t>(&prices[k]);
                simd_t factor_batch = stdx::simd_load<simd_t>(&factors[j]);
                sum += price_batch * factor_batch * multiplier;
            }
        }
        stdx::simd_store(sum, &ema[i]);
    }
    
    return ema;
}

// Свертка изображения с ядром фильтра
void convolution(const std::vector<float>& input, std::vector<float>& output,
                int width, int height, const std::vector<float>& kernel, int kernel_size) {
    int k_radius = kernel_size / 2;
    
    for (int y = k_radius; y < height - k_radius; ++y) {
        for (int x = k_radius; x < width - k_radius; ++x) {
            float sum = 0.0f;
            for (int ky = 0; ky < kernel_size; ++ky) {
                for (int kx = 0; kx < kernel_size; ++kx) {
                    int input_idx = (y + ky - k_radius) * width + (x + kx - k_radius);
                    int kernel_idx = ky * kernel_size + kx;
                    sum += input[input_idx] * kernel[kernel_idx];
                }
            }
            output[y * width + x] = sum;
        }
    }
}

void convolutionSIMD(std::mdspan<const float, std::extents<size_t, std::dynamic_extent, std::dynamic_extent>> input,
                    std::mdspan<float, std::extents<size_t, std::dynamic_extent, std::dynamic_extent>> output,
                    std::mdspan<const float, std::extents<size_t, std::dynamic_extent, std::dynamic_extent>> kernel) {
    int height = input.extent(0);
    int width = input.extent(1);
    int k_height = kernel.extent(0);
    int k_width = kernel.extent(1);
    int k_radius_y = k_height / 2;
    int k_radius_x = k_width / 2;
    
    using simd_t = stdx::native_simd<float>;
    const int simd_width = simd_t::size();
    
    for (int y = k_radius_y; y < height - k_radius_y; ++y) {
        // Обрабатываем по simd_width пикселей за раз
        for (int x = k_radius_x; x < width - k_radius_x; x += simd_width) {
            std::array<simd_t, 1> sum = {simd_t(0.0f)};
            
            for (int ky = 0; ky < k_height; ++ky) {
                for (int kx = 0; kx < k_width; ++kx) {
                    float k_val = kernel(ky, kx);
                    simd_t k_simd(k_val);
                    
                    // Загружаем блок пикселей
                    simd_t input_simd;
                    for (int sx = 0; sx < simd_width; ++sx) {
                        input_simd[sx] = input(y + ky - k_radius_y, x + sx + kx - k_radius_x);
                    }
                    
                    sum[0] += input_simd * k_simd;
                }
            }
            
            // Сохраняем результат
            for (int sx = 0; sx < simd_width; ++sx) {
                if (x + sx < width - k_radius_x)
                    output(y, x + sx) = sum[0][sx];
            }
        }
    }
}

// Наивная попытка ускорения - может оказаться медленнее!
stdx::native_simd<float> process(stdx::native_simd<float> data) {
  // Операции, которые плохо векторизуются
  return data * stdx::sin(data) / (stdx::cos(data) + 0.01f);
}

#ifdef DEBUG_SIMD
  using debug_simd_t = stdx::fixed_size_simd<float, 4, stdx::simd_abi::scalar>;
#else
  using debug_simd_t = stdx::native_simd<float>;
#endif

// Проверка, что операция безопасна для параллельного выполнения
#ifdef PARALLEL_CHECK
  std::atomic<bool> accessing{false};
  bool was_accessing = accessing.exchange(true);
  assert(!was_accessing && "Parallel access detected!");
  // Операция с SIMD-данными
  accessing = false;
#else
  // Обычное выполнение операции
#endif

template<typename T>
bool is_aligned(const T* ptr, size_t alignment) {
  return reinterpret_cast<uintptr_t>(ptr) % alignment == 0;
}
 
// Использование
if (!is_aligned(data.data(), 32)) {
  std::cerr << "Warning: Unaligned data may degrade SIMD performance" << std::endl;
}

void debug_mask(const stdx::native_simd_mask<float>& mask) {
  std::cout << "Mask: ";
  for (size_t i = 0; i < mask.size(); ++i) {
    std::cout << (mask[i] ? '1' : '0');
  }
  std::cout << std::endl;
}