Skip to main content

pounce_solve_report/
console.rs

1//! Ipopt-style banner / problem-stats / final-summary printing,
2//! shared by the `pounce` CLI (which re-exports this module as
3//! `print`) and the Python bindings (`pounce.print_banner`,
4//! `Problem.print_problem_stats`, `Solver.print_summary`). Output is
5//! structured to match upstream Ipopt's console layout closely enough
6//! that anyone familiar with `ipopt` can spot at a glance whether
7//! POUNCE is converging similarly.
8
9use pounce_common::types::Number;
10use pounce_nlp::return_codes::ApplicationReturnStatus;
11use pounce_nlp::solve_statistics::SolveStatistics;
12use pounce_nlp::tnlp::{IndexStyle, NlpInfo, SparsityRequest, TNLP};
13use pounce_nlp::tnlp_adapter::{FixedVarTreatment, TNLPAdapter};
14use std::cell::RefCell;
15use std::rc::Rc;
16
17/// Same sentinel Ipopt uses for "no bound": ±1e19. Only referenced by the
18/// unit tests now that `collect_stats` derives bounds from the adapter
19/// classification rather than re-thresholding raw bounds itself.
20#[cfg(test)]
21const BOUND_INF: f64 = 1.0e19;
22
23/// `PartialEq` is load-bearing, not a convenience: `IpoptApplication`
24/// compares each attempt's block against the last one it printed so a retry
25/// does not reprint an identical header. See `emit_problem_stats`.
26#[derive(Debug, Clone, Copy, PartialEq, Eq)]
27pub struct ProblemStats {
28    pub n: i32,
29    pub m: i32,
30    pub nnz_jac_eq: i32,
31    pub nnz_jac_ineq: i32,
32    pub nnz_h_lag: i32,
33    pub var_lower_only: i32,
34    pub var_upper_only: i32,
35    pub var_both: i32,
36    pub var_free: i32,
37    pub n_eq: i32,
38    pub n_ineq: i32,
39    pub ineq_lower_only: i32,
40    pub ineq_upper_only: i32,
41    pub ineq_both: i32,
42}
43
44/// Gather everything the banner block needs, reported over the **reduced**
45/// problem that the algorithm actually solves — i.e. after
46/// `fixed_variable_treatment` removes fixed (`x_l == x_u`) variables under
47/// `make_parameter`. This mirrors Ipopt, whose banner is computed from the
48/// post-`IpTNLPAdapter` problem; computing it from the raw TNLP instead made
49/// pounce over-report variables and bucket fixed vars as "lower and upper
50/// bounds" (#140).
51///
52/// To stay byte-for-byte consistent with the solve, the counts are taken from
53/// a throwaway [`TNLPAdapter`] built with the same options — reusing the exact
54/// production classification (including the `make_parameter → relax_bounds`
55/// auto-switch). The Jacobian / Hessian nnz are read from the raw structure and
56/// filtered to drop entries in fixed-variable columns (the columns Ipopt
57/// removes). Returns `None` if any required TNLP call fails.
58pub fn collect_stats(
59    tnlp: &Rc<RefCell<dyn TNLP>>,
60    lo_inf: Number,
61    up_inf: Number,
62    fixed_treatment: FixedVarTreatment,
63) -> Option<ProblemStats> {
64    let adapter =
65        TNLPAdapter::new_with_options(Rc::clone(tnlp), lo_inf, up_inf, fixed_treatment).ok()?;
66    let cls = adapter.classification();
67    let info: NlpInfo = *adapter.nlp_info();
68    let n_full_x = cls.n_full_x as usize;
69    let m = info.m as usize;
70    let one_based = matches!(info.index_style, IndexStyle::Fortran);
71
72    // --- Variable bound buckets over the reduced (non-fixed) variable set.
73    // `x_l_map` / `x_u_map` hold positions in `x_var` that carry a finite
74    // lower / upper bound. A fixed var under `make_parameter` is absent from
75    // both (it was dropped from `x_var`); under `relax_bounds` it lands in
76    // both — matching Ipopt's banner in either mode.
77    let nv = cls.n_x_var() as usize;
78    let mut has_l = vec![false; nv];
79    let mut has_u = vec![false; nv];
80    for &p in &cls.x_l_map {
81        has_l[p as usize] = true;
82    }
83    for &p in &cls.x_u_map {
84        has_u[p as usize] = true;
85    }
86    let (mut var_lower_only, mut var_upper_only, mut var_both, mut var_free) = (0, 0, 0, 0);
87    for k in 0..nv {
88        match (has_l[k], has_u[k]) {
89            (true, true) => var_both += 1,
90            (true, false) => var_lower_only += 1,
91            (false, true) => var_upper_only += 1,
92            (false, false) => var_free += 1,
93        }
94    }
95
96    // --- Constraint counts / inequality bound buckets, straight from the
97    // classification (equality = `c_map`, inequality = `d_map`).
98    let n_eq = cls.n_c;
99    let n_ineq = cls.n_d;
100    let nd = cls.n_d as usize;
101    let mut has_dl = vec![false; nd];
102    let mut has_du = vec![false; nd];
103    for &p in &cls.d_l_map {
104        has_dl[p as usize] = true;
105    }
106    for &p in &cls.d_u_map {
107        has_du[p as usize] = true;
108    }
109    let (mut ineq_lower_only, mut ineq_upper_only, mut ineq_both) = (0, 0, 0);
110    for k in 0..nd {
111        match (has_dl[k], has_du[k]) {
112            (true, true) => ineq_both += 1,
113            (true, false) => ineq_lower_only += 1,
114            (false, true) => ineq_upper_only += 1,
115            // A "free" inequality has no finite bound on either side (e.g. an
116            // `.nl` range row left fully open). It is still counted in
117            // `n_ineq`, so bucket it under "both" or the printed breakdown
118            // won't sum to the total.
119            (false, false) => ineq_both += 1,
120        }
121    }
122
123    // Which raw rows are equality rows (for the Jacobian split).
124    let mut row_is_eq = vec![false; m];
125    for &r in &cls.c_map {
126        row_is_eq[r as usize] = true;
127    }
128
129    // --- Jacobian split: read the raw structure once, drop fixed-variable
130    // columns (the ones `make_parameter` removes), tally per-row.
131    let nnz_total = info.nnz_jac_g as usize;
132    let (mut nnz_jac_eq, mut nnz_jac_ineq) = (0, 0);
133    if nnz_total > 0 && m > 0 {
134        let mut irow = vec![0_i32; nnz_total];
135        let mut jcol = vec![0_i32; nnz_total];
136        let mut t = tnlp.borrow_mut();
137        if t.eval_jac_g(
138            None,
139            true,
140            SparsityRequest::Structure {
141                irow: &mut irow,
142                jcol: &mut jcol,
143            },
144        ) {
145            for k in 0..nnz_total {
146                let col = if one_based {
147                    (jcol[k] - 1) as usize
148                } else {
149                    jcol[k] as usize
150                };
151                // Skip nonzeros in fixed-variable columns: `full_to_var[col]`
152                // is `-1` for a dropped fixed var.
153                if col >= n_full_x || cls.full_to_var[col] < 0 {
154                    continue;
155                }
156                let row = if one_based {
157                    (irow[k] - 1) as usize
158                } else {
159                    irow[k] as usize
160                };
161                if row < m && row_is_eq[row] {
162                    nnz_jac_eq += 1;
163                } else {
164                    nnz_jac_ineq += 1;
165                }
166            }
167        }
168    }
169
170    // --- Hessian nnz over the reduced problem: drop any entry touching a
171    // fixed-variable row or column. If the TNLP supplies no Hessian
172    // (`eval_h` → false), fall back to the raw count.
173    let mut nnz_h_lag = info.nnz_h_lag;
174    let nnz_h = info.nnz_h_lag as usize;
175    if nnz_h > 0 {
176        let mut irow = vec![0_i32; nnz_h];
177        let mut jcol = vec![0_i32; nnz_h];
178        let mut t = tnlp.borrow_mut();
179        if t.eval_h(
180            None,
181            true,
182            1.0,
183            None,
184            true,
185            SparsityRequest::Structure {
186                irow: &mut irow,
187                jcol: &mut jcol,
188            },
189        ) {
190            let mut kept = 0_i32;
191            for k in 0..nnz_h {
192                let r = if one_based {
193                    (irow[k] - 1) as usize
194                } else {
195                    irow[k] as usize
196                };
197                let c = if one_based {
198                    (jcol[k] - 1) as usize
199                } else {
200                    jcol[k] as usize
201                };
202                if r < n_full_x
203                    && c < n_full_x
204                    && cls.full_to_var[r] >= 0
205                    && cls.full_to_var[c] >= 0
206                {
207                    kept += 1;
208                }
209            }
210            nnz_h_lag = kept;
211        }
212    }
213
214    Some(ProblemStats {
215        n: cls.n_x_var(),
216        m: info.m,
217        nnz_jac_eq,
218        nnz_jac_ineq,
219        nnz_h_lag,
220        var_lower_only,
221        var_upper_only,
222        var_both,
223        var_free,
224        n_eq,
225        n_ineq,
226        ineq_lower_only,
227        ineq_upper_only,
228        ineq_both,
229    })
230}
231
232/// POUNCE wordmark in block letters, printed above the copyright banner.
233const LOGO: [&str; 5] = [
234    "####    ###   #   #  #   #   ####  #####",
235    "#   #  #   #  #   #  ##  #  #      #    ",
236    "####   #   #  #   #  # # #  #      #### ",
237    "#      #   #  #   #  #  ##  #      #    ",
238    "#       ###    ###   #   #   ####  #####",
239];
240
241/// Width of the copyright banner's asterisk rules — wide enough to span
242/// the longest banner text line. The wordmark is centered against this,
243/// and a matching rule is printed above it.
244const BANNER_WIDTH: usize = 80;
245
246/// Print the branded POUNCE ASCII wordmark, mimicking the project logo.
247///
248/// Block letters get a top-lit **steel** sheen (light silver → dark
249/// steel down the rows); three diagonal **molten claw** slashes rake
250/// upper-right → lower-left, glowing bright gold at the top into deep
251/// red at the bottom — the brand logo's look. Emitted through
252/// `anstream::stdout()`, which strips the ANSI when stdout is redirected
253/// or `NO_COLOR` is set (non-TTY sinks get the plain text), with a
254/// 256-color downgrade on non-truecolor terminals. The metallic letters
255/// are tuned for a dark terminal background.
256pub fn print_logo() {
257    use std::io::Write as _;
258    let width = LOGO
259        .iter()
260        .map(|l| l.chars().count())
261        .max()
262        .unwrap_or(1)
263        .max(2);
264    let mut out = anstream::stdout();
265    // Leading rule matching the copyright banner's width, then a blank
266    // line, then the centered wordmark. The rule is left in the terminal's
267    // default color (like the banner's own rules) so it stays distinct on
268    // any background. `anstream` strips the styling when stdout isn't a TTY.
269    let _ = writeln!(out, "{}", "*".repeat(BANNER_WIDTH));
270    let _ = writeln!(out);
271    let pad = " ".repeat(BANNER_WIDTH.saturating_sub(width) / 2);
272    for row in logo_rows(true) {
273        let _ = writeln!(out, "{pad}{row}");
274    }
275    let _ = writeln!(out);
276}
277
278/// Render the POUNCE wordmark as styled rows (one `String` per line):
279/// steel-sheen letters with three molten claw slashes, in the project
280/// palette. Emits ANSI styling only when `color`; otherwise plain
281/// `#`/`/` block characters. Shared by the solve header ([`print_logo`])
282/// and the interactive debugger's open banner (rendered to stderr).
283pub fn logo_rows(color: bool) -> Vec<String> {
284    use pounce_common::style::{ALPHA_HOT, BRIGHT_YEL, TIGER_ORANGE, downgrade, truecolor_enabled};
285
286    fn lerp(a: u8, b: u8, t: f64) -> u8 {
287        (a as f64 + (b as f64 - a as f64) * t)
288            .round()
289            .clamp(0.0, 255.0) as u8
290    }
291    fn mix(a: anstyle::RgbColor, b: anstyle::RgbColor, t: f64) -> anstyle::RgbColor {
292        anstyle::RgbColor(lerp(a.0, b.0, t), lerp(a.1, b.1, t), lerp(a.2, b.2, t))
293    }
294    // Steel sheen (top-lit): light silver at the top row → dark steel at
295    // the bottom. Molten ramp: gold → tiger-orange → deep red top-to-bottom.
296    const STEEL_HI: anstyle::RgbColor = anstyle::RgbColor(0xd2, 0xd6, 0xdc);
297    const STEEL_LO: anstyle::RgbColor = anstyle::RgbColor(0x5c, 0x60, 0x68);
298
299    let rows = LOGO.len();
300    let width = LOGO
301        .iter()
302        .map(|l| l.chars().count())
303        .max()
304        .unwrap_or(1)
305        .max(2);
306    let vfrac = |r: usize| {
307        if rows <= 1 {
308            0.0
309        } else {
310            r as f64 / (rows - 1) as f64
311        }
312    };
313    // Molten color for a claw cell at row `r` (0 = top, hottest).
314    let molten = |r: usize| {
315        let t = vfrac(r);
316        if t < 0.5 {
317            mix(BRIGHT_YEL, TIGER_ORANGE, t / 0.5)
318        } else {
319            mix(TIGER_ORANGE, ALPHA_HOT, (t - 0.5) / 0.5)
320        }
321    };
322
323    let mut grid: Vec<Vec<Option<(char, anstyle::RgbColor)>>> = vec![vec![None; width]; rows];
324    for (r, line) in LOGO.iter().enumerate() {
325        let steel = mix(STEEL_HI, STEEL_LO, vfrac(r));
326        for (c, ch) in line.chars().enumerate() {
327            if ch != ' ' {
328                grid[r][c] = Some((ch, steel));
329            }
330        }
331    }
332    // Three parallel molten claw slashes, upper-right → lower-left (`/`).
333    for &start in &[width / 4, width / 4 + 6, width / 4 + 12] {
334        for r in 0..rows {
335            let c = start + (rows - 1 - r);
336            if c < width {
337                grid[r][c] = Some(('/', molten(r)));
338            }
339        }
340    }
341
342    let truecolor = truecolor_enabled();
343    grid.iter()
344        .map(|row| {
345            let mut rendered = String::new();
346            for cell in row {
347                match cell {
348                    Some((ch, rgb)) if color => {
349                        let style = anstyle::Style::new()
350                            .bold()
351                            .fg_color(Some(downgrade(*rgb, truecolor)));
352                        rendered.push_str(&format!(
353                            "{}{}{}",
354                            style.render(),
355                            ch,
356                            style.render_reset()
357                        ));
358                    }
359                    Some((ch, _)) => rendered.push(*ch),
360                    None => rendered.push(' '),
361                }
362            }
363            rendered.trim_end().to_string()
364        })
365        .collect()
366}
367
368pub fn print_banner(linear_solver: &str) {
369    use std::io::IsTerminal as _;
370
371    // OSC 8 hyperlink so supporting terminals make the URL clickable;
372    // only emitted to a TTY so redirected output stays plain text.
373    const URL: &str = "https://github.com/jkitchin/pounce";
374    let link = if std::io::stdout().is_terminal() {
375        format!("\x1b]8;;{URL}\x1b\\{URL}\x1b]8;;\x1b\\")
376    } else {
377        URL.to_string()
378    };
379
380    let rule = "*".repeat(BANNER_WIDTH);
381    println!("{rule}");
382    println!("This program contains POUNCE, a pure-Rust interior-point optimization solver");
383    println!("for nonlinear, conic, and global problems (its NLP core is ported from Ipopt).");
384    println!("Released under the Eclipse Public License (EPL) — drop-in compatible with Ipopt.");
385    println!("         For more information visit {link}");
386    println!("{rule}");
387    println!();
388    println!(
389        "This is POUNCE version {}, running with linear solver {}.",
390        env!("CARGO_PKG_VERSION"),
391        linear_solver
392    );
393    println!();
394}
395
396pub fn print_problem_stats(s: &ProblemStats) {
397    println!(
398        "Number of nonzeros in equality constraint Jacobian...: {:>8}",
399        s.nnz_jac_eq
400    );
401    println!(
402        "Number of nonzeros in inequality constraint Jacobian.: {:>8}",
403        s.nnz_jac_ineq
404    );
405    println!(
406        "Number of nonzeros in Lagrangian Hessian.............: {:>8}",
407        s.nnz_h_lag
408    );
409    println!();
410    println!(
411        "Total number of variables............................: {:>8}",
412        s.n
413    );
414    println!(
415        "                     variables with only lower bounds: {:>8}",
416        s.var_lower_only
417    );
418    println!(
419        "                variables with lower and upper bounds: {:>8}",
420        s.var_both
421    );
422    println!(
423        "                     variables with only upper bounds: {:>8}",
424        s.var_upper_only
425    );
426    println!(
427        "Total number of equality constraints.................: {:>8}",
428        s.n_eq
429    );
430    println!(
431        "Total number of inequality constraints...............: {:>8}",
432        s.n_ineq
433    );
434    println!(
435        "        inequality constraints with only lower bounds: {:>8}",
436        s.ineq_lower_only
437    );
438    println!(
439        "   inequality constraints with lower and upper bounds: {:>8}",
440        s.ineq_both
441    );
442    println!(
443        "        inequality constraints with only upper bounds: {:>8}",
444        s.ineq_upper_only
445    );
446    println!();
447}
448
449/// Evaluation-callback tallies for the end-of-run summary, decoupled
450/// from any particular TNLP wrapper so both the CLI's `CountingTnlp`
451/// and the Python bindings' callback problem can supply them.
452#[derive(Debug, Clone, Copy, Default)]
453pub struct EvalCounts {
454    pub n_obj: u64,
455    pub n_grad_f: u64,
456    pub n_g: u64,
457    pub n_jac_g: u64,
458    pub n_h: u64,
459}
460
461/// Emit the "the point is outside the model you wrote" warning, in bold red.
462///
463/// This is the one line in the summary block that reports a *defect in the
464/// answer* rather than a residual of the solve, and it is easy to read past in
465/// a wall of sixteen-digit residuals — which is the whole failure mode it
466/// exists to catch: on the LISWET/YAO family a `1e-8` widening moves the
467/// optimum by 25% under `EXIT: Optimal Solution Found`, and the only trace in
468/// the log is a number nobody looks at.
469///
470/// Written through `anstream::stdout()`, which strips the ANSI when stdout is
471/// not a TTY or `NO_COLOR` is set — so redirected logs, the benchmark
472/// harness's scrapes and every `assert!(stdout.contains(...))` in the test
473/// suite see the plain text, byte-for-byte as before. It is deliberately kept
474/// *out* of the upstream-compatible residual table above: that block is
475/// diffed against `ipopt`'s own output, and this line has no counterpart
476/// there to diff against.
477fn print_declared_violation(value: Number) {
478    use std::io::Write as _;
479    let style = anstyle::Style::new()
480        .bold()
481        .fg_color(Some(anstyle::AnsiColor::Red.into()));
482    let mut out = anstream::stdout();
483    let _ = writeln!(out);
484    let _ = writeln!(
485        out,
486        "{}Violation of the model as declared (before the \
487         bound_relax_factor widening): {}{}",
488        style.render(),
489        fmt_ipopt(value),
490        style.render_reset()
491    );
492}
493
494pub fn print_summary(
495    status: ApplicationReturnStatus,
496    stats: &SolveStatistics,
497    counters: &EvalCounts,
498) {
499    println!();
500    println!();
501    println!("Number of Iterations....: {}", stats.iteration_count);
502    println!();
503    println!("                                   (scaled)                 (unscaled)");
504    let row = |label: &str, scaled: f64, unscaled: f64| {
505        println!(
506            "{label}:   {}    {}",
507            fmt_ipopt(scaled),
508            fmt_ipopt(unscaled)
509        );
510    };
511    row(
512        "Objective...............",
513        stats.final_scaled_objective,
514        stats.final_objective,
515    );
516    // The residual rows carry genuinely different values in the two columns
517    // whenever `nlp_scaling` is active, and printing the scaled number twice
518    // hid exactly the discrepancy that matters: on gh #200's `quartc` the
519    // objective row correctly showed `2.49e-06` / `2.49e+02` while dual
520    // infeasibility read `8.38e-09` in *both* columns — when the unscaled value
521    // is `0.84`, eight orders above it. A user auditing a suspicious
522    // certificate was shown a report that agreed with the certificate. The
523    // unscaled statistics were already computed and already surfaced through
524    // the Python bindings; only the console dropped them.
525    row(
526        "Dual infeasibility......",
527        stats.final_dual_inf,
528        stats.final_unscaled_dual_inf,
529    );
530    row(
531        "Constraint violation....",
532        stats.final_constr_viol,
533        stats.final_unscaled_constr_viol,
534    );
535    // Was a hardcoded `0.0` until gh#900: the row asserted a measurement the
536    // solver never took, and it read `0.00e+00` on exactly the solves where it
537    // matters — a `bound_relax_factor`-widened run whose point sits outside
538    // the box the caller wrote. Variable bounds carry no scaling (POUNCE
539    // scales the objective and the constraint rows only), so the one number is
540    // correct in both columns rather than being duplicated for want of a
541    // second.
542    row(
543        "Variable bound violation",
544        stats.final_declared_box_viol,
545        stats.final_declared_box_viol,
546    );
547    row(
548        "Complementarity.........",
549        stats.final_compl,
550        stats.final_unscaled_compl,
551    );
552    row(
553        "Overall NLP error.......",
554        stats.final_kkt_error,
555        stats.final_unscaled_kkt_error,
556    );
557    // gh #528. The strict gate judges the primal term against what each row
558    // can actually represent in floating point, so on a model whose constraint
559    // values run to `~1e8` the raw error above can sit above `tol` beside an
560    // `EXIT: Optimal Solution Found`. Print the number that was tested, but
561    // only when it differs — every `O(1)` model, and every run with the floor
562    // switched off, keeps the summary block byte-identical to upstream's.
563    if stats.final_kkt_error_above_noise.is_finite()
564        && stats.final_kkt_error.is_finite()
565        && stats.final_kkt_error_above_noise != stats.final_kkt_error
566    {
567        println!(
568            "  ...above the per-row floating-point noise floor:   {}",
569            fmt_ipopt(stats.final_kkt_error_above_noise),
570        );
571        println!(
572            "  (the strict convergence test judges this value; the residual \
573             below it is finer than the row's own arithmetic can resolve. \
574             Set primal_noise_floor_kappa = 0 to disable.)"
575        );
576    }
577    // How far outside the model AS DECLARED the returned point sits. This arm
578    // applies the `bound_relax_factor` widening by design -- a feasible-iterate
579    // log-barrier needs `x` strictly inside its bounds -- so the residuals
580    // above are honest residuals of the WIDENED model, and on a row-degenerate
581    // one the two differ by orders: netlib `wood1p` prints `1.71e-14` above at
582    // a point `9.84e-09` outside the declared model. Printed only when a
583    // widening was applied and it actually moved the number, so the block
584    // stays byte-identical to upstream's on everything else.
585    if stats.final_declared_constr_viol.is_finite()
586        && stats.final_declared_constr_viol > stats.final_constr_viol * 10.0
587        && stats.final_declared_constr_viol > 0.0
588    {
589        print_declared_violation(stats.final_declared_constr_viol);
590    }
591    println!();
592    println!();
593    println!(
594        "Number of objective function evaluations             = {}",
595        counters.n_obj
596    );
597    println!(
598        "Number of objective gradient evaluations             = {}",
599        counters.n_grad_f
600    );
601    println!(
602        "Number of equality constraint evaluations            = {}",
603        counters.n_g
604    );
605    println!(
606        "Number of inequality constraint evaluations          = {}",
607        counters.n_g
608    );
609    println!(
610        "Number of equality constraint Jacobian evaluations   = {}",
611        counters.n_jac_g
612    );
613    println!(
614        "Number of inequality constraint Jacobian evaluations = {}",
615        counters.n_jac_g
616    );
617    println!(
618        "Number of Lagrangian Hessian evaluations             = {}",
619        counters.n_h
620    );
621    // gh #819. `Number of Iterations` above is the index of the last row the
622    // iteration table printed, `r` rows included — Ipopt's rule, and now
623    // POUNCE's on every exit path. That total is the right headline, but it
624    // hides the split, and the split is what tells a reader whether an
625    // eleven-second solve was eleven seconds of Newton steps or a restoration
626    // grind. Printed only when restoration actually ran, so a summary from a
627    // solve that never entered it stays byte-identical to upstream's.
628    if stats.restoration_inner_iters > 0 || stats.restoration_calls > 0 {
629        println!(
630            "Number of restoration iterations                     = {} (in {} call{})",
631            stats.restoration_inner_iters,
632            stats.restoration_calls,
633            if stats.restoration_calls == 1 {
634                ""
635            } else {
636                "s"
637            },
638        );
639    }
640    // gh#857. Printed only when the linear solver actually escalated, so a
641    // summary from a solve that never did stays byte-identical to
642    // upstream's — the same rule the restoration split above follows.
643    // Worth a line of its own because an escalation is not visible
644    // anywhere else in a default-verbosity run: it reroutes the rest of
645    // the solve (FERAL's ladder changes which pivots are taken and never
646    // steps back down), and until this line a reader comparing two runs
647    // had no way to tell that apart from an ordinary trajectory
648    // difference.
649    if stats.quality_escalations > 0 {
650        println!(
651            "Number of linear solver quality escalations          = {}",
652            stats.quality_escalations,
653        );
654    }
655    // gh#884. Same rule: printed only when it happened, so a summary from
656    // a solve that never saw the signature stays byte-identical to
657    // upstream's. Worth its own line for the same reason as the
658    // escalation count — a promoted retry reports the *retry's* iteration
659    // count and residuals, so without this the second solve leaves no
660    // trace at all, and a reader comparing two runs would see a
661    // trajectory that came from nowhere.
662    //
663    // This line reports the **signature**, never the retry's verdict, and
664    // that is forced rather than chosen: this block runs once per attempt,
665    // from inside the solve, while the promotion is decided after the last
666    // attempt returns. `dual_divergence_retry_promoted` therefore reads
667    // `false` here even on the run that promotes, and an earlier draft that
668    // printed it produced a summary contradicting the JSON report beside
669    // it. `run_with_dual_divergence_retry` prints the verdict, being the
670    // only place that knows it.
671    //
672    // Two further narrowings. No number is printed: the flag is sticky
673    // across the attempts of one solve, so on a promoted run it is `true`
674    // while *this* attempt's residuals are clean, and a residual beside the
675    // word "detected" would read as its evidence while belonging to the
676    // other attempt. And it is suppressed on `Solve_Succeeded`, because
677    // passing *through* a biactive runaway and recovering is routine on an
678    // MPCC lowering — `mpcc_qpec_small_biactive.nl` sets the flag at
679    // default options and then converges in 29 iterations. A warning over a
680    // correct answer is noise.
681    if stats.dual_divergence_signature && !matches!(status, ApplicationReturnStatus::SolveSucceeded)
682    {
683        println!("Biactive dual divergence (gh#884)                    = detected");
684    }
685    println!(
686        "Total seconds in POUNCE                              = {:.3}",
687        stats.total_wallclock_time_secs
688    );
689}
690
691/// The two lines that end a run: `EXIT: <verdict>` and
692/// `POUNCE <version>: <verdict>`.
693///
694/// Split out of [`print_summary`] because they are the only part of the block
695/// that is about the RUN rather than about the attempt. Every retry driver
696/// re-enters the solve routine, so printing them per attempt produced a
697/// mid-run `POUNCE 0.11.0: Solved To Acceptable Level.` that reads as the
698/// final answer on a run that goes on to report `Optimal Solution Found` —
699/// and a consumer keeping the last `EXIT:` line (p3_control.py does; see
700/// `issue_508_infeasibility_gap_status`) had to hope the last attempt was the
701/// one that shipped. The statistics block above still prints per attempt,
702/// because it describes what that attempt achieved.
703///
704/// The caller decides when the run is over; see
705/// `IpoptApplication::defer_end_verdict`.
706pub fn print_exit_verdict(status: ApplicationReturnStatus) {
707    println!();
708    println!("EXIT: {}", status_message(status));
709    println!();
710    println!(
711        "POUNCE {}: {}",
712        env!("CARGO_PKG_VERSION"),
713        status_message(status)
714    );
715}
716
717/// Emit an Ipopt-style end-of-run summary for the dedicated convex
718/// (LP / QP / conic) IPM path. That path otherwise prints only a compact
719/// one-line result, so the `Number of Iterations....:` and
720/// `Objective...............:` lines the general NLP path emits are missing.
721/// Downstream consumers that parse Ipopt's summary block — notably the
722/// benchmark harness's `extract_obj`/`extract_iters` in
723/// `benchmarks/scripts/run_nl_bench.sh` — then see a null objective and zero
724/// iterations even though the solve succeeded. This prints the same labelled
725/// lines (objective + KKT residual rows) so those consumers capture the real
726/// values. The convex solver reports a single (unscaled, user-sense) objective
727/// and residuals, so the "(scaled)"/"(unscaled)" columns carry the same value.
728pub fn print_convex_summary(
729    iterations: usize,
730    objective: f64,
731    primal_inf: f64,
732    dual_inf: f64,
733    complementarity: f64,
734    kkt_error: f64,
735    // How far the returned `x` sits outside the **declared** variable box —
736    // `QpResiduals::bound_violation` measured against a re-extraction at
737    // `BoundRelax::NONE`, falling back to the solved model's own box when no
738    // widening was applied (the two coincide there). Ipopt's `Variable bound
739    // violation` line, which this arm printed as a hardcoded `0.0` until
740    // gh#900.
741    bound_violation: f64,
742    // How far outside the model AS DECLARED the point sits, when the solve
743    // applied the `bound_relax_factor` widening and the two differ. The
744    // residuals above measure the widened model the solver was handed — the
745    // model its convergence test is about — so without this line a reader
746    // takes `Constraint violation....: 8.68e-13` for their own model's
747    // feasibility when the point is `4.99e-06` outside a declared row
748    // (netlib `afiro`, gh #744/#745). `None` when no widening applied.
749    declared_primal_inf: Option<f64>,
750) {
751    println!();
752    println!();
753    println!("Number of Iterations....: {iterations}");
754    println!();
755    println!("                                   (scaled)                 (unscaled)");
756    let row = |label: &str, v: f64| {
757        println!("{label}:   {}    {}", fmt_ipopt(v), fmt_ipopt(v));
758    };
759    row("Objective...............", objective);
760    row("Dual infeasibility......", dual_inf);
761    row("Constraint violation....", primal_inf);
762    row("Variable bound violation", bound_violation);
763    row("Complementarity.........", complementarity);
764    row("Overall NLP error.......", kkt_error);
765    // Only when the widening actually moved the number: on the vast majority
766    // of models the point satisfies the declared rows to the same order and
767    // an extra line would be noise.
768    if let Some(d) = declared_primal_inf {
769        if d > primal_inf * 10.0 && d > 0.0 {
770            print_declared_violation(d);
771        }
772    }
773    println!();
774}
775
776/// The end-of-run verdict block for the convex path: wall-clock total,
777/// `EXIT:` banner, and the `POUNCE <version>:` line — the same three lines,
778/// in the same order and spelling, that [`print_summary`] ends the NLP
779/// path's log with.
780///
781/// gh #767: the convex log stopped after the residual block, with no `EXIT:`
782/// banner and no terminal status of any kind. `benchmarks/scripts/run_nl_bench.sh`
783/// already compensates with a ladder of convex-specific stdout scrapes, so
784/// every *other* consumer of the CLI had to reimplement that ladder to learn
785/// how a convex-routed solve ended. The phrases come from [`status_message`],
786/// the table the NLP path prints from, so a consumer that recognises Ipopt's
787/// end-of-run vocabulary needs no convex-specific case at all.
788///
789/// The caller supplies the status already mapped onto the NLP-side
790/// [`ApplicationReturnStatus`] (the CLI's `qp_status_to_ars`), for the same
791/// reason: one status vocabulary across both engines.
792pub fn print_convex_end(status: ApplicationReturnStatus, total_seconds: f64) {
793    println!();
794    println!("Total seconds in POUNCE                              = {total_seconds:.3}");
795    println!();
796    println!("EXIT: {}", status_message(status));
797    println!();
798    println!(
799        "POUNCE {}: {}",
800        env!("CARGO_PKG_VERSION"),
801        status_message(status)
802    );
803}
804
805/// Format a number in Ipopt's scientific notation: 16-digit mantissa,
806/// signed 2-digit exponent (e.g. `3.7952009505566139e+03`). Rust's
807/// `{:.16e}` is close but emits a 1-digit exponent without leading
808/// sign, which makes side-by-side diffs against `ipopt` output messy.
809pub fn fmt_ipopt(v: f64) -> String {
810    if v.is_nan() {
811        return "nan".to_string();
812    }
813    if v.is_infinite() {
814        return if v > 0.0 { "inf".into() } else { "-inf".into() };
815    }
816    let s = format!("{:.16e}", v);
817    let Some(e_pos) = s.rfind('e') else {
818        return s;
819    };
820    let (mantissa, exp_part) = s.split_at(e_pos);
821    let exp_str = &exp_part[1..];
822    let (sign, digits) = if let Some(rest) = exp_str.strip_prefix('-') {
823        ('-', rest)
824    } else if let Some(rest) = exp_str.strip_prefix('+') {
825        ('+', rest)
826    } else {
827        ('+', exp_str)
828    };
829    let padded = if digits.len() < 2 {
830        format!("0{digits}")
831    } else {
832        digits.to_string()
833    };
834    format!("{mantissa}e{sign}{padded}")
835}
836
837pub fn status_message(s: ApplicationReturnStatus) -> &'static str {
838    match s {
839        ApplicationReturnStatus::SolveSucceeded => "Optimal Solution Found.",
840        ApplicationReturnStatus::SolvedToAcceptableLevel => "Solved To Acceptable Level.",
841        ApplicationReturnStatus::InfeasibleProblemDetected => {
842            "Converged to a point of local infeasibility. Problem may be infeasible."
843        }
844        ApplicationReturnStatus::SearchDirectionBecomesTooSmall => {
845            "Search Direction is becoming Too Small."
846        }
847        ApplicationReturnStatus::DivergingIterates => {
848            "Iterates diverging; problem might be unbounded."
849        }
850        ApplicationReturnStatus::UserRequestedStop => "Stopping optimization at user request.",
851        ApplicationReturnStatus::FeasiblePointFound => "Feasible Point Found.",
852        ApplicationReturnStatus::MaximumIterationsExceeded => {
853            "Maximum Number of Iterations Exceeded."
854        }
855        ApplicationReturnStatus::RestorationFailed => "Restoration Failed!",
856        ApplicationReturnStatus::ErrorInStepComputation => "Error in step computation.",
857        ApplicationReturnStatus::MaximumCpuTimeExceeded => "Maximum CPU time exceeded.",
858        ApplicationReturnStatus::MaximumWallTimeExceeded => "Maximum wallclock time exceeded.",
859        ApplicationReturnStatus::NotEnoughDegreesOfFreedom => "Not Enough Degrees of Freedom.",
860        ApplicationReturnStatus::InvalidProblemDefinition => "Invalid Problem Definition.",
861        ApplicationReturnStatus::InvalidOption => "Invalid Option.",
862        ApplicationReturnStatus::InvalidNumberDetected => {
863            "Invalid number in NLP function or derivative detected."
864        }
865        ApplicationReturnStatus::UnrecoverableException => "Unrecoverable Exception.",
866        ApplicationReturnStatus::NonIpoptExceptionThrown => "Exception of type generic.",
867        ApplicationReturnStatus::InsufficientMemory => "Insufficient memory.",
868        ApplicationReturnStatus::InternalError => "INTERNAL ERROR: Unknown SolverReturn value.",
869    }
870}
871
872#[cfg(test)]
873mod inequality_tally_tests {
874    //! Regression test for code review L26: the inequality bound-type
875    //! breakdown (`lower_only` / `both` / `upper_only`) must always sum to
876    //! `n_ineq`. A "free" inequality row (no finite bound on either side)
877    //! previously fell through to a no-op arm, so the breakdown summed to
878    //! *less* than the total whenever such a row was present.
879    use super::*;
880    use pounce_common::types::{Index, Number};
881    use pounce_nlp::tnlp::{BoundsInfo, IndexStyle, IpoptCq, IpoptData, Solution, StartingPoint};
882
883    /// Two free variables, three inequality rows of distinct bound types:
884    /// row 0 lower-only, row 1 both, row 2 *free* (the bug trigger). No
885    /// equality rows. The breakdown must sum to `n_ineq == 3`.
886    struct FreeIneqRow;
887    impl TNLP for FreeIneqRow {
888        fn get_nlp_info(&mut self) -> Option<NlpInfo> {
889            Some(NlpInfo {
890                n: 2,
891                m: 3,
892                nnz_jac_g: 3,
893                nnz_h_lag: 0,
894                index_style: IndexStyle::C,
895            })
896        }
897        fn get_bounds_info(&mut self, b: BoundsInfo<'_>) -> bool {
898            b.x_l.iter_mut().for_each(|v| *v = -BOUND_INF);
899            b.x_u.iter_mut().for_each(|v| *v = BOUND_INF);
900            // row 0: lower-only  [0, +inf)
901            b.g_l[0] = 0.0;
902            b.g_u[0] = BOUND_INF;
903            // row 1: both        [0, 1]
904            b.g_l[1] = 0.0;
905            b.g_u[1] = 1.0;
906            // row 2: free        (-inf, +inf) — the regression trigger
907            b.g_l[2] = -BOUND_INF;
908            b.g_u[2] = BOUND_INF;
909            true
910        }
911        fn get_starting_point(&mut self, sp: StartingPoint<'_>) -> bool {
912            sp.x.iter_mut().for_each(|v| *v = 0.0);
913            true
914        }
915        fn eval_f(&mut self, _x: &[Number], _new_x: bool) -> Option<Number> {
916            Some(0.0)
917        }
918        fn eval_grad_f(&mut self, _x: &[Number], _new_x: bool, grad_f: &mut [Number]) -> bool {
919            grad_f.iter_mut().for_each(|v| *v = 0.0);
920            true
921        }
922        fn eval_g(&mut self, _x: &[Number], _new_x: bool, g: &mut [Number]) -> bool {
923            g.iter_mut().for_each(|v| *v = 0.0);
924            true
925        }
926        fn eval_jac_g(
927            &mut self,
928            _x: Option<&[Number]>,
929            _new_x: bool,
930            mode: SparsityRequest<'_>,
931        ) -> bool {
932            match mode {
933                SparsityRequest::Structure { irow, jcol } => {
934                    // one entry per row so the eq/ineq Jacobian split also
935                    // visits each row.
936                    irow.copy_from_slice(&[0, 1, 2]);
937                    jcol.copy_from_slice(&[0, 0, 0]);
938                }
939                SparsityRequest::Values { values } => {
940                    values.copy_from_slice(&[1.0, 1.0, 1.0]);
941                }
942            }
943            true
944        }
945        fn finalize_solution(&mut self, _sol: Solution<'_>, _d: &IpoptData, _q: &IpoptCq) {}
946    }
947
948    #[test]
949    fn free_inequality_row_keeps_breakdown_summing_to_total() {
950        let tnlp: Rc<RefCell<dyn TNLP>> = Rc::new(RefCell::new(FreeIneqRow));
951        let s = collect_stats(
952            &tnlp,
953            -BOUND_INF,
954            BOUND_INF,
955            FixedVarTreatment::MakeParameter,
956        )
957        .expect("collect_stats succeeds");
958
959        assert_eq!(s.n_eq, 0, "no equality rows");
960        assert_eq!(s.n_ineq, 3, "all three rows are inequalities");
961        // The headline invariant L26 flagged: the three printed buckets must
962        // account for every inequality row.
963        let bucket_sum: Index = s.ineq_lower_only + s.ineq_both + s.ineq_upper_only;
964        assert_eq!(
965            bucket_sum, s.n_ineq,
966            "ineq bound-type breakdown ({} lower + {} both + {} upper) must sum to n_ineq={}",
967            s.ineq_lower_only, s.ineq_both, s.ineq_upper_only, s.n_ineq
968        );
969        // The free row is bucketed under "both" alongside the genuine
970        // both-bounded row 1.
971        assert_eq!(s.ineq_lower_only, 1);
972        assert_eq!(s.ineq_upper_only, 0);
973        assert_eq!(s.ineq_both, 2);
974    }
975
976    /// #140 regression. Three variables, the middle one fixed
977    /// (`x_l == x_u`). Under the default `make_parameter` the banner must
978    /// report the *reduced* problem: the fixed var is dropped from the total,
979    /// is NOT bucketed as "lower and upper bounds", and its Jacobian column is
980    /// excluded from the nnz tally — matching Ipopt (and the problem the
981    /// algorithm actually solves). Previously the banner walked the raw TNLP
982    /// and over-reported all three.
983    struct OneFixedVar;
984    impl TNLP for OneFixedVar {
985        fn get_nlp_info(&mut self) -> Option<NlpInfo> {
986            Some(NlpInfo {
987                n: 3,
988                m: 1,
989                nnz_jac_g: 3,
990                nnz_h_lag: 0,
991                index_style: IndexStyle::C,
992            })
993        }
994        fn get_bounds_info(&mut self, b: BoundsInfo<'_>) -> bool {
995            // var 0: lower-only [0, +inf)   var 1: FIXED at 2   var 2: free
996            b.x_l[0] = 0.0;
997            b.x_u[0] = BOUND_INF;
998            b.x_l[1] = 2.0;
999            b.x_u[1] = 2.0;
1000            b.x_l[2] = -BOUND_INF;
1001            b.x_u[2] = BOUND_INF;
1002            // one equality row (keeps n_x_var=2 >= n_c=1, so no relax switch)
1003            b.g_l[0] = 0.0;
1004            b.g_u[0] = 0.0;
1005            true
1006        }
1007        fn get_starting_point(&mut self, sp: StartingPoint<'_>) -> bool {
1008            sp.x.iter_mut().for_each(|v| *v = 0.0);
1009            true
1010        }
1011        fn eval_f(&mut self, _x: &[Number], _new_x: bool) -> Option<Number> {
1012            Some(0.0)
1013        }
1014        fn eval_grad_f(&mut self, _x: &[Number], _new_x: bool, grad_f: &mut [Number]) -> bool {
1015            grad_f.iter_mut().for_each(|v| *v = 0.0);
1016            true
1017        }
1018        fn eval_g(&mut self, _x: &[Number], _new_x: bool, g: &mut [Number]) -> bool {
1019            g.iter_mut().for_each(|v| *v = 0.0);
1020            true
1021        }
1022        fn eval_jac_g(
1023            &mut self,
1024            _x: Option<&[Number]>,
1025            _new_x: bool,
1026            mode: SparsityRequest<'_>,
1027        ) -> bool {
1028            match mode {
1029                // The equality row touches all three columns, including the
1030                // fixed var (col 1) that must be filtered out.
1031                SparsityRequest::Structure { irow, jcol } => {
1032                    irow.copy_from_slice(&[0, 0, 0]);
1033                    jcol.copy_from_slice(&[0, 1, 2]);
1034                }
1035                SparsityRequest::Values { values } => {
1036                    values.copy_from_slice(&[1.0, 1.0, 1.0]);
1037                }
1038            }
1039            true
1040        }
1041        fn finalize_solution(&mut self, _sol: Solution<'_>, _d: &IpoptData, _q: &IpoptCq) {}
1042    }
1043
1044    #[test]
1045    fn make_parameter_banner_reports_reduced_problem() {
1046        let tnlp: Rc<RefCell<dyn TNLP>> = Rc::new(RefCell::new(OneFixedVar));
1047        let s = collect_stats(
1048            &tnlp,
1049            -BOUND_INF,
1050            BOUND_INF,
1051            FixedVarTreatment::MakeParameter,
1052        )
1053        .expect("collect_stats succeeds");
1054
1055        // Fixed var removed: 3 raw vars → 2 optimized.
1056        assert_eq!(s.n, 2, "fixed variable must be dropped from the total");
1057        assert_eq!(s.var_both, 0, "fixed var must NOT count as lower-and-upper");
1058        assert_eq!(s.var_lower_only, 1, "var 0 is lower-only");
1059        assert_eq!(s.var_free, 1, "var 2 is free");
1060        assert_eq!(s.var_upper_only, 0);
1061        // The fixed column is excluded from the Jacobian nnz.
1062        assert_eq!(s.n_eq, 1);
1063        assert_eq!(
1064            s.nnz_jac_eq, 2,
1065            "fixed-var column dropped from the Jacobian"
1066        );
1067        assert_eq!(s.nnz_jac_ineq, 0);
1068    }
1069
1070    #[test]
1071    fn relax_bounds_banner_keeps_fixed_variable() {
1072        // Under relax_bounds the fixed var stays in the optimization and is
1073        // reported as a lower-and-upper-bounded variable — matching Ipopt.
1074        let tnlp: Rc<RefCell<dyn TNLP>> = Rc::new(RefCell::new(OneFixedVar));
1075        let s = collect_stats(&tnlp, -BOUND_INF, BOUND_INF, FixedVarTreatment::RelaxBounds)
1076            .expect("collect_stats succeeds");
1077
1078        assert_eq!(s.n, 3, "relax_bounds keeps the fixed variable");
1079        assert_eq!(
1080            s.var_both, 1,
1081            "fixed var reported as lower-and-upper bounded"
1082        );
1083        assert_eq!(s.var_lower_only, 1);
1084        assert_eq!(s.var_free, 1);
1085        assert_eq!(s.nnz_jac_eq, 3, "all columns retained under relax_bounds");
1086    }
1087}