Эксперименты с hypot III

В посте Эксперименты с hypot я допустил существенную неточность: указал ссылку на реализацию std::hypot (double) без проверки реализации std::hypotf (float), в то время как для сравнения использовалась как раз версия одинарной точности. Итак, hypotf.

Оказывается, это, по сути, та же hypot_v3! Я было думал, что hypotf использует ту же методику, что и hypot. В этой мысли меня поддерживало и её низкое быстродействие (в 5 раз медленнее hypot_v3). Видимо, дополнительные проверки оказались столь дороги. Таким образом, с одной стороны, очевидно, что результаты hypot_v3 и hypotf (для чисел) совпадут. С другой стороны, разработчики руководствовались теми же соображениями, что и я. Функция hypot_v3 даёт корректно округлённый результат в одинарной точности, и её можно использовать для экспериментальной проверки прочих вариантов.

Вариант в двойной точности не мог быть создан с помощью того же приёма, так как четверная точность не поддерживается аппаратно на большинстве процессоров, а программная реализация очень медленная и требует подключения специальной библиотеки (хотя она и есть в комплекте gcc). И вот тут возник вопрос: float vs double даёт возможность проверить качество того или иного алгоритма, а если мы его применяем только в double, то мы, по сути, ограничены теоретическим анализом.

Я, наконец, прочёл комментарий в e_hypot.c, который, вроде бы (с кодом не сравнивал), описывает применённый там метод. А метод основан на соображении следующего вида:

«If (assume round-to-nearest) z=x*x+y*y has error less than sqrt(2)/2 ulp, than sqrt(z) has error less than 1 ulp (exercise).»

Т.е. если вычислить правильно округлённую сумму квадратов, то и корень будет правильно округлён (докажите, ага). Конечно, чисто формально можно разложить в ряд:

\sqrt{z+\varepsilon}-\sqrt{z} = \dfrac{\varepsilon}{2\sqrt{z}} + O(\varepsilon^2).

Если привести абсолютную погрешность к относительной, то получим

\dfrac{\sqrt{z(1+\varepsilon)} - \sqrt{z}}{\sqrt{z}} = \sqrt{1+\varepsilon}-1 = \dfrac{\varepsilon}{2} + O(\varepsilon^2).

В первом приближении извлечение корня уполовинивает относительную погрешность! С такой точки зрения, у нас должен получаться правильно округлённый результат даже в том случае, если последняя цифра z неверна. Впрочем, следует учесть погрешность вычисления самого корня. Если считать, что эпсилон = 0.5 ULP, то у нас появляется дополнительная относительная погрешность \varepsilon:

\left(1 + \dfrac{\varepsilon}{2} + O(\varepsilon^2)\right)(1 + \varepsilon) = 1 + \dfrac{3}{2}\varepsilon + O(\varepsilon^2).

Т.е. всё же ошибка в 0.75 ULP. Я выше написал «правильно округлён», но это, увы, не то же самое, что «ошибка меньше 1 ULP». Ошибку меньше 1 ULP мы получили, но для «правильно округлён» требуется ошибка не более 0.5 ULP, чего мы не достигаем. Эксперимент показывает, что, действительно, вычисление вида

float hypot_v2c(float x, float y)
{
  x = std::fabs(x);
  y = std::fabs(y);

  float
    min_ = x < y? x: y,
    max_ = x < y? y: x,
    mul, invmul;

  if (max_ + min_ == max_)
    return max_;

  // max_ != 0
  // min_ >= 0x1p-24 * max_

  if (max_ <= 0x1p-126f)
  {
    invmul = 0x1p-126f;  // 0x1p-23 <= max_ <= 1
    mul = 0x1p+126f;     // 0x1p-23 <= min_
  }
  else if (max_ <= 0x1p-63f)
  {
    invmul = 0x1p-63f;   // 0x1p-63 <= max_ <= 1
    mul = 0x1p+63f;      // 0x1p-87 <= min_
  }
  else if (max_ <= 1.f)
  {
    invmul = 1.f;        // 0x1p-63 <= max_ <= 1
    mul = 1.f;           // 0x1p-87 <= min_
  }
  else if (max_ <= 0x1p+63f)
  {
    invmul = 0x1p+63f;   // 0x1p-63 <= max_ <= 1
    mul = 0x1p-63f;      // 0x1p-87 <= min_
  }
  else
  {
    invmul = 0x1p+127f;  // 0x1p-64 <= max_ <  2
    mul = 0x1p-127f;     // 0x1p-88 <= min_
  }

  max_ *= mul;
  min_ *= mul;

  const double a = max_, b = min_, s = a*a + b*b;
  return invmul * std::sqrt((float)s);
}

даёт ошибку в последней двоичной цифре, хотя и немного реже, чем вариант hypot_v2b. Следовательно, то же поведение можно ожидать от std::hypot: ошибка в последней цифре. Т.е. она немногим лучше прямого вычисления в духе hypot_v2b. В связи с этим возникает вопрос: можно ли более-менее эффективно вычислить действительно правильно округлённую гипотенузу, не прибегая к арифметике удвоенной разрядности. И ещё: как быть с длиной n-мерного вектора? Метод hypot_v2b можно применить к n-мерному вектору: найти максимум по абсолютной величине, выбрать порядок, вычислить сумму промасштабированных квадратов с компенсацией и потом уже взять корень, и можно надеяться на погрешность в одну-две последних цифры.

Реклама

Эксперименты с hypot II

Небольшое дополнение к предыдущему посту. Ввиду низкого быстродействия frexp/ldexp может быть разумно заменить их на выбор множителя (степени двойки) банальными if-else.

