params <-
list(family = "lapis", preset = "homage")

## ----setup-opts, include = FALSE----------------------------------------------
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  message = FALSE,
  warning = FALSE,
  fig.width = 7,
  fig.height = 4.1,
  fig.align = "center",
  out.width = "92%",
  dpi = 100
)

## ----albers-classes, echo=FALSE, results='asis'-------------------------------
cat(sprintf(
  paste0(
    '<script>document.addEventListener("DOMContentLoaded",function(){',
    'document.body.classList.remove("palette-red","palette-lapis","palette-ochre","palette-teal","palette-green","palette-violet","preset-homage","preset-interaction","preset-study","preset-structural","preset-adobe","preset-midnight");',
    'document.body.classList.add("palette-%s","preset-%s");',
    '});</script>'
  ),
  params$family,
  params$preset
))

## ----helpers, include = FALSE-------------------------------------------------
# Shared plotting helpers. Kept out of the reader's view: the reader sees
# the API call and the picture, never the plotting boilerplate.
ec_blue <- "#2166ac"; ec_grey <- "grey78"

# A computed slice (blue) drawn on top of the full dense spectrum (grey).
plot_spectrum <- function(all_vals, computed, ylab = "value", main = NULL) {
  s <- sort(all_vals, decreasing = TRUE)
  k <- length(computed)
  op <- par(mar = c(4, 4.6, if (is.null(main)) 1 else 2.4, 1)); on.exit(par(op))
  plot(seq_along(s), s, pch = 19, cex = 0.4, col = ec_grey, bty = "n",
       xlab = "rank (largest to smallest)", ylab = ylab, main = main)
  points(seq_len(k), sort(computed, decreasing = TRUE),
         pch = 19, cex = 1.2, col = ec_blue)
  legend("topright", c("full spectrum (dense reference)", "computed by eigencore"),
         pch = 19, pt.cex = c(0.7, 1.2), col = c(ec_grey, ec_blue),
         bty = "n", cex = 0.9)
}

## ----setup--------------------------------------------------------------------
library(eigencore)

## ----make-A-------------------------------------------------------------------
set.seed(1)
n <- 200
A <- crossprod(matrix(rnorm(n * n), n, n)) / n + diag(n)
fit <- eig_partial(A, k = 5, target = largest())
fit

## ----first-certificate--------------------------------------------------------
fit$certificate$passed

## ----spectrum, echo = FALSE, fig.cap = "The five largest eigenvalues (blue) located within the full spectrum of A (grey). eigencore computes only the requested slice, then certifies it.", fig.alt = "Scatter plot of all 200 eigenvalues sorted from largest to smallest in grey, with the five largest highlighted in blue at the top-left."----
all_vals <- eigen(A, symmetric = TRUE, only.values = TRUE)$values
plot_spectrum(all_vals, fit$values, ylab = "eigenvalue")

## ----generalized--------------------------------------------------------------
B <- diag(seq(1, 5, length.out = n))
fit_gen <- eig_partial(A, k = 5, target = largest(), B = B,
                       method = lobpcg(maxit = 200))
fit_gen

## ----svd----------------------------------------------------------------------
M <- matrix(rnorm(400 * 50), 400, 50)
svd_fit <- svd_partial(M, rank = 5, target = largest())
svd_fit

## ----svd-scree, echo = FALSE, fig.cap = "The five leading singular values (blue) computed by eigencore, shown against the full singular-value spectrum of M (grey).", fig.alt = "Scatter plot of all 50 singular values of M sorted descending in grey, with the top five highlighted in blue."----
all_sv <- svd(M, nu = 0, nv = 0)$d
plot_spectrum(all_sv, svd_fit$d, ylab = "singular value")

## ----rspectra-----------------------------------------------------------------
res <- eigs_sym(A, k = 5, which = "LA")
str(res, max.level = 1)

## ----rspectra-cert------------------------------------------------------------
res$certificate

## ----as-operator--------------------------------------------------------------
Aop <- as_operator(A)
Aop

## ----problem------------------------------------------------------------------
P <- eigen_problem(A, structure = hermitian(), target = largest())

## ----plan---------------------------------------------------------------------
plan <- plan_solver(P, k = 5)
plan

## ----solve--------------------------------------------------------------------
same_fit <- solve(P, k = 5)
isTRUE(all.equal(same_fit$values, fit$values))