code/hyperloglog/hyperloglog_optimized.cpp - master - AMaftuna/Cardinality-Estimation-Algorithm

code/hyperloglog/hyperloglog_optimized.cpp
328 строк10 KB
Abdulloeva Maftuna
improved
22 июн 2025, 02:46
22 июн 2025, 02:46f7fca03
#include "hyperloglog_optimized.h"
#include <cstring>
#include <cmath>
#include <algorithm>

#include <array>

uint64_t MurmurHash64A(const void* key, int len, uint64_t seed) {
    const uint64_t m = 0xc6a4a7935bd1e995ULL;
    const int r = 47;
    uint64_t h = seed ^ (len * m);

    const uint64_t* data = (const uint64_t*)key;
    const uint64_t* end = data + (len / 8);

    while(data != end) {
        uint64_t k = *data++;

        k *= m;
        k ^= k >> r;
        k *= m;

        h ^= k;
        h *= m;
    }

    const unsigned char* data2 = (const unsigned char*)data;
    switch(len & 7) {
        case 7: h ^= (uint64_t)(data2[6]) << 48;
        case 6: h ^= (uint64_t)(data2[5]) << 40;
        case 5: h ^= (uint64_t)(data2[4]) << 32;
        case 4: h ^= (uint64_t)(data2[3]) << 24;
        case 3: h ^= (uint64_t)(data2[2]) << 16;
        case 2: h ^= (uint64_t)(data2[1]) << 8;
        case 1: h ^= (uint64_t)(data2[0]);
                h *= m;
    };

    h ^= h >> r;
    h *= m;
    h ^= h >> r;
    return h;
}

BitPackedRegisters::BitPackedRegisters(int numRegisters)
    : m(numRegisters)
{
    int totalBits = m * 6;
    int totalBytes = (totalBits + 7) / 8;
    data.resize(totalBytes, 0);
}

uint8_t BitPackedRegisters::get(int index) const {
    int bitPos = index * 6;
    size_t bytePos = bitPos / 8;
    int offset = bitPos % 8;
    uint16_t value = data[bytePos];
    if (bytePos + 1 < data.size()) {
        value |= (uint16_t)data[bytePos + 1] << 8;
    }
    uint8_t reg = (value >> offset) & 0x3F;
    return reg;
}

void BitPackedRegisters::set(int index, uint8_t newVal) {
    int bitPos = index * 6;
    size_t bytePos = bitPos / 8;
    int offset = bitPos % 8;
    uint16_t value = data[bytePos];
    if (bytePos + 1 < data.size()) {
        value |= (uint16_t)data[bytePos + 1] << 8;
    }
    uint8_t current = (value >> offset) & 0x3F;
    if(newVal <= current) return;
    uint16_t mask = ~(0x3F << offset);
    value = (value & mask) | ((newVal & 0x3F) << offset);
    data[bytePos] = value & 0xFF;
    if (bytePos + 1 < data.size()) {
        data[bytePos + 1] = (value >> 8) & 0xFF;
    }
}

void BitPackedRegisters::merge(const BitPackedRegisters& other) {
    if (other.m != m) throw std::invalid_argument("Error");
    for (int i = 0; i < m; i++) {
        uint8_t a = get(i);
        uint8_t b = other.get(i);
        if (b > a) {
            set(i, b);
        }
    }
}

int BitPackedRegisters::getSize() const {
    return m;
}

const std::vector<uint8_t>& BitPackedRegisters::RawData() const { 
    return data;
}

int BitPackedRegisters::ByteSize() const { 
    return static_cast<int>(data.size());
}
void BitPackedRegisters::LoadRaw(const std::vector<uint8_t>& buf) {
    if ((int)buf.size() != ByteSize()) throw std::invalid_argument("size mismatch");
    data = buf;
}

using TLookup = std::array<double, 256>;

struct TCorrection {
    TLookup Estimations;
    TLookup Biases;

