7

what is the most efficient way to compute a logarithm in solidity? Is there a library that implements it or a built in function?

6

Though there's no current implementation (and I couldn't see one in the dapp-bin either), you could implement your own using a Taylor Series, as suggested by Vitalik in this old Reddit thread.

For the example in the thread, it effectively comes down to something like the following:

 x = msg.data[0]
 LOG = 0
 while x >= 1500000:
     LOG = LOG + 405465
     x = x * 2 / 3
 x = x - 1000000
 y = x
 i = 1
 while i < 10:
      LOG = LOG + (y / i)
      i = i + 1
      y = y * x / 1000000
      LOG = LOG - (y / i)
      i = i + 1
      y = y * x / 1000000
return(LOG)

Bear in mind the number of steps and the associated gas cost...

12

Here's a very efficient (< 700 gas) way to calculate the ceiling of log_2:

function log2(uint x) returns (uint y){
   assembly {
        let arg := x
        x := sub(x,1)
        x := or(x, div(x, 0x02))
        x := or(x, div(x, 0x04))
        x := or(x, div(x, 0x10))
        x := or(x, div(x, 0x100))
        x := or(x, div(x, 0x10000))
        x := or(x, div(x, 0x100000000))
        x := or(x, div(x, 0x10000000000000000))
        x := or(x, div(x, 0x100000000000000000000000000000000))
        x := add(x, 1)
        let m := mload(0x40)
        mstore(m,           0xf8f9cbfae6cc78fbefe7cdc3a1793dfcf4f0e8bbd8cec470b6a28a7a5a3e1efd)
        mstore(add(m,0x20), 0xf5ecf1b3e9debc68e1d9cfabc5997135bfb7a7a3938b7b606b5b4b3f2f1f0ffe)
        mstore(add(m,0x40), 0xf6e4ed9ff2d6b458eadcdf97bd91692de2d4da8fd2d0ac50c6ae9a8272523616)
        mstore(add(m,0x60), 0xc8c0b887b0a8a4489c948c7f847c6125746c645c544c444038302820181008ff)
        mstore(add(m,0x80), 0xf7cae577eec2a03cf3bad76fb589591debb2dd67e0aa9834bea6925f6a4a2e0e)
        mstore(add(m,0xa0), 0xe39ed557db96902cd38ed14fad815115c786af479b7e83247363534337271707)
        mstore(add(m,0xc0), 0xc976c13bb96e881cb166a933a55e490d9d56952b8d4e801485467d2362422606)
        mstore(add(m,0xe0), 0x753a6d1b65325d0c552a4d1345224105391a310b29122104190a110309020100)
        mstore(0x40, add(m, 0x100))
        let magic := 0x818283848586878898a8b8c8d8e8f929395969799a9b9d9e9faaeb6bedeeff
        let shift := 0x100000000000000000000000000000000000000000000000000000000000000
        let a := div(mul(x, magic), shift)
        y := div(mload(add(m,sub(255,a))), shift)
        y := add(y, mul(256, gt(arg, 0x8000000000000000000000000000000000000000000000000000000000000000)))
    }  
}
4

Here is high-precision ln(x) implementation for 128.128 fixed point numbers:

/**
 * 2^127.
 */
uint128 private constant TWO127 = 0x80000000000000000000000000000000;

/**
 * 2^128 - 1.
 */
uint128 private constant TWO128_1 = 0xFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF;

/**
 * ln(2) * 2^128.
 */
uint128 private constant LN2 = 0xb17217f7d1cf79abc9e3b39803f2f6af;

/**
 * Return index of most significant non-zero bit in given non-zero 256-bit
 * unsigned integer value.
 *
 * @param x value to get index of most significant non-zero bit in
 * @return index of most significant non-zero bit in given number
 */
function mostSignificantBit (uint256 x) pure internal returns (uint8 r) {
  require (x > 0);

  if (x >= 0x100000000000000000000000000000000) {x >>= 128; r += 128;}
  if (x >= 0x10000000000000000) {x >>= 64; r += 64;}
  if (x >= 0x100000000) {x >>= 32; r += 32;}
  if (x >= 0x10000) {x >>= 16; r += 16;}
  if (x >= 0x100) {x >>= 8; r += 8;}
  if (x >= 0x10) {x >>= 4; r += 4;}
  if (x >= 0x4) {x >>= 2; r += 2;}
  if (x >= 0x2) r += 1; // No need to shift x anymore
}
/*
function mostSignificantBit (uint256 x) pure internal returns (uint8) {
  require (x > 0);

  uint8 l = 0;
  uint8 h = 255;

  while (h > l) {
    uint8 m = uint8 ((uint16 (l) + uint16 (h)) >> 1);
    uint256 t = x >> m;
    if (t == 0) h = m - 1;
    else if (t > 1) l = m + 1;
    else return m;
  }

  return h;
}
*/

