hypercore-kb / opencog-kb / language-learning / mvp / mi-kernel.cl

/*
 * mi-kernel.cl -- GPU-accelerated MI computation
 *
 * Computes MI = log2[count * N / (left_marg * right_marg)]
 *            = log2(count) + log2(N) - log2(left_marg) - log2(right_marg)
 *
 * Input vectors (all same length = number of pairs):
 *   counts      - pair counts count(x,y)
 *   left_margs  - left marginals count(x,*)
 *   right_margs - right marginals count(*,y)
 *   n_value     - single element vector containing N (total count)
 *
 * Output:
 *   mi_values   - computed MI for each pair
 *
 * Note: atomese-simd automatically adds a size parameter (sz) to all kernels.
 * Use double precision to match atomese-simd's FloatValue representation.
 */

#if defined(cl_khr_fp64)
#  pragma OPENCL EXTENSION cl_khr_fp64: enable
#elif defined(cl_amd_fp64)
#  pragma OPENCL EXTENSION cl_amd_fp64: enable
#else
#  error double precision is not supported
#endif

__kernel void compute_mi(
    __global double* mi_values,
    __global const double* counts,
    __global const double* left_margs,
    __global const double* right_margs,
    __global const double* n_value,
    const unsigned long sz)
{
    size_t i = get_global_id(0);
    if (i >= sz) return;

    double count = counts[i];
    double left = left_margs[i];
    double right = right_margs[i];
    double n = n_value[0];

    // Avoid log(0) - use small epsilon
    double eps = 1e-10;
    count = fmax(count, eps);
    left = fmax(left, eps);
    right = fmax(right, eps);
    n = fmax(n, eps);

    // MI = log2(count) + log2(N) - log2(left) - log2(right)
    // log2(x) = log(x) / log(2) = log(x) * 1.4426950408889634
    double log2_factor = 1.4426950408889634;

    double mi = (log(count) + log(n) - log(left) - log(right)) * log2_factor;

    mi_values[i] = mi;
}

/*
 * Compute PMI_raw (N-independent part)
 * PMI_raw = log2(count) - log2(left_marg) - log2(right_marg)
 */
__kernel void compute_pmi_raw(
    __global double* pmi_values,
    __global const double* counts,
    __global const double* left_margs,
    __global const double* right_margs,
    const unsigned long sz)
{
    size_t i = get_global_id(0);
    if (i >= sz) return;

    double count = counts[i];
    double left = left_margs[i];
    double right = right_margs[i];

    double eps = 1e-10;
    count = fmax(count, eps);
    left = fmax(left, eps);
    right = fmax(right, eps);

    double log2_factor = 1.4426950408889634;

    double pmi = (log(count) - log(left) - log(right)) * log2_factor;

    pmi_values[i] = pmi;
}

/*
 * Add log2(N) to PMI_raw values to get MI
 * MI = PMI_raw + log2(N)
 */
__kernel void adjust_for_n(
    __global double* mi_values,
    __global const double* pmi_raw_values,
    __global const double* n_value,
    const unsigned long sz)
{
    size_t i = get_global_id(0);
    if (i >= sz) return;

    double n = n_value[0];
    double log2_n = log(fmax(n, 1e-10)) * 1.4426950408889634;

    mi_values[i] = pmi_raw_values[i] + log2_n;
}