Как вычислить трехмерное число Мортона (чередовать биты по 3 дюйма)

Я ищу быстрый способ вычисления числа Мортона 3D. На этом сайте есть хитрость, основанная на магических числах, для двумерных чисел Мортона, но не представляется очевидным, как распространить ее на трехмерные.

Таким образом, в основном у меня есть 3 10-битных числа, которые я хочу чередовать в одно 30-битное число с минимальным количеством операций.

10 ответов

Решение

Вы можете использовать ту же технику. Я предполагаю, что переменные содержат 32-битные целые числа с максимальными 22 битами, установленными в 0 (что немного более ограничительно, чем необходимо). Для каждой переменной x содержащий одно из трех 10-битных целых чисел, мы делаем следующее:

x = (x | (x << 16)) & 0x030000FF;
x = (x | (x <<  8)) & 0x0300F00F;
x = (x | (x <<  4)) & 0x030C30C3;
x = (x | (x <<  2)) & 0x09249249;

Затем с x,y а также z три манипулированных 10-битных целых числа, которые мы получаем, взяв:

x | (y << 1) | (z << 2)

Способ работы этой техники заключается в следующем. Каждый из x = ... строки выше "разбивают" группы битов пополам так, чтобы между ними было достаточно места для битов других целых чисел. Например, если мы рассмотрим три 4-битных целых числа, мы разделим одно с битами 1234 на 000012000034, где нули зарезервированы для других целых чисел. На следующем шаге мы разделим 12 и 34 таким же образом, чтобы получить 001002003004. Даже если 10 битов не дают хорошего повторного деления на две группы, вы можете просто считать это 16 битами, где вы теряете самые высокие в конце,

Как видно из первой строки, вам нужно только это для каждого входного целого числа x он считает, что x & 0x03000000 == 0,

Вот мое решение с помощью скрипта Python:

Я понял намек на его комментарий: Фабиан Риг Гизен
Прочитайте длинный комментарий ниже! Нам нужно отслеживать, какие биты должны идти, как далеко!
Затем на каждом шаге мы выбираем эти биты, перемещаем их и применяем битовую маску (см. Комментарий в последних строках), чтобы маскировать их!

Bit Distances: [0, 2, 4, 6, 8, 10, 12, 14, 16, 18]
Bit Distances (binary): ['0', '10', '100', '110', '1000', '1010', '1100', '1110', '10000', '10010']
Shifting bits by 1   for bits idx: []
Shifting bits by 2   for bits idx: [1, 3, 5, 7, 9]
Shifting bits by 4   for bits idx: [2, 3, 6, 7]
Shifting bits by 8   for bits idx: [4, 5, 6, 7]
Shifting bits by 16  for bits idx: [8, 9]
BitPositions: [0, 1, 2, 3, 4, 5, 6, 7, 8, 9]
Shifted bef.:   0000 0000 0000 0000 0000 0011 0000 0000 hex: 0x300
Shifted:        0000 0011 0000 0000 0000 0000 0000 0000 hex: 0x3000000
NonShifted:     0000 0000 0000 0000 0000 0000 1111 1111 hex: 0xff
Bitmask is now: 0000 0011 0000 0000 0000 0000 1111 1111 hex: 0x30000ff

Shifted bef.:   0000 0000 0000 0000 0000 0000 1111 0000 hex: 0xf0
Shifted:        0000 0000 0000 0000 1111 0000 0000 0000 hex: 0xf000
NonShifted:     0000 0011 0000 0000 0000 0000 0000 1111 hex: 0x300000f
Bitmask is now: 0000 0011 0000 0000 1111 0000 0000 1111 hex: 0x300f00f

Shifted bef.:   0000 0000 0000 0000 1100 0000 0000 1100 hex: 0xc00c
Shifted:        0000 0000 0000 1100 0000 0000 1100 0000 hex: 0xc00c0
NonShifted:     0000 0011 0000 0000 0011 0000 0000 0011 hex: 0x3003003
Bitmask is now: 0000 0011 0000 1100 0011 0000 1100 0011 hex: 0x30c30c3

