Can anyone come up with a nice way of computing a solution to the linear diophantine equation $ax + by = c$where $a,b,c \in \mathbb{Z}$ and $\gcd(a,b) \mid c$, such that all the calculations are carried out without any of the intermediate results exceeding $\max\{|a|, |b|, |c|\}$ in absolute value?I.e. if $a,b,c$ are 64-bit integers, all the calculations should be done using 64-bit integers only.The standard method is of course to use Euclidean algorithm to first find a solution to $ax + by = \gcd(a,b)$ and then multiplying $x$and $y$by $\frac{c}{\gcd(a,b)}$. The Euclidean algorithm is fine, but this last multiplication might go out of range.
Avoiding overflow in finding a solution to$ax + by = c$
ds.algorithms
Here is a modified version of the Euclidean algorithm that is stable.When considering the equation$$ax + by = c \quad (1)$$we may w.l.o.g. assume that $a, b, c \ge 0$, $b \ge 1$ and $b \ge a$. From now on we will only consider such equations. Moreover, we will set $d = \gcd(a,b)$.We call a solution $(x,y)$ good if $0 \le x < \frac{b}{d}$. Such a solution always exists, since all solutions are parametrized by $(x + t \cdot \frac{b}{d}, y - t \cdot \frac{a}{d})$.If $(x,y)$ is a good solution, then$$\frac{c}{b} - \frac{a}{d} < y \le \frac{c}{b} \quad (2)$$This follows simply from the bounds on $x$ and the fact that $y = \frac{c}{b} - \frac{ax}{b}$. Notice that in particular $|y| \le \max(a,c)$.We will now solve the problem by induction. Assume first that $a=0$. Then $x=0$, $y = \frac{c}{b}$ is a good solution. Let us then suppose by induction that we can obtain good solutions for all coefficients smaller than $a,b,c$. We let $(x',y')$ be a good solution to the equation$$(b - \left\lfloor \frac{b}{a} \right\rfloor a) x' + a y' = c \quad (3)$$Then $(\tilde{x}, \tilde{y}) = (y' - \left\lfloor \frac{b}{a} \right\rfloor x',  x')$ is a solution to (1). First of all we can calculate $\left\lfloor \frac{b}{a} \right\rfloor x'$ without overflow since $0 \le x' < \frac{a}{d}$. We note that (2) applied to the equation (3) gives us$$\frac{c}{a} - \frac{b - \left\lfloor \frac{b}{a}\right\rfloor a}{d} < y' \le \frac{c}{a}.$$Hence it follows that$$\frac{c}{a} - \frac{b}{d} < y' - \left\lfloor \frac{b}{a} \right\rfloor x' \le \frac{c}{a} \quad (4)$$Thus $|\tilde{x}| \le \max(b,c)$ and $\tilde{x}$ can be computed without overflow.Finally it remains to normalize $(\tilde{x}, \tilde{y})$ to a good solution. This involves computing $x = \tilde{x} + t \cdot \frac{b}{d}$ and $y = \tilde{y} - t \cdot \frac{a}{d}$ where $t = \left\lceil \frac{-\tilde{x}}{b/d} \right\rceil$. Assume first that $\tilde{x} < 0$. From (4) we see that then $t=1$ and everything can be done without overflow. Assume then that $\tilde{x}\ge 0$. Then $|t| \le \frac{\tilde{x}}{b/d}$. Because $b \ge a$, both products $t \cdot \frac{b}{d}$ and $t \cdot \frac{a}{d}$ can be calculated without overflow. We are done.