I don’t see how that formula could predict the input power. vqiq is output power (to the rotor, not at the shaft, since the rotor inertia, air resistance, and bearing friction will eat up some torque before it makes it to the shaft). Input power is that plus the mosfet+winding resistance multiplied by total current squared, plus eddy current loss which is probably small enough to ignore in most cases, or at least ignore its scaling and lump it in with a general safety margin. So if I were to hazard a guess at a formula for input power, it would be vqiq + (id²+iq²)*(mosfet_resistance + phase_resistance) + safety_margin
EDIT: Oops, I forgot it’s current squared times resistance, so no square root in the current sum. Which means that even if we assume id=0, we need quadratic formula to do the inverse calculation.
r = mosfet_resistance + phase_resistance
iq²*r + vqiq - power = 0
iq = (-vq ± sqrt(vq² + 4*r*power)) / 2*r
Not sure offhand whether the + or - solution is the one we want or if it will need a branch.