lambda d,P=2**32:[s:=min(range(32),key=lambda n:abs(ceil((P<<n)/d)*d-(P<<n))),m:=ceil((P<<s)/d),0]*(1<m<=P and m*d>>(32+s))or[s:=min([n for n in range(32)if 1<2**(33+n)/d-P<=P],key=lambda n:abs(ceil(2**(33+n)/d)*d-(2*P<<n))),m:=floor(2**(33+s)/d+1)-P,1]
from math import*
Try it online!
Probably severely golfable.
Returns a list of three elements in the order [s, m, f].
Here is a calculator that uses my function. You input a Python division expression (dividend/divisor), and it uses my function to return a quotient without directly dividing the dividend by the divisor.
Thanks @pxeger for suggesting the use of bitwise operators (2**(32+s) -> P<<s) in some places, which saves ten bytes.
Explanation
The problem can be solved by comparing the "magic quotient" expression to the actual quotient. Since f is either zero or one, we can split the problem into two parts.
If f is zero, then the "magic quotient" expression becomes:
$$\left\lfloor \frac{y}{2^s} \right\rfloor = \left\lfloor \frac{\left\lfloor \frac{mx}{2^{32}} \right\rfloor}{2^s} \right\rfloor$$
We can ignore the floor functions so that the problem becomes one of approximation rather than complete accuracy. Then, the magic quotient becomes:
$$\frac{mx}{2^{32+s}}$$
What we want to do is find \$m\$ and \$s\$ such that the value above most closely approximates \$\frac{x}{d}\$. We want
$$\frac{mx}{2^{32+s}} \approx \frac{x}{d} \implies mdx \approx x2^{32+s} \implies md \approx 2^{32+s}$$
This problem is a lot easier to solve than the daunting task of dealing with floors and correcting errors, or, even worse, dealing with those pesky fractions. In fact, this part is very simple. We can just loop through the 32 possible values of \$s\$ and find the one with the lowest
$$\left\lceil \frac{2^{32+s}}{d} \right\rceil \times d - 2^{32+s}$$
We use the ceiling, and not the floor, because later on we use the floor to compute the quotient, and the ceiling allows us to get a more precise value.
The above expression is evaluated for each \$s\$ from 0 to 31, and the \$s\$ for which that expression is the lowest is the one that is to be outputted (or is at least considered for output). Next, we evaluate \$m\$, which is the ceiling expression above.
Then \$f\$ is zero (assumed above).
Now, of course, \$f\$ cannot always be zero, so some checks are done to determine whether or not \$f = 0\$ in fact works. \$f = 0\$ works if:
- \$m\$ is within the required range; and
- the magic quotient of \$d\$ by \$d\$ is equal to one (if it isn't, then we need some additional value in the numerator, which is where \$f = 1\$ comes in).
The second condition is sufficient because, considering some nonnegative integer \$ad+b\$ where \$a\$ and \$b\$ are nonnegative integers, \$b\$ being less than \$d\$, the value of
$$\left\lfloor \frac{\left\lfloor \frac{mad+mb}{2^{32}} \right\rfloor}{2^s} \right\rfloor = \left\lfloor \frac{\left\lfloor \frac{amd}{2^{32}} \right\rfloor}{2^s} \right\rfloor$$
This is because \$mb\$ is less than \$md\$, so \$\frac{mb}{2^{32}}\$ is less than one.
If the magic quotient of \$d\$ by \$d\$ is one, then
$$\left\lfloor \frac{\left\lfloor \frac{md}{2^{32}} \right\rfloor}{2^s} \right\rfloor = 1$$
Consequently,
$$\left\lfloor \frac{\left\lfloor \frac{amd}{2^{32}} \right\rfloor}{2^s} \right\rfloor = a \left\lfloor \frac{\left\lfloor \frac{md}{2^{32}} \right\rfloor}{2^s} \right\rfloor = a$$
which is indeed the integer part of \$\frac{ad+b}{d}\$. Since this works for every \$ad+b\$, the \$s\$, \$m\$ and \$f\$ values are satisfactory.
Otherwise, we need to consider \$f = 1\$. This value of \$f\$ will definitely work.
The problem is more complicated if \$f = 1\$, but it is still a problem of approximation. I'll leave the work of finding the simple equation to you, but here it is for those who don't have the time / are too lazy to do the work themselves:
$$d(m+2^{32}) \approx 2^{33+s}$$
From here, the work is very similar to before, except that since we must get a good value, we restrict the search space to the values of \$s\$ such that the \$m\$ corresponding to \$s\$ is within the required range. Another difference is that \$m\$ is found using a Python equivalent of -~ rather than ceil, the difference between them being that -~ increments integers whereas ceil leaves them unchanged. I'm still unsure as to why this is needed, but it works so I'm not going to question it.
I don't do these long mathematical explanations very often, so feel free to critique my mathematics. If the LaTEX looks weird, please do tell me so I can fix it. My reasoning might also be wrong, so if it is, please tell me.
