|
| 1 | +""" pyplots.ai |
| 2 | +root-locus-basic: Root Locus Plot for Control Systems |
| 3 | +Library: bokeh 3.9.0 | Python 3.14.3 |
| 4 | +Quality: 90/100 | Created: 2026-03-20 |
| 5 | +""" |
| 6 | + |
| 7 | +import numpy as np |
| 8 | +from bokeh.io import export_png, save |
| 9 | +from bokeh.models import ColumnDataSource, HoverTool, Label, Range1d, Span |
| 10 | +from bokeh.plotting import figure |
| 11 | +from bokeh.resources import Resources |
| 12 | + |
| 13 | + |
| 14 | +# Data - Transfer function: G(s) = 1 / (s(s+1)(s+3)) |
| 15 | +# Open-loop poles at s = 0, -1, -3; no zeros |
| 16 | +open_loop_poles = np.array([0.0, -1.0, -3.0]) |
| 17 | + |
| 18 | +# Characteristic equation: s^3 + 4s^2 + 3s + K = 0 |
| 19 | +# Compute root locus by solving for poles at each gain K |
| 20 | +gains = np.concatenate( |
| 21 | + [np.linspace(0, 0.5, 200), np.linspace(0.5, 5, 600), np.linspace(5, 20, 600), np.linspace(20, 80, 600)] |
| 22 | +) |
| 23 | + |
| 24 | +# Coefficients: s^3 + 4s^2 + 3s + K |
| 25 | +all_roots = np.zeros((len(gains), 3), dtype=complex) |
| 26 | +for i, k in enumerate(gains): |
| 27 | + coeffs = [1, 4, 3, k] |
| 28 | + all_roots[i] = np.sort_complex(np.roots(coeffs)) |
| 29 | + |
| 30 | +# Organize into branches by tracking continuity |
| 31 | +n_branches = 3 |
| 32 | +branch_real = [] |
| 33 | +branch_imag = [] |
| 34 | +branch_gain = [] |
| 35 | +branch_id = [] |
| 36 | + |
| 37 | +for b in range(n_branches): |
| 38 | + reals = all_roots[:, b].real |
| 39 | + imags = all_roots[:, b].imag |
| 40 | + branch_real.extend(reals) |
| 41 | + branch_imag.extend(imags) |
| 42 | + branch_gain.extend(gains) |
| 43 | + branch_id.extend([f"Branch {b + 1}"] * len(gains)) |
| 44 | + |
| 45 | +branch_real = np.array(branch_real) |
| 46 | +branch_imag = np.array(branch_imag) |
| 47 | +branch_gain = np.array(branch_gain) |
| 48 | + |
| 49 | +# Colors |
| 50 | +branch_colors = ["#306998", "#E07B39", "#E07B39"] |
| 51 | +colors = [] |
| 52 | +for b in range(n_branches): |
| 53 | + colors.extend([branch_colors[b]] * len(gains)) |
| 54 | + |
| 55 | +source = ColumnDataSource( |
| 56 | + data={"real": branch_real, "imag": branch_imag, "gain": branch_gain, "branch": branch_id, "color": colors} |
| 57 | +) |
| 58 | + |
| 59 | +# Open-loop poles source |
| 60 | +pole_source = ColumnDataSource(data={"real": open_loop_poles, "imag": np.zeros_like(open_loop_poles)}) |
| 61 | + |
| 62 | +# Imaginary axis crossings (stability boundary) |
| 63 | +# For s^3 + 4s^2 + 3s + K = 0 via Routh criterion: K_critical = 12, w = sqrt(3) |
| 64 | +w_crit = np.sqrt(3) |
| 65 | +crossing_source = ColumnDataSource(data={"real": [0.0, 0.0], "imag": [w_crit, -w_crit]}) |
| 66 | + |
| 67 | +# Breakaway point: derivative of K = -(s^3 + 4s^2 + 3s) → 3s^2 + 8s + 3 = 0 |
| 68 | +# s = (-8 ± sqrt(64-36)) / 6 → s ≈ -0.451 (on real axis between -1 and 0) |
| 69 | +breakaway_s = (-8 + np.sqrt(28)) / 6 |
| 70 | +breakaway_K = -(breakaway_s**3 + 4 * breakaway_s**2 + 3 * breakaway_s) |
| 71 | + |
| 72 | +# Plot |
| 73 | +y_max = max(abs(branch_imag)) * 1.15 |
| 74 | +p = figure( |
| 75 | + width=4800, |
| 76 | + height=2700, |
| 77 | + title="root-locus-basic · bokeh · pyplots.ai", |
| 78 | + x_axis_label="Real Axis (σ)", |
| 79 | + y_axis_label="Imaginary Axis (jω)", |
| 80 | + x_range=Range1d(-5.5, 1.5), |
| 81 | + y_range=Range1d(-y_max, y_max), |
| 82 | + tools="pan,wheel_zoom,box_zoom,reset,save", |
| 83 | + active_scroll="wheel_zoom", |
| 84 | + match_aspect=True, |
| 85 | +) |
| 86 | + |
| 87 | +# Constant damping ratio lines (radial lines from origin into left half-plane) |
| 88 | +for zeta in [0.2, 0.4, 0.6, 0.8]: |
| 89 | + r_max = 5.0 |
| 90 | + x_end = -r_max * zeta |
| 91 | + y_end = r_max * np.sqrt(1 - zeta**2) |
| 92 | + p.line(x=[0, x_end], y=[0, y_end], line_color="#B0B0B0", line_width=1.5, line_alpha=0.3, line_dash="dotted") |
| 93 | + p.line(x=[0, x_end], y=[0, -y_end], line_color="#B0B0B0", line_width=1.5, line_alpha=0.3, line_dash="dotted") |
| 94 | + label_r = 4.5 |
| 95 | + lx = -label_r * zeta |
| 96 | + ly = label_r * np.sqrt(1 - zeta**2) |
| 97 | + p.add_layout( |
| 98 | + Label( |
| 99 | + x=lx, |
| 100 | + y=ly + 0.15, |
| 101 | + text=f"ζ={zeta}", |
| 102 | + text_font_size="24pt", |
| 103 | + text_color="#999999", |
| 104 | + text_alpha=0.6, |
| 105 | + text_font="Helvetica", |
| 106 | + ) |
| 107 | + ) |
| 108 | + |
| 109 | +# Constant natural frequency arcs (semicircles) |
| 110 | +for wn in [1, 2, 3, 4]: |
| 111 | + theta = np.linspace(np.pi / 2, np.pi, 60) |
| 112 | + arc_x = wn * np.cos(theta) |
| 113 | + arc_y = wn * np.sin(theta) |
| 114 | + p.line( |
| 115 | + x=arc_x.tolist(), y=arc_y.tolist(), line_color="#B0B0B0", line_width=1.5, line_alpha=0.25, line_dash="dotted" |
| 116 | + ) |
| 117 | + p.line( |
| 118 | + x=arc_x.tolist(), y=(-arc_y).tolist(), line_color="#B0B0B0", line_width=1.5, line_alpha=0.25, line_dash="dotted" |
| 119 | + ) |
| 120 | + p.add_layout( |
| 121 | + Label( |
| 122 | + x=-0.15, |
| 123 | + y=wn + 0.15, |
| 124 | + text=f"ωn={wn}", |
| 125 | + text_font_size="22pt", |
| 126 | + text_color="#999999", |
| 127 | + text_alpha=0.6, |
| 128 | + text_font="Helvetica", |
| 129 | + ) |
| 130 | + ) |
| 131 | + |
| 132 | +# Real axis segments of root locus |
| 133 | +# Poles at 0, -1, -3 → segments: [-1, 0] and [-∞, -3] |
| 134 | +p.segment( |
| 135 | + x0=[-1], y0=[0], x1=[0], y1=[0], line_color="#306998", line_width=6, line_alpha=0.35, legend_label="Real-axis locus" |
| 136 | +) |
| 137 | +p.segment(x0=[-5.5], y0=[0], x1=[-3], y1=[0], line_color="#306998", line_width=6, line_alpha=0.35) |
| 138 | + |
| 139 | +# Locus branches |
| 140 | +scatter = p.scatter( |
| 141 | + x="real", |
| 142 | + y="imag", |
| 143 | + source=source, |
| 144 | + size=6, |
| 145 | + color="color", |
| 146 | + alpha=0.7, |
| 147 | + line_color=None, |
| 148 | + legend_label="Complex branches", |
| 149 | +) |
| 150 | + |
| 151 | +# Open-loop poles (× markers) |
| 152 | +p.scatter( |
| 153 | + x="real", |
| 154 | + y="imag", |
| 155 | + source=pole_source, |
| 156 | + size=45, |
| 157 | + marker="x", |
| 158 | + color="#2B2B2B", |
| 159 | + line_width=5, |
| 160 | + legend_label="Open-loop poles", |
| 161 | +) |
| 162 | + |
| 163 | +# Imaginary axis crossings (stability boundary markers) |
| 164 | +p.scatter( |
| 165 | + x="real", |
| 166 | + y="imag", |
| 167 | + source=crossing_source, |
| 168 | + size=35, |
| 169 | + marker="diamond", |
| 170 | + color="#E74C3C", |
| 171 | + line_color="white", |
| 172 | + line_width=3, |
| 173 | + legend_label="Stability crossing (K=12)", |
| 174 | +) |
| 175 | + |
| 176 | +# Breakaway point marker |
| 177 | +p.scatter( |
| 178 | + x=[breakaway_s], |
| 179 | + y=[0], |
| 180 | + size=30, |
| 181 | + marker="square", |
| 182 | + color="#306998", |
| 183 | + line_color="white", |
| 184 | + line_width=3, |
| 185 | + legend_label=f"Breakaway (K={breakaway_K:.2f})", |
| 186 | +) |
| 187 | + |
| 188 | +# Stability boundary - imaginary axis |
| 189 | +stability_line = Span( |
| 190 | + location=0, dimension="height", line_color="#E74C3C", line_width=2.5, line_alpha=0.2, line_dash="dashed" |
| 191 | +) |
| 192 | +p.add_layout(stability_line) |
| 193 | + |
| 194 | +# Direction arrows on branches indicating increasing gain |
| 195 | +for b in range(n_branches): |
| 196 | + arrow_idx = len(gains) * 2 // 3 |
| 197 | + r = all_roots[arrow_idx, b] |
| 198 | + r_next = all_roots[min(arrow_idx + 20, len(gains) - 1), b] |
| 199 | + dx = r_next.real - r.real |
| 200 | + dy = r_next.imag - r.imag |
| 201 | + length = np.sqrt(dx**2 + dy**2) |
| 202 | + if length > 0.01: |
| 203 | + p.scatter( |
| 204 | + x=[r.real], |
| 205 | + y=[r.imag], |
| 206 | + size=28, |
| 207 | + marker="triangle", |
| 208 | + color=branch_colors[b], |
| 209 | + angle=[np.arctan2(dy, dx) - np.pi / 2], |
| 210 | + alpha=0.9, |
| 211 | + ) |
| 212 | + |
| 213 | +# Annotations at key engineering points |
| 214 | +p.add_layout( |
| 215 | + Label( |
| 216 | + x=0.2, |
| 217 | + y=w_crit + 0.25, |
| 218 | + text=f"jω-crossing: K={12}, ω=√3≈{w_crit:.2f}", |
| 219 | + text_font_size="28pt", |
| 220 | + text_color="#E74C3C", |
| 221 | + text_font="Helvetica", |
| 222 | + text_font_style="bold", |
| 223 | + ) |
| 224 | +) |
| 225 | + |
| 226 | +p.add_layout( |
| 227 | + Label( |
| 228 | + x=breakaway_s + 0.05, |
| 229 | + y=-0.45, |
| 230 | + text=f"Breakaway: σ={breakaway_s:.3f}, K={breakaway_K:.2f}", |
| 231 | + text_font_size="26pt", |
| 232 | + text_color="#306998", |
| 233 | + text_font="Helvetica", |
| 234 | + text_font_style="bold", |
| 235 | + ) |
| 236 | +) |
| 237 | + |
| 238 | +# Asymptote annotation |
| 239 | +# Asymptote angles: (2k+1)*180/n for n=3 poles, 0 zeros → 60°, 180°, 300° |
| 240 | +# Centroid: (-0 -1 -3) / 3 = -4/3 |
| 241 | +centroid = sum(open_loop_poles) / len(open_loop_poles) |
| 242 | +p.add_layout( |
| 243 | + Label( |
| 244 | + x=centroid - 0.05, |
| 245 | + y=0.35, |
| 246 | + text=f"Centroid σ={centroid:.2f}", |
| 247 | + text_font_size="24pt", |
| 248 | + text_color="#555555", |
| 249 | + text_font="Helvetica", |
| 250 | + text_font_style="italic", |
| 251 | + ) |
| 252 | +) |
| 253 | +p.scatter(x=[centroid], y=[0], size=18, marker="circle", color="#555555", line_color="white", line_width=2, alpha=0.7) |
| 254 | + |
| 255 | +# HoverTool |
| 256 | +hover = HoverTool( |
| 257 | + renderers=[scatter], |
| 258 | + tooltips=[("Pole", "@real{0.00} + @imag{0.00}j"), ("Gain K", "@gain{0.00}"), ("Branch", "@branch")], |
| 259 | + point_policy="snap_to_data", |
| 260 | + mode="mouse", |
| 261 | +) |
| 262 | +p.add_tools(hover) |
| 263 | + |
| 264 | +# Style |
| 265 | +p.title.text_font_size = "72pt" |
| 266 | +p.title.text_color = "#2B2B2B" |
| 267 | +p.title.text_font = "Helvetica" |
| 268 | + |
| 269 | +p.xaxis.axis_label_text_font_size = "48pt" |
| 270 | +p.yaxis.axis_label_text_font_size = "48pt" |
| 271 | +p.xaxis.axis_label_text_font = "Helvetica" |
| 272 | +p.yaxis.axis_label_text_font = "Helvetica" |
| 273 | +p.xaxis.major_label_text_font_size = "36pt" |
| 274 | +p.yaxis.major_label_text_font_size = "36pt" |
| 275 | +p.xaxis.axis_label_text_color = "#3A3A3A" |
| 276 | +p.yaxis.axis_label_text_color = "#3A3A3A" |
| 277 | +p.xaxis.major_label_text_color = "#555555" |
| 278 | +p.yaxis.major_label_text_color = "#555555" |
| 279 | + |
| 280 | +p.xaxis.axis_line_color = None |
| 281 | +p.yaxis.axis_line_color = None |
| 282 | +p.xaxis.major_tick_line_color = None |
| 283 | +p.yaxis.major_tick_line_color = None |
| 284 | +p.xaxis.minor_tick_line_color = None |
| 285 | +p.yaxis.minor_tick_line_color = None |
| 286 | + |
| 287 | +p.grid.grid_line_alpha = 0.12 |
| 288 | +p.grid.grid_line_width = 1.5 |
| 289 | +p.grid.grid_line_color = "#AAAAAA" |
| 290 | + |
| 291 | +p.background_fill_color = "#FAFAFA" |
| 292 | +p.border_fill_color = "white" |
| 293 | +p.outline_line_color = None |
| 294 | + |
| 295 | +p.xaxis.ticker.desired_num_ticks = 10 |
| 296 | +p.yaxis.ticker.desired_num_ticks = 10 |
| 297 | + |
| 298 | +# Legend styling |
| 299 | +p.legend.location = "top_right" |
| 300 | +p.legend.label_text_font_size = "28pt" |
| 301 | +p.legend.label_text_font = "Helvetica" |
| 302 | +p.legend.label_text_color = "#3A3A3A" |
| 303 | +p.legend.background_fill_alpha = 0.85 |
| 304 | +p.legend.background_fill_color = "white" |
| 305 | +p.legend.border_line_color = "#CCCCCC" |
| 306 | +p.legend.border_line_width = 2 |
| 307 | +p.legend.glyph_width = 40 |
| 308 | +p.legend.glyph_height = 40 |
| 309 | +p.legend.spacing = 8 |
| 310 | +p.legend.padding = 15 |
| 311 | +p.legend.margin = 20 |
| 312 | + |
| 313 | +# Save |
| 314 | +export_png(p, filename="plot.png") |
| 315 | +save(p, filename="plot.html", resources=Resources(mode="cdn"), title="Root Locus Plot") |
0 commit comments