inline float hypot_v2b(float x, float y)
{
  x = std::fabs(x);
  y = std::fabs(y);

  float
    min_ = x < y? x: y,
    max_ = x < y? y: x,
    mul, invmul;

  if (max_ <= 0x1p-126f)
  {
    invmul = 0x1p-126f;
    mul = 0x1p+126f;
  }
  else if (max_ <= 0x1p-63f)
  {
    invmul = 0x1p-63f;
    mul = 0x1p+63f;
  }
  else if (max_ <= 1.f)
  {
    invmul = 1.f;
    mul = 1.f;
  }
  else if (max_ <= 0x1p+63f)
  {
    invmul = 0x1p+63f;
    mul = 0x1p-63f;
  }
  else
  {
    invmul = 0x1p+126f;
    mul = 0x1p-126f;
  }

  max_ *= mul;
  min_ *= mul;

  return invmul * std::sqrt(max_* max_ + min_* min_);
}

Данный код предотвращает некорректное переполнение и исчезновение и обеспечивает погрешность 1 ULP. Конечно, гребёнка if-ов сказывается на скорости работы — у меня получилось, что hypot_v2b примерно втрое медленнее hypot_v1 (но всё же в десять раз быстрее hypot_v2a). Эта оценка довольно условна, так как скорость исполнения ветвлений зависит от множества факторов.

Аналогичный способ можно применить для замены frexp/ldexp при вычислении произведений. Например, если произведение получилось меньше по абсолютной величине наименьшего положительного нормального числа, то пересчитать его, предварительно домножив на 0x1p+63f; если произведение стало бесконечным (переполнение) или приблизилось по величине к наибольшему представимому числу, то пересчитать его, домножив на 0x1p-63f.

Эксперименты с hypot

В стандартной библиотеке математических функций C++ (унаследованной от C) имеется функция hypot(x, y), возвращающая евклидову длину вектора (x, y).

Зачем особая функция, если можно просто определить «в лоб»?

inline float hypot_v1(float x, float y)
{
  return std::sqrt(x * x + y * y);
}

(Здесь я использую float и только float по причине, о которой будет сказано ниже.)

Причин может быть несколько:

  1. Гарантировать отсутствие переполнения (из-за возведения в квадрат большого числа или сложения больших чисел) в случае, когда конечный результат (длина) представим.
  2. Гарантировать отсутствие влияния исчезновения (из-за возведения в квадрат слишком маленького числа) на результат. В частности, в случае, когда один из аргументов x, y есть нуль, а другой — число, hypot(x, y) должен возвращать модуль другого аргумента, даже если это субнормальное число.
  3. Обеспечивать погрешность лучше 1 ULP.

Легко показать, что hypot_v1 не удовлетворяет всем этим пунктам. Экспериментально установлено, что погрешность hypot_v1 на основной части диапазона, где не возникает ситуация вредоносного переполнения или исчезновения составляет 1 ULP (проверка для y из набора { 0.f, 1e-40f, 1e-30f, 1e-20f, 1e-15f, 1e-6f, 1e-2f, 1.f, 1e+2f, 1e+6f, 1e+15f, 1e+20f, 1e+30f } и для всех возможных значений x).

В активе hypot_v1 сравнительно большая скорость вычисления (у меня получилось, что она примерно в 10 раз быстрее, чем стандартная hypot на g++ -O3).

Можно ли получить какой-то иной компромисс в плане скорости и качества? В Википедии можно наблюдать следующий вариант (на момент написания). Ниже дан исправленный код:

inline float hypot_v2(float x, float y)
{
  x = std::fabs(x);
  y = std::fabs(y);

  const float
    min_ = x < y? x: y,
    max_ = x < y? y: x;

  if (max_ == 0.f || std::isinf(max_))
    return max_;

  const float q = min_ / max_;
  return max_ * std::sqrt(1.f + q*q);
}

Хорош ли этот код? Он страхует нас от патологического переполнения и исчезновения. Что можно сказать о его скорости? В моём эксперименте получилось, что hypot_v2 примерно вдвое медленнее hypot_v1, что можно считать приемлемым. Что можно сказать о погрешности? Экспериментально получено максимальное значение 2 ULP.

В стремлении убрать операции, вносящие погрешность, я написал иной вариант, использующий вынесение степеней двойки:

inline float hypot_v2a(float x, float y)
{
  x = std::fabs(x);
  y = std::fabs(y);

  if (x == 0.f)
    return y;
  if (y == 0.f)
    return x;

  float
    min_ = x < y? x: y,
    max_ = x < y? y: x;

  int min_e, max_e;
  min_ = std::frexp(min_, &min_e);
  max_ = std::frexp(max_, &max_e);

  min_ *= min_;
  max_ *= max_;

  return std::ldexp(
    std::sqrt(std::ldexp(min_, 2*(min_e - max_e)) + max_),
    max_e);
}

Этот вариант действительно обеспечивает погрешность <= 1 ULP (экспериментально). К сожалению, этот код в моём случае работал втрое медленнее стандартной hypot. (Что может говорить о потенциально низком быстродействии кода предложенного мной для перемножения массивов чисел.)

Теперь следует прояснить вопрос об экспериментальном вычислении погрешности. Так как используются числа float (IEEE-754 binary32, 23 двоичных разряда после запятой), и доступны числа double (IEEE-754 binary64, 52 двоичных разряда после запятой), то нетрудно организовать вычисление ближайшего к истинной величине значения в формате float, попросту выполнив вычисление с точностью double и округлив затем результат до float:

inline float hypot_v3(float x, float y)
{
  const double x_ = x, y_ = y;
  return (float)std::sqrt(x_ * x_ + y_ * y_);
}

