Хруев Антон. Технология SEQ|MPI.Многошаговая схема решения двумерных задач глобальной оптимизации. Распараллеливание по характеристикам. Вариант 13 (#294)

<!--
Требования к названию pull request:

"<Фамилия> <Имя>. Технология <TECHNOLOGY_NAME:SEQ|OMP|TBB|STL|MPI>.
<Полное название задачи>. Вариант <Номер>"
-->

## Описание
<!--
Пожалуйста, предоставьте подробное описание вашей реализации, включая:
 - основные детали решения (описание выбранного алгоритма)
 - применение технологии параллелизма (если применимо)
-->

- **Задача**: Многошаговая схема решения двумерных задач глобальной
оптимизации. Распараллеливание по характеристикам.
- **Вариант**: 13
- **Технология**:  SEQ,  MPI
- **Описание**:  
**SEQ:** Основная идея алгоритма: Двумерная задач оптимизации сводится к
одномерной с помощью развертки Гильберта. Алгоритм состоит из:
инициализации, где в список испытаний добавляются точки ```t=0``` и
```t=1```, для которых вычисляются значения целевой функции; оценки
константы Липшица; вычисления характеристик интервалов и выбора
интервала с максимальной характеристикой; построение новой точки.
**MPI:** Распараллеливание выполняется по интервалам с наибольшими
характеристиками. Алгоритм включает инициализацию, вычисление параметра
М (константа Липшица), вычисление характеристик, распределение работы,
где каждый процесс получает интервал для вычисления новой точки. После
проводится обмен результатами через ```MPI_Allgather```.
 Более подробное описание находится в отчёте.

---

## Чек-лист
<!--
Пожалуйста, убедитесь, что следующие пункты выполнены **до** отправки
pull request'а и запроса его ревью:
-->

- [x] **Статус CI**: Все CI-задачи (сборка, тесты, генерация отчёта)
успешно проходят на моей ветке в моем форке
- [x] **Директория и именование задачи**: Я создал директорию с именем
`<фамилия>_<первая_буква_имени>_<короткое_название_задачи>`
- [x] **Полное описание задачи**: Я предоставил полное описание задачи в
теле pull request
- [x] **clang-format**: Мои изменения успешно проходят `clang-format`
локально в моем форке (нет ошибок форматирования)
- [x] **clang-tidy**: Мои изменения успешно проходят `clang-tidy`
локально в моем форке (нет предупреждений/ошибок)
- [x] **Функциональные тесты**: Все функциональные тесты успешно
проходят локально на моей машине
- [x] **Тесты производительности**: Все тесты производительности успешно
проходят локально на моей машине
- [x] **Ветка**: Я работаю в ветке, названной точно так же, как
директория моей задачи (например, `nesterov_a_vector_sum`), а не в
`master`
- [x] **Правдивое содержание**: Я подтверждаю, что все сведения,
указанные в этом pull request, являются точными и достоверными

<!--
ПРИМЕЧАНИЕ: Ложные сведения в этом чек-листе могут привести к отклонению
PR и получению нулевого балла за соответствующую задачу.
-->

---------

Co-authored-by: brb484ul <pppphp2800@gmail.com>

Author: pt11241

tasks/khruev_a_global_opt/common/include/common.hpp

@@ -0,0 +1,137 @@
+#pragma once
+
+#include <algorithm>
+#include <cmath>
+#include <cstdint>
+#include <string>
+#include <tuple>
+
+#include "task/include/task.hpp"
+
+namespace khruev_a_global_opt {
+
+struct SearchData {
+  int func_id;
+  double ax, bx;
+  double ay, by;
+  double epsilon;
+  int max_iter;
+  double r;
+};
+
+struct SearchResult {
+  double x;
+  double y;
+  double value;
+  int iter_count;
+};
+
+using InType = SearchData;
+using OutType = SearchResult;
+using BaseTask = ppc::task::Task<InType, OutType>;
+using TestType = std::tuple<std::string, int, double, double, double, double, double>;
+
+// 1. Кривая Гильберта: переводит t [0,1] -> (x, y) [0,1]x[0,1]
+// Вспомогательная функция для поворота/отражения квадранта
+// n: размер текущего квадрата (степень двойки)
+// x, y: координаты
+// rx, ry: биты, определяющие квадрант
+inline void Rotate(uint64_t n, uint64_t &x, uint64_t &y, uint64_t rx, uint64_t ry) {
+  if (ry == 0) {
+    if (rx == 1) {
+      x = n - 1 - x;
+      y = n - 1 - y;
+    }
+    // Swap x и y
+    uint64_t t = x;
+    x = y;
+    y = t;
+  }
+}
+
+// 1. Кривая Гильберта: переводит t [0,1] -> (x, y) [0,1]x[0,1]
+// Используем 32-битный порядок кривой (сетка 2^32 x 2^32), что покрывает точность double
+inline void D2xy(double t, double &x, double &y) {
+  // Ограничиваем t диапазоном [0, 1]
+  t = std::max(t, 0.0);
+  t = std::min(t, 1.0);
+
+  // Порядок кривой N = 32. Всего точек 2^64 (влезает в uint64_t).
+  // Масштабируем t [0, 1] в целое число s [0, 2^64 - 1]
+  // Используем (2^64 - 1) как множитель.
+  const auto max_dist = static_cast<uint64_t>(-1);
+  const double two_to_64 = ldexp(1.0, 64);
+  double val = t * two_to_64;
+  uint64_t s = 0;
+  if (val >= two_to_64) {
+    s = max_dist;
+  } else {
+    s = static_cast<uint64_t>(val);
+  }
+
+  uint64_t ix = 0;
+  uint64_t iy = 0;
+
+  // Итеративный алгоритм от младших битов к старшим (bottom-up)
+  // n - это размер под-квадрата на текущем уровне
+  for (uint64_t nn = 1; nn < (1ULL << 32); nn <<= 1) {
+    uint64_t rx = 1 & (s / 2);
+    uint64_t ry = 1 & (s ^ rx);
+
+    Rotate(nn, ix, iy, rx, ry);
+
+    ix += nn * rx;
+    iy += nn * ry;
+
+    s /= 4;
+  }
+
+  // Нормализуем обратно в [0, 1]
+  // Делим на 2^32 (размер сетки)
+  const double scale = 1.0 / 4294967296.0;  // 1.0 / 2^32
+  x = static_cast<double>(ix) * scale;
+  y = static_cast<double>(iy) * scale;
+}
+
+inline double TargetFunction(int id, double x, double y) {
+  if (id == 1) {
+    // квадратичная
+    return ((x - 0.5) * (x - 0.5)) + ((y - 0.5) * (y - 0.5));
+  }
+  if (id == 2) {
+    // Rastrigin-like
+    constexpr double kA = 10.0;
+    constexpr double kTwoPi = 6.2831853071795864769;
+    return (2.0 * kA) + ((x * x) - (kA * std::cos(kTwoPi * x))) + ((y * y) - (kA * std::cos(kTwoPi * y)));
+  }
+  if (id == 3) {
+    // BoothFunc
+    const double t1 = x + (2.0 * y) - 7.0;
+    const double t2 = (2.0 * x) + y - 5.0;
+    return (t1 * t1) + (t2 * t2);
+  }
+  if (id == 4) {
+    // MatyasFunc
+    return (0.26 * ((x * x) + (y * y))) - (0.48 * x * y);
+  }
+  if (id == 5) {  // Himme
+    const double t1 = ((x * x) + y - 11.0);
+    const double t2 = (x + (y * y) - 7.0);
+    return (t1 * t1) + (t2 * t2);
+  }
+
+  return 0.0;
+}
+
+// Структура для хранения испытаний
+struct Trial {
+  double x;  // координата на отрезке [0, 1]
+  double z;  // значение функции
+  int id;    // исходный индекс
+
+  bool operator<(const Trial &other) const {
+    return x < other.x;
+  }
+};
+
+}  // namespace khruev_a_global_opt

tasks/khruev_a_global_opt/info.json

@@ -0,0 +1,9 @@
+{
+  "student": {
+    "first_name": "Антон",
+    "last_name": "Хруев",
+    "middle_name": "Олегович",
+    "group_number": "3823Б1ФИ1",
+    "task_number": "3"
+  }
+}

tasks/khruev_a_global_opt/mpi/include/ops_mpi.hpp

@@ -0,0 +1,49 @@
+#pragma once
+
+#include <vector>
+
+#include "khruev_a_global_opt/common/include/common.hpp"
+#include "task/include/task.hpp"
+
+namespace khruev_a_global_opt {
+
+class KhruevAGlobalOptMPI : public BaseTask {
+ public:
+  static constexpr ppc::task::TypeOfTask GetStaticTypeOfTask() {
+    return ppc::task::TypeOfTask::kMPI;
+  }
+
+  explicit KhruevAGlobalOptMPI(const InType &in);
+  bool ValidationImpl() override;
+  bool PreProcessingImpl() override;
+  bool RunImpl() override;
+  bool PostProcessingImpl() override;
+
+ private:
+  struct IntervalInfo {
+    double R;
+    int index;  // индекс в массиве trials_ (указывает на правый конец интервала)
+
+    bool operator>(const IntervalInfo &other) const {
+      return R > other.R;  // Для сортировки по убыванию
+    }
+  };
+
+  struct Point {
+    double x, z;
+  };
+
+  std::vector<Trial> trials_;
+  OutType result_{};
+
+  double CalculateFunction(double t);
+  void AddTrialUnsorted(double t, double z);
+
+  double ComputeM();
+  [[nodiscard]] std::vector<IntervalInfo> ComputeIntervals(double m) const;
+  bool LocalShouldStop(const std::vector<IntervalInfo> &intervals, int num_to_check);
+  [[nodiscard]] double GenerateNewX(int idx, double m) const;
+  void CollectAndAddPoints(const std::vector<Point> &global_res, int &k);
+};
+
+}  // namespace khruev_a_global_opt

tasks/khruev_a_global_opt/mpi/src/ops_mpi.cpp

@@ -0,0 +1,186 @@
+#include "khruev_a_global_opt/mpi/include/ops_mpi.hpp"
+
+#include <mpi.h>
+
+#include <algorithm>
+#include <cmath>
+#include <cstddef>
+#include <functional>
+#include <utility>
+#include <vector>
+
+#include "khruev_a_global_opt/common/include/common.hpp"
+
+namespace khruev_a_global_opt {
+
+KhruevAGlobalOptMPI::KhruevAGlobalOptMPI(const InType &in) {
+  SetTypeOfTask(GetStaticTypeOfTask());
+  GetInput() = in;
+}
+
+bool KhruevAGlobalOptMPI::ValidationImpl() {
+  return GetInput().max_iter > 0 && GetInput().epsilon > 0 && GetInput().r > 1.0;
+}
+
+bool KhruevAGlobalOptMPI::PreProcessingImpl() {
+  trials_.clear();
+  return true;
+}
+
+double KhruevAGlobalOptMPI::CalculateFunction(double t) {
+  double u = 0;
+  double v = 0;
+  D2xy(t, u, v);
+  double real_x = GetInput().ax + (u * (GetInput().bx - GetInput().ax));
+  double real_y = GetInput().ay + (v * (GetInput().by - GetInput().ay));
+  return TargetFunction(GetInput().func_id, real_x, real_y);
+}
+
+void KhruevAGlobalOptMPI::AddTrialUnsorted(double t, double z) {
+  Trial tr{};
+  tr.x = t;
+  tr.z = z;
+  trials_.push_back(tr);
+}
+
+double KhruevAGlobalOptMPI::ComputeM() {
+  double max_slope = 0.0;
+  for (size_t i = 1; i < trials_.size(); ++i) {
+    double dx = trials_[i].x - trials_[i - 1].x;
+    double dz = std::abs(trials_[i].z - trials_[i - 1].z);
+    if (dx > 1e-12) {
+      max_slope = std::max(max_slope, dz / dx);
+    }
+  }
+  return (max_slope > 0) ? GetInput().r * max_slope : 1.0;
+}
+
+std::vector<KhruevAGlobalOptMPI::IntervalInfo> KhruevAGlobalOptMPI::ComputeIntervals(double m) const {
+  std::vector<IntervalInfo> intervals;
+  intervals.reserve(trials_.size() - 1);
+  for (size_t i = 1; i < trials_.size(); ++i) {
+    double dx = trials_[i].x - trials_[i - 1].x;
+    double z_r = trials_[i].z;
+    double z_l = trials_[i - 1].z;
+
+    // Классическая формула Стронгина
+    double r = (m * dx) + (((z_r - z_l) * (z_r - z_l)) / (m * dx)) - (2.0 * (z_r + z_l));
+    intervals.push_back({r, static_cast<int>(i)});
+  }
+  return intervals;
+}
+
+bool KhruevAGlobalOptMPI::LocalShouldStop(const std::vector<IntervalInfo> &intervals, int num_to_check) {
+  bool local_stop = true;
+  for (int i = 0; i < num_to_check; ++i) {
+    int idx = intervals[i].index;
+    if (trials_[idx].x - trials_[idx - 1].x > GetInput().epsilon) {
+      local_stop = false;
+      break;
+    }
+  }
+  return local_stop;
+}
+
+double KhruevAGlobalOptMPI::GenerateNewX(int idx, double m) const {
+  double x_r = trials_[idx].x;
+  double x_l = trials_[idx - 1].x;
+  double z_r = trials_[idx].z;
+  double z_l = trials_[idx - 1].z;
+
+  // Формула вычисления новой точки испытания
+  double new_x = (0.5 * (x_r + x_l)) - ((z_r - z_l) / (2.0 * m));
+
+  // Защита от выхода за границы интервала
+  if (new_x <= x_l || new_x >= x_r) {
+    new_x = 0.5 * (x_r + x_l);
+  }
+  return new_x;
+}
+
+void KhruevAGlobalOptMPI::CollectAndAddPoints(const std::vector<Point> &global_res, int &k) {
+  for (auto global_re : global_res) {
+    if (global_re.x >= 0.0) {
+      AddTrialUnsorted(global_re.x, global_re.z);
+      k++;
+    }
+  }
+}
+
+bool KhruevAGlobalOptMPI::RunImpl() {
+  int rank = 0;
+  int size = 0;
+  MPI_Comm_rank(MPI_COMM_WORLD, &rank);
+  MPI_Comm_size(MPI_COMM_WORLD, &size);
+
+  // Начальные испытания
+  if (trials_.empty()) {
+    AddTrialUnsorted(0.0, CalculateFunction(0.0));
+    AddTrialUnsorted(1.0, CalculateFunction(1.0));
+  }
+  std::ranges::sort(trials_, [](const Trial &a, const Trial &b) { return a.x < b.x; });
+
+  int k = 2;  // Уже провели 2 испытания
+
+  while (k < GetInput().max_iter) {
+    // 1. Вычисляем M (одинаково на всех процессах)
+    double m = ComputeM();
+
+    // 2. Вычисляем характеристики R
+    auto intervals = ComputeIntervals(m);
+
+    // 3. Сортируем интервалы по убыванию характеристики
+    std::ranges::sort(intervals, std::greater<>());
+
+    // 4. Проверка критерия остановки (по самому широкому из лучших интервалов)
+    int num_to_check = std::min(static_cast<int>(intervals.size()), size);
+    bool local_stop = LocalShouldStop(intervals, num_to_check);
+
+    // Синхронная остановка
+    int stop_signal = 0;
+    int local_stop_signal = local_stop ? 1 : 0;
+    MPI_Allreduce(&local_stop_signal, &stop_signal, 1, MPI_INT, MPI_LAND, MPI_COMM_WORLD);
+    if (stop_signal != 0) {
+      break;
+    }
+
+    // 5. Генерация новых точек
+    Point my_point{.x = -1.0, .z = 0.0};
+
+    if (std::cmp_less(rank, static_cast<int>(intervals.size()))) {
+      int idx = intervals[rank].index;
+      double new_x = GenerateNewX(idx, m);
+
+      my_point.x = new_x;
+      my_point.z = CalculateFunction(new_x);
+    }
+
+    // 6. Сбор результатов со всех процессов
+    std::vector<Point> global_res(size);
+    MPI_Allgather(&my_point, 2, MPI_DOUBLE, global_res.data(), 2, MPI_DOUBLE, MPI_COMM_WORLD);
+
+    // 7. Обновление списка
+    CollectAndAddPoints(global_res, k);
+    std::ranges::sort(trials_, [](const Trial &a, const Trial &b) { return a.x < b.x; });
+  }
+
+  // Поиск финального минимума среди всех испытаний
+  auto it = std::ranges::min_element(trials_, {}, &Trial::z);
+
+  double u = 0;
+  double v = 0;
+  D2xy(it->x, u, v);
+  result_.x = GetInput().ax + (u * (GetInput().bx - GetInput().ax));
+  result_.y = GetInput().ay + (v * (GetInput().by - GetInput().ay));
+  result_.value = it->z;
+  result_.iter_count = k;
+
+  GetOutput() = result_;
+  return true;
+}
+
+bool KhruevAGlobalOptMPI::PostProcessingImpl() {
+  return true;
+}
+
+}  // namespace khruev_a_global_opt