Sqrt.opcodes

// Gas efficient SQRT method, computes floor(sqrt(x)).
// Constant gas cost of 355
// 
// For testing use this tool: https://www.evm.codes/playground?fork=cancun
// Author: Lohann Ferreira
PUSH0
CALLDATALOAD
PUSH1 1
DUP2

// r = floor(log2(x))
PUSH16 0xFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF
DUP2
GT
PUSH1 7
SHL
SWAP1
DUP2
SHR

PUSH8 0xFFFFFFFFFFFFFFFF
DUP2
GT
PUSH1 6
SHL
SWAP1
DUP2
SHR

PUSH4 0xFFFFFFFF
DUP2
GT
PUSH1 5
SHL
SWAP1
DUP2
SHR

PUSH2 0xFFFF
DUP2
GT
PUSH1 4
SHL
SWAP1
DUP2
SHR

PUSH1 0xFF
DUP2
GT
PUSH1 3
SHL
SWAP1
DUP2
SHR

PUSH1 0xF
DUP2
GT
PUSH1 2
SHL
SWAP1
DUP2
SHR

PUSH1 0x3
SWAP1
GT
DUP8 // PUSH1 1
SHL

OR
OR
OR
OR
OR
OR

// SQRT
DUP2
SHR
SHL

// Netwon's Method
// r = (r + x/r) / 2
DUP1
DUP3
DIV
ADD
PUSH1 1
SHR

DUP1
DUP3
DIV
ADD
PUSH1 1
SHR

DUP1
DUP3
DIV
ADD
PUSH1 1
SHR

DUP1
DUP3
DIV
ADD
PUSH1 1
SHR

DUP1
DUP3
DIV
ADD
PUSH1 1
SHR

DUP1
DUP3
DIV
ADD
PUSH1 1
SHR

DUP1
DUP3
DIV
ADD
PUSH1 1
SHR

// y = x/r
DUP1
SWAP2
DIV

// min(r, x/r)
DUP2
GT
SWAP1
SUB

// Return result
PUSH0
MSTORE
PUSH1 32
PUSH0
RETURN

Sqrt.sol

// SPDX-License-Identifier: GPL-2.0
// Gas efficient SQRT method, computes floor(sqrt(x)).
// Constant gas cost of 407 + solidity overhead.
// 
// For testing use this tool: https://www.evm.codes/playground?fork=cancun
// Author: Lohann Ferreira
pragma solidity >=0.7.0 <0.9.0;

contract Sqrt {
    function sqrt(uint256 x) public pure returns (uint256 r) {
        assembly ("memory-safe") {
            // r = floor(log2(x))
            r := shl(7, gt(x, 0xFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF))
            let xx := shr(r, x)

            let rr := shl(6, gt(x, 0xFFFFFFFFFFFFFFFF))
            xx := shr(rr, xx)
            r := or(r, rr)

            rr := shl(5, gt(xx, 0xFFFFFFFF))
            xx := shr(rr, xx)
            r := or(r, rr)

            rr := shl(4, gt(xx, 0xFFFF))
            xx := shr(rr, xx)
            r := or(r, rr)

            rr := shl(3, gt(xx, 0xFF))
            xx := shr(rr, xx)
            r := or(r, rr)

            rr := shl(2, gt(xx, 0x0F))
            xx := shr(rr, xx)
            r := or(r, rr)

            rr := shl(1, gt(xx, 0x03))
            xx := shr(rr, xx)
            r := or(r, rr)
            
            r := shl(shr(1, r), 1)
            
            // Newton's method
            r := shr(1, add(r, div(x, r)))
            r := shr(1, add(r, div(x, r)))
            r := shr(1, add(r, div(x, r)))
            r := shr(1, add(r, div(x, r)))
            r := shr(1, add(r, div(x, r)))
            r := shr(1, add(r, div(x, r)))
            r := shr(1, add(r, div(x, r)))

            // r = min(r, x/r)
            r := sub(r, gt(r, div(x, r)))
        }
    }
}