# Efficient Algorithm to Find the Closest Integer Representation, in the Form $A\times\frac{N}{D}$ for a Value

## The Problem

I am working on a problem that boils down to finding the closest representation of an arbitrary number ($$x$$) in the form:

$$x = A\times\frac{N}{D}$$

Where $$A$$ is a 32-bit integer, and $$N$$ and $$D$$ are 16-bit integers.

Typical inputs are real numbers between 1 and 1 billion.

## My Attempts

I first attempted to use $$N$$ and $$D$$ to represent the fractional part of the number, but I can only correct for errors less than one part in 65 thousand, not very useful.

I then looked at naively testing all the possible combinations of $$N$$ and $$D$$ to find one with the lowest remainder, even when skipping like fractions (i.e., 2/4, 3/6, etc.), and skipping values of $$N$$ that would result in a fraction that is too large ($$N \leq \frac{x}{D}$$). I still need to search more than a billion combinations.

Some test cases I have found challenging; clearly, it's pretty easy to come up with these.

• $$1000\times \frac{255}{256}$$
• $$123456 \times \frac{123}{456}$$

## Backround

This problem comes up as I am trying to configure a PLL to convert arbitrary input frequencies to a uniform output frequency with the lowest error possible.

## Question

Is there a much more efficient way to solve this problem? Does this kind of problem have a name that I can look up?

• How is the real number given as input in your problem? An input to an algorithm should be a finite string. Jul 1, 2021 at 21:28
• @Vincenzo Well, in my application, it would be a double-precision float. Here I use the term "Algorithm" in its colloquial sense. I also reviewed the tag description to see if this question was a good fit, "An algorithm is a sequence of well-defined steps that defines an abstract solution to a problem. Use this tag when your issue is related to design and analysis of algorithms". Based on that, I thought the tag was a good fit. Jul 1, 2021 at 21:35
• I am looking for an algorithm that is better (ideally substantially) than brute force for finding the optimal set of integers $A$, $N$, $D$ that gives the lowest error. Jul 1, 2021 at 21:49

Generally speaking, this problem is called diophantine approximation.

However, your inputs are not just real numbers but double-precision floating point (aka double) numbers, which are in fact rational. Any double number has the form:

$$(-1)^{sign} ( 1 . b_{51} b_{50} ... b_{0}) \times 2^{e - 1023}$$ where $$b$$ encodes the mantissa (52 bits) and $$e'=e-1023$$ encodes the exponent (11 bits). You mention that your typical inputs are at least 1, so your exponent $$e'$$ is always nonnegative. The mantissa is always an integer multiple of $$2^{-52}$$, and will be an integer multiple of $$2^{-\ell}$$ whenever the least significant nonzero bit of the mantissa is the $$\ell$$-th one.
The bottom line is that to achieve an exact representation you could do the following:

• look at the binary representation of your floating point number $$x$$, isolating the mantissa bit sequence $$b$$ and the exponent $$e'=e-1023$$;
• let $$\ell$$ be the index of the least significant nonzero bit of the mantissa $$b$$;
• your number $$x$$ is an integer multiple of $$1/D = 2^{e' - \ell}$$, and no fraction with a smaller denominator than $$D$$ can represent it exactly;
• divide $$x$$ by $$1/D = 2^{e'-\ell}$$ to find $$N \times A$$.

In the worst case, this may however not guarantee that $$N$$ and $$A$$ satisfy your size requirement, since $$\ell- e'$$ could be up to $$52$$, while you only have $$32 + 16 = 48$$ bits available for $$N \times A$$. A practical solution (but perhaps non-optimal) would be to round the least significant $$52-48 = 4$$ bits in the mantissa of $$x$$.

If you insist on finding the best approximation given a constraint on the largest denominator $$D$$, we need to turn to diophantine approximation techniques and continued fractions.

A fraction $$A/B$$ is a best rational approximation to real number $$x$$ if $$|A - Bx| \le |C - Dx|$$ for all integers $$C$$ and $$D$$ where $$0 < D \le B$$.

Then it is known that the best rational approximation of $$x$$ is a so-called convergent to $$x$$, which you can recover from the continued fraction expansion of $$x$$: $$x = a_0 + \frac{1}{a_1 + \frac{1}{a_2 + \ldots}}$$ The $$m$$-th convergent $$C_m$$ is defined by the fraction $$A_m/B_m$$ where $$A_k$$ and $$B_k$$ are defined recursively by

$$A_{-1}=1$$, $$A_0=a_0$$, $$A_{k} = a_{k} A_{k-1} + A_{k-2}$$

$$B_{-1}=0$$, $$B_0=1$$, $$B_{k} = a_{k} B_{k-1} + B_{k-2}$$.

In other words, for each $$k=0,1,2,\ldots$$:

• generate $$a_k$$ by setting $$a_k = \lfloor x_k \rfloor$$, $$x_{k+1} = 1/(x_k - \lfloor x_k \rfloor)$$, where $$x_0 = x$$;
• compute $$A_k$$ and $$B_k$$ with the formulas above
• the $$k$$-th convergent is $$A_k/B_k$$.

This will give you a sequence of closer and closer approximations of $$x$$, each of which is a best rational approximation. In your case, you want to stop when $$B_k$$ does not fit anymore into $$16$$ bits.

In your examples, we get the following:

First example.

$$x = 1000 \times 255 / 256 \approx 996.09375$$

$$A_0 / B_0 = 996 / 1 \approx 996.0$$

$$A_1 / B_1 = 9961 / 10 \approx 996.1$$

$$A_2 / B_2 = 10957 / 11 \approx 996.0909090909091$$

$$A_3 / B_3 = 20918 / 21 \approx 996.0952380952381$$

$$A_4 / B_4 = 31875 / 32 \approx 996.09375$$

Second example.

$$x = 123456 \times 123 / 456 \approx 33300.63157894737$$

$$A_0 / B_0 = 33300 / 1 \approx 33300.0$$

$$A_1 / B_1 = 33301 / 1 \approx 33301.0$$

$$A_2 / B_2 = 66601 / 2 \approx 33300.5$$

$$A_3 / B_3 = 99902 / 3 \approx 33300.666666666664$$

$$A_4 / B_4 = 266405 / 8 \approx 33300.625$$

$$A_5 / B_5 = 632712 / 19 \approx 33300.63157894737$$

• Thanks! This seems to work well for my needs in practice. However, for around 0.08% of inputs, it does not find an optimal solution. For example, with the input $795082588 \times (1763/65225)$, the result is $988571985/46$, with an error of $1.5 \times 10^{-14}$. You can see my implementation here. Jul 3, 2021 at 0:39

The "best" rational approximations to a real number are given by it's continuous fraction, check that out.