    double GetBias(double e) const {
        for (size_t idx = 0;; ++idx) {
            const auto estr = Estimations[idx];
            if (estr >= e) {
                if (idx == 0) {
                    return Biases[0];
                }
                const auto estl = Estimations[idx - 1];
                const auto biasl = Biases[idx - 1];
                const auto biasr = Biases[idx];
                const auto de = estr - estl;
                const auto db = biasr - biasl;
                const auto scale = e - estl;
                return biasl + scale * db / de;
            } else if (std::fabs(estr) < 1e-4) {
                //limiter
                return Biases[idx - 1];
            }
        }
    }
};

double EstimateBias(double e, unsigned precision) {
    static const TCorrection CORRECTIONS[1 + ImprovedHyperLogLog::PRECISION_MAX - ImprovedHyperLogLog::PRECISION_MIN] = {
#include "hyperloglog_corrections.inc"
    };
    if (precision < ImprovedHyperLogLog::PRECISION_MIN || precision > ImprovedHyperLogLog::PRECISION_MAX) {
        return 0.;
    }

    return CORRECTIONS[precision - ImprovedHyperLogLog::PRECISION_MIN].GetBias(e);
}


ImprovedHyperLogLog::ImprovedHyperLogLog(unsigned precision)
    : precision_(precision),
      num_buckets_(1 << precision),
      registers_(num_buckets_),
      sparseMode_(true),
      sparseThreshold_(num_buckets_ / 4)
{}

 ImprovedHyperLogLog::ImprovedHyperLogLog(const ImprovedHyperLogLog& o)
     : precision_(o.precision_)
     , num_buckets_(o.num_buckets_)
     , registers_(o.registers_)
     , sparseMode_(o.sparseMode_)
     , sparseRegs_(o.sparseRegs_)
     , sparseThreshold_(o.sparseThreshold_)
 {}

 ImprovedHyperLogLog& ImprovedHyperLogLog::operator=(const ImprovedHyperLogLog& o) {
     if (this != &o) {
         precision_       = o.precision_;
         num_buckets_     = o.num_buckets_;
         registers_       = o.registers_;
         sparseMode_      = o.sparseMode_;
         sparseRegs_      = o.sparseRegs_;
         sparseThreshold_ = o.sparseThreshold_;
     }
     return *this;
 }

void ImprovedHyperLogLog::Add(const std::string& value) {
    uint64_t hash = HashString(value);
    uint32_t idx = hash >> (64 - precision_);
    uint64_t remaining = hash << precision_;
    uint8_t rank = CountLeadingZeros(remaining) + 1;
    std::lock_guard<std::mutex> lock(mtx_);
    if (sparseMode_) {
        AddSparse(idx, rank);
        if (sparseRegs_.size() > sparseThreshold_) ToDense();
    } else {
        registers_.set(idx, rank);
    }
}

void ImprovedHyperLogLog::AddSparse(uint32_t idx, uint8_t rank) {
    for (auto &p : sparseRegs_) {
        if (p.first == idx) {
            p.second = std::max(p.second, rank);
            return;
        }
    }
    sparseRegs_.emplace_back(idx, rank);
}

void ImprovedHyperLogLog::ToDense() {
    for (auto &p : sparseRegs_) {
        registers_.set(p.first, p.second);
    }
    sparseRegs_.clear();
    sparseMode_ = false;
}

double ImprovedHyperLogLog::EstimateRaw() const {
    double sum = 0.0;
    if (sparseMode_) {
        sum = double(num_buckets_);
        for (auto &p : sparseRegs_) {
            sum -= 1;
            sum += std::exp2(-p.second);
        }
    } else {
        for (int i = 0; i < num_buckets_; ++i) sum += std::exp2(-registers_.get(i));
    }
    double alpha = GetAlphaConstant(num_buckets_);
    return alpha * num_buckets_ * num_buckets_ / sum;
}

double ImprovedHyperLogLog::BiasCorrected(double E) const {
    return E - EstimateBias(E, GetPrecision());
}