Кстати, этот вариант также не страдает от переполнения и исчезновения и вычисляется лишь вдвое медленнее hypot_v1, обеспечивая погрешность <= 0.5 ULP. Он плох тем, что не масштабируется на double (хотя можно задействовать софтовое вычисление с точностью IEEE-754 binary128, заплатив за это резким падением производительности), зато очень удобен для определения качества других вариантов, включая стандартную hypot.

В случае g++ результат hypot везде совпал с результатом hypot_v3 (т.е. точность <= 0.5 ULP стандартной hypot обеспечена). Интересно, что для её вычисления используется древний (1993г.) код от Sun Microsystems, довольно зубодробительного вида.

В принципе, альтернативой hypot и hypot_v1 может стать их гибрид, использующий hypot_v1 там, где она работает хорошо (с погрешностью 1 ULP), а в остальных случаях (предположительно, редких) вызывающий hypot:

inline float hypot_v4(float x, float y)
{
  x = std::fabs(x);
  y = std::fabs(y);

  if (1e-19f <= x && x <= 1e+19f &&
      1e-19f <= y && y <= 1e+19f)
    return std::sqrt(x*x + y*y);

  return std::hypot(x, y);
}

Такой код способен работать со скоростью, в среднем близкой к скорости hypot_v1.

Вычисление чисел Фибоначчи во время компиляции С++14

Шесть лет назад я привёл пример вычисления целочисленных констант в рамках стандарта C++03. Последний стандарт C++14 ввёл новый вид шаблона: шаблоны переменных. Данная конструкция не только позволяет писать то же самое намного более коротко, но и использовать практически произвольные типы данных.

Пример чисел Фибоначчи, вычисляемых во время компиляции в плавающей точке:

template <unsigned N>
const double fib = fib<N-1> + fib<N-2>;

template <>
const double fib<0> = 0.;

template <>
const double fib<1> = 1.;

Конечно, если константы требуется использовать как параметры шаблонов (для чего, собственно, они и вычисляются во время компиляции), то они должны быть целочисленными и, кроме того, использовать ключевое слово constexpr:

template <unsigned N>
constexpr unsigned long long fib = fib<N-1> + fib<N-2>;

template <>
constexpr unsigned long long fib<0> = 0;

template <>
constexpr unsigned long long fib<1> = 1;

Борьба с переполнением и исчезновением при вычислении произведений

Пусть дан массив чисел с плавающей запятой, требуется найти произведение.

double product(const double x[], size_t n)
{
  double p = 1.;
  for (size_t i = 0; i < n; ++i)
    p *= x[i];
  return p;
}

Какие недостатки у приведённого выше очевидного решения?

Из трёх компиляторов (g++ 7.1, clang 4.0, icc 17.0) векторизацию данного цикла (параметры -O3 -mavx2) выполнил только icc, что лишь подтверждает известный тезис о том, что icc порой действует чересчур агрессивно: операция умножения в плавающей точке не ассоциативна.

Для любителей обобщённого кода можно записать то же самое, но для «любых чисел»:

#include <iterator>
#include <functional>
#include <numeric>

template <class It>
auto product(It from, It to)
{
  return std::accumulate(from, to,
    typename std::iterator_traits<It>::value_type(1),
    std::multiplies<>{});
}

Потенциальная проблема прямого «наивного» перемножения заключается в возможном переполнении или исчезновении промежуточного результата, в то время как реальный окончательный результат был бы представим в диапазоне используемых чисел (например, перемножаем много больших чисел, затем много маленьких, или наоборот).

Одно из возможных решений — складывать логарифмы, затем сумму возвести в степень:

\prod\limits_i{x_i} = \exp\sum\limits_i{\log{x_i}}.

Данный способ, вероятно, лучше подходит для работы с числами с фиксированной запятой (на практике, естественно, лучше использовать логарифм и возведение в степень по основанию 2).

В случае плавающей запятой на основе средств стандартной библиотеки C доступен способ получше: комбинация frexp и ldexp.

Лучше он как с точки зрения скорости (по сути это просто битовые операции над представлением числа), так и с точки зрения точности (не привносит дополнительной погрешности по сравнению с обычным произведением).

#include <cmath>
#include <iterator>

template <class It, class NT>
NT product(It from, It to, NT first)
{
  using std::frexp;
  using std::ldexp;
  int e = 0;
  NT p = frexp(first, &e);
  for (int e2 = 0; from != to; ++from)
  {
    p = frexp(p * *from, &e2);
    e += e2;
  }
  return ldexp(p, e);
}

template <class It> inline
auto product(It from, It to)
{
  return product(from, to,
    typename std::iterator_traits<It>::value_type(1));
}

Или то же самое в виде функтора:

template <class NT, class ET = long>
class Frexp_ldexp_product
{
  NT p;
  ET e;

public:
  using Number = NT;
  using Exponent = ET;

  explicit Frexp_ldexp_product(Number first = Number(1))
  {
    using std::frexp;
    int first_e = 0;
    p = frexp(first, &first_e);
    e = first_e;
  }

  Number operator()() const
  {
    using std::ldexp;
    return ldexp(p, e);
  }

  void operator()(Number next)
  {
    using std::frexp;
    int next_e = 0;
    p = frexp(p * next, &next_e);
    e += next_e;
  }

  Number significand() const
  {
    return p;
  }

  Exponent exponent() const
  {
    return e;
  }
};

// Пример использования: при обработке диапазона.
template <class It, class NT>
NT product(It from, It to, NT first)
{
  return std::for_each(from, to,
    Frexp_ldexp_product<NT>{first})();
}

Функтор удобнее, если значения в произведении генерируются на ходу, а не лежат в памяти. Кроме того, есть возможность получить отдельно мантиссу ∈ (0.5, 1] и порядок, даже если результирующее произведение не помещается в диапазон используемого типа чисел с плавающей точкой.