/**
 * Calculate log_2 (x / 2^128) * 2^128.
 *
 * @param x parameter value
 * @return log_2 (x / 2^128) * 2^128
 */
function log_2 (uint256 x) pure internal returns (int256) {
  require (x > 0);

  uint8 msb = mostSignificantBit (x);

  if (msb > 128) x >>= msb - 128;
  else if (msb < 128) x <<= 128 - msb;

  x &= TWO128_1;

  int256 result = (int256 (msb) - 128) << 128; // Integer part of log_2

  int256 bit = TWO127;
  for (uint8 i = 0; i < 128 && x > 0; i++) {
    x = (x << 1) + ((x * x + TWO127) >> 128);
    if (x > TWO128_1) {
      result |= bit;
      x = (x >> 1) - TWO127;
    }
    bit >>= 1;
  }

  return result;
}

/**
 * Calculate ln (x / 2^128) * 2^128.
 *
 * @param x parameter value
 * @return ln (x / 2^128) * 2^128
 */
function ln (uint256 x) pure internal returns (int256) {
  require (x > 0);

  int256 l2 = log_2 (x);
  if (l2 == 0) return 0;
  else {
    uint256 al2 = uint256 (l2 > 0 ? l2 : -l2);
    uint8 msb = mostSignificantBit (al2);
    if (msb > 127) al2 >>= msb - 127;
    al2 = (al2 * LN2 + TWO127) >> 128;
    if (msb > 127) al2 <<= msb - 127;

    return int256 (l2 >= 0 ? al2 : -al2);
  }
}
0

I was interested in a 2log: So the first one does do practically the same thing as above. But the other does binary search in the logarithmic value in [0,256]. The gas cost of the latter is very constant. The gas cost of the first seems pretty linear to the output. The best would be to do binary search up to a certain precision and then finishing with the simple searching probably.

    function findLogSimple(uint b) returns (uint){
        for(uint i=0;2**i<=b;i++){}
        return i-1;
    }
    function findLogBinarySearch(uint b) returns (uint){
        uint up = 256;
        uint down = 0;
        uint attempt = (up+down)/2;
        while (up>down+1){
            if(b>=(2**attempt)){
                down=attempt;
            }else{
                up=attempt;
            }
            attempt=(up+down)/2;
        }
        return attempt;
    }

And even slightly better is a mixture of the two:

function findLogMix(uint b) returns (uint){
    //if(b==0){return 0;}
    uint up = 256;
    uint down = 0;
    uint attempt = (up+down)/2;
    while (up>down+4){
        if(b>=(2**attempt)){
            down=attempt;
        }else{
            up=attempt;
        }
        attempt=(up+down)/2;
    }
uint temp = 2**down;
while(temp<=b){
    down++;
    temp=temp*2;
}
return down;
}
0

Compute the largest integer smaller than or equal to the binary logarithm of a given input:

uint256 constant ONE = 1;

function floorLog(uint256 n) pure returns (uint8) {
    uint8 res = 0;

    if (n < 256) {
        // At most 8 iterations
        while (n > 1) {
            n >>= 1;
            res += 1;
        }
    }
    else {
        // Exactly 8 iterations
        for (uint8 s = 128; s > 0; s >>= 1) {
            if (n >= (ONE << s)) {
                n >>= s;
                res |= s;
            }
        }
    }

    return res;
}

Compute the binary logarithm of a given input scaled up by PRECISION bits:

uint8 constant PRECISION = 127;
uint256 constant FIXED_1 = 0x080000000000000000000000000000000; // (1<<(PRECISION)
uint256 constant FIXED_2 = 0x100000000000000000000000000000000; // (2<<(PRECISION)
uint256 constant MAX_NUM = 0x1ffffffffffffffffffffffffffffffff; // (1<<(256-PRECISION))-1

function log(uint256 numerator, uint256 denominator) pure returns (uint256) {
    uint256 res = 0;

    assert(numerator <= MAX_NUM);
    uint256 x = numerator * FIXED_1 / denominator;

    // If x >= 2, then we compute the integer part of log2(x), which is larger than 0.
    if (x >= FIXED_2) {
        uint8 count = floorLog(x / FIXED_1);
        x >>= count; // now x < 2
        res = count * FIXED_1;
    }

    // If x > 1, then we compute the fraction part of log2(x), which is larger than 0.
    if (x > FIXED_1) {
        for (uint8 i = PRECISION; i > 0; --i) {
            x = (x * x) / FIXED_1; // now 1 < x < 4
            if (x >= FIXED_2) {
                x >>= 1; // now 1 < x < 2
                res += ONE << (i - 1);
            }
        }
    }

    return res;
}

