Post History

Standard Disclaimer: As always for this kind of thing, for production use, you should use an existing library written by experts to do as much of this as possible.

Meta Note: I also highly encourage looking into ECC (and particularly [Curve25519](https://cr.yp.to/papers.html#curve25519)) as it has a lot of usability and performance benefits. Curve25519 keys are only 32-bytes long, which makes using the public key *directly* as an identifier much more practical, and can be generated by simply taking 256 bits of random data (though you'd usually normalize with some trivial bit masking). This is *much* simpler and less error-prone than finding some large primes as RSA requires. Being able to pick a convenient modulus, like $2^{255}-19$, can allow significant simplifications and speed-ups in modular operations. The paper discusses Bernstein's rationale on this choice.

To your explicit question, I don't know if someone has implemented a constant-time multiplication using Java's `BigInteger` library. I suspect not. Really, the question is whether Java's `BigInteger` library even provides enough timing guarantees for that to be possible, either by specification or in practice. Obviously, two `BigInteger`s can take up different amounts of memory depending on their magnitudes. This, by itself, will almost certainly lead to timing differences. You could attempt to finesse this by storing $x$ as $2^N + x$ for some fixed, suitably large $N$ such that $2^N$ is larger than any *intermediate* involving computation with the $x$s would be. This (along with reductions mod $2^N$ when necessary) would *probably* guarantee a consistent memory size but, I suspect, would still not guarantee consistent times.

In practice, for a production system, you wouldn't want to use `BigInteger` anyway. For any particular algorithm, the number of bits you need to represent the numbers is fixed, and it makes more sense to use a representation that represents these with a fixed amount of memory independent of magnitude.

Here's a route you could take that would be fairly easy to implement. First, I'm going to assume the modulus is public, so we don't need to worry about timing for operations that only involve the modulus (and other public information). Let $W$ be the number of bits in the primitive words you use, e.g. 32 or 64, and $N = BW$ be the number of bits you'll use to represent numbers, e.g. 2048, 4096 for RSA or 256 for something like Curve25519. Here I'm assuming its a multiple, $B$, of the number of bits in the primitive words. In other words, we're representing numbers as degree < $B$ polynomials with coefficients in $\{0,\dots, 2^W-1\}$, i.e. $x = \sum_{k=0}^{B-1} a_k 2^{kW}$. This polynomial is simply represented by an array of length $B$ of $W$-bit words corresponding to the $a_k$. We can view these arrays as "vectors" (technically elements of a $\mathbb Z$-[module](https://en.wikipedia.org/wiki/Module_(mathematics))). (This polynomial-oriented view is very powerful and I highly recommend [Multidigit Multiplication for Mathematicians](https://cr.yp.to/papers.html#m3) by Dan Bernstein for more discussion on this.)

Let $n$ be our public modulus which we'll assume is coprime to $2^N$ though that will not matter for most of the following but is generally true in these contexts. Multiplication by $N$ mod $n$ is a linear transformation. Let $M_N$ be the matrix that represents multiplying by $2^N$ mod $n$. You can pre-compute this using the `BigInteger` library. The $i$th column of the matrix will simply $2^{N+iW}$ mod $n$ as computed by `BigInteger` and converted into the "vector" representation. Let $r_N = 2^N \pmod n$, and also $r_j = 2^j \pmod n$ for various $j < N$ such that $2^j > n$.

To multiply $x$ and $y$ mod $n$, you can first just use the [grade-school algorithm](https://en.wikipedia.org/wiki/Multiplication_algorithm) which is essentially the naive convolution algorithm with carrying to multiply $x$ and $y$ as (multi-digit) integers. The result will require a $2B$ word long array. We can view this result as $a2^N + b$. We then simply compute $M_N a + b$ which will produce a number less than $2^{N+1}$. Unfortunately, this number may be greater than or equal to $2^N$ so we need to do a step of carry. To avoid data-dependent timing, we add in $cr_N$ where $c$ is the carry bit produced by the addition, i.e. it's either $0$ or $1$, while dropping the carry. In symbols, $M_N a + b = c2^N + y = cr_N + y \pmod n$ for $y < 2^N$. Unfortunately again, this addition could carry as well, so we do the same thing. Luckily, it can't carry again. (Why not?)

**NOTE**, the result of this is a number (represented in vector form) less than $2^N$ but ***NOT*** necessarily less than $n$. The result is correct mod $n$, but isn't necessarily fully reduced. Generally, this is fine for intermediate calculations. To perform a full reduction, we use $x = a + b2^{j_0} + c2^{j_1}$ where $2^{j_0} < n < 2^{j_0 + 1}$ and $j_0 < j_1 < N$. $a$ will be a $j_0$-bit number, $b$ a $(j_1 - j_0)$-bit number, and $c$ a $(N-j_1)$-bit number. We have $x = a + b2^{j_0} + cr_{j_1} \pmod n$, but since $r_{j_1} < n$, $x < a + (b+2c)2^{j_0}$. $b+2c$ has at most $\max(j_1 - j_0, N-j_1+1)+1$ bits. If $j_1$ was chosen as roughly halfway between $N$ and $j_0$, this means it takes roughly half the bits to represent the part beyond bit $j_0$. We can iterate this a fixed number of times in a fixed pattern (dependent on $N$ and $j_0$) roughly halving the number of bits required beyond $j_0$ each time. Eventually we can finish with a fixed few rounds of subtracting by $n$ with borrow where we *add* back an $n$ if the borrow flag is set. Again, we do this in a data-independent way by *always* adding back an $n$ multiplied by the borrow flag.

You can see why doing a full reduction is not desirable for every intermediate operation. In the approach above, modular multiplication is the usual grade-school algorithm, followed by adding a matrix multiplied version of the upper part to the lower part, followed by adding a constant multiplied by the carry flag. You'd need to implement conversion to/from `BigInteger`, addition, subtraction, and scaling by a "small" scalar all with carries, but as just normal multi-digit integer arithmetic. For convenience, you could use an array of 64-bit numbers (say) but use them the represent 32-bit chunks. That is, you have $x = \sum_{k=0}^{B-1} a_k 2^{32k}$ but now $a_k$ is (inefficiently) represented as a 64-bit number. The benefit of doing this is you don't need to perform carries in your addition, subtraction, and multiplication, algorithms, you can instead have a separate `propagate_carries` (or `reduce`) function that normalizes all the coefficients to $\{0,\dots,2^{32}\}$ and performs carries then. Worst-case, you do an explicit `propagate_carries` after every operation anyway, and it just lets you write a bit simpler code. In many cases, it's easy to see that you can delay a carry until later saving some computation. Of course, this is at the cost of storing and performing 64-bit operations everywhere.

As a final warning, and another reason you shouldn't do this yourself for production code without a lot of care, the above algorithms use a fixed pattern of constant-time operations (given $N$ and $n$) so are not data-dependent. However, this doesn't stop the compiler from doing some "optimization" that *does* lead to data-dependent timings. For example, maybe the compiler realizes that the carry will only ever be $0$ or $1$ and decides to "optimize" as `if (carry == 0) return result; else return add(x, r_N);` Short of writing assembly code, it's hard to be certain something like this isn't happening.