{"spec_id":"root-locus-basic","library":"ggplot2","language":"r","code":"#' anyplot.ai\n#' root-locus-basic: Root Locus Plot for Control Systems\n#' Library: ggplot2 3.5.1 | R 4.4.1\n#' Quality: 86/100 | Created: 2026-06-18\n\nlibrary(ggplot2)\nlibrary(dplyr)\nlibrary(ragg)\n\nset.seed(42)\n\n# Theme tokens\nTHEME       <- Sys.getenv(\"ANYPLOT_THEME\", \"light\")\nPAGE_BG     <- if (THEME == \"light\") \"#FAF8F1\" else \"#1A1A17\"\nELEVATED_BG <- if (THEME == \"light\") \"#FFFDF6\" else \"#242420\"\nINK         <- if (THEME == \"light\") \"#1A1A17\" else \"#F0EFE8\"\nINK_SOFT    <- if (THEME == \"light\") \"#4A4A44\" else \"#B8B7B0\"\nINK_MUTED   <- if (THEME == \"light\") \"#6B6A63\" else \"#A8A79F\"\n\nIMPRINT_PALETTE <- c(\n  \"#009E73\",  # 1 — brand green\n  \"#C475FD\",  # 2 — lavender\n  \"#4467A3\",  # 3 — blue\n  \"#BD8233\",  # 4 — ochre\n  \"#AE3030\"   # 5 — matte red\n)\n\n# Data: third-order system  G(s) = K / (s(s+2)(s+4))\n# Open-loop poles: 0, -2, -4. No finite zeros.\n# Root locus computed numerically by sweeping K.\ngain_vals <- c(seq(0, 2, length.out = 300),\n               seq(2, 10, length.out = 400),\n               seq(10, 40, length.out = 300))\ngain_vals <- sort(unique(gain_vals))\n\n# Characteristic polynomial: s^3 + 6s^2 + 8s + K = 0\n# Roots via polyroot for each K\nroot_locus_rows <- lapply(gain_vals, function(K) {\n  # coefficients: s^3 + 6s^2 + 8s + K\n  r <- polyroot(c(K, 8, 6, 1))\n  # Sort by imaginary part to keep branch assignment consistent\n  r <- r[order(Im(r))]\n  list(\n    gain   = rep(K, 3),\n    branch = c(\"Branch 1\", \"Branch 2\", \"Branch 3\"),\n    re     = Re(r),\n    im     = Im(r)\n  )\n})\n\ndf <- bind_rows(lapply(root_locus_rows, function(x) {\n  data.frame(gain = x$gain, branch = x$branch, re = x$re, im = x$im,\n             stringsAsFactors = FALSE)\n}))\n\n# Keep only points within a reasonable window\ndf <- df[abs(df$re) <= 8 & abs(df$im) <= 8, ]\n\n# Open-loop poles (K=0)\npoles <- data.frame(re = c(0, -2, -4), im = c(0, 0, 0))\n\n# Imaginary axis crossing: K at which Routh criterion gives jw roots\n# For s^3+6s^2+8s+K: Routh row 1=[1,8], row 2=[6,K], pivot=(48-K)/6\n# Crossing when K=48, auxiliary eq 6s^2+48=0 -> s=+/-j*sqrt(8)\ncrossing <- data.frame(re = c(0, 0), im = c(sqrt(8), -sqrt(8)))\n\n# Arrow data: show direction of increasing gain along each branch\narrow_k_vals <- c(6, 18, 6)  # one per branch at a distinctive location\nbranch_names <- c(\"Branch 1\", \"Branch 2\", \"Branch 3\")\n\narrow_df <- bind_rows(lapply(seq_along(branch_names), function(i) {\n  k0 <- arrow_k_vals[i]\n  dk <- 0.6\n  seg_start <- df[abs(df$gain - k0) < 0.15 & df$branch == branch_names[i], ][1, ]\n  seg_end   <- df[abs(df$gain - (k0 + dk)) < 0.15 & df$branch == branch_names[i], ][1, ]\n  if (!is.na(seg_start$re) && !is.na(seg_end$re)) {\n    data.frame(x = seg_start$re, y = seg_start$im,\n               xend = seg_end$re, yend = seg_end$im,\n               branch = branch_names[i])\n  }\n}))\n\n# Damping ratio reference lines: radial rays from origin at angle arccos(zeta)\n# from the negative real axis. A pole with damping ratio zeta lies at angle\n# (pi - arccos(zeta)) from the positive real axis => direction (-zeta, sqrt(1-zeta^2)).\nzeta_vals <- c(0.2, 0.4, 0.6, 0.8)\nplot_re_max <- 7\nplot_im_max <- 6\n\nzeta_df <- bind_rows(lapply(zeta_vals, function(z) {\n  re_end <- -plot_re_max * z\n  im_end <- plot_re_max * sqrt(1 - z^2)\n  # Clip to plot bounds\n  if (abs(im_end) > plot_im_max) {\n    sc <- plot_im_max / abs(im_end)\n    re_end <- re_end * sc\n    im_end <- im_end * sc\n  }\n  rbind(\n    data.frame(x = 0, y = 0, xend = re_end, yend =  im_end),\n    data.frame(x = 0, y = 0, xend = re_end, yend = -im_end)\n  )\n}))\n\n# Natural frequency circles: concentric circles at omega_n = 1..5\nomega_df <- bind_rows(lapply(1:5, function(w) {\n  theta <- seq(0, 2 * pi, length.out = 200)\n  data.frame(re = w * cos(theta), im = w * sin(theta), grp = paste0(\"wn\", w))\n}))\n\n# Title\nplot_title <- \"root-locus-basic · r · ggplot2 · anyplot.ai\"\ntitle_size <- 12\n\n# Plot\np <- ggplot(df, aes(x = re, y = im, color = branch, group = branch)) +\n  # Natural frequency circles (dashed, muted — behind locus)\n  geom_path(\n    data = omega_df,\n    aes(x = re, y = im, group = grp),\n    color = INK_MUTED, linewidth = 0.3, linetype = \"dashed\",\n    alpha = 0.55, inherit.aes = FALSE\n  ) +\n  # Damping ratio rays (dashed, muted — behind locus)\n  geom_segment(\n    data = zeta_df,\n    aes(x = x, y = y, xend = xend, yend = yend),\n    color = INK_MUTED, linewidth = 0.3, linetype = \"dashed\",\n    alpha = 0.55, inherit.aes = FALSE\n  ) +\n  # Axis reference lines\n  geom_hline(yintercept = 0, color = INK_SOFT, linewidth = 0.4) +\n  geom_vline(xintercept = 0, color = INK_SOFT, linewidth = 0.4) +\n  # Root locus branches\n  geom_path(linewidth = 1.0, alpha = 0.85) +\n  # Direction arrows\n  geom_segment(\n    data = arrow_df,\n    aes(x = x, y = y, xend = xend, yend = yend, color = branch),\n    arrow = arrow(length = unit(0.18, \"cm\"), type = \"closed\"),\n    linewidth = 1.2,\n    inherit.aes = FALSE\n  ) +\n  # Imaginary axis crossings (stability boundary)\n  geom_point(\n    data = crossing,\n    aes(x = re, y = im),\n    shape = 18, size = 4, color = IMPRINT_PALETTE[5],\n    inherit.aes = FALSE\n  ) +\n  # Open-loop poles\n  geom_point(\n    data = poles,\n    aes(x = re, y = im),\n    shape = 4, size = 4, stroke = 1.5, color = INK,\n    inherit.aes = FALSE\n  ) +\n  # Color scale (3 branches)\n  scale_color_manual(\n    values = IMPRINT_PALETTE[1:3],\n    name   = \"Branch\"\n  ) +\n  # Equal axes to preserve geometry\n  coord_equal(xlim = c(-7, 2), ylim = c(-6, 6)) +\n  labs(\n    title = plot_title,\n    x     = \"Real Axis\",\n    y     = \"Imaginary Axis\"\n  ) +\n  theme_minimal(base_size = 8) +\n  theme(\n    plot.background   = element_rect(fill = PAGE_BG,     color = PAGE_BG),\n    panel.background  = element_rect(fill = PAGE_BG,     color = NA),\n    panel.grid.major  = element_line(color = INK_MUTED,  linewidth = 0.15),\n    panel.grid.minor  = element_line(color = INK_MUTED,  linewidth = 0.08),\n    panel.border      = element_blank(),\n    axis.line         = element_line(color = INK_SOFT,   linewidth = 0.4),\n    axis.title        = element_text(color = INK,        size = 10),\n    axis.text         = element_text(color = INK_SOFT,   size = 8),\n    plot.title        = element_text(color = INK,        size = title_size,\n                                     margin = margin(b = 8)),\n    legend.background = element_rect(fill = ELEVATED_BG, color = INK_SOFT,\n                                     linewidth = 0.3),\n    legend.text       = element_text(color = INK_SOFT,  size = 8),\n    legend.title      = element_text(color = INK,       size = 9),\n    legend.position   = \"right\",\n    plot.margin       = margin(12, 14, 10, 10)\n  )\n\n# Save\nggsave(\n  filename = sprintf(\"plot-%s.png\", THEME),\n  plot     = p,\n  device   = ragg::agg_png,\n  width    = 8,\n  height   = 4.5,\n  units    = \"in\",\n  dpi      = 400\n)\n"}