Calculating roots in solidity

I’ve been looking for a way to calculate roots of ERC20 caliber numbers, it’s definitely a challenge.

I found Bancor has a solution but that rounds down and yields no decimals.

I found another that yields all the decimals you’d want for small numbers, but fails completely for ERC20 caliber numbers (18 decimals).

Looks like I may have to implement my own, but would love to hear thoughts from the community here.

1 Like

Hi @realisation,

What roots do you need, is it just square roots or is it roots of higher degree (https://en.wikipedia.org/wiki/Nth_root)?

I am not sure I understand what you mean by ERC20 caliber numbers? Do you mean 18 decimals?

In your smart contract I assumed that your calculations are done in the base token units and don’t use the decimals.

Tagging community members who may have dealt with this type of maths: @bugduino, @alsco77

Yeah, exactly, numbers of 18 decimal places – and yes, Nth roots.

I have been using numbers 18 decimal and using wmul and wdiv from DSMath, which is a nice function they have which handles the decimals of precision.

It seems like I’ll need to create an nth power function to enable nth roots with these types of numbers.

1 Like

I’ve never dealt with roots for the moment in Solidity, but when searching for solutions on other math problems I found this issue https://github.com/andytudhope/Recollections/issues/4 which may be related to your problem.
As last resort you could always use an iterative approach using the Bisection method or Newton-Raphson.

1 Like

Here’s an implementation of Newton’s method -

    function nthRoot(uint _a, uint _n, uint _dp, uint _maxIts) public returns(uint) {
        assert (_n > 1);

        // The scale factor is a crude way to turn everything into integer calcs.
        // Actually do (a * (10 ^ ((dp + 1) * n))) ^ (1/n)
        // We calculate to one extra dp and round at the end
        uint one = 10 ** (1 + _dp);
        uint a0 = one ** _n * _a;

        // Initial guess: 1.0
        uint xNew = one;
        uint x;

        uint iter = 0;
        while (xNew != x && iter < _maxIts) {
            uint x = xNew;
            uint t0 = x ** (_n - 1);
            if (x * t0 > a0) {
                xNew = x - (x - a0 / t0) / _n;

            } else {
                xNew = x + (a0 / t0 - x) / _n;
            }
            ++iter;
        }

        // Round to nearest in the last dp.
        return (xNew + 5) / 10;
    }

It is fairly gas intensive, and doesn’t work a lot of the time. If you put in a different value for the decimal places, then it can never converges, even if it is a smaller value. Or if you ask for a smaller root, or a larger one, then it will converge quickly but to the wrong answer. And then for other inputs, mostly ones that are small in general (not that many decimal places, only 4th or 5th roots, not a billion for the base) it’s better. But it’s hard to have confidence in this one. Maybe it is worth trying to improve - depending on how other stuff goes.

1 Like

The Bisection method is the only one that surely converges, but it’s slow, for the Newthon method you have to start with good initial guess otherwise it could diverge.

1 Like