Shifted bef.:   0000 0010 0000 1000 0010 0000 1000 0010 hex: 0x2082082
Shifted:        0000 1000 0010 0000 1000 0010 0000 1000 hex: 0x8208208
NonShifted:     0000 0001 0000 0100 0001 0000 0100 0001 hex: 0x1041041
Bitmask is now: 0000 1001 0010 0100 1001 0010 0100 1001 hex: 0x9249249

x &= 0x3ff
x = (x | x << 16) & 0x30000ff   <<< THIS IS THE MASK for shifting 16 (for bit 8 and 9)
x = (x | x << 8) & 0x300f00f
x = (x | x << 4) & 0x30c30c3
x = (x | x << 2) & 0x9249249

Таким образом, для 10-битного числа и 2 чередующихся битов (для 32-битных) вам нужно сделать следующее!:

x &= 0x3ff
x = (x | x << 16) & 0x30000ff   #<<< THIS IS THE MASK for shifting 16 (for bit 8 and 9)
x = (x | x << 8) & 0x300f00f
x = (x | x << 4) & 0x30c30c3
x = (x | x << 2) & 0x9249249

А для 21-битного числа и 2 чередующихся битов (для 64-битных) вам нужно сделать следующее!:

x &= 0x1fffff
x = (x | x << 32) & 0x1f00000000ffff
x = (x | x << 16) & 0x1f0000ff0000ff
x = (x | x << 8) & 0x100f00f00f00f00f
x = (x | x << 4) & 0x10c30c30c30c30c3
x = (x | x << 2) & 0x1249249249249249

А для 42-битного числа и 2 чередующихся битов (для 128-битного) вам нужно сделать следующее (на случай, если вам это нужно;-)):

x &= 0x3ffffffffff
x = (x | x << 64) & 0x3ff0000000000000000ffffffffL
x = (x | x << 32) & 0x3ff00000000ffff00000000ffffL
x = (x | x << 16) & 0x30000ff0000ff0000ff0000ff0000ffL
x = (x | x << 8) & 0x300f00f00f00f00f00f00f00f00f00fL
x = (x | x << 4) & 0x30c30c30c30c30c30c30c30c30c30c3L
x = (x | x << 2) & 0x9249249249249249249249249249249L

Скрипт Python для создания и проверки шаблонов чередования!!!

def prettyBinString(x,d=32,steps=4,sep=".",emptyChar="0"):
    b = bin(x)[2:]
    zeros = d - len(b)


    if zeros <= 0: 
        zeros = 0
        k = steps - (len(b) % steps)
    else:
        k = steps - (d % steps)

    s = ""
    #print("zeros" , zeros)
    #print("k" , k)
    for i in range(zeros): 
        #print("k:",k)
        if(k%steps==0 and i!= 0):
            s+=sep
        s += emptyChar
        k+=1

    for i in range(len(b)):
        if( (k%steps==0 and i!=0 and zeros == 0) or  (k%steps==0 and zeros != 0) ):
            s+=sep
        s += b[i]
        k+=1
    return s    

def binStr(x): return prettyBinString(x,32,4," ","0")


def computeBitMaskPatternAndCode(numberOfBits, numberOfEmptyBits):
    bitDistances=[ i*numberOfEmptyBits for i in range(numberOfBits) ]
    print("Bit Distances: " + str(bitDistances))
    bitDistancesB = [bin(dist)[2:] for dist in  bitDistances]
    print("Bit Distances (binary): " + str(bitDistancesB))
    moveBits=[] #Liste mit allen Bits welche aufsteigend um 2, 4,8,16,32,64,128 stellen geschoben werden müssen

    maxLength = len(max(bitDistancesB, key=len))
    abort = False
    for i in range(maxLength):
        moveBits.append([])
        for idx,bits in enumerate(bitDistancesB):
            if not len(bits) - 1 < i:
                if(bits[len(bits)-i-1] == "1"):
                    moveBits[i].append(idx)

    for i in range(len(moveBits)):
        print("Shifting bits by " + str(2**i) + "\t for bits idx: " + str(moveBits[i]))

    bitPositions = range(numberOfBits);
    print("BitPositions: " + str(bitPositions))
    maskOld = (1 << numberOfBits) -1

    codeString = "x &= " + hex(maskOld) + "\n"
    for idx in xrange(len(moveBits)-1, -1, -1):
        if len(moveBits[idx]):


           shifted = 0
           for bitIdxToMove in moveBits[idx]:
                shifted |= 1<<bitPositions[bitIdxToMove];
                bitPositions[bitIdxToMove] += 2**idx; # keep track where the actual bit stands! might get moved several times

           # Get the non shifted part!     
           nonshifted = ~shifted & maskOld

           print("Shifted bef.:\t" + binStr(shifted) + " hex: " + hex(shifted))
           shifted = shifted << 2**idx
           print("Shifted:\t" + binStr(shifted)+ " hex: " + hex(shifted))

           print("NonShifted:\t" + binStr(nonshifted) + " hex: " + hex(nonshifted))
           maskNew =  shifted | nonshifted
           print("Bitmask is now:\t" + binStr(maskNew) + " hex: " + hex(maskNew) +"\n")
           #print("Code: " + "x = x | x << " +str(2**idx)+ " & " +hex(maskNew))

           codeString += "x = (x | x << " +str(2**idx)+ ") & " +hex(maskNew) + "\n"
           maskOld = maskNew
    return codeString


