Корень из целого числа в compile-time часть 2

Не так давно в этом блоге мелькала статья на тему вычисления квадратного корня из целого числа в compile-time. Сейчас будет небольшое продолжение.

Некоторое время назад постоянный читатель этого блога прислал мне ссылку на задачку. Сайт Insidecpp иногда очень интересно почитать, и даже пафос и надменность автора как-то можно перетерпеть, когда видишь сколько он там всего понаписал. Прореагировал я примерно так: "чтобы моё решение подошло к этой задаче, нужно исправить всего 1 строчку". А потом задумался.

Автор предлагает реализовать вычисление квадратного корня в compile-time, но только с формулировкой "найти целое, ближайшее к правильному ответу", а не "наибольшее целое, квадрат которого не превышает":

Напишите вычисление квадратного корня в compile-time-е. Чтобы усложнить задачу я предлагаю следующее. Если квадратный корень от заданного числа не может быть вычислен в целых, то вашим результатом должно быть ближайшее целое к реальному результату. Обратите внимание, что речь идет именно о ближайшем целом. То есть не ближайшее меньшее результата, а просто ближайшее. Оно может быть как больше, так и меньше реального результата.

Т.е., если в моей задаче корень из 3 - это решительно 1, то в задаче автора - это должно быть 2. Решение автора меня очень порадовало своей краткостью, но к сожалению это решение при ближайшем рассмотрении оказалось неверным. Сейчас я расскажу как поправить код из предыдущей статьи, чтобы получить решение для задачки с insidecpp.ru.

Пусть дано:
y * y = x
0 < a < y < b
a и b - это целые числа, ближайшие к y. Нам нужно из a и b выбрать число, которое отличается от y минимально.
Т.е. d1=y-a, d2=b-y
Если d1 < d2, тогда нужно выбрать a, если d1 > d2, тогда - b.

т.е. если буквой z обозначить знак, который ставится между d1 и d2 при сравнении, то получим:
d1 z d2
y-a z b-y
2*y z a+b
4*y*y z (a+b)*(a+b)
как известно, y*y это x, т.е.:
4*x z (a+b)*(a+b)
т.к. в нашем случае вместо a и b используются соответственно lo и lo+1, результат:
4*x z 4*lo*lo + 4*lo + 1
или, чтобы избежать переполнения,
4*(x-lo*lo) z 4*lo + 1
если 4*(x-lo*lo) > 4*lo + 1, выбераем - lo+1. Иначе - lo.

Единственное место в коде, которое необходимо поправить - это:

template<unsigned x, unsigned lo, unsigned y>
struct sqrt<x, lo, 1, y, less>
{  
  enum { value = lo * lo + 2 * lo <= x - 1 ? (lo + 1) : lo };  
};

Меняем логику на:

template<unsigned x, unsigned lo, unsigned y>
struct sqrt<x, lo, 1, y, less>
{  
  enum { value = 4*(x-lo*lo) > 4*lo + 1 ? lo+1:lo };  
};

Задача решена.