Напоследок хочется отметить, что если проблема переполнения/исчезновения снята, то рост накопленной погрешности остановить не так просто. Хотя результат произведения зависит от порядка, оценка относительной погрешности от порядка выполнения операций не зависит:

(x(1+n\varepsilon))(y(1+m\varepsilon)) \rightarrow xy(1+n\varepsilon)(1+m\varepsilon)(1+\varepsilon) \rightarrow xy(1+(n+m+1)\varepsilon).

В каком порядке не перемножай, всё равно получится относительная погрешность ne, где n — число сомножителей, e — половина «эпсилона» (ULP).

ChaiScript и mingw32

ChaiScript — сравнительно простой C++-подобный скриптовый язык, предназначенный для встраивания в C++-приложения. Недавно я попробовал его задействовать с одной из сборок mingw32 (gcc-6.3.0-64).

Итак, берём исходник helloworld с головной страницы сайта ChaiScript:

#include <chaiscript/chaiscript.hpp>

std::string helloWorld(const std::string &t_name) {
  return "Hello " + t_name + "!";
}

int main() {
  chaiscript::ChaiScript chai;
  chai.add(chaiscript::fun(&helloWorld), "helloWorld");

  chai.eval(R"(
    puts(helloWorld("Bob"));
  )");
}

и вперёд, с песней… а откуда столько ошибок компиляции?

error: 'mutex' in namespace 'std' does not name a type

что, правда что ли? А как же полная поддержка C++14 в gcc 6.3 и всё такое? Открываем стандартный заголовочный файл mutex и видим