numberOfBits = 10;
numberOfEmptyBits = 2;
codeString = computeBitMaskPatternAndCode(numberOfBits,numberOfEmptyBits);
print(codeString)

def partitionBy2(x):
    exec(codeString)
    return x

def checkPartition(x):
    print("Check partition for: \t" + binStr(x))
    part = partitionBy2(x);
    print("Partition is : \t\t" + binStr(part))
    #make the pattern manualy
    partC = long(0);
    for bitIdx in range(numberOfBits):
        partC  = partC | (x & (1<<bitIdx)) << numberOfEmptyBits*bitIdx
    print("Partition check is :\t" + binStr(partC))
    if(partC == part):
        return True
    else:
        return False

checkError = False        
for i in range(20):
    x = random.getrandbits(numberOfBits);
    if(checkPartition(x) == False):
        checkError = True
        break
if not checkError:
    print("CHECK PARTITION SUCCESSFUL!!!!!!!!!!!!!!!!...")
else:
    print("checkPartition has ERROR!!!!")

Я также добавлю код декодирования через некоторое время!

Самым простым, вероятно, является справочная таблица, если у вас есть 4K свободного места:

static uint32_t t [ 1024 ] = { 0, 0x1, 0x8, ... };

uint32_t m ( int a, int b, int c )
{
    return t[a] | ( t[b] << 1 ) | ( t[c] << 2 );
}

Хакерство битов использует сдвиги и маски для распределения битов, поэтому каждый раз, когда он сдвигает значение и / или упорядочивает его, копируя некоторые биты в пустые места, затем маскируя комбинации, чтобы остались только оригинальные биты.

например:

x = 0xabcd;
  = 0000_0000_0000_0000_1010_1011_1100_1101    

x = (x | (x << S[3])) & B[3]; 

  = ( 0x00abcd00 | 0x0000abcd ) & 0xff00ff 
  = 0x00ab__cd & 0xff00ff 
  = 0x00ab00cd
  = 0000_0000_1010_1011_0000_0000_1100_1101
x = (x | (x << S[2])) & B[2]; 
  = ( 0x0ab00cd0 | 0x00ab00cd) & 0x0f0f0f0f 
  =   0x0a_b_c_d & 0x0f0f0f0f 
  = 0x0a0b0c0d 
  = 0000_1010_0000_1011_0000_1100_0000_1101
x = (x | (x << S[1])) & B[1]; 
  = ( 0000_1010_0000_1011_0000_1100_0000_1101 | 
      0010_1000_0010_1100_0011_0000_0011_0100 ) &
      0011_0011_0011_0011_0011_0011_0011_0011
  =   0010_0010_0010_0011_0011_0000_0011_0001
x = (x | (x << S[0])) & B[0]; 
  = ( 0010_0010_0010_0011_0011_0000_0011_0001 | 
      0100_0100_0100_0110_0110_0000_0110_0010 ) &
      0101_0101_0101_0101_0101_0101_0101_0101
  =   0100_0010_0100_0101_0101_0000_0101_0001

