This week I implemented the velusqrt algorithms for the rust library im working on: https://github.com/GiacomoPope/isogeny_rs

Good news is that it works generally and is pretty fast (~100s of microseconds for an isogeny in CSIDH compared to 10s of milliseconds in Sage) but it's definitely not standing up to the state of the art. There's a few reasons for this:

1. I am not using the scaled remainder tree for multipoint evaluation. A proper implementation of this needs a handful of specialised polynomial arithmetic functions to compute the low, inverse and middle elements of the product of two polynomials f and g.

2. I probably have made some silly mistakes in Rust itself and I'm wasting time with clones and allocation when I could be borrowing more or using scratch buffers

3. My estimates about inversion costs and multiplications through sets is slightly off and I should invert less and use projective things more. This is a little harder to guesstimate in a general setting without picking a precise prime.

If anyone wants to collaborate on optimising this please let me know. I’m sure there’s loads I could be doing better!