double ImprovedHyperLogLog::Estimate() const {
    double E = EstimateRaw();
    double threshold = 2.5 * num_buckets_;
    double thresholdLarge = (1.0 / 30.0) * std::exp2(64);
    if (E <= threshold) {
        size_t V = sparseMode_ ? (num_buckets_ - sparseRegs_.size()) : 0;
        if (!sparseMode_) {
            for (int i = 0; i < num_buckets_; ++i)
                if (registers_.get(i) == 0) ++V;
        }
        if (V > 0) {
            double lc = num_buckets_ * std::log(double(num_buckets_) / V);
            return std::floor(lc + 0.5);
        }
    } 
    if (E <= thresholdLarge) {
        return std::floor(BiasCorrected(E) + 0.5);
    }
    double n = -std::exp2(64) * std::log(1 - E / std::exp2(64));
    return std::floor(n + 0.5);
}

void ImprovedHyperLogLog::Merge(const ImprovedHyperLogLog& other) {
     if (other.precision_ != precision_) throw std::invalid_argument("Precision mismatch");
    std::lock_guard<std::mutex> lock(mtx_);
    if (sparseMode_ && other.sparseMode_) {
        for (auto &p : other.sparseRegs_) AddSparse(p.first, p.second);
        if (sparseRegs_.size() > sparseThreshold_) ToDense();
    } else {
        if (sparseMode_) ToDense();
        for (int i = 0; i < num_buckets_; ++i) {
            uint8_t a = registers_.get(i);
            uint8_t b = other.sparseMode_ ? (
                [&](){ for (auto &p: other.sparseRegs_) if(p.first==unsigned(i)) return p.second; return other.registers_.get(i);} )()
                : other.registers_.get(i);
            registers_.set(i, std::max(a, b));
        }
    }
}

double ImprovedHyperLogLog::EstimateIntersection(const ImprovedHyperLogLog& a, const ImprovedHyperLogLog& b) {
    ImprovedHyperLogLog uni(a);
    uni.Merge(b);
    return a.Estimate() + b.Estimate() - uni.Estimate();
}

double ImprovedHyperLogLog::EstimateDifference(const ImprovedHyperLogLog& a, const ImprovedHyperLogLog& b) {
    ImprovedHyperLogLog uni(a);
    uni.Merge(b);
    return uni.Estimate() - b.Estimate();
}

unsigned ImprovedHyperLogLog::GetPrecision() const {
    return precision_;
}

void ImprovedHyperLogLog::Serialize(IOutputStream& out) const {
    uint32_t p = precision_;
    out.Write(&p, sizeof(p));
    uint32_t m = static_cast<uint32_t>(num_buckets_);
    out.Write(&m, sizeof(m));

    const auto& buf = registers_.RawData();
    uint32_t byteCount = static_cast<uint32_t>(buf.size());
    out.Write(&byteCount, sizeof(byteCount));
    out.Write(buf.data(), byteCount);
}

ImprovedHyperLogLog ImprovedHyperLogLog::Deserialize(IInputStream& in) {
    uint32_t p;
    in.Read(&p, sizeof(p));
    uint32_t m;
    in.Read(&m, sizeof(m));

    uint32_t byteCount;
    in.Read(&byteCount, sizeof(byteCount));
    std::vector<uint8_t> buf(byteCount);
    in.Read(buf.data(), byteCount);

    ImprovedHyperLogLog hll(p);
    hll.registers_.LoadRaw(buf);
    hll.ToDense();
    return hll;
}

uint8_t ImprovedHyperLogLog::CountLeadingZeros(uint64_t value) const {
    return value ? __builtin_clzll(value) : 64;
}

double ImprovedHyperLogLog::GetAlphaConstant(int m) const {
    if (m == 16)
        return 0.673;
    else if (m == 32)
        return 0.697;
    else if (m == 64)
        return 0.709;
    else
        return 0.7213 / (1 + 1.079 / m);
}