#ifdef _GLIBCXX_HAS_GTHREADS

  // Common base class for std::recursive_mutex and std::recursive_timed_mutex
  class __recursive_mutex_base
  {
...

Легко убедиться, что макрос _GLIBCXX_HAS_GTHREADS, оказывается, не определён.
Что делать? Ответы со StackOverflow:

  1. на Windows использовать MSVC (ага, уже);
  2. использовать Boost.Thread или ещё какую-нибудь стороннюю библиотеку (ну да, и перелицовывать на неё весь чужой код);
  3. найти другую сборку mingw32 с каким-то образом прикрученными C++threads, например, TDM GCC (да, только вот TDM GCC застрял на версии компилятора 5.1… вообще странно — OpenMP работает, а C++threads не шмогла?..);
  4. сделать свою сборку mingw32 с каким-то образом прикрученными C++threads (спасибо, может быть когда-нибудь потом…);
  5. использовать drop-in замену на основе Win32, например, mingw-std-threads (о, это может сработать!).

Итак, копирую файлы mingw.thread.h, mingw.mutex.h, mingw.condition_variable.h. Добавляю инклюды перед чайскриптом:

#include <mingw.thread.h>
#include <mingw.mutex.h>
#include <mingw.condition_variable.h>
#include <chaiscript/chaiscript.hpp>

Компилируем… Опять куча ошибок, например,

mingw.condition_variable.h 192 error: 'unique_lock' has not been declared

Ну ладно, это всё тот же mutex, добавим и его инклюд после mingw.mutex.h. И полюбуемся на новую кучу ошибок, например:

std_mutex.h 132 error: redefinition of 'struct std::defer_lock_t'

Ок, поменяем местами инклюды мьютексов. Компилирум… И наблюдаем ошибки снова в коде ChaiScript, например:

chaiscript_stdlib.hpp 57 error: invalid use of incomplete type 'class std::future'

Srsly? Incomplete type? Хммм… А где вообще этот future определён? Во future? Лезем туда и любуемся на старый добрый _GLIBCXX_HAS_GTHREADS, закрывающий основную часть кода (включая async, кстати). А не обмануть ли тебя? Попробуем.

#ifndef _GLIBCXX_HAS_GTHREADS
#include <mingw.thread.h>
#include <mutex>
#include <mingw.mutex.h>
#include <mingw.condition_variable.h>

#define _GLIBCXX_HAS_GTHREADS
#include <future>
#undef _GLIBCXX_HAS_GTHREADS
#endif // _GLIBCXX_HAS_GTHREADS

И получаем опять переопределения.

thread 61 error: redefinition of 'class std::thread'

Добавим ещё инклюдов

#ifndef _GLIBCXX_HAS_GTHREADS
#include <thread>
#include <mingw.thread.h>
#include <mutex>
#include <mingw.mutex.h>
#include <condition_variable>
#include <mingw.condition_variable.h>

#define _GLIBCXX_HAS_GTHREADS
#include <future>
#undef _GLIBCXX_HAS_GTHREADS
#endif // _GLIBCXX_HAS_GTHREADS

Число ошибок уменьшилось… Оказывается, нет стандартных определений once_flag и call_once:

future 312 error: 'once_flag' does not name a type

Т.е. в mingw-std-threads кое-чего не хватает. Стандартные определения (в mutex) использовать нельзя, так как они опираются на отсутствующую библиотеку pthreads (кстати, на всякий случай упомяну, что попытки её подключить через параметры компилятора не увенчались успехом).

Естественная мысль: а может сделать это уже руками?.. Я попробовал. Результат — одна ошибка

future 548 error: expected class-name before '{' token

фактически там встретился неопределённый идентификатор __at_thread_exit_elt (имя класса)… В общем, реализовывать ещё и замены подобных потрохов мне стало как-то резко влом.

В принципе, ChaiScript позволяет попросту отключить C++threads (и внутреннюю синхронизацию) через макрос CHAISCRIPT_NO_THREADS. Итак, выкинем всю требуху:

#define CHAISCRIPT_NO_THREADS
#include <chaiscript/chaiscript.hpp>

Компилируем… Долго… И фейл!

Fatal error: can't close obj\Debug\main.o: File too big

WTF?! StackOverflow подсказывает, что проблема в ограничениях формата экзешников PE/COFF (от «coffin», что ли?). К счастью, в этом конкретном случае она решается включением оптимизации (хотя бы -Og или -O2), благодаря чему размер объектов получается меньше. Собираем. Три минуты (helloworld? да, здесь у меня не шибко быстрый процессор, но всё же…)? Exe в 40Мб? Круто, чё. Но оно наконец-то работает — сообщение выводится.

Справедливости ради: на -O3 сборка занимает около одной минуты и exe имеет размер 2.7Мб.

Для «надёжности» можно добавить ещё один параметр компилятора: -Wa,-mbig-obj (именно так, через запятую без пробелов). Но как-то всё это не воодушевляет, да.

Conar 6: Parser — окончание

Начало

Предыдущая часть

Последнее, что осталось «дописать» — собственно организация процесса разбора, а именно: поиск ключей и выборка параметров. При поиске ключей проверяем четыре варианта: «ключ» «…» «другой ключ», «ключ=значение», «-ключключ…ключ» и «ключзначение» (именно в таком порядке). Третий вариант применяется только для однобуквенных ключей.

Все нижеприведённые функции являются членами класса Parser.

Функция, проверяющая параметр на удовлетворение первому варианту («ключ»):

bool is_direct_key(const std::string &param) const
{
  return mapper.count(param) != 0;
}

Функция, проверяющая параметр на удовлетворение второму варианту («ключ=значение»), возвращает позицию знака «=»:

size_t is_key_equals(const std::string &param) const
{
  const auto eq_pos = param.find('=');
  if (eq_pos == std::string::npos)
    return 0;

  return is_direct_key(param.substr(0, eq_pos))? eq_pos: 0;
}

Функция, проверяющая параметр на удовлетворение третьему варианту («-ключключ…ключ»):

bool is_key_sequence(const std::string &param) const
{
  const auto psz = param.size();
  if (psz < 3 || param[0] != '-')
    return false;

  char key[3] {'-'};
  for (size_t i = 1; i < psz; ++i)
  {
    key[1] = param[i];
    if (mapper.count(key) == 0)
      return false;
  }

  return true;
}

Наконец, функция, проверяющая параметр на удовлетворение четвёртому варианту («ключзначение»), возвращает позицию первого символа значения:

size_t is_key_value(std::string param) const
{
  while (!param.empty())
  {
    param.pop_back();
    if (is_direct_key(param))
      return param.size();
  }

  return 0;
}

Данная функция проверяет «с конца», пытаясь подобрать самый длинный ключ, что может быть нежелаемым поведением (например, мы хотим использовать в такой манере только однобуквенные ключи — такую опцию я хочу добавить позже).

Используя определённые выше функции можно набросать скелет цикла обработки последовательности параметров:

  Unparsed_parameters unparsed;
  Parameters params(from, to);
  for (auto key = begin(params), params_end = end(params),
            prev_key = params_end;
    key != params_end;)
  {
    // Check all possible variants one by one.
    if (is_direct_key(*key))
    {
      // The only situation in which we may have a sequence
      // of argument parameters after the key.

      prev_key = key;
    }
    else if (const auto eq_pos = is_key_equals(*key))
    {    }
    else if (is_key_sequence(*key))
    {    }
    else if (const auto val_pos = is_key_value(*key))
    {    }
    else
    {
      // Nope, this is not a valid parameter of any form.
      unparsed.emplace_back(move(*key));
      ++key; // next...
    }
  }
  return unparsed;

Если сработало условие is_direct_key(*key), то надо найти следующий ключ, а всё что между передать диапазоном обработчику опции, соответствующей *key. Первая мысль: запоминать последнюю такую позицию (итератор prev_key) и вызывать обработчик, когда был найден правый конец. Недостаток: это надо делать под каждым if’ом, хотелось бы избежать такого дублирования кода…

Подумав, как можно разобраться с данной ситуацией, я остановился на варианте «умного» кода (здесь «умный» — не комплимент), а именно, использовать приём, известный как «хэширование»:

    size_t pos;
    if
    (
      int situation =
        is_direct_key(*key)?          1
      : (pos = is_key_equals(*key))?  2
      : is_key_sequence(*key)?        3
      : (pos = is_key_value(*key))?   4
      : 0 // --> else
    )
    {   }
    else
    {
      // Nope, this is not a valid parameter of any form.
      unparsed.emplace_back(move(*key));
      ++key; // next...
    }

Теперь я могу избежать дублирования кода, а необходимая информация о случае хранится в situation и pos. Хотя нет, есть ещё случай, когда ключ — последний, дальше до конца идут только его аргументы… Продублировать вызов обработчика после цикла? Думаю, лучше ситуацию «конец параметров» тоже занести в situation:

  for (auto key = begin(params), params_end = end(params),
            prev_key = params_end;;)
  {
    size_t pos;
    if
    (
      int situation =
        key == params_end?            1
      : is_direct_key(*key)?          2
      : (pos = is_key_equals(*key))?  3
      : is_key_sequence(*key)?        4
      : (pos = is_key_value(*key))?   5
      : 0 // --> else
    )

Интересно, как можно было бы записать то же самое чище и без дублирования кода…

Определим для удобства записи функцию, «вытаскивающая» ссылку обработчик по ключу:

template <class Key>
Handler& handler_of(const Key &key)
{
  return std::get<0>(options[mapper.at(key)]);
}

Теперь я могу, наконец, записать, что делается под if’ом:

  if (prev_key != params_end)
  {
    for (auto last = handler_of(*prev_key)
                     (prev_key + 1, key);
         last != key; ++last)
    { unparsed.emplace_back(move(*last)); }
    prev_key = params_end;
  }

  switch (situation)
  {
  case 1:
    return unparsed;

  case 2:
    prev_key = key;
    break;

  case 3: case 5:
    {
      const int has_equals = situation == 3;
      auto k = key->substr(0, pos);
      key->erase(0, pos + has_equals);
      if (handler_of(k)(key, key + 1) == key)
        unparsed.emplace_back(move(*key));
    }
    break;

  case 4:
    for (size_t i = 1, sz = key->size(); i < sz; ++i)
    {
      char k[3] { '-', (*key)[i] };
      handler_of(k)(params_end, params_end);
    }
  }

Красотой этот код не отличается, увы. И надо написать для него тесты.

В результате тестирования была исправлена ветка else:

    else
    {
      // Nope, this is not a valid parameter of any form.
      if (prev_key == params_end)
        unparsed.emplace_back(move(*key));
    }

При тестировании со строковыми параметрами обнаружилась старая ошибка: оператор >> читает «словами», однако в случае строковых параметров их надо не «читать» (с пробельными символами-разделителями, да), а просто «отдавать» целиком. Для этого я заменю явное чтение с помощью >> на вызов новой функции impl::read_value:

/// Tries to read a generic val from string str.
template <class S, class T>
inline bool read_value
  (std::istringstream &reader, S &&str, T &val)
{
  reader.str(std::forward<S>(str));
  return reader >> val;
}

/// Forwards a string instead of reading it.
inline bool read_value
  (std::istringstream&, std::string &&str, std::string &val)
{
  val = std::move(str);
  return true;
}

Забавно, что после этой замены сработал тест value_out (код ошибки 7) — из-за того, что теперь «чтение» забирает строки с помощью move, а не копирует их.

Вместе с тестами таки вылезло за пределы 1000 строк. В итоге, я не решился использовать C++17 (только static_assert без сообщения, а можно было декомпозиционное определение и string_view, например).

Код данного варианта целиком.

О затенении перегрузок функций

Недавно вновь столкнулся с ошибкой компиляции, о причинах которой к тому времени я уже успел забыть. Пришлось вспомнить :)

Рассмотрим простенький код: определим оператор вывода вектора в поток.

#include <iostream>
#include <iterator>
#include <vector>
using namespace std;

ostream& operator<<(ostream &os, const vector<int> &v)
{
	for (auto x: v)
		os << x << ' ';
	return os;
}

int main()
{
	cout << vector<int>{ 1, 2, 3, 4 } << endl;
	return 0;
}

Он работает так, как и предполагается, никаких проблем нет. Range-for — удобная штука: до его появления обход контейнера с помощью итераторов выглядел страшненько. Вместо явного цикла можно было воспользоваться возможностями STL: алгоритмом copy и адаптером потока ostream_iterator.

#include <iostream>
#include <iterator>
#include <vector>
#include <algorithm>
using namespace std;

ostream& operator<<(ostream &os, const vector<int> &v)
{
	ostream_iterator<int> oi(os, " ");
	copy(begin(v), end(v), oi);
	return os;
}

int main()
{
	cout << vector<int>{ 1, 2, 3, 4 } << endl;
	return 0;
}

Этот вариант работает точно так же, как и предыдущий. Всё в порядке.

Попробуем обобщить его, например, на вектор векторов…

#include <iostream>
#include <iterator>
#include <vector>
#include <algorithm>
using namespace std;

ostream& operator<<(ostream &os, const vector<int> &v)
{
	ostream_iterator<int> oi(os, " ");
	copy(begin(v), end(v), oi);
	return os;
}

ostream& operator<<
	(ostream &os, const vector<vector<int>> &m)
{
	ostream_iterator<vector<int>> oi(os, "\n");
	copy(begin(m), end(m), oi);
	return os;
}

int main()
{
	using VI = vector<int>;
	vector<VI> m
	{
		VI{ 1, 2, 3, 4 },
		VI{ 5, 6, 7, 8 },
		VI{ 9, 0, 1, 2 },
	};
	
	cout << m << endl;
	return 0;
}

И вот она, ошибка… Этот код уже не компилируется. Как часто бывает в таких случаях, компилятор может высыпать ворох маловразумительных ошибок. Проблема в операторе <<. Точнее, в том факте, что в том же пространстве имён std, где определены ostream_iterator и vector, определён и набор стандартных перегрузок operator<<, который «затеняет» operator<< из глобального пространства имён.

Простейший способ избавиться от ошибки — поместить свою версию operator<< также в пространство имён std. Например, такой код уже благополучно компилируется и работает правильно:

#include <iostream>
#include <iterator>
#include <vector>
#include <algorithm>
using namespace std;

namespace std
{
	ostream& operator<<(ostream &os, const vector<int> &v)
	{
		ostream_iterator<int> oi(os, " ");
		copy(begin(v), end(v), oi);
		return os;
	}
}

ostream& operator<<
	(ostream &os, const vector<vector<int>> &m)
{
	ostream_iterator<vector<int>> oi(os, "\n");
	copy(begin(m), end(m), oi);
	return os;
}

int main()
{
	using VI = vector<int>;
	vector<VI> m
	{
		VI{ 1, 2, 3, 4 },
		VI{ 5, 6, 7, 8 },
		VI{ 9, 0, 1, 2 },
	};
	
	cout << m << endl;
	return 0;
}

Я намеренно не стал помещать второй operator<< внутрь namespace std, чтобы показать, что здесь проблема была с первым operator<<, который должен вызываться некоей функцией, определённой в std.

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

Следующий пример определяет рекурсивную версию operator<<, которую тоже приходится поместить в namespace std:

#include <iostream>
#include <iterator>
#include <vector>
#include <algorithm>
using namespace std;

namespace std
{
	template <class T>
	ostream& operator<<(ostream &os, const vector<T> &m)
	{
		ostream_iterator<T> oi(os.put('{'), " ");
		copy(begin(m), end(m), oi);
		return os.put('}');
	}
}

int main()
{
	using VI = vector<int>;
	vector<VI> m
	{
		VI{ 1, 2, 3, 4 },
		VI{ 5, 6, 7, 8 },
		VI{ 9, 0, 1, 2 },
	};
	
	cout << m << endl;
	return 0;
}

Вообще, «по умолчанию» считается плохим тоном вносить собственные определения в пространство имён std (есть исключения, например, частные специализации hash).

Для операций с собственными типами, определёнными в собственном пространстве имён, возможно положиться на ADL, поместив перегрузки операторов в то же самое пространство имён. В вышеприведённых примерах использовался std::vector, а определение операций для стандартных классов приводит к опасности коллизии разных определений для одних и тех же шаблонных классов, появившихся в разных местах проекта. Лично мне бы хотелось иметь в C++ «strong typedef», тем более, что для него даже новых ключевых слов вводить не надо:

namespace my
{
	using Data = new vector<int>;
	// Теперь Data != vector<int> и
	// определён в пространстве имён my.
	ostream& operator<<(ostream &os, const Data &d);
}

Мечты-мечты.

Conar 5: структура данных Parser

Начало

Предыдущая часть

Займёмся классом Parser.

При добавлении опции требуется строить отображение ключ -> (обработчик, инфо). Чтобы не дублировать пары (обработчик, инфо) (что важно для «обработчика», который может хранить изменяемое состояние), можно или оборачивать их в shared_ptr и отображать ключ в shared_ptr, или складывать в отдельный контейнер и отображать ключ в индекс или итератор. Второй вариант мне нравится больше.

Обработчик будем приводить к типу Handler, объявленному следующим образом:

using Parameters = std::vector<std::string>;
using PIt = Parameters::iterator;
using Handler = std::function<PIt(PIt, PIt)>;

Теперь определим собственно структуры данных, содержащиеся в объекте Parser. Mapper отображает ключи в индексы вектора Options.

using Option = std::tuple<Handler, Option_info>;
using Options = std::vector<Option>;
using Mapper = std::map<std::string, std::size_t>;
Options options;
Mapper mapper;

Напишем код, добавляющий опцию (тройку «обработчик, ключи, инфо«).

/// Add an option description.
/// Use functions flag, value and seq to make an option.
template <class Opt>
Parser& operator()(Opt option)
{
  using namespace std;
  static_assert(tuple_size<Opt>::value == 3);

  // Add the option to options.
  const auto index = options.size();
  options.emplace_back(move(get<0>(option)),
                       move(get<2>(option)));

  // Register keys of the option.
  for (auto &key: get<1>(option))
    mapper.emplace(key_prep(move(key)), index);

  return *this;
}

Функция key_prep добавляет знаки «-» в начало строки:

static std::string key_prep(std::string &&key)
{
  if (!key.empty() && std::isalnum(key[0]))
  {
    if (key.size() == 1)
      key.insert(0, 1, '-');
    else
      key.insert(0, "--");
  }

  return move(key);
}

Теперь можно заполнить тело функции help, которую можно использовать и в целях тестирования. Подумав, я решил убрать параметр «ширина строки», поскольку это не очень нужное усложнение. Реальную ширину строки консоли нельзя получить стандартными методами, а насильная вставка перевода строки может приводить к появлению избыточных пустых строк. Функция получилась длинная, но простая:

/// Construct a help message.
std::string help(
  std::size_t key_column_width = 20,
  const char *line_sep = "\n\n") const
{
  using namespace std;
  ostringstream write;
  write.setf(ios::left);

  // Option index -> first key (iterator).
  vector<Mapper::const_iterator> key_by_index(options.size());
  for (auto p = mapper.end(), b = mapper.begin();;)
  {
    --p;
    key_by_index.at(p->second) = p;
    if (p == b)
      break;
  }

  // Make the text.
  for (auto p = mapper.begin(), e = mapper.end(); p != e; ++p)
  {
    if (p->first.size() < key_column_width)
    {
      write.width(key_column_width);
      write << p->first;
    }
    else
    {
      write.width(0);
      write << p->first << '\n';
    }

    const auto pivot = key_by_index[p->second];
    if (pivot == p)
      write << get<1>(options[p->second]);
    else
      write << ("See " + pivot->first + ".");
    write << line_sep;
  }

  return write.str();
}

Самую сложную часть парсера — собственно разбор списка параметров — из-за недостатка времени в данный момент я оставил для следующей части.

Код данного варианта целиком.

Продолжение следует…

Conar 4: больше SFINAE

Начало

Предыдущая часть

Займёмся функциями, обеспечивающими чтение последовательностей. Первой такой функцией у меня будет не seq, а вариант value, который записывает каждое распознанное значение по итератору вывода.

template <class T, class OutIt> inline
auto value(OutIt dest,
  Possible_keys keys, Option_info info = {},
  std::enable_if_t<
    std::is_base_of<
      std::output_iterator_tag,
    typename std::iterator_traits<OutIt>::iterator_category>
      ::value, int> = 0)
{
  using std::move;
  return value_parser<T>(
    [dest](T val) mutable { *dest++ = std::move(val); },
    move(keys), move(info));
}

Мелкая неприятность: пользователь должен явно указывать тип считываемых значений T. В принципе, можно сделать два варианта: аналогичный приведённому выше и выводящий тип значения с помощью iterator_traits::value_type. К сожалению, стандартные итераторы вывода, реализующие вставку в контейнеры, определяют value_type как void, делая тип элементов контейнера невыводимым для функции value. Кроме того, заметим, что forward_iterator_tag наследует от input_iterator_tag, но не от output_iterator_tag.

Введём простенькие вспомогательные определения.

/// Check if iterator It is compatible with category Cat.
template <class It, class Cat>
using Has_iterator_category =
  is_base_of<Cat,
    typename iterator_traits<It>::iterator_category>;

/// Check if iterator It is an output iterator.
template <class It>
using Is_output_iterator =
  Has_iterator_category<It,
    output_iterator_tag>;

/// Check if iterator It is a forward iterator.
template <class It>
using Is_forward_iterator =
  Has_iterator_category<It,
    forward_iterator_tag>;

/// Check if iterator It declares void as its value type.
template <class It>
using Has_void_value =
  is_same<void,
    typename iterator_traits<It>::value_type>;

Теперь два варианта value, пишущих через итераторы можно определить следующим образом:

/// Construct a value option, which writes values through an iterator.
/// Explicitly typed variant when the type is not derivable from OutIt.
template <class T, class OutIt> inline
auto value(OutIt dest,
  Possible_keys keys, Option_info info = {},
  std::enable_if_t<
    impl::Is_output_iterator<OutIt>::value &&
    impl::Has_void_value<OutIt>::value,
    int> = 0)
{
  using std::move;
  return value_parser<T>(
    [dest](T val) mutable { *dest++ = std::move(val); },
    move(keys), move(info));
}

/// Construct a value option, which writes values through an iterator.
template <class OutIt> inline
auto value(OutIt dest,
  Possible_keys keys, Option_info info = {},
  std::enable_if_t<
    (impl::Is_output_iterator<OutIt>::value ||
     impl::Is_forward_iterator<OutIt>::value) &&
    !impl::Has_void_value<OutIt>::value,
    int> = 0)
{
  using std::move;
  using T = typename std::iterator_traits<OutIt>::value_type;
  return value_parser<T>(
    [dest](T val) mutable { *dest++ = std::move(val); },
    move(keys), move(info));
}

Тестирование этих определений вскрыло небольшое упущение: value_parser также должна создавать mutable-замыкание, если мы хотим упаковать в него своё mutable-замыкание, что требуется, в свою очередь, из-за необходимости передвигать итератор на каждом вызове и сохранять его изменённую версию между вызовами. Соответственно, я также добавил mutable в flag_parser и seq_parser.

Теперь можно перейти к seq.

Вариант seq, заполняющий условно бесконечную последовательность по итератору, написать не составит большого труда. Полная аналогия с вариантами value, указанными выше. Как и в случае с value, предполагается, что «места хватит» (интересно, сколько параметров командной строки можно скормить в распространённых ОС/интерпретаторах?).

Приведу одну из них. Вторая абсолютно аналогична.

template <class T, class OutIt> inline
auto seq(OutIt dest,
  Possible_keys keys, Option_info info = {},
  std::enable_if_t<
    impl::Is_output_iterator<OutIt>::value &&
    impl::Has_void_value<OutIt>::value,
    int> = 0)
{
  using std::move;
  return seq_parser<T>(
    [dest](T *from, T *to) mutable
    {
      while (from != to)
        *dest++ = std::move(*from++);
    },
    move(keys), move(info));
}

Для удобства пользователя можно добавить ещё три варианта seq, которые принимают непосредственно контейнеры, извлекают тип элементов через вложенное объявление value_type и выбирают способ вставки (push_back, insert или push) в зависимости от того, какой подходит. Для различения между push_back и insert (линейными и ассоциативными контейнерами) нам в качестве заготовки пригодится SFINAE-костыль Is_incrementable. Последний вариант seq подразумевает вставку в адаптер вроде std::stack с помощью вызова push.

/// Check if container C supports push_back(T&&).
template <class C, class T>
struct Supports_push_back
{
  template <class Y>
  static true_type check(decltype((void)(
    declval<Y>().push_back(declval<T>())
  ), 1));

  template <class Y>
  static false_type check(unsigned);

  using type = decltype(check<C&>(1));
  static constexpr bool value = type::value;
};

Два других детектора (Supports_insert и Supports_push) полностью аналогичны (да, макрос сюда так и просится). Вспомогательные определения для краткости записи сигнатур функций seq:

template <class C>
using Use_push_back =
  Supports_push_back<C, typename C::value_type>;

template <class C>
using Use_insert = bool_constant<
  Supports_insert<C, typename C::value_type>::value &&
  !Use_push_back<C>::value>;

template <class C>
using Use_push = bool_constant<
  Supports_push<C, typename C::value_type>::value &&
  !Use_push_back<C>::value &&
  !Use_insert<C>::value>;

Теперь вариант seq, который заполняет контейнер, предоставляющий функцию push_back, может быть записан следующим образом:

/// Construct a seq option,
/// which appends a "push_back" container.
template <class Container> inline
auto seq(Container &cont,
  Possible_keys keys, Option_info info = {},
  std::enable_if_t<
    impl::Use_push_back<Container>::value,
    int> = 0)
{
  using std::move;
  using T = typename Container::value_type;
  return seq_parser<T>(
    [&cont](T *from, T *to) mutable
    {
      while (from != to)
        cont.push_back(std::move(*from++));
    },
    move(keys), move(info));
}

Ещё два варианта seq полностью аналогичны.

Код данного варианта целиком.

Продолжение следует.