lkb-pipeline/core/eigensolid.R

source("boot.R")
box_use(./laurent2)
box_use(./lkb)

CROSSING_PAIRS <- list(c(0,1), c(2,3), c(4,5), c(6,7))

cross_step <- function(strands, gens, gens_inv) {
  n <- length(strands)
  if (n != 8) stop("eigensolid: only B8 (8 strands) supported, got B", n)
  new_strands <- vector("list", n)
  for (pair in CROSSING_PAIRS) {
    i <- pair[[1]] + 1L
    j <- pair[[2]] + 1L
    si <- strands[[i]]
    sj <- strands[[j]]
    m <- laurent2$l2_mat_mul(si, sj)
    new_strands[[i]] <- m
    new_strands[[j]] <- laurent2$l2_mat_identity(1)
  }
  new_strands
}

extract_crossing_matrix <- function(strands) {
  m <- matrix(0L, nrow = 8, ncol = 8)
  for (i in seq_along(strands)) {
    for (j in seq_along(strands)) {
      if (i == j) next
      m[i, j] <- laurent2$l2_equal(strands[[i]], strands[[j]])
    }
  }
  m
}

residual_from_strand <- function(strand_mat) {
  entry <- laurent2$l2_str(strand_mat[[1, 3]])
  nterms <- length(strand_mat$coef)
  q16_16_resid <- min(nterms, 65535L)
  q16_16_resid
}

eigensolid_converged <- function(strands, step) {
  if (step == 0) return(FALSE)
  paired <- crossStep(strands)
  old_m <- extract_crossing_matrix(strands)
  new_m <- extract_crossing_matrix(paired)
  identical(old_m, new_m)
}

run_eigensolid <- function(con, max_jobs = Inf) {
  cat(sprintf("\n=== Eigensolid convergence check ===\n"))

  done <- DBI::dbGetQuery(con, "
    SELECT equation_id, result_entry
    FROM lkb.compute_queue
    WHERE status = 'done' AND result_entry IS NOT NULL
  ")
  cat(sprintf("  %d completed LKB entries to process\n", nrow(done)))

  processed <- 0L
  for (i in seq_len(nrow(done))) {
    if (processed >= max_jobs) break

    eq_id <- done$equation_id[[i]]
    entry_str <- done$result_entry[[i]]

    already <- DBI::dbGetQuery(con,
      "SELECT COUNT(*) AS n FROM lkb.eigensolid_receipts WHERE equation_id = $1",
      params = list(eq_id))$n[[1]]
    if (already > 0) {
      next
    }

    if (entry_str == "0" || entry_str == "" || is.na(entry_str)) {
      crossing_m <- matrix(0L, 8, 8)
      sidon_slack <- 128L
      step_count <- 0L
      residuals <- list(0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L)
      scar_absent <- TRUE
    } else {
      row <- DBI::dbGetQuery(con,
        "SELECT sidon_set, braid_word, steps_done FROM lkb.compute_queue WHERE equation_id = $1",
        params = list(eq_id))
      sidon_set <- as.integer(strsplit(row$sidon_set[[1]], ",")[[1]])
      steps_done <- row$steps_done[[1]]
      n <- row$steps_done[[1]]

      gens <- lkb$lkb_generators(8)
      gens_inv <- lkb$lkb_generators_inverse(8)
      N <- 28L
      I_mat <- laurent2$l2_mat_identity(N)

      word_str <- row$braid_word[[1]]
      if (word_str == "ID" || is.na(word_str) || word_str == "") {
        step_count <- 0L
        strands <- replicate(8, I_mat, simplify = FALSE)
      } else {
        steps <- as.integer(strsplit(word_str, ",")[[1]])
        mat <- I_mat
        for (g in steps[seq_len(min(steps_done, length(steps)))]) {
          if (g > 0) mat <- laurent2$l2_mat_mul(mat, gens[[g]])
          else mat <- laurent2$l2_mat_mul(mat, gens_inv[[-g]])
        }
        step_count <- steps_done
        strands <- replicate(8, mat, simplify = FALSE)
      }

      crossing_m <- extract_crossing_matrix(strands)
      sidon_slack <- if (length(sidon_set) > 0) 128L - max(sidon_set, na.rm = TRUE) else 0L
      residuals <- lapply(strands, function(s) residual_from_strand(s))
      scar_absent <- TRUE

      cat(sprintf("  [%s] step %d, slack %d\n", eq_id, step_count, sidon_slack))
    }

    crossing_json <- jsonlite::toJSON(crossing_m, auto_unbox = TRUE)
    residuals_json <- jsonlite::toJSON(unlist(residuals), auto_unbox = TRUE)

    DBI::dbExecute(con,
      "INSERT INTO lkb.eigensolid_receipts
         (equation_id, crossing_matrix, sidon_slack, step_count, residuals, scar_absent)
       VALUES ($1, $2::JSONB, $3, $4, $5::JSONB, $6)
       ON CONFLICT (equation_id) DO UPDATE SET
         crossing_matrix = EXCLUDED.crossing_matrix,
         sidon_slack = EXCLUDED.sidon_slack,
         step_count = EXCLUDED.step_count,
         residuals = EXCLUDED.residuals,
         scar_absent = EXCLUDED.scar_absent",
      params = list(eq_id, crossing_json, sidon_slack, step_count, residuals_json, scar_absent))

    processed <- processed + 1L
    if (processed %% 10 == 0) {
      cat(sprintf("  ... %d receipts generated\n", processed))
    }
  }

  total <- DBI::dbGetQuery(con, "SELECT COUNT(*) AS n FROM lkb.eigensolid_receipts")$n[[1]]
  cat(sprintf("  Done: %d/%d receipts generated (total: %d)\n", processed, nrow(done), total))
  invisible(processed)
}

if (interactive() && requireNamespace("DBI", quietly = TRUE) && requireNamespace("RPostgres", quietly = TRUE)) {
  con <- DBI::dbConnect(RPostgres::Postgres(),
    host = Sys.getenv("PGHOST", "100.79.14.103"),
    dbname = "research_stack",
    user = Sys.getenv("PGUSER", "kaleidoscope"),
    password = Sys.getenv("PGPASSWORD", "kaleidoscope"))
  on.exit(DBI::dbDisconnect(con))
  run_eigensolid(con, max_jobs = Inf)
}