В каждой итерации каждый блок делится на два, самый правый бит самой левой половины блока перемещается в свою конечную позицию, и применяется маска, поэтому остаются только требуемые биты.

После того как вы разложили входы, сместить их так, чтобы значения одного попадали в нули другого, легко.

Чтобы расширить эту технику более чем на два бита между значениями в конечном результате, вы должны увеличить сдвиги между тем, где эти биты заканчиваются. Это становится немного сложнее, так как размер начального блока не является степенью 2, поэтому вы можете разделить его на середину или на границу степени 2.

Так что такая эволюция может сработать:

0000_0000_0000_0000_0000_0011_1111_1111    
0000_0011_0000_0000_0000_0000_1111_1111    
0000_0011_0000_0000_1111_0000_0000_1111    
0000_0011_0000_1100_0011_0000_1100_0011    
0000_1001_0010_0100_1001_0010_0100_1001    

// 0000_0000_0000_0000_0000_0011_1111_1111    
x = ( x | ( x << 16 ) ) & 0x030000ff;
// 0000_0011_0000_0000_0000_0000_1111_1111    
x = ( x | ( x << 8 ) ) & 0x0300f00f;
// 0000_0011_0000_0000_1111_0000_0000_1111    
x = ( x | ( x << 4 ) ) & 0x030c30c3;
// 0000_0011_0000_1100_0011_0000_1100_0011    
x = ( x | ( x << 2 ) ) & 0x09249249;
// 0000_1001_0010_0100_1001_0010_0100_1001    

Выполните одно и то же преобразование для входов, сдвиньте один за другим и два на два или их вместе, и все готово.

Я взял вышеизложенное и изменил его, чтобы объединить 3 16-разрядных числа в 48- (действительно 64-разрядное) число. Возможно, это спасет кого-то от размышлений, чтобы добраться туда.

#include <inttypes.h>
#include <assert.h>
uint64_t zorder3d(uint64_t x, uint64_t y, uint64_t z){
     static const uint64_t B[] = {0x00000000FF0000FF, 0x000000F00F00F00F,
                                    0x00000C30C30C30C3, 0X0000249249249249};           
     static const int S[] =  {16, 8, 4, 2}; 
     static const uint64_t MAXINPUT = 65536;

     assert( ( (x < MAXINPUT) ) && 
      ( (y < MAXINPUT) ) && 
      ( (z < MAXINPUT) )
     );

     x = (x | (x << S[0])) & B[0];
     x = (x | (x << S[1])) & B[1];
     x = (x | (x << S[2])) & B[2];
     x = (x | (x << S[3])) & B[3];

     y = (y | (y << S[0])) & B[0];
     y = (y | (y << S[1])) & B[1];
     y = (y | (y << S[2])) & B[2];
     y = (y | (y << S[3])) & B[3];

     z = (z | (z <<  S[0])) & B[0];
     z = (z | (z <<  S[1])) & B[1];
     z = (z | (z <<  S[2])) & B[2];
     z = (z | (z <<  S[3])) & B[3];

     return ( x | (y << 1) | (z << 2) );
    }

Хорошее время, я только что сделал это в прошлом месяце!

Ключ должен был сделать две функции. Один разносит биты на каждый третий бит. Затем мы можем объединить три из них вместе (со сдвигом для последних двух), чтобы получить окончательное значение Мортона с чередованием.

Этот код чередуется, начиная с битов HIGH (что более логично для значений с фиксированной точкой). Если ваше приложение имеет только 10 битов на компонент, просто сдвиньте каждое значение влево на 22, чтобы оно началось с старших битов.

/* Takes a value and "spreads" the HIGH bits to lower slots to seperate them.
   ie, bit 31 stays at bit 31, bit 30 goes to bit 28, bit 29 goes to bit 25, etc.
   Anything below bit 21 just disappears. Useful for interleaving values
   for Morton codes. */
inline unsigned long spread3(unsigned long x)
{
  x=(0xF0000000&x) | ((0x0F000000&x)>>8) | (x>>16); // spread top 3 nibbles
  x=(0xC00C00C0&x) | ((0x30030030&x)>>4);
  x=(0x82082082&x) | ((0x41041041&x)>>2);
  return x;
}

