Нормализуем деление вместо двоичного поиска

big_int.cpp

@@ -220,8 +220,9 @@ static vector<Digit> div_by_base(const vector<Digit>& magnitude)
     return vector<Digit>(magnitude.cbegin() + 1, magnitude.cend());
 }
 
-/// Делит кусок числителя на знаменатель (подбором).
-/// Числитель и знаменатель таковы, что в результате всегда одна цифра
+/// Делит кусок числителя на знаменатель.
+/// Числитель и знаменатель таковы, что в результате всегда одна цифра.
+/// К тому же аргументы должны быть нормализованы (старшая цифра знаменателя >= base / 2)
 static Digit div_chunk(const vector<Digit>& chunk, const vector<Digit>& denominator)
 {
     if (first_is_less(chunk, denominator))
@@ -231,27 +232,29 @@ static Digit div_chunk(const vector<Digit>& chunk, const vector<Digit>& denomina
            || (chunk.size() == denominator.size() + 1 // Или числитель длиннее знаменателя на 1 цифру
                && first_is_less(div_by_base(chunk), denominator))); // И эта цифра не лишняя
 
-    Digit left = 0;
-    Digit right = base - 1;
-    Digit best_digit = 0;
+    assert(denominator.back() >= base / 2);
 
-    while (left <= right)
-    {
-        Digit middle = left + (right - left) / 2;
-        vector<Digit> digit{middle};
+    // Две старшие цифры куска.
+    // Если длина куска равна длине знаменателя, то старшая цифра будет 0
+    DDigit chunk_2_digits = chunk.back();
+    if (chunk.size() != denominator.size())
+        chunk_2_digits = chunk_2_digits * base + chunk[chunk.size() - 2];
 
-        if (first_is_less(chunk, mul_magnitudes(digit, denominator)))
-        {
-            right = middle - 1;
-        }
-        else
-        {
-            best_digit = middle;
-            left = middle + 1;
-        }
+    // Делим две старшие цифры куска на старшую цифру знаменателя
+    Digit digit = (Digit)min(chunk_2_digits / denominator.back(), DDigit(base - 1));
+
+#ifndef NDEBUG
+    size_t misses_count = 0; // Число промахов
+#endif
+
+    // while (chunk < digit * denominator) // Пока цифра слишком большая
+    while (first_is_less(chunk, mul_magnitudes(vector<Digit>{digit}, denominator)))
+    {
+        --digit;
+        assert(++misses_count <= 2);
     }
 
-    return best_digit;
+    return digit;
 }
 
 /// Возвращает неполное частное и остаток (делит столбиком).
@@ -270,32 +273,46 @@ static pair<vector<Digit>, vector<Digit>> div_mod_magnitudes(const vector<Digit>
     if (first_is_less(numerator, denominator))
         return {vector<Digit>{0}, numerator};
 
+    // Нормализующий множитель
+    Digit scale = base / (denominator.back() + 1);
+    assert(scale <= base / 2); // Старшая цифра denominator не может быть 0, т.е. минимум 1
+
+    // Нормализованный числитель
+    vector<Digit> norm_numerator = (scale == 1)
+        ? numerator
+        : mul_magnitudes(vector<Digit>{scale}, numerator);
+
+    // Нормализованный знаменатель
+    vector<Digit> norm_denominator = (scale == 1)
+        ? denominator
+        : mul_magnitudes(vector<Digit>{scale}, denominator);
+
     // Неполное частное
     vector<Digit> quotient;
-    quotient.reserve(numerator.size() - denominator.size() + 1);
+    quotient.reserve(norm_numerator.size() - norm_denominator.size() + 1);
 
     // Текущий кусок числителя (и одновременно остаток)
     vector<Digit> chunk;
-    quotient.reserve(denominator.size() + 1);
+    quotient.reserve(norm_denominator.size() + 1);
 
     // Копируем цифры числителя, начиная с конца (со старших разрядов)
-    for (size_t i = numerator.size() - 1; i != size_t(-1); --i)
+    for (size_t i = norm_numerator.size() - 1; i != size_t(-1); --i)
     {
         // Копируем очередную цифру числителя в начало куска (так как цифры хранятся в обратном порядке)
-        chunk.insert(chunk.begin(), numerator[i]); // chunk = chunk * base + numerator[i]
+        chunk.insert(chunk.begin(), norm_numerator[i]); // chunk = chunk * base + norm_numerator[i]
 
         // Убираем ведущие нули (когда chunk поделился без остатка и стал равен 0, а потом добавили 0)
         while (chunk.size() > 1 && chunk.back() == 0)
             chunk.pop_back();
 
         // Делим кусок на знаменатель и получаем очередную цифру результата
-        Digit digit = div_chunk(chunk, denominator);
+        Digit digit = div_chunk(chunk, norm_denominator);
 
         // Добавляем цифру к результату
         quotient.push_back(digit);
 
         // кусок -= цифра * знаменатель ⇔ кусок %= знаменатель
-        chunk = sub_magnitudes(chunk, mul_magnitudes(vector<Digit>{digit}, denominator));
+        chunk = sub_magnitudes(chunk, mul_magnitudes(vector<Digit>{digit}, norm_denominator));
     }
 
     // Мы добавляли цифры в конец вектора, а не в начало, поэтому реверсируем
@@ -305,6 +322,14 @@ static pair<vector<Digit>, vector<Digit>> div_mod_magnitudes(const vector<Digit>
     while (quotient.size() > 1 && quotient.back() == 0)
         quotient.pop_back();
 
+    // Денормализуем остаток
+    if (scale != 1)
+    {
+        pair<vector<Digit>, vector<Digit>> denorm_chunk = div_by_digit(chunk, scale);
+        assert(denorm_chunk.second == vector<Digit>{0}); // Должно поделиться без остатка
+        chunk = denorm_chunk.first;
+    }
+
     return {quotient, chunk}; // В chunk находится остаток
 }