/*
* gpu-connector-cosine.cl — Brute-force pairwise cosine on dense connector matrix
*
* For L0 word clustering: individual connectors ("the-", "sat+") as dimensions.
* Dense matrix: N_WORDS rows × N_CONN columns (float, row-major).
* N_WORDS ≤ 2048, N_CONN ≤ 4096 → ~32MB matrix, 2M pairs — trivial for GPU.
*
* Kernel 1: compute_connector_norms — one thread per word
* Kernel 2: connector_cosine_brute — one thread per pair (i,j) where i < j
*/
/* ─── Kernel 1: Compute L2 norms per word ─── */
__kernel void compute_connector_norms(
__global const float* conn_matrix, /* [n_words × n_conn] row-major */
__global float* word_norms, /* [n_words] output: sqrt(sum of squares) */
const unsigned int n_words,
const unsigned int n_conn)
{
unsigned int w = get_global_id(0);
if (w >= n_words) return;
float sum_sq = 0.0f;
unsigned int base = w * n_conn;
for (unsigned int c = 0; c < n_conn; c++) {
float v = conn_matrix[base + c];
sum_sq += v * v;
}
word_norms[w] = sqrt(sum_sq);
}
/* ─── Kernel 2: Brute-force pairwise cosine ─── */
__kernel void connector_cosine_brute(
__global const float* conn_matrix, /* [n_words × n_conn] row-major */
__global const float* word_norms, /* [n_words] pre-computed norms */
const unsigned int n_words,
const unsigned int n_conn,
const float threshold, /* minimum cosine to output */
__global unsigned int* out_wa, /* output: word index a */
__global unsigned int* out_wb, /* output: word index b */
__global float* out_cos, /* output: cosine value */
__global unsigned int* out_count, /* atomic counter */
const unsigned int max_output)
{
/* Linearized pair index: tid maps to (i, j) where i < j */
unsigned int tid = get_global_id(0);
unsigned int total_pairs = (n_words * (n_words - 1)) / 2;
if (tid >= total_pairs) return;
/* Decode linear index to (i, j) */
/* Using the formula: i = n - 2 - floor(sqrt(-8*tid + 4*n*(n-1) - 7) / 2 - 0.5)
j = tid + i + 1 - n*(n-1)/2 + (n-i)*((n-i)-1)/2 */
unsigned int i = n_words - 2 - (unsigned int)floor(sqrt((float)(2 * total_pairs - 2 * tid) - 0.75f) - 0.5f);
unsigned int j = tid - i * (n_words - 1) + (i * (i + 1)) / 2 + 1;
/* Safety check */
if (i >= n_words || j >= n_words || i >= j) return;
/* Skip if either norm is zero */
float norm_i = word_norms[i];
float norm_j = word_norms[j];
if (norm_i < 1e-6f || norm_j < 1e-6f) return;
/* Compute dot product */
float dot = 0.0f;
unsigned int base_i = i * n_conn;
unsigned int base_j = j * n_conn;
for (unsigned int c = 0; c < n_conn; c++) {
dot += conn_matrix[base_i + c] * conn_matrix[base_j + c];
}
float cosine = dot / (norm_i * norm_j);
/* Clamp to [0, 1] */
if (cosine < 0.0f) cosine = 0.0f;
if (cosine > 1.0f) cosine = 1.0f;
/* Output if above threshold */
if (cosine >= threshold) {
unsigned int idx = atomic_add(out_count, 1);
if (idx < max_output) {
out_wa[idx] = i;
out_wb[idx] = j;
out_cos[idx] = cosine;
}
}
}