inline unsigned long morton(unsigned long x, unsigned long y, unsigned long z)
{
  return spread3(x) | (spread3(y)>>1) | (spread3(z)>>2);
}

Следующий код находит число Мортона трех 10-битных входных чисел. Он использует idee из вашей ссылки и выполняет распределение битов в шагах 5-5, 3-2-3-2, 2-1-1-1-2-1-1-1 и 1-1-1- 1-1-1-1-1-1-1, потому что 10 не является степенью двойки.

......................9876543210
............98765..........43210
........987....56......432....10
......98..7..5..6....43..2..1..0
....9..8..7..5..6..4..3..2..1..0

Выше вы можете увидеть расположение каждого бита до первого и после каждого из четырех шагов.

public static Int32 GetMortonNumber(Int32 x, Int32 y, Int32 z)
{
    return SpreadBits(x, 0) | SpreadBits(y, 1) | SpreadBits(z, 2);
}

public static Int32 SpreadBits(Int32 x, Int32 offset)
{
    if ((x < 0) || (x > 1023))
    {
        throw new ArgumentOutOfRangeException();
    }

    if ((offset < 0) || (offset > 2))
    {
        throw new ArgumentOutOfRangeException();
    }

    x = (x | (x << 10)) & 0x000F801F;
    x = (x | (x <<  4)) & 0x00E181C3;
    x = (x | (x <<  2)) & 0x03248649;
    x = (x | (x <<  2)) & 0x09249249;

    return x << offset;
}

У меня была похожая проблема сегодня, но вместо 3 чисел я должен объединить произвольное число чисел любой битовой длины. Я использовал свой собственный алгоритм распределения и маскирования битов и применил его к C# BigIntegers. Вот код, который я написал. На этапе компиляции он вычисляет магические числа и маску для заданного количества измерений и битовой глубины. Затем вы можете повторно использовать объект для нескольких преобразований.

/// <summary>
/// Convert an array of integers into a Morton code by interleaving the bits.
/// Create one Morton object for a given pair of Dimension and BitDepth and reuse if when encoding multiple 
/// Morton numbers.
/// </summary>  
public class Morton
{
    /// <summary>
    /// Number of bits to use to represent each number being interleaved.
    /// </summary>
    public int BitDepth { get; private set; }

    /// <summary>
    /// Count of separate numbers to interleave into a Morton number.
    /// </summary>
    public int Dimensions { get; private set; }

    /// <summary>
    /// The MagicNumbers spread the bits out to the right position.
    /// Each must must be applied and masked, because the bits would overlap if we only used one magic number.
    /// </summary>
    public BigInteger LargeMagicNumber { get; private set; }
    public BigInteger SmallMagicNumber { get; private set; }

    /// <summary>
    /// The mask removes extraneous bits that were spread into positions needed by the other dimensions.
    /// </summary>
    public BigInteger Mask { get; private set; }

    public Morton(int dimensions, int bitDepth)
    {
        BitDepth = bitDepth;
        Dimensions = dimensions;
        BigInteger magicNumberUnit = new BigInteger(1UL << (int)(Dimensions - 1));
        LargeMagicNumber = magicNumberUnit;
        BigInteger maskUnit = new BigInteger(1UL << (int)(Dimensions - 1));
        Mask = maskUnit;
        for (var i = 0; i < bitDepth - 1; i++)
        {
            LargeMagicNumber = (LargeMagicNumber << (Dimensions - 1)) | (i % 2 == 1 ? magicNumberUnit : BigInteger.Zero);
            Mask = (Mask << Dimensions) | maskUnit;       
        }
        SmallMagicNumber = (LargeMagicNumber >> BitDepth) << 1; // Need to trim off pesky ones place bit.
    }

