Is there a widely available implementation of an adaptive, recursive, numerically stable IIR filter?
While I may not be the final authority, I have looked into this. In short, no. 

From my own research, the problem is stability. When the filter is extended to have a wider bandwidth, which of course means that the y[n] values and lags 'move quickly', and then is quenched to a lower bandwidth, where the y[n] values and lags 'move slowly', the y[n] buffer may be loaded with fast oscillated entries that appear as initial conditions to the low bandwidth filter. The filter can easily go unstable. If the bandwidth is altered in an adiabatic manner this instability may be avoided (I'm supposing) but then there's an additional constraint on your adaptivity. 

Now, you don't have to take my word for it. In "Adaptive Filters: Theory and Applications" (2nd) page 5 (Amazon), Boroujeny writes, "However, as we shall see in the later chapters, because of the many difficulties involved in adaptation of IIR filters, their application in the area of adaptive filters is rather limited." The author states that the adaptation process can place (digital) poles outside of the unit circle, even if they started within the unit circle. The rest of his book is mostly dedicated to FIR filters.

I hope this information helps.