Finding minimum residue of a linear function in O(log M) time
Разница между en1 и en2, 88 символ(ов) изменены
I saw [a blog (that I can no longer find](https://codeforces.me/blog/entry/90585) (thanks to [user:bugman,2021-05-14] for finding it!) earlier today asking how to solve the following problem:↵

_Given a linear function $ax + b$, find the minimum $x \geq 0$ such that $ax + b \in [L, R]\ (\text{mod } M)$_.↵

To solve that problem, we can make the following reduction: If $gcd(a, M) > 1$, we divide everything by the GCD. The $x$ for which $ax + b = L\ (\text{mod } M)$ is $(L - b) a^{-1}\ (\text{mod } M)$. Denote this value by $b_0$. Then, the minimum $x$ to get to $L + y$ is $a^{-1} y + b_0\ (\text{mod } M)$. This gives us a reduced problem:↵

_Find the minimum value of $ay + b \text{ mod } M$ over $y \leq k$_.↵

This seemed pretty hard, but surprisingly I figured out how to do it in $\mathcal{O}(\log M)$ time! The algorithm is as follows:↵

In every step, we reduce the modulo from $M$ to $\min(a, M - a) \leq M / 2$. This guarantees we do at most $\mathcal{O}(\log M)$ steps.↵

To reduce it to $a$, we consider the first value $s$ among $[0, a)$ achieved by $ay + b$. If $b < a$, it is $b$. Otherwise, it is $b - M\text{ mod } a$. We check in $\mathcal{O}(1)$ if we reach $s$ for some $y \leq k$. If we do, set $M$ to $a$, $a$ to $-M \text{ mod } a$ and $b$ to $s$. Otherwise, output $b$ as the first $a$ values are the only values such that the previous value was larger, thus if we never attain them, the first value we attain is the largest we do.↵

To reduce it to $M - a$, we consider the first value $s = b \text{ mod } M - a$ among $[0, M - a)$ achieved by $ay + b$. If we reach $s$ for some $y \leq k$, set $M$ to $M - a$, $a$ to $a \text{ mod } M - a$ and $b$ to $s$. Otherwise, binary search the smallest value $b - t(M - a)$ that we attain for some $y \leq k$ and return it.↵

<spoiler summary="code">↵

~~~~~↵

#include <bits/stdc++.h>↵
using namespace std;↵
using ll = long long;↵

pair<ll, ll> extEucMod(ll a, ll b, ll p) {↵
if (b == 0) return {1, 0};↵
ll m = a / b;↵
auto sub = extEucMod(b, a - b * m, p);↵
return {sub.second, (sub.first - m*sub.second) % p};↵
}↵
ll modInv(ll a, ll p) {↵
ll res = extEucMod(p, a, p).second;↵
return (res < 0 ? res + p : res);↵
}↵
ll gcd(ll a, ll b) { return (b == 0 ? a : gcd(b, a % b)); }↵
ll mSub(ll a, ll b, ll m) { return (a >= b ? a - b : a - b + m); }↵
ll posMod(ll a, ll m) { ll res = a % m; return (res < 0 ? res + m : res); }↵
ll getSteps(ll t, ll ia, ll b, ll m) { return mSub(t, b, m) * ia % m; }↵

// Returns minimum value of ax + b (mod m) for x \in [0, k]. O(log m) time↵
ll minRem(ll a0, ll b0, ll m0, ll k) {↵
ll g = gcd(a0, m0);↵
if (g > 1) return g * minRem(a0 / g, b0 / g, m0 / g, k) + (b0 % g);↵
for (ll a = a0, b = b0, m = m0, ia0 = modInv(a0, m0), s; a; b = s) {↵
if (a > m-a) {↵
m = m - a;↵
s = b % m;↵
if (getSteps(s, ia0, b0, m0) > k) {↵
ll low = 1, high = (b - s) / m;↵
while(low != high) {↵
ll mid = (low + high) >> 1;↵
if (getSteps(s + m * mid, ia0, b0, m0) <= k) high = mid;↵
else low = mid + 1;↵
}↵
return s + m * low;↵
}↵
a %= m;↵
} else {↵
s = (b < a ? b : posMod(b - m, a));↵
if (getSteps(s, ia0, b0, m0) > k) return b;↵

ll tmp = posMod(-m, a);↵
m = a;↵
a = tmp;↵
}↵
}↵
return 0;↵
}↵

// Returns minimum x such that ax + b (mod m) \in [le, ri] or -1 if there is no such x. O(log m) time↵
ll firstInRange(ll a, ll b, ll m, ll le, ll ri) {↵
ll g = gcd(a, m);↵
if (g > 1) {↵
le = le / g + (le % g > b % g);↵
ri = ri / g - (ri % g < b % g);↵
a /= g; b /= g; m /= g;↵
}↵
if (le > ri) return -1; // impossible↵

ll ia = modInv(a, m);↵
ll x0 = mSub(le, b, m) * ia % m;↵
return minRem(ia, x0, m, ri - le);↵
}↵

int main() {↵
ios_base::sync_with_stdio(false);↵
cin.tie(0);↵

ll a0, b0, m0, k;↵
cin >> a0 >> b0 >> m0 >> k;↵
cout << minRem(a0, b0, m0, k) << '\n';↵

/*↵
ll a, b, m, le, ri;↵
cin >> a >> b >> m >> le >> ri;↵
cout << firstInRange(a, b, m, le, ri) << '\n';↵
*/↵
}↵
~~~~~↵

</spoiler>↵

What I wanted to ask is, is this a known problem, and is there a simpler or perhaps even faster solution to it? The problem seems fairly simple, so I doubt nobody has thought about it before.

История

 
 
 
 
Правки
 
 
  Rev. Язык Кто Когда Δ Комментарий
en5 Английский mango_lassi 2021-05-15 13:18:17 829 Remove binary search
en4 Английский mango_lassi 2021-05-15 12:45:36 1404 Reverted to en2
en3 Английский mango_lassi 2021-05-15 02:10:10 1408 Turns out the other case handled in my original solution was pointless >_>. I removed it, much simplifying the explanation and code :)
en2 Английский mango_lassi 2021-05-14 06:13:03 88 blog found
en1 Английский mango_lassi 2021-05-14 04:28:14 4162 Initial revision (published)