    /// <summary>
    /// Interleave the bits from several integers into a single BigInteger.
    /// The high-order bit from the first number becomes the high-order bit of the Morton number.
    /// The high-order bit of the second number becomes the second highest-ordered bit in the Morton number.
    /// 
    /// How it works.
    /// 
    /// When you multupliy by the magic numbers you make multiple copies of the the number they are multplying, 
    /// each shifted by a different amount.
    /// As it turns out, the high order bit of the highest order copy of a number is N bits to the left of the 
    /// second bit of the second copy, and so forth. 
    /// This is because each copy is shifted one bit less than N times the copy number.
    /// After that, you apply the AND-mask to unset all bits that are not in position.
    /// 
    /// Two magic numbers are needed because since each copy is shifted one less than the bitDepth, consecutive
    /// copies would overlap and ruin the algorithm. Thus one magic number (LargeMagicNumber) handles copies 1, 3, 5, etc, while the 
    /// second (SmallMagicNumber) handles copies 2, 4, 6, etc.
    /// </summary>
    /// <param name="vector">Integers to combine.</param>
    /// <returns>A Morton number composed of Dimensions * BitDepth bits.</returns>
    public BigInteger Interleave(int[] vector)
    {
        if (vector == null || vector.Length != Dimensions)
            throw new ArgumentException("Interleave expects an array of length " + Dimensions, "vector");
        var morton = BigInteger.Zero;
        for (var i = 0; i < Dimensions; i++)
        {
            morton |= (((LargeMagicNumber * vector[i]) & Mask) | ((SmallMagicNumber * vector[i]) & Mask)) >> i;
        }
        return morton;
    }


    public override string ToString()
    {
        return "Morton(Dimension: " + Dimensions + ", BitDepth: " + BitDepth 
            + ", MagicNumbers: " + Convert.ToString((long)LargeMagicNumber, 2) + ", " + Convert.ToString((long)SmallMagicNumber, 2)
            + ", Mask: " + Convert.ToString((long)Mask, 2) + ")";
    }
}

Ниже приведен фрагмент кода для генерации ключа Мортона размером 64 бита для трехмерной точки.

using namespace std;

unsigned long long spreadBits(unsigned long long x)
{
    x=(x|(x<<20))&0x000001FFC00003FF;
    x=(x|(x<<10))&0x0007E007C00F801F;
    x=(x|(x<<4))&0x00786070C0E181C3;
    x=(x|(x<<2))&0x0199219243248649;
    x=(x|(x<<2))&0x0649249249249249;
    x=(x|(x<<2))&0x1249249249249249;
    return x;
}

int main()
{
    unsigned long long x,y,z,con=1;
    con=con<<63;
    printf("%#llx\n",(spreadBits(x)|(spreadBits(y)<<1)|(spreadBits(z)<<2))|con);    
}

Вот генератор, который я сделал в Ruby для создания методов кодирования произвольной длины:

def morton_code_for(bits)
  method = ''

  limit_mask = (1 << (bits * 3)) - 1
  split = (2 ** ((Math.log(bits) / Math.log(2)).to_i + 1)).to_i
  level = 1

  puts "// Coding for 3 #{bits}-bit values"

  loop do
    shift = split
    split /= 2
    level *= 2

    mask = ([ '1' * split ] * level).join('0' * split * 2).to_i(2) & limit_mask

    expression = "v = (v | (v << %2d)) & 0x%016x;" % [ shift, mask ]

    method << expression

    puts "%s // 0b%064b" % [ expression, mask ]

    break if (split <= 1)
  end

  puts
  print "// Test of method results: "
  v = (1 << bits) - 1
  puts eval(method).to_s(2)
end

morton_code_for(21)

Выходные данные являются достаточно общими и могут быть адаптированы по мере необходимости. Образец вывода:

// Coding for 3 21-bit values
v = (v | (v << 32)) & 0x7fff00000000ffff; // 0b0111111111111111000000000000000000000000000000001111111111111111
v = (v | (v << 16)) & 0x00ff0000ff0000ff; // 0b0000000011111111000000000000000011111111000000000000000011111111
v = (v | (v <<  8)) & 0x700f00f00f00f00f; // 0b0111000000001111000000001111000000001111000000001111000000001111
v = (v | (v <<  4)) & 0x30c30c30c30c30c3; // 0b0011000011000011000011000011000011000011000011000011000011000011
v = (v | (v <<  2)) & 0x1249249249249249; // 0b0001001001001001001001001001001001001001001001001001001001001001

// Test of method results: 1001001001001001001001001001001001001001001001001001001001001