A few notes on the function above:

  • It returns floor(log(numerator/denominator)*2^PRECISION)
  • Both input values must range between 1 and 2^(256-PRECISION)-1
  • The output value ranges between 0 and floor(log(2^(256-PRECISION)-1)*2^PRECISION)
  • It assumes numerator >= denominator, because the output would be negative otherwise

Compute the natural logarithm of a given input scaled up to PRECISION bits:

uint256 private constant LOG_NUMERATOR   = 0x3f80fe03f80fe03f80fe03f80fe03f8;
uint256 private constant LOG_DENOMINATOR = 0x5b9de1d10bf4103d647b0955897ba80;
function ln(uint256 numerator, uint256 denominator) pure returns (uint256) {
    return log(numerator, denominator) * LOG_NUMERATOR / LOG_DENOMINATOR;
}

Note that all the constant values above should be computed according to the value of PRECISION.

If you right-shift the result by PRECISION bits, then you get a rough (integer) approximation.

If you divide it by 2 ^ PRECISION on a calculator, then you get a much better approximation.

The maximum value of PRECISION safe for use is 127.

0

A fast and branch-free implementation to compute ceiling of log2. Translated from https://stackoverflow.com/questions/3272424/compute-fast-log-base-2-ceiling

pragma solidity ^0.4.23;


library MathUtil {
  function ceilLog2(uint _x) pure internal returns(uint) {
    require(_x > 0);

    uint x = _x;
    uint y = (((x & (x - 1)) == 0) ? 0 : 1);
    uint j = 128;
    uint k = 0;

    k = (((x & 0xFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF00000000000000000000000000000000) == 0) ? 0 : j);
    y += k;
    x >>= k;
    j >>= 1;

    k = (((x & 0xFFFFFFFFFFFFFFFF0000000000000000) == 0) ? 0 : j);
    y += k;
    x >>= k;
    j >>= 1;

    k = (((x & 0xFFFFFFFF00000000) == 0) ? 0 : j);
    y += k;
    x >>= k;
    j >>= 1;

    k = (((x & 0x00000000FFFF0000) == 0) ? 0 : j);
    y += k;
    x >>= k;
    j >>= 1;

    k = (((x & 0x000000000000FF00) == 0) ? 0 : j);
    y += k;
    x >>= k;
    j >>= 1;

    k = (((x & 0x00000000000000F0) == 0) ? 0 : j);
    y += k;
    x >>= k;
    j >>= 1;

    k = (((x & 0x000000000000000C) == 0) ? 0 : j);
    y += k;
    x >>= k;
    j >>= 1;

    k = (((x & 0x0000000000000002) == 0) ? 0 : j);
    y += k;
    x >>= k;
    j >>= 1;

    return y;
  }
}
  • Gives 0 instead of 255 for 0x80b20f6c0271cf198d45a6b3b25e59e60fd67d663aba08a56c2c88fa8c5f2427 – k06a Aug 7 at 8:39
0

I've started working on a fixed-point math library for solidity. It's open source (apache 2) and you're invited to use and contribute. This library can solve your problem by providing a set of basic math functions such as log, power, root, etc., that use fixed point decimal numbers of the kind customarily used for ERC20 tokens and automatic market makers.

Please find the code here: https://github.com/extraterrestrial-tech/fixidity

It's low on documentation right now, so please contact me if you have any questions.

0

Improved answer of Mikhail Vladimirov by optimizing mostSignificantBit method. Now log2 costs 25K gas instead of 30K gas:

function mostSignificantBit(uint256 x) pure public returns (uint8 r) {
    uint t;
    if ((t = (x >> 128)) > 0) { x = t; r += 128; }
    if ((t = (x >> 64)) > 0) { x = t; r += 64; }
    if ((t = (x >> 32)) > 0) { x = t; r += 32; }
    if ((t = (x >> 16)) > 0) { x = t; r += 16; }
    if ((t = (x >> 8)) > 0) { x = t; r += 8; }
    if ((t = (x >> 4)) > 0) { x = t; r += 4; }
    if ((t = (x >> 2)) > 0) { x = t; r += 2; }
    if ((t = (x >> 1)) > 0) { x = t; r += 1; }
}
  • This sounds doubtful to me. Function log_2 calls mostSignificantBit only once, and my variant of mostSignificantBit consumes about 2K gas in worst case (Solidity 0.5.10 with optimizations). Your variant consumes 830 gas (after adding require (x > 0)). So, improvement should be about 1.2K gas, not 75K. BTW, I optimized mostSignificantBit and now it costs 742 gas in worst case. – Mikhail Vladimirov Aug 9 at 7:08
  • @MikhailVladimirov sorry, fixed estimation of original code to be 30K gas. – k06a Aug 9 at 8:37

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.