Для платформ с достаточно быстрым целочисленным умножением и суперскалярным выполнением следующий код обычно работает хорошо. Основная идея здесь заключается в том, что генерация n -мерного кода Мортона может быть сгенерирована с помощью строительного блока, который распределяет исходные биты на n битов друг от друга. Этот процесс расширения включает сдвиг битов влево и объединение промежуточных результатов с помощью ИЛИ. Пока нет вредных коллизий, можно использовать ADD вместо OR, а последовательность зависимых SHIFT-ADD по сути является умножением. Операнды-источники и операнды-приемники должны быть соответствующим образом замаскированы, чтобы свести к минимуму вероятность коллизий в определенных битовых позициях и потенциального выноса из таких битовых позиций.

При компиляции для x86-64 с использованием текущих компиляторов gcc, clang и Intel classic при максимальной оптимизации (-O3),spread_ten_bits_by_3()компилируется в последовательность из 11 команд (безRET).

      #include <stdio.h>
#include <stdlib.h>
#include <stdint.h>

/* Spread least significant ten bits of argument to every third bit of result.
   For i in [0,9]: r<3*i+0> = a<i>, r<3*i+1> = 0, r<3*i+2> = 0. Bits a<31:10> 
   are assumed to be zero! 
*/
uint32_t spread_ten_bits_by_3 (uint32_t a)
{
    const uint32_t sm1 = 0x333ul;      // source mask: select bits 0,1,4,5,8,9
    const uint32_t sm2 = 0x0ccul;      // source mask: select bits 2,3,6,7
    const uint32_t ml1 = 0x50005ul;    // multiplier: spread to bits 0-11,16-27
    const uint32_t ml2 = 0x00500ul;    // multiplier: spread to bits 8-19
    const uint32_t ml3 = 0x05050ul;    // multiplier: spread to bits 6-13,14-21
    const uint32_t rm1 = 0x09009009ul; // result mask: bits 0,3,12,15,24,27
    const uint32_t rm2 = 0x00240240ul; // result mask: bits 6,9,18,21
    uint32_t t = a & sm1;
    return ((((t * ml1) | (t * ml2)) & rm1) | (((a & sm2) * ml3) & rm2));
}

/* Compute 30-bit 3D Morton code from the least significant 10 bits of each of
   the source operands. All other source operand bits are assumed to be zero!
*/
uint32_t morton3d_ten_bits (uint32_t x, uint32_t y, uint32_t z)
{
    return ((spread_ten_bits_by_3 (x) << 0) |
            (spread_ten_bits_by_3 (y) << 1) |
            (spread_ten_bits_by_3 (z) << 2));    
}

uint32_t spread_ten_bits_by_3_ref (uint32_t a)
{
    uint32_t r = 0;
    r |= (((a >> 0) & 1) <<  0);
    r |= (((a >> 1) & 1) <<  3);
    r |= (((a >> 2) & 1) <<  6);
    r |= (((a >> 3) & 1) <<  9);
    r |= (((a >> 4) & 1) << 12);
    r |= (((a >> 5) & 1) << 15);
    r |= (((a >> 6) & 1) << 18);
    r |= (((a >> 7) & 1) << 21);
    r |= (((a >> 8) & 1) << 24);
    r |= (((a >> 9) & 1) << 27);
    return r;
}

int main (void)
{
    for (uint32_t i = 0; i < 1024; i++) {
        uint32_t res = spread_ten_bits_by_3 (i);
        uint32_t ref = spread_ten_bits_by_3_ref (i);
        if (res != ref) {
            printf ("error @ %08lx: res=%08lx ref=%08lx\n", 
                    (unsigned long int)i, 
                    (unsigned long int)res, 
                    (unsigned long int)ref);
            return EXIT_FAILURE;
        }
    }
    return EXIT_SUCCESS;
}

Пример кода, сгенерированного clang 17.0.1 для x86-64, предоставлен Compiler Explorer:

      spread_ten_bits_by_3:                   # @spread_ten_bits_by_3
        mov     eax, edi
        and     eax, 819
        imul    ecx, eax, 327685
        shl     eax, 8
        lea     edx, [rax + 4*rax]
        or      edx, ecx
        and     edx, 151031817
        and     edi, 204
        imul    eax, edi, 20560
        and     eax, 2359872
        or      eax, edx
        ret
Другие вопросы по тегам