Data: Solver Sensitivity report (#501)

* Data: Solver Sensitivity report

Builds on the Answer and Limits reports. reports: [sensitivity] (on
data.solver and data.solverReport) writes a Sensitivity report for
models solved by Simplex LP or GRG Nonlinear without integer
constraints and with a Max or Min objective.

Simplex LP: the model is linearised at the solution and solved again,
and the report is read from the final tableau: each variable's reduced
cost and the range of its objective coefficient over which the basis
stays optimal, and each constraint's shadow price and the range of its
right side over which that price holds (1E+30 for no limit).

GRG Nonlinear: difference gradients at the solution give each binding
constraint's Lagrange multiplier (a least-squares fit on the variables
inside their bounds, dropping constraints whose multiplier has the
wrong sign) and each variable's reduced gradient.

The Solver Results dialog offers the report.

* Data: Solver Results offers the Sensitivity report only where it applies

Integer and binary models, the Evolutionary method and Value Of goals
have no Sensitivity report, and choosing it with other reports made the
whole report request fail. Solver Results now disables the Sensitivity
choice for such models (with the reason on hover) and leaves it out of
the request, so the other chosen reports are still written.

sensitivity_applies (in cmd::solver) reads the model from data.solver
params, as the dialog holds them; a test checks that it agrees with the
report for each kind of model.

* Data: Solver bounds the nonlinear Sensitivity report's work

The Lagrange multipliers of a GRG Nonlinear model are found by least
squares over the binding constraints; a constraint whose multiplier has
the wrong sign is dropped and the rest are solved again. Each pass costs
about 2·k²·f operations for k binding constraints and f free variables,
and each pass drops at most one constraint, so a model at the
1,000-constraint limit could take minutes and hold the window.

The passes now share a limit of about 10⁹ operations (about a second).
Past it the report is refused with a message, before anything changes:
the workbook and the Undo history stay as they were, as for the other
report refusals. The limit is lowered in tests so the refusal can be
reached with a small model.

* Data: Solver refuses a Simplex Sensitivity report for a point that is not optimal

The Simplex LP Sensitivity report solves the linearised model again to
read its final tableau, but shows the final values it was given.
data.solverReport accepts any values, so a report could pair values that
are not optimal (or not feasible) with the shadow prices, reduced costs
and ranges of the true optimum. The report now checks that the given
point satisfies every constraint and bound and reaches the optimal
objective (both within the constraint precision; another optimal point
with the same objective is accepted) and is refused otherwise, before
anything changes.

* Data: Solver compares the Sensitivity report's objective in both directions

The check that a point is optimal compared its objective with the
optimum in one direction only, as for an inequality. A point just
outside a bound (within the constraint precision) can still have an
objective far below the minimum when its costs are large, and that
point was accepted. The objective is now compared by the size of the
difference, so such a point is refused.

* Data: Solver refuses a nonlinear Sensitivity report for a point that is not optimal

data.solverReport accepts any finite variable values. The Simplex LP
report already refuses values that are not optimal, but the GRG
Nonlinear report described whatever point it was given: for example, a
minimisation of A1^2 given A1 = 1 found no binding rows and reported a
reduced gradient of 2 instead of requiring the optimum at A1 = 0.

The nonlinear report now checks, from the gradients it already
computes, that the values are a Karush-Kuhn-Tucker point before
anything changes: every constraint holds within the tolerance the
nonlinear method accepts a solution within, every bound holds within
the constraint precision, each variable inside its bounds has a reduced
gradient near zero, and a variable at a bound has a reduced gradient
that points only out of that bound for the objective's direction. A
binding row whose multiplier has the wrong sign is still dropped, and
the gradient it would have explained now causes the refusal. The
stationarity tolerance is 1E-3 of the largest objective-gradient
component (at least 1), well above the reduced gradients left at optima
the method solves.

* Data: Solver writes the other reports when the Sensitivity report does not apply to the solution

Since the nonlinear Sensitivity report refuses values that are not
optimal, data.solver with reports that include sensitivity failed as a
whole when the run stopped before an optimum (an evaluation limit) or
the values were not one, so the solution and the Answer and Limits
reports were lost. Solver Results also offered Sensitivity after any
solve, and ticking it with Answer after such a run wrote no report.

data.solver now writes the Sensitivity report only when the run reached
an optimal or converged solution that the report accepts. Otherwise the
solve still applies, the other reports are written, and the result
lists the report in skippedReports with the reason. Models that have no
Sensitivity report (integer, Evolutionary, Value Of) are refused as
before. data.solverReport still refuses values that are not optimal;
its new skipNotOptimal parameter leaves the report out instead. Solver
Results offers Sensitivity only after an optimal or converged status,
with the reason on hover, passes skipNotOptimal, and shows the reason
when the report is left out.

* Data: Solver skips a Sensitivity report in a protected workbook without refusing the call

Structure protection is now checked only when a report sheet will be written. A data.solverReport request with skipNotOptimal whose only report is skipped changes nothing, so it succeeds in a workbook whose structure is protected and returns skippedReports; a request that would write a sheet is still refused.

* Data: Solver refuses report sheets in a protected workbook before building them

A request that must write a report sheet (any Answer or Limits report, or a Sensitivity report that is not allowed to be skipped) is refused before any report is built, so a protected workbook does not run the Limits report's optimisations first. The check after the reports are built remains for a Sensitivity report that may be skipped.

* Data: Solver reports Simplex sensitivity only at the vertex it describes

The Simplex LP Sensitivity report solves the linearised model again and
takes its prices and ranges from that final tableau, but it discarded
the vertex the solve returned and showed the values it was given. When
the optimum is not unique, another optimal point passed the optimality
check, so the report paired those values with the ranges of a different
vertex. The report now also requires the given values to match the
re-solved vertex within the constraint precision, and refuses another
optimum with a message that says so. A Sensitivity report that data.solver
or Solver Results may skip is left out instead, as for other points the
report does not support. Values from a Simplex LP solve still match, so
the usual Solver Results flow is unchanged.

* Data: Solver gives a redundant constraint's right side no allowable change

The Simplex LP Sensitivity report ranges each right side over the basic
values of the real columns only. A row that the other rows make
redundant (max x with x = 2 and 2x = 4) keeps its artificial column
basic at zero after phase one, and changing either right side alone
moves that artificial, so no change has a solution. The report showed
1E+30 and 2 for the first row instead. When a basic artificial has a
nonzero entry in the column of the right side being changed, both
allowable changes are now 0. The objective-coefficient ranges are
unaffected: such a row has no cost and no entries in the real columns,
and a cost change does not affect feasibility.

* Data: Solver builds the Simplex Sensitivity report's program at the solve's start values

The Simplex LP solve linearises the model at the start values, but the
Sensitivity report linearised it again at the solution. The difference
coefficients then differ in their last bits, and when the optimum is
not unique that can break a tie between pivots the other way: the
report's solve returned the other vertex and refused the solve's own
values. With max 0.1x + 0.1y over a tied edge, more than half of the
sampled starting points were refused.

The solve and the report now build the program with one function,
simplex_program, and the report builds it from the start values it is
given (data.solver's start values, which Solver Results passes back as
original), so its solve repeats the solve's pivots and returns the same
vertex. Start values at which the model cannot be evaluated fall back
to the solution, and values that are not the solve's are still refused
by the vertex check.

* Data: Solver's Sensitivity report sets its noise levels from the coefficients, not the objective's value

The figures a Sensitivity report prints as 0 were those below 1E-9 of the objective's value (1E-6 for the nonlinear figures), so a constant in the objective hid real prices: max 1E9 + x with x <= 10 wrote the row's shadow price, exactly 1, as 0, although the constant changes nothing about the optimisation. A shadow price and a reduced cost are combinations of the objective coefficients, so Simplex LP now takes the noise level from the largest coefficient; a reduced gradient is a part of the objective's gradient and a Lagrange multiplier objective per unit of its row's gradient, so GRG Nonlinear takes it from the gradient's size, and the multipliers' sign test compares lambda times the row's gradient size with it. Allowable changes keep the tableau's pivot tolerance as before.

* Data: Solver Sensitivity report tests for variables at a bound and for both goals

Hand-derived cases pin the report's conventions. Simplex LP: a maximum's variable at its upper bound (reduced cost 1, Allowable Increase 1E+30, Allowable Decrease 1), a minimum with a >= row and a variable at its upper bound (reduced cost -1, Allowable Increase 1, Allowable Decrease 1E+30), and Winston's Dakota Furniture problem, whose variable at zero has reduced cost -5 with Allowable Increase 5 and Allowable Decrease 1E+30, with every shadow price and allowable range worked out from the final basis. GRG Nonlinear: the Lagrange multiplier is the derivative of the optimal objective with respect to the right side and the reduced gradient that with respect to the variable's bound for a minimum (-1) and the negated maximum (+1), and the mirror image of the existing minimum with a variable at its bound (multiplier +2, reduced gradient -4).

* Docs: ROADMAP.md records the Solver Sensitivity report

The detail documents (docs/gaps.md, docs/target-app-parity.md) record the report, but the roadmap summary did not: its header change note, progress log and revision history now carry it.

* Data: Solver's Sensitivity report bounds each figure's rounding noise by its own terms, not the model's largest coefficient or gradient

The noise level a figure had to exceed to print was 1E-9 of the largest objective coefficient (Simplex LP) or 1E-6 of the largest gradient component (GRG Nonlinear) for every figure, so an unrelated large coefficient hid real values: max 1E9*x + y with x <= 1 and y <= 1 as rows wrote the second row's shadow price, exactly 1, as 0, and the same model's Lagrange multiplier 1 and a third variable's reduced cost or gradient -1 went the same way. The rounding error of a sum scales with the sum of its terms' sizes, so each figure now has its own bound from the terms it is made of. Simplex LP: a shadow price y_i = sum over the basis rows of c_B(r) * Binv(r, i) is bounded by 1E-9 of those terms' sizes (the row's unit column of the final tableau holds Binv's column i); a reduced cost cf_j - sum of a_ij * y_i by 1E-9 of |cf_j| plus |a_ij| times each price's term sizes, so a price that is a small difference of large terms carries their error into the reduced costs it is used in. GRG Nonlinear: a multiplier solved from its row of the normal equations is bounded by 1E-6 of its own terms (the row's gradient times the objective's, and its coupling terms with the other active rows times their multipliers) over the row's squared norm; the multipliers' sign test uses the same bound, so a row leaves the active set on the criterion the report rounds with; a reduced gradient is bounded by 1E-6 of its own component and the multipliers' terms plus the multipliers' bounds times |g_ij|. Nearly parallel active gradients amplify every multiplier's error beyond a bound local to a row; such multipliers print as computed. A constant in the objective still changes no bound, since none depends on the objective's value. The nonlinear method stops on the 1E9 model from (0, 0) with "unbounded", so its two regression tests ask for the report at the optimum through data.solverReport.

* Data: Solver's Sensitivity report tests a nonlinear point's optimality by each variable's and row's own scale, not the model's largest

The Karush-Kuhn-Tucker check that refuses a Sensitivity report for a point that is not optimal accepted a reduced gradient within 1E-3 of the largest gradient component, and a row violation within ten times the constraint precision of the largest constant right side. Both scales were the model's, not the figure's: on max 1E9*x + y with x <= 1 and y <= 1 as rows the tolerance was 1E6, so (1, 0.5), where raising y gains 1 per unit, passed as optimal and the report was written for it; on max y - (x - 2)^2 with x <= 1E9 and y <= 1, (2, 1.0001) passed although the second row is over by 1E-4 against its own precision of 2E-5. Each reduced gradient is now tested against 1E-3 of the terms it is made of (the objective's component and the multipliers' parts, at least 1, the same floor as before), and each row against ten times the constraint precision of its own right side. The nonlinear method itself accepts a solution relative to the largest constant right side; a row over by more than its own precision is not satisfied whatever the other rows' sizes, so the report refuses it. The term sizes are shared with the reduced gradients' rounding bound. Two regression tests, each proving the workbook is left as it was and that the optimum is still reported.

* Data: Solver refuses a Sensitivity report the model has none of before any report is built or the solve runs

The refusals that apply to the model rather than the solution (integer or binary constraints, the Evolutionary method, a Value Of goal) were made inside the report loop, when it reached the Sensitivity report, so a request listing the Limits or Answer report first built those reports for an error: a Value Of model asking for limits and sensitivity ran the Limits report's two optimisations per variable and received only the refusal, and data.solver ran its whole solve before it. prepare_reports now checks the model's refusal before the loop, in both modes (a refusal of the model was already an error in both; only the solution's not-optimal refusal is skipped under skipNotOptimal), and data.solver checks it before the solve. The regression test asks for limits and sensitivity on an integer model, which now gets the Sensitivity refusal instead of the Limits report's, and on a Value Of model, which gets the Sensitivity refusal with no sheet written.

---------

Co-authored-by: FredApps <FredApps@users.noreply.github.com>

ROADMAP.md

@@ -2,7 +2,7 @@
 
 **Stage: alpha** · next: beta, ~25 points of "ready for real work" and ~160–240 h away
 
-> **Last reviewed:** 2026-11-01 · **Last updated:** 2026-11-02 · **Change:** minor (floating objects: one selection, eight resize handles, keyboard delete; earlier: Chart data tables: live source values beneath column, bar, line, area and compatible combo charts; legend keys and borders; XLSX round-trip and print/PDF; earlier: Formulas › Function Library as labelled buttons with category menus; Recently Used and Insert Function categories; XLSX workbook-view selection, bounded tab ratios and sheet-edit regressions; readiness estimates unchanged; earlier: CSV and text open detects each column's date order, so day-first dates are not silently swapped) · **Target:** Microsoft Excel (Microsoft 365)
+> **Last reviewed:** 2026-11-01 · **Last updated:** 2026-11-02 · **Change:** minor (What-if and analysis: Solver Sensitivity report for Simplex LP and GRG Nonlinear; earlier: floating objects: one selection, eight resize handles, keyboard delete; earlier: Chart data tables: live source values beneath column, bar, line, area and compatible combo charts; legend keys and borders; XLSX round-trip and print/PDF; earlier: Formulas › Function Library as labelled buttons with category menus; Recently Used and Insert Function categories; XLSX workbook-view selection, bounded tab ratios and sheet-edit regressions; readiness estimates unchanged; earlier: CSV and text open detects each column's date order, so day-first dates are not silently swapped) · **Target:** Microsoft Excel (Microsoft 365)
 
 GridCraft aims at full Microsoft Excel parity, and at being better in places: faster, open (XLSX
 is the native format), fully scriptable by agents (MCP, control channel, CLI), and on the web.
@@ -135,6 +135,7 @@ Ranked; estimates in Opus 5.5 agent-hours. Detail and the full beta-gate list in
 
 | Date | What landed |
 |---|---|
+| 2026-11-02 | Solver writes a Sensitivity report: Simplex LP shadow prices, reduced costs and allowable ranges from the final tableau; GRG Nonlinear reduced gradients and Lagrange multipliers from difference gradients, at a Karush-Kuhn-Tucker point only; Solver Results offers it where it applies (no score change) |
 | 2026-11-02 | Floating objects: charts, pictures, shapes and slicers share one selection with a frame of eight resize handles (each moves only its own edges, also on right-to-left sheets); Delete or Backspace deletes the selected object and Escape drops the selection; slicers follow row and column edits and sorts like the other objects; chart series follow cut-moves and sheet renames and deletes; a copied sheet gives its objects fresh ids; `ui.inspect` reports `selectedObject` and `ui.object.select` selects any object (#190, #213) |
 | 2026-11-01 | Chart data tables are available from Chart Design › Add Chart Element and Format Chart Area. Tables follow live source values and number formats, keep blank cells blank, offer legend keys and border settings, and survive XLSX save/reopen; the shared renderer includes them in print and PDF. Readiness estimates unchanged. |
 | 2026-11-01 | XLSX keeps the first visible workbook window's active sheet and display settings together, falling back to the first view when all are hidden. Tab ratios are normalized to 0–1000 on import and export; invalid numeric text uses 600. Synthetic round trips cover isolated attributes, missing views, and active tabs after sheet deletion/reordering; firstSheet remains separate #2264 work (#3057). |
@@ -252,6 +253,7 @@ Ranked; estimates in Opus 5.5 agent-hours. Detail and the full beta-gate list in
 
 | Date | Change | Summary |
 |---|---|---|
+| 2026-11-02 | minor | What-if and analysis: Solver Sensitivity report for Simplex LP and GRG Nonlinear (no score change) |
 | 2026-11-02 | minor | Progress log: one selection and eight resize handles for every floating object |
 | 2026-11-01 | minor | Chart data tables are available from Chart Design › Add Chart Element and Format Chart Area. Tables follow live source values and number formats, keep blank cells blank, offer legend keys and border settings, and survive XLSX save/reopen; the shared renderer includes them in print and PDF. Readiness estimates unchanged. |
 | 2026-11-01 | trivial | Workbook views: selected window settings stay together, tab ratios are bounded, and regressions cover isolated settings and sheet edits; readiness estimates unchanged (#3057) |

crates/engine/src/cmd/solver.rs

@@ -40,15 +40,16 @@ mod model_range;
 mod nlp;
 mod problem;
 mod reports;
+mod sensitivity;
 
 use evolutionary::{EvoEnd, evolutionary};
-use lp::{Outcome, Sense, branch_and_bound, linearize, solve_lp};
+use lp::{Outcome, SimplexProgram, branch_and_bound, simplex_program, solve_lp};
-pub use model::saved_model;
 use model::{Goal, Method, Model, model_params, parse_model, save_model};
+pub use model::{saved_model, sensitivity_applies, sensitivity_status};
 use model_range::{load_model_cells, reset_model, save_model_cells};
 use nlp::{nlp, project};
 use problem::{Evaluator, as_num, build_problem};
-use reports::{ReportInput, insert_report_selections, prepare_reports, report, report_kinds, write_reports};
+use reports::{NotOptimal, ReportInput, insert_report_selections, prepare_reports, report, report_kinds, write_reports};
 
 const CMD: &str = "data.solver";
 const REPORT_CMD: &str = "data.solverReport";
@@ -78,7 +79,7 @@ pub fn specs() -> Vec<CommandSpec> {
             "Solver…",
             ["Data", "Analysis"],
             None,
-            "{objective: \"B10\" (formula cell), goal?: max|min|value (default max), value?: number (for goal value), variables: \"B2:B4,D1\" (≤200 constant cells, same sheet), constraints?: [\"B5 <= 10\", \"B2:B4 >= C2:C4\", \"B2 int\", \"B3 bin\", {lhs: \"B6\", op: \"=\", rhs: 3 | \"C6\"}], nonNegative?: bool (default true: variables without a lower bound get 0), method?: grg|simplex|evolutionary (default grg), populationSize?: 100, mutationRate?: 0.075, randomSeed?: 0 (a fixed seed, so runs repeat), convergence?: 0.0001, generationsWithoutImprovement?: 50, requireBounds?: bool (default true: every variable needs both bounds; false searches within ±100 of the start) (Evolutionary options), keep?: bool (default true; false restores the original values), saveModel?: bool (default true: hidden names, solver_* as Excel names them, with the options; the evaluation limit and generations without improvement, which Excel lacks, as gridcraft_solver_evals and gridcraft_solver_stall), maxEvaluations?: 20000 (at most 200000), integerOptimality?: percent (default 0), saved?: bool (solve the sheet's saved model; other params override), reports?: [answer|limits] (written on new sheets when a solution is found)} → {status: optimal|converged|cannotImprove|iterationLimit|infeasible|unbounded|notLinear, message, objective, values: {cell: value}, original: {cell: value}, model: the model as params, evaluations, changed, reports?: [sheet names]}. No `objective` → {model: the saved model} (the menu opens the dialog).",
+            "{objective: \"B10\" (formula cell), goal?: max|min|value (default max), value?: number (for goal value), variables: \"B2:B4,D1\" (≤200 constant cells, same sheet), constraints?: [\"B5 <= 10\", \"B2:B4 >= C2:C4\", \"B2 int\", \"B3 bin\", {lhs: \"B6\", op: \"=\", rhs: 3 | \"C6\"}], nonNegative?: bool (default true: variables without a lower bound get 0), method?: grg|simplex|evolutionary (default grg), populationSize?: 100, mutationRate?: 0.075, randomSeed?: 0 (a fixed seed, so runs repeat), convergence?: 0.0001, generationsWithoutImprovement?: 50, requireBounds?: bool (default true: every variable needs both bounds; false searches within ±100 of the start) (Evolutionary options), keep?: bool (default true; false restores the original values), saveModel?: bool (default true: hidden names, solver_* as Excel names them, with the options; the evaluation limit and generations without improvement, which Excel lacks, as gridcraft_solver_evals and gridcraft_solver_stall), maxEvaluations?: 20000 (at most 200000), integerOptimality?: percent (default 0), saved?: bool (solve the sheet's saved model; other params override), reports?: [answer|sensitivity|limits] (written on new sheets when a solution is found)} → {status: optimal|converged|cannotImprove|iterationLimit|infeasible|unbounded|notLinear, message, objective, values: {cell: value}, original: {cell: value}, model: the model as params, evaluations, changed, reports?: [sheet names], skippedReports?: [{report: \"sensitivity\", reason}] (a Sensitivity report left out because the run did not reach an optimum or the values are not one; the solve and the other reports still apply)}. No `objective` → {model: the saved model} (the menu opens the dialog).",
             has_doc,
             run
         ),
@@ -87,7 +88,7 @@ pub fn specs() -> Vec<CommandSpec> {
             "Solver Reports",
             [],
             None,
-            "{reports: [answer|limits], values?: {cell: number} (the solution; default the variable cells' current values), original?: {cell: number} (the values before solving; default `values`), message?: text (the result line), model?: data.solver params (default the sheet's saved model, with any data.solver params given here over it)} → {sheets: [names]}. Writes Answer and Limits reports on new sheets before the model's sheet; Limits needs Simplex LP or GRG Nonlinear and no integer constraints. data.solver's result can be passed as is.",
+            "{reports: [answer|sensitivity|limits], values?: {cell: number} (the solution; default the variable cells' current values), original?: {cell: number} (the values before solving; default `values`), message?: text (the result line), skipNotOptimal?: bool (default false: a Sensitivity report for values that are not optimal refuses the call; true leaves it out and writes the others), model?: data.solver params (default the sheet's saved model, with any data.solver params given here over it)} → {sheets: [names], skippedReports?: [{report, reason}]}. Writes Answer, Sensitivity and Limits reports on new sheets before the model's sheet; Sensitivity and Limits need Simplex LP or GRG Nonlinear and no integer constraints, and Sensitivity a Max or Min objective. data.solver's result can be passed as is.",
             has_doc,
             report
         ),
@@ -146,24 +147,13 @@ fn optimise(ev: &mut Evaluator, m: &Model, goal: Goal, x0: &[f64], bounds: Bound
         Outcome::Infeasible
     } else {
         match m.method {
-            Method::Simplex => match linearize(ev, x0) {
+            Method::Simplex => match simplex_program(ev, goal, x0) {
                 Err("nonlinear") | Err("error") => {
                     not_linear = true;
                     Outcome::Infeasible
                 }
                 Err(_) => Outcome::Limit,
-                Ok(((cf, f_const), lin)) => {
+                Ok(SimplexProgram { c, rows, mut work, .. }) => {
-                    let mut rows: Vec<(Vec<f64>, Sense, f64)> = lin.iter().zip(&ev.rows).map(|((a, k), r)| (a.clone(), r.sense, -k)).collect();
-                    let c: Vec<f64> = match goal {
-                        Goal::Max => cf.iter().map(|v| -v).collect(),
-                        Goal::Min => cf.clone(),
-                        Goal::Value(v) => {
-                            rows.push((cf.clone(), Sense::Eq, v - f_const));
-                            vec![0.0; cf.len()]
-                        }
-                    };
-                    let size = (rows.len() as u64 + 2 * x0.len() as u64 + 1) * (4 * x0.len() as u64 + 2 * rows.len() as u64 + 1);
-                    let mut work = LP_WORK.max(size);
                     let mut relax = |l: &[f64], h: &[f64]| solve_lp(&c, &rows, l, h, &mut work);
                     if any_int { branch_and_bound(lo, hi, integer, MAX_NODES_LP, gap, &mut relax).0 } else { relax(lo, hi) }
                 }
@@ -218,6 +208,14 @@ fn run(s: &mut Session, p: &Json) -> Result<Json> {
     let save = bool_param(&params, "saveModel").unwrap_or(true);
     let reports = report_kinds(CMD, &params)?;
     let prob = build_problem(&model);
+    // A Sensitivity report the model has none of (integer, Evolutionary, Value Of) fails the
+    // call, so it is refused before the solve runs rather than after it.
+    if let Some(kinds) = &reports
+        && kinds.contains(&reports::ReportKind::Sensitivity)
+        && let Some(why) = sensitivity::sensitivity_refusal(&model, &prob)
+    {
+        return Err(EngineError::Other(why.into()));
+    }
     let d = s.doc()?;
     let sh = d.wb.sheet(model.sheet).ok_or(EngineError::NoDocument)?;
     let x0: Vec<f64> = model.vars.iter().map(|c| as_num(&sh.value(*c))).map(|v| if v.is_finite() { v } else { 0.0 }).collect();
@@ -277,12 +275,14 @@ fn run(s: &mut Session, p: &Json) -> Result<Json> {
     }
     let original: serde_json::Map<String, Json> = model.vars.iter().zip(&x0).map(|(c, v)| (c.a1(), json!(v))).collect();
     let changed = keep && solution.is_some();
-    let built = match (&reports, &solution) {
+    // A Sensitivity report the solution doesn't support is left out (with its reason), so the
+    // solve and the other reports still apply.
+    let (built, skipped) = match (&reports, &solution) {
         (Some(kinds), Some(x)) => {
             let inp = ReportInput { model: &model, prob: &prob, original: x0, solution: x.clone(), message: message.to_string(), budget };
-            prepare_reports(s.doc()?, kinds, &inp)?
+            prepare_reports(s.doc()?, kinds, &inp, NotOptimal::Skip { optimum: sensitivity_status(status) })?
         }
-        _ => Vec::new(),
+        _ => (Vec::new(), Vec::new()),
     };
     let count = built.len();
     let before = s.doc()?.wb.clone();
@@ -325,6 +325,9 @@ fn run(s: &mut Session, p: &Json) -> Result<Json> {
     if reports.is_some() && solution.is_some() {
         out["reports"] = json!(names);
     }
+    if !skipped.is_empty() {
+        out["skippedReports"] = json!(skipped);
+    }
     Ok(out)
 }
 

crates/engine/src/cmd/solver/lp.rs

@@ -1,10 +1,14 @@
+use super::model::Goal;
 use super::problem::Evaluator;
-use super::{CellRef, EPS};
+use super::{CellRef, EPS, LP_WORK};
 
 // ---------------------------------------------------------------- evaluation
 
 type BoundRow = (Vec<(usize, f64)>, Sense, f64);
 
+/// A row of an LP: coefficients, sense and right side.
+pub(super) type LpRow = (Vec<f64>, Sense, f64);
+
 /// Objective (gradient, constant) and per-row (gradient, constant).
 pub(super) type Linearized = ((Vec<f64>, f64), Vec<(Vec<f64>, f64)>);
 
@@ -34,9 +38,10 @@ pub(super) enum Outcome {
     Limit,
 }
 
-struct Map {
+/// A column of the simplex tableau standing for a variable: `x_var = offset + Σ sign·y`.
-    var: usize,
+pub(super) struct Map {
-    sign: f64,
+    pub(super) var: usize,
+    pub(super) sign: f64,
 }
 
 struct VariableMap {
@@ -45,6 +50,23 @@ struct VariableMap {
     bound_rows: Vec<BoundRow>,
 }
 
+/// The final tableau of an optimal LP, for the Sensitivity report.
+pub(super) struct FinalTableau {
+    /// Rows `B⁻¹A | B⁻¹b`.
+    pub(super) t: Vec<Vec<f64>>,
+    pub(super) basis: Vec<usize>,
+    /// The phase-2 reduced costs.
+    pub(super) obj: Vec<f64>,
+    pub(super) width: usize,
+    /// The first artificial column.
+    pub(super) first_art: usize,
+    pub(super) cols: Vec<Map>,
+    /// Per input row: the column that started as its unit column, and whether the row was
+    /// negated to make its right side non-negative.
+    pub(super) ident: Vec<usize>,
+    pub(super) flipped: Vec<bool>,
+}
+
 struct Tableau {
     rows: Vec<Vec<f64>>,
     basis: Vec<usize>,
@@ -82,7 +104,9 @@ fn map_variables(n: usize, lo: &[f64], hi: &[f64]) -> VariableMap {
     VariableMap { columns, offset, bound_rows }
 }
 
-fn normalize_rows(rows: &[(Vec<f64>, Sense, f64)], bound_rows: Vec<BoundRow>, columns: &[Map], offset: &[f64]) -> Vec<(Vec<f64>, Sense, f64)> {
+/// The rows over the tableau's columns, each with a non-negative right side, and per row
+/// whether it was negated to get there.
+fn normalize_rows(rows: &[(Vec<f64>, Sense, f64)], bound_rows: Vec<BoundRow>, columns: &[Map], offset: &[f64]) -> (Vec<LpRow>, Vec<bool>) {
     let width = columns.len();
     let mut transformed = Vec::with_capacity(rows.len() + bound_rows.len());
     for (coefficients, sense, rhs) in rows {
@@ -105,7 +129,8 @@ fn normalize_rows(rows: &[(Vec<f64>, Sense, f64)], bound_rows: Vec<BoundRow>, co
         transformed.push((row, sense, rhs));
     }
     // Rows with a negative right side are negated so every right side is >= 0.
-    transformed
+    let flipped = transformed.iter().map(|row| row.2 < 0.0).collect();
+    let normalized = transformed
         .into_iter()
         .map(|(row, sense, rhs)| {
             if rhs < 0.0 {
@@ -119,7 +144,8 @@ fn normalize_rows(rows: &[(Vec<f64>, Sense, f64)], bound_rows: Vec<BoundRow>, co
                 (row, sense, rhs)
             }
         })
-        .collect()
+        .collect();
+    (normalized, flipped)
 }
 
 fn build_tableau(normalized: &[(Vec<f64>, Sense, f64)], decision_cols: usize, slack_cols: usize, artificial_cols: usize) -> Tableau {
@@ -388,20 +414,33 @@ fn simplex(t: &mut [Vec<f64>], basis: &mut [usize], obj: &mut [f64], allowed: &[
 
 /// Minimises `c·x` subject to `rows` (coefficients, sense, rhs) and `lo ≤ x ≤ hi`.
 pub(super) fn solve_lp(c: &[f64], rows: &[(Vec<f64>, Sense, f64)], lo: &[f64], hi: &[f64], budget: &mut u64) -> Outcome {
+    solve_lp_tableau(c, rows, lo, hi, budget).0
+}
+
+/// [`solve_lp`], with the final tableau when it is optimal.
+pub(super) fn solve_lp_tableau(
+    c: &[f64],
+    rows: &[(Vec<f64>, Sense, f64)],
+    lo: &[f64],
+    hi: &[f64],
+    budget: &mut u64,
+) -> (Outcome, Option<FinalTableau>) {
     if lo.iter().zip(hi).any(|(l, h)| l > h) {
-        return Outcome::Infeasible;
+        return (Outcome::Infeasible, None);
     }
     // x_j = offset + Σ sign·y_k with y ≥ 0.
     let mapping = map_variables(c.len(), lo, hi);
     let VariableMap { columns: cols, offset, bound_rows: extra } = mapping;
     let ny = cols.len();
     // Rows over y: (coefficients, sense, rhs).
-    let normalized = normalize_rows(rows, extra, &cols, &offset);
+    let (normalized, flipped) = normalize_rows(rows, extra, &cols, &offset);
     let n_slack = normalized.iter().filter(|row| row.1 != Sense::Eq).count();
     let n_art = normalized.iter().filter(|row| row.1 != Sense::Le).count();
     let mut tableau = build_tableau(&normalized, ny, n_slack, n_art);
+    // Each row's starting basic column is its unit column.
+    let ident = tableau.basis.clone();
     if let Err(outcome) = phase_one(&mut tableau, &normalized, budget) {
-        return outcome;
+        return (outcome, None);
     }
     let Tableau { rows: mut t, mut basis, width, decision_cols: ny, slack_cols: n_slack, .. } = tableau;
     let is_art = |column: usize| column >= ny + n_slack && column < width;
@@ -423,8 +462,8 @@ pub(super) fn solve_lp(c: &[f64], rows: &[(Vec<f64>, Sense, f64)], lo: &[f64], h
     }
     let allowed: Vec<bool> = (0..width).map(|j| !is_art(j)).collect();
     match simplex(&mut t, &mut basis, &mut obj, &allowed, budget) {
-        Step::Limit => return Outcome::Limit,
+        Step::Limit => return (Outcome::Limit, None),
-        Step::Unbounded => return Outcome::Unbounded,
+        Step::Unbounded => return (Outcome::Unbounded, None),
         Step::Optimal => {}
     }
     let mut y = vec![0.0; width];
@@ -442,7 +481,43 @@ pub(super) fn solve_lp(c: &[f64], rows: &[(Vec<f64>, Sense, f64)], lo: &[f64], h
         }
     }
     let val: f64 = c.iter().zip(&x).map(|(a, b)| a * b).sum();
-    Outcome::Solved(x, val, true)
+    let tableau = FinalTableau { t, basis, obj, width, first_art: ny + n_slack, cols, ident, flipped };
+    (Outcome::Solved(x, val, true), Some(tableau))
+}
+
+/// The LP to minimise for `goal`, from the linearised model: the costs, and the rows
+/// (coefficients, sense, right side) in the order of `model_rows`.
+pub(super) fn linear_program(goal: Goal, cf: &[f64], f_const: f64, lin: &[(Vec<f64>, f64)], model_rows: &[Row]) -> (Vec<f64>, Vec<LpRow>) {
+    let mut rows: Vec<LpRow> = lin.iter().zip(model_rows).map(|((a, k), r)| (a.clone(), r.sense, -k)).collect();
+    let c = match goal {
+        Goal::Max => cf.iter().map(|v| -v).collect(),
+        Goal::Min => cf.to_vec(),
+        Goal::Value(v) => {
+            rows.push((cf.to_vec(), Sense::Eq, v - f_const));
+            vec![0.0; cf.len()]
+        }
+    };
+    (c, rows)
+}
+
+/// The linear program Simplex LP solves: the model linearised at `x0`, its objective
+/// coefficients, the costs and rows for the goal, and the work limit for its size.
+pub(super) struct SimplexProgram {
+    pub(super) cf: Vec<f64>,
+    pub(super) c: Vec<f64>,
+    pub(super) rows: Vec<LpRow>,
+    pub(super) work: u64,
+}
+
+/// Builds the [`SimplexProgram`] at `x0`. The solve and the Sensitivity report both build it
+/// here from the start values, so the report's solve repeats the solve's pivots and returns the
+/// same vertex: linearising at another point changes the coefficients in their last bits, and
+/// that can break a tie between pivots the other way.
+pub(super) fn simplex_program(ev: &mut Evaluator, goal: Goal, x0: &[f64]) -> std::result::Result<SimplexProgram, &'static str> {
+    let ((cf, f_const), lin) = linearize(ev, x0)?;
+    let (c, rows) = linear_program(goal, &cf, f_const, &lin, &ev.rows);
+    let size = (rows.len() as u64 + 2 * x0.len() as u64 + 1) * (4 * x0.len() as u64 + 2 * rows.len() as u64 + 1);
+    Ok(SimplexProgram { cf, c, rows, work: LP_WORK.max(size) })
 }
 
 /// Linearises the model at `x0` (n + 1 evaluations) and checks it at one more point. Returns the

crates/engine/src/cmd/solver/model.rs

@@ -305,6 +305,32 @@ pub(super) fn parse_model(s: &Session, p: &Json) -> Result<Model> {
     Ok(Model { sheet, objective, goal, vars, vars_text, cons, non_negative, method, evo: evo_options(p), max_evals, gap })
 }
 
+/// Whether a `data.solver` result status is a solution the Sensitivity report can describe: an
+/// optimum (Simplex LP) or a converged local optimum (GRG Nonlinear). A run that stopped at its
+/// evaluation limit, or an Evolutionary run, has none.
+pub fn sensitivity_status(status: &str) -> bool {
+    matches!(status, "optimal" | "converged")
+}
+
+/// Whether `data.solver` params (such as the `model` of its result) describe a model that can
+/// have a Sensitivity report: Simplex LP or GRG Nonlinear, a Max or Min objective and no
+/// integer or binary constraints.
+pub fn sensitivity_applies(p: &Json) -> bool {
+    let method = str_param(p, "method").unwrap_or("grg").to_ascii_lowercase().replace([' ', '_', '-'], "");
+    let goal = str_param(p, "goal").or_else(|| str_param(p, "to")).unwrap_or("max").to_ascii_lowercase();
+    let integer = p.get("constraints").and_then(Json::as_array).is_some_and(|list| {
+        list.iter().any(|c| {
+            let rel = match c {
+                Json::String(t) => split_constraint(t).map(|(_, rel, _)| rel),
+                Json::Object(_) => str_param(c, "op").or_else(|| str_param(c, "relation")).and_then(Rel::parse),
+                _ => None,
+            };
+            matches!(rel, Some(Rel::Int | Rel::Bin))
+        })
+    });
+    !integer && !matches!(method.as_str(), "evolutionary" | "evo" | "genetic") && !matches!(goal.as_str(), "value" | "valueof" | "value of")
+}
+
 // ---------------------------------------------------------------- saved model
 
 /// `"B2:B4,D1"` as absolute references, each qualified by `sheet_name` when given.

crates/engine/src/cmd/solver/reports.rs

@@ -1,9 +1,10 @@
 use super::lp::Outcome;
 use super::model::{Constraint, Goal, Method, Model, Rel, Rhs, parse_model, saved_model};
 use super::problem::{Evaluator, Problem, as_num, build_problem};
+use super::sensitivity::{sensitivity_refusal, sensitivity_report};
 use super::{
-    Cell, CellRef, Ctx, DocState, EngineError, Json, PRECISION, REPORT_CMD, Result, Session, Sheet, Trial, Value, bad, edit, json, optimise,
+    Cell, CellRef, Ctx, DocState, EngineError, Json, PRECISION, REPORT_CMD, Result, Session, Sheet, Trial, Value, bad, bool_param, edit, json,
-    str_param,
+    optimise, str_param,
 };
 
 // ---------------------------------------------------------------- reports
@@ -11,6 +12,7 @@ use super::{
 #[derive(Clone, Copy, Debug, PartialEq)]
 pub(super) enum ReportKind {
     Answer,
+    Sensitivity,
     Limits,
 }
 
@@ -18,6 +20,7 @@ impl ReportKind {
     fn title(self) -> &'static str {
         match self {
             ReportKind::Answer => "Answer",
+            ReportKind::Sensitivity => "Sensitivity",
             ReportKind::Limits => "Limits",
         }
     }
@@ -26,13 +29,14 @@ impl ReportKind {
 /// The `reports` list (`None` when absent).
 pub(super) fn report_kinds(cmd: &'static str, p: &Json) -> Result<Option<Vec<ReportKind>>> {
     let Some(list) = p.get("reports") else { return Ok(None) };
-    let list = list.as_array().ok_or_else(|| bad(cmd, "`reports` must be an array of answer and limits"))?;
+    let list = list.as_array().ok_or_else(|| bad(cmd, "`reports` must be an array of answer, sensitivity and limits"))?;
     let mut kinds = Vec::new();
     for k in list {
         let kind = match k.as_str().map(|t| t.trim().to_ascii_lowercase()).as_deref() {
             Some("answer") => ReportKind::Answer,
+            Some("sensitivity") => ReportKind::Sensitivity,
             Some("limits") => ReportKind::Limits,
-            _ => return Err(bad(cmd, format!("unknown report {k} (answer or limits)"))),
+            _ => return Err(bad(cmd, format!("unknown report {k} (answer, sensitivity or limits)"))),
         };
         if !kinds.contains(&kind) {
             kinds.push(kind);
@@ -42,14 +46,14 @@ pub(super) fn report_kinds(cmd: &'static str, p: &Json) -> Result<Option<Vec<Rep
 }
 
 /// A report's cells: (row, column, value, bold).
-type Lines = Vec<(u32, u32, Value, bool)>;
+pub(super) type Lines = Vec<(u32, u32, Value, bool)>;
 
-fn abs_ref(c: CellRef) -> String {
+pub(super) fn abs_ref(c: CellRef) -> String {
     format!("${}${}", gridcraft_core::col_to_letters(c.col), u64::from(c.row) + 1)
 }
 
 /// A cell's name in the reports: the nearest text to its left and the nearest text above it.
-fn cell_name(sh: &Sheet, c: CellRef) -> String {
+pub(super) fn cell_name(sh: &Sheet, c: CellRef) -> String {
     let text = |c: CellRef| match sh.value(c) {
         Value::Text(t) if !t.trim().is_empty() => Some(t.to_string()),
         _ => None,
@@ -68,7 +72,7 @@ fn method_name(m: Method) -> &'static str {
 }
 
 /// Header labels in bold from column 1, then the rows below them; returns the next free row.
-fn table(out: &mut Lines, row: u32, headers: &[&str], rows: Vec<Vec<Value>>) -> u32 {
+pub(super) fn table(out: &mut Lines, row: u32, headers: &[&str], rows: Vec<Vec<Value>>) -> u32 {
     for (k, h) in headers.iter().enumerate() {
         out.push((row, 1 + k as u32, Value::text(*h), true));
     }
@@ -92,7 +96,7 @@ pub(super) struct ReportInput<'a> {
     pub(super) budget: u64,
 }
 
-fn title(out: &mut Lines, kind: ReportKind, sheet: &str) -> u32 {
+pub(super) fn title(out: &mut Lines, kind: ReportKind, sheet: &str) -> u32 {
     out.push((0, 0, Value::text(format!("GridCraft {} Report", kind.title())), true));
     out.push((1, 0, Value::text(format!("Worksheet: {sheet}")), false));
     3
@@ -255,18 +259,67 @@ fn limits_report(d: &DocState, sh: &Sheet, inp: &ReportInput) -> Result<Lines> {
 }
 
 /// Builds every report before the solution or saved model can change.
-pub(super) fn prepare_reports(d: &DocState, kinds: &[ReportKind], inp: &ReportInput) -> Result<Vec<(String, Lines)>> {
+/// What happens to a Sensitivity report the model allows but the solution doesn't: `Refuse`
+/// fails the call before anything changes; `Skip` leaves the report out, with its reason, and
+/// writes the others. `optimum` is false when the solve stopped before an optimum, so the
+/// report isn't attempted.
+#[derive(Clone, Copy, Debug, PartialEq)]
+pub(super) enum NotOptimal {
+    Refuse,
+    Skip { optimum: bool },
+}
+
+/// The reports left out under [`NotOptimal::Skip`], as `{report, reason}`.
+pub(super) type Skipped = Vec<Json>;
+
+pub(super) fn prepare_reports(
+    d: &DocState,
+    kinds: &[ReportKind],
+    inp: &ReportInput,
+    not_optimal: NotOptimal,
+) -> Result<(Vec<(String, Lines)>, Skipped)> {
     if kinds.is_empty() {
-        return Ok(Vec::new());
+        return Ok((Vec::new(), Vec::new()));
     }
-    if d.wb.protected_structure {
+    // A request that must write a sheet is refused before any report is built (a Limits report
-        return Err(EngineError::Other("The workbook is protected and cannot be changed.".into()));
+    // runs two optimisations per variable); only a Sensitivity report that may be skipped waits
+    // for the check after the loop.
+    let must_write = kinds.iter().any(|k| *k != ReportKind::Sensitivity || not_optimal == NotOptimal::Refuse);
+    if must_write && d.wb.protected_structure {
+        return Err(protected());
+    }
+    // A model that has no Sensitivity report (integer, Evolutionary, Value Of) is refused in
+    // both modes before any report is built: a Limits report listed first would otherwise run
+    // its two optimisations per variable for an error. Only the solution's own refusal (not
+    // optimal) waits for the loop, where `Skip` may leave the report out.
+    if kinds.contains(&ReportKind::Sensitivity)
+        && let Some(why) = sensitivity_refusal(inp.model, inp.prob)
+    {
+        return Err(EngineError::Other(why.into()));
     }
     let sh = d.wb.sheet(inp.model.sheet).ok_or(EngineError::NoDocument)?;
     let mut built = Vec::new();
+    let mut skipped = Vec::new();
     for k in kinds {
         let lines = match k {
             ReportKind::Answer => answer_report(d, sh, inp),
+            ReportKind::Sensitivity => match not_optimal {
+                NotOptimal::Skip { optimum } => {
+                    let got = if optimum {
+                        sensitivity_report(d, sh, inp).map_err(|e| e.to_string())
+                    } else {
+                        Err("Solver did not reach an optimal solution, so the Sensitivity report was not written.".to_string())
+                    };
+                    match got {
+                        Ok(lines) => lines,
+                        Err(reason) => {
+                            skipped.push(json!({"report": "sensitivity", "reason": reason}));
+                            continue;
+                        }
+                    }
+                }
+                NotOptimal::Refuse => sensitivity_report(d, sh, inp)?,
+            },
             ReportKind::Limits => limits_report(d, sh, inp)?,
         };
         let name = (1..=10_000u32)
@@ -275,7 +328,16 @@ pub(super) fn prepare_reports(d: &DocState, kinds: &[ReportKind], inp: &ReportIn
             .ok_or_else(|| EngineError::Other("Too many report sheets.".into()))?;
         built.push((name, lines));
     }
-    Ok(built)
+    // Only a report sheet changes the workbook's structure, so a request whose reports were all
+    // skipped still succeeds in a protected workbook.
+    if !built.is_empty() && d.wb.protected_structure {
+        return Err(protected());
+    }
+    Ok((built, skipped))
+}
+
+fn protected() -> EngineError {
+    EngineError::Other("The workbook is protected and cannot be changed.".into())
 }
 
 /// Inserts prepared reports before the model's sheet, which stays active.
@@ -292,7 +354,7 @@ pub(super) fn write_reports(cx: &mut Ctx, at: usize, built: Vec<(String, Lines)>
             let style = if b { bold } else { gridcraft_model::StyleId::DEFAULT };
             sheet.set_cell(CellRef::new(r, c), Cell { value: v, formula: None, style });
         }
-        for (c, w) in [(0u32, 14.0f32), (1, 110.0), (2, 140.0), (3, 120.0), (4, 120.0), (5, 120.0), (6, 120.0), (7, 120.0)] {
+        for (c, w) in [(0u32, 14.0f32), (1, 110.0), (2, 140.0), (3, 120.0), (4, 120.0), (5, 155.0), (6, 140.0), (7, 140.0)] {
             sheet.cols.insert(c, gridcraft_model::LineInfo { size: Some(w), ..Default::default() });
         }
         cx.wb.sheets.insert(at + i, std::sync::Arc::new(sheet));
@@ -321,9 +383,9 @@ pub(super) fn insert_report_selections(s: &mut Session, at: usize, count: usize)
 
 /// `data.solverReport`: reports for a solution of the sheet's saved model.
 pub(super) fn report(s: &mut Session, p: &Json) -> Result<Json> {
-    let kinds = report_kinds(REPORT_CMD, p)?.ok_or_else(|| bad(REPORT_CMD, "missing `reports` (answer, limits)"))?;
+    let kinds = report_kinds(REPORT_CMD, p)?.ok_or_else(|| bad(REPORT_CMD, "missing `reports` (answer, sensitivity, limits)"))?;
     if kinds.is_empty() {
-        return Err(bad(REPORT_CMD, "`reports` is empty (answer, limits)"));
+        return Err(bad(REPORT_CMD, "`reports` is empty (answer, sensitivity, limits)"));
     }
     let d = s.doc()?;
     // The given model as is, else the saved model with any params given here over it.
@@ -365,8 +427,14 @@ pub(super) fn report(s: &mut Session, p: &Json) -> Result<Json> {
     let budget = model.max_evals;
     let message = str_param(p, "message").unwrap_or("").chars().take(500).collect();
     let inp = ReportInput { model: &model, prob: &prob, original, solution, message, budget };
-    let built = prepare_reports(s.doc()?, &kinds, &inp)?;
+    let not_optimal = if bool_param(p, "skipNotOptimal").unwrap_or(false) { NotOptimal::Skip { optimum: true } } else { NotOptimal::Refuse };
-    let sheets = edit(s, |cx| Ok(write_reports(cx, model.sheet, built)))?;
+    let (built, skipped) = prepare_reports(s.doc()?, &kinds, &inp, not_optimal)?;
+    // Nothing left to write (every report skipped): no edit, so no Undo step.
+    let sheets = if built.is_empty() { Vec::new() } else { edit(s, |cx| Ok(write_reports(cx, model.sheet, built)))? };
     insert_report_selections(s, model.sheet, sheets.len());
-    Ok(json!({"sheets": sheets}))
+    let mut out = json!({"sheets": sheets});
+    if !skipped.is_empty() {
+        out["skippedReports"] = json!(skipped);
+    }
+    Ok(out)
 }

crates/engine/src/cmd/solver/sensitivity.rs

@@ -0,0 +1,614 @@
+//! The Sensitivity report. Simplex LP: each variable's reduced cost and the range of its
+//! objective coefficient over which the basis stays optimal, and each constraint's shadow price
+//! and the range of its right side over which the shadow price holds, from the final tableau of
+//! the LP linearised at the solution. GRG Nonlinear: the reduced gradients and the Lagrange
+//! multipliers of the binding constraints at the solution, from difference gradients.
+
+use super::lp::{FinalTableau, LpRow, Outcome, Row, Sense, SimplexProgram, simplex_program, solve_lp_tableau};
+use super::model::{Goal, Method, Model};
+use super::problem::{Evaluator, Point, Problem, as_num};
+use super::reports::{Lines, ReportInput, ReportKind, abs_ref, cell_name, table, title};
+use super::{CellRef, DocState, EPS, EngineError, Result, Sheet, Trial, Value};
+
+/// Excel prints an unlimited allowable change as 1E+30.
+const NO_LIMIT: f64 = 1e30;
+
+/// Most arithmetic the Lagrange multipliers of a nonlinear model may take (about a second): a
+/// pass over k binding rows and f free variables costs about 3·k²·f (the normal equations and
+/// each multiplier's noise), and each pass drops at most one row, so a model at the constraint
+/// limit could otherwise take minutes.
+#[cfg(not(test))]
+const MULTIPLIER_WORK: u64 = 1_000_000_000;
+/// A low limit in tests, so the refusal is reachable with a small model.
+#[cfg(test)]
+const MULTIPLIER_WORK: u64 = 20_000;
+
+/// Why a model has no Sensitivity report (`None` when it has one).
+pub(super) fn sensitivity_refusal(m: &Model, prob: &Problem) -> Option<&'static str> {
+    if prob.integer.iter().any(|b| *b) {
+        Some("The Sensitivity report doesn't apply to models with integer or binary constraints.")
+    } else if m.method == Method::Evolutionary {
+        Some("The Sensitivity report is available for Simplex LP and GRG Nonlinear.")
+    } else if matches!(m.goal, Goal::Value(_)) {
+        Some("The Sensitivity report needs an objective to maximise or minimise.")
+    } else {
+        None
+    }
+}
+
+/// A sensitivity figure: values within `noise` of zero (rounding) print as 0, unlimited ones
+/// as 1E+30. The noise level is local to the figure: it comes from the terms the figure is a
+/// sum of (the rounding error of a sum scales with the sum of its terms' sizes), never from
+/// the model's largest coefficient or gradient component, so a small price or gradient keeps
+/// its value beside a large unrelated one, and never from the objective's value, so adding a
+/// constant to the objective changes nothing the report shows.
+fn sens_num(v: f64, noise: f64) -> Value {
+    if !v.is_finite() || v >= NO_LIMIT {
+        Value::number(NO_LIMIT)
+    } else if v.abs() <= noise {
+        Value::number(0.0)
+    } else {
+        Value::number(v)
+    }
+}
+
+fn num(v: f64) -> Value {
+    if v.is_finite() { Value::number(v) } else { Value::Error(gridcraft_core::CellError::Value) }
+}
+
+/// The smallest ratio `num / den` over the `den` above zero (the step to the first bound).
+fn min_ratio(pairs: impl Iterator<Item = (f64, f64)>) -> f64 {
+    pairs.filter(|(_, den)| *den > EPS).map(|(num, den)| num.max(0.0) / den).fold(f64::INFINITY, f64::min)
+}
+
+/// Solves the square system `a·x = b` by Gaussian elimination with partial pivoting; unknowns
+/// of a singular system that have no pivot are 0.
+fn solve_dense(mut a: Vec<Vec<f64>>, mut b: Vec<f64>) -> Vec<f64> {
+    let n = b.len();
+    let scale = a.iter().flatten().fold(0.0f64, |m, v| m.max(v.abs())).max(f64::MIN_POSITIVE);
+    let mut pivots: Vec<Option<usize>> = vec![None; n];
+    let mut row = 0;
+    for col in 0..n {
+        let Some(p) = (row..n).max_by(|i, k| {
+            let v = |r: &usize| a.get(*r).and_then(|x| x.get(col)).map_or(0.0, |v| v.abs());
+            v(i).total_cmp(&v(k))
+        }) else {
+            break;
+        };
+        if a.get(p).and_then(|x| x.get(col)).is_none_or(|v| v.abs() <= 1e-12 * scale) {
+            continue;
+        }
+        a.swap(row, p);
+        b.swap(row, p);
+        let (Some(prow), Some(pb)) = (a.get(row).cloned(), b.get(row).copied()) else { break };
+        let pv = prow.get(col).copied().unwrap_or(1.0);
+        for (i, (r, bi)) in a.iter_mut().zip(b.iter_mut()).enumerate() {
+            if i == row {
+                continue;
+            }
+            let f = r.get(col).copied().unwrap_or(0.0) / pv;
+            if f != 0.0 {
+                for (v, p) in r.iter_mut().zip(&prow) {
+                    *v -= f * p;
+                }
+                *bi -= f * pb;
+            }
+        }
+        if let Some(e) = pivots.get_mut(col) {
+            *e = Some(row);
+        }
+        row += 1;
+    }
+    pivots
+        .iter()
+        .enumerate()
+        .map(|(col, p)| match p {
+            Some(r) => b.get(*r).copied().unwrap_or(0.0) / a.get(*r).and_then(|x| x.get(col)).copied().unwrap_or(1.0),
+            None => 0.0,
+        })
+        .collect()
+}
+
+fn out_of_evaluations() -> EngineError {
+    EngineError::Other("Solver ran out of evaluations.".into())
+}
+
+fn not_optimal() -> EngineError {
+    EngineError::Other("The Sensitivity report needs an optimal solution: these values are not optimal for the model. Solve it first.".into())
+}
+
+/// What both report bodies share: the model, the solution, its objective value and the rows'
+/// final values and right sides.
+struct At<'a> {
+    inp: &'a ReportInput<'a>,
+    sh: &'a Sheet,
+    point: Point,
+    finals: Vec<(f64, f64)>,
+}
+
+impl At<'_> {
+    fn var_head(&self, j: usize, c: CellRef) -> Vec<Value> {
+        vec![Value::text(abs_ref(c)), Value::text(cell_name(self.sh, c)), num(self.inp.solution.get(j).copied().unwrap_or(f64::NAN))]
+    }
+    fn row_head(&self, i: usize, r: &Row) -> Vec<Value> {
+        let fin = self.finals.get(i).map_or(f64::NAN, |f| f.0);
+        vec![Value::text(abs_ref(r.lhs)), Value::text(cell_name(self.sh, r.lhs)), num(fin)]
+    }
+}
+
+/// The variable rows, the constraint rows and their headers.
+type Body = (Vec<Vec<Value>>, Vec<Vec<Value>>, Vec<&'static str>, Vec<&'static str>);
+
+pub(super) fn sensitivity_report(d: &DocState, sh: &Sheet, inp: &ReportInput) -> Result<Lines> {
+    let m = inp.model;
+    if let Some(why) = sensitivity_refusal(m, inp.prob) {
+        return Err(EngineError::Other(why.into()));
+    }
+    let n = m.vars.len();
+    let x = &inp.solution;
+    let mut ev = Evaluator {
+        trial: Trial::new(d),
+        sheet: m.sheet,
+        vars: m.vars.clone(),
+        objective: m.objective,
+        rows: inp.prob.rows.clone(),
+        evals: 0,
+        // Room for the gradients (GRG) or two linearisations (Simplex LP) at the least.
+        budget: inp.budget.max(3 * n as u64 + 4),
+        target: None,
+    };
+    let point = ev.eval(x).ok_or_else(out_of_evaluations)?;
+    // The rows' final values and right sides at the solution.
+    let mut cells = Vec::new();
+    for r in &inp.prob.rows {
+        cells.push(r.lhs);
+        cells.extend(r.rhs);
+    }
+    let mut trial = Trial::new(d);
+    let sets: Vec<(usize, CellRef, Value)> = m.vars.iter().zip(x).map(|(c, v)| (m.sheet, *c, Value::number(*v))).collect();
+    let reads: Vec<(usize, CellRef)> = cells.iter().map(|c| (m.sheet, *c)).collect();
+    let mut vals = trial.run(&sets, &reads).into_iter().map(|v| as_num(&v));
+    let finals: Vec<(f64, f64)> = inp
+        .prob
+        .rows
+        .iter()
+        .map(|r| {
+            let l = vals.next().unwrap_or(f64::NAN);
+            (l, if r.rhs.is_some() { vals.next().unwrap_or(f64::NAN) } else { r.rhs_const })
+        })
+        .collect();
+    let at = At { inp, sh, point, finals };
+    let (var_rows, con_rows, var_headers, con_headers) =
+        if m.method == Method::Simplex { simplex_body(&mut ev, &at)? } else { nonlinear_body(&mut ev, &at)? };
+    let mut out = Lines::new();
+    let r = title(&mut out, ReportKind::Sensitivity, &sh.name);
+    out.push((r, 0, Value::text("Variable Cells"), true));
+    let r = table(&mut out, r + 1, &var_headers, var_rows);
+    out.push((r, 0, Value::text("Constraints"), true));
+    table(&mut out, r + 1, &con_headers, con_rows);
+    Ok(out)
+}
+
+/// Simplex LP: reduced costs, shadow prices and their ranges from the final tableau.
+fn simplex_body(ev: &mut Evaluator, at: &At) -> Result<Body> {
+    let (inp, m) = (at.inp, at.inp.model);
+    // The program the solve built, from the same start values (data.solver's x0, which Solver
+    // Results passes back as `original`), so the solve below repeats its pivots. Start values
+    // the model can't be evaluated at (given only to the report) fall back to the solution.
+    let start: Vec<f64> = inp.original.iter().map(|v| if v.is_finite() { *v } else { 0.0 }).collect();
+    let built = match simplex_program(ev, m.goal, &start) {
+        Err("nonlinear" | "error") if start != inp.solution => simplex_program(ev, m.goal, &inp.solution),
+        other => other,
+    };
+    let SimplexProgram { cf, c, rows, mut work } = built.map_err(|e| {
+        if e == "budget" {
+            out_of_evaluations()
+        } else {
+            EngineError::Other("The model isn't linear in the variable cells, so Simplex LP can't report its sensitivity.".into())
+        }
+    })?;
+    let (outcome, tab) = solve_lp_tableau(&c, &rows, &inp.prob.lo, &inp.prob.hi, &mut work);
+    let (Outcome::Solved(vertex, best, _), Some(tab)) = (outcome, tab) else {
+        return Err(EngineError::Other("Solver could not find an optimal basis for the Sensitivity report.".into()));
+    };
+    // The report describes the optimum, so the values it shows must be one.
+    if !optimal_point(&c, &rows, (&inp.prob.lo, &inp.prob.hi), &inp.solution, best) {
+        return Err(not_optimal());
+    }
+    // The prices and ranges come from the tableau of `vertex`, so the values shown must be that
+    // vertex. Another optimum with the same objective (when the optimum is not unique) has other
+    // ranges, so it is refused too, although it is optimal.
+    let near = |a: f64, b: f64| (a - b).abs() <= super::PRECISION * (1.0 + a.abs().max(b.abs()));
+    if vertex.len() != inp.solution.len() || !vertex.iter().zip(&inp.solution).all(|(a, b)| near(*a, *b)) {
+        return Err(EngineError::Other(
+            "The Sensitivity report needs the values at the Simplex LP solution: these values are another optimal solution, whose ranges differ. Solve the model with Simplex LP first."
+                .into(),
+        ));
+    }
+    // Minimising `sign · objective`.
+    let sign = if m.goal == Goal::Max { -1.0 } else { 1.0 };
+    // A shadow price is the dual yᵢ = Σᵣ c_B(r)·B⁻¹(r, i): the row's unit column of the final
+    // tableau holds B⁻¹'s column i, and c_B(r) is the cost of the column basic in row r (a
+    // slack or artificial costs nothing). The rounding error of a sum scales with the sum of
+    // its terms' sizes, so a price below 1E-9 of Σᵣ |c_B(r)·B⁻¹(r, i)| is noise. The bound is
+    // the row's own: a small price keeps its value beside a large coefficient that plays no
+    // part in it, and a constant in the objective (no coefficient) changes no bound. Allowable
+    // changes to a coefficient or a right side are ratios the tableau's own pivot tolerance
+    // (`EPS`) already bounds.
+    let cost_size = |b: usize| tab.cols.get(b).map_or(0.0, |mp| cf.get(mp.var).copied().unwrap_or(0.0).abs());
+    let price_terms = |col: usize| (0..tab.t.len()).map(|r| tab.basis.get(r).map_or(0.0, |b| cost_size(*b)) * entry(&tab, r, col).abs()).sum::<f64>();
+    let (shadow, terms): (Vec<f64>, Vec<f64>) = (0..rows.len())
+        .map(|i| {
+            let col = tab.ident.get(i).copied().unwrap_or(usize::MAX);
+            (sign * flip(&tab, i) * -tab.obj.get(col).copied().unwrap_or(0.0), price_terms(col))
+        })
+        .unzip();
+    let noise = |terms: f64| 1e-9 * terms.max(f64::MIN_POSITIVE);
+    let con_rows = ev
+        .rows
+        .iter()
+        .enumerate()
+        .map(|(i, r)| {
+            // The basic values move along the row's unit column as its right side changes.
+            let col = tab.ident.get(i).copied().unwrap_or(usize::MAX);
+            let dir: Vec<(f64, f64)> =
+                (0..tab.t.len()).filter(|r| real_basic(&tab, *r)).map(|r| (entry(&tab, r, tab.width), flip(&tab, i) * entry(&tab, r, col))).collect();
+            // An artificial still basic after phase one stands for a row the others make
+            // redundant, and it is zero there. When this right side moves it, any change alone
+            // leaves no solution, so both ranges are 0.
+            let pinned = (0..tab.t.len()).any(|r| !real_basic(&tab, r) && entry(&tab, r, col).abs() > EPS);
+            let (increase, decrease) =
+                if pinned { (0.0, 0.0) } else { (min_ratio(dir.iter().map(|(v, d)| (*v, -d))), min_ratio(dir.iter().map(|(v, d)| (*v, *d)))) };
+            let b = at.finals.get(i).map_or(f64::NAN, |f| f.1);
+            let mut row = at.row_head(i, r);
+            let price = sens_num(shadow.get(i).copied().unwrap_or(0.0), noise(terms.get(i).copied().unwrap_or(0.0)));
+            row.extend([price, num(b), sens_num(increase, EPS), sens_num(decrease, EPS)]);
+            row
+        })
+        .collect();
+    let nonbasic: Vec<usize> = (0..tab.first_art).filter(|l| !tab.basis.contains(l)).collect();
+    let var_rows = m
+        .vars
+        .iter()
+        .enumerate()
+        .map(|(j, c)| {
+            let coef = cf.get(j).copied().unwrap_or(0.0);
+            let used: f64 = rows.iter().zip(&shadow).map(|((a, _, _), p)| a.get(j).copied().unwrap_or(0.0) * p).sum();
+            // The reduced cost cf_j − Σᵢ a_ij·yᵢ is a sum of cf_j and the prices' terms times
+            // a_ij: the prices' terms rather than the prices, since a price that is a small
+            // difference of large terms carries their rounding error into every reduced cost
+            // it is used in.
+            let used_terms: f64 = rows.iter().zip(&terms).map(|((a, _, _), t)| a.get(j).copied().unwrap_or(0.0).abs() * t).sum();
+            // A change Δ in the coefficient changes each nonbasic reduced cost by Δ·e.
+            let own = |col: usize| tab.cols.get(col).filter(|mp| mp.var == j).map_or(0.0, |mp| mp.sign);
+            let moves: Vec<(f64, f64)> = nonbasic
+                .iter()
+                .map(|&l| {
+                    let through_basis: f64 = (0..tab.t.len()).map(|r| tab.basis.get(r).map_or(0.0, |b| own(*b)) * entry(&tab, r, l)).sum();
+                    (tab.obj.get(l).copied().unwrap_or(0.0), sign * (own(l) - through_basis))
+                })
+                .collect();
+            let increase = min_ratio(moves.iter().map(|(d, e)| (*d, -e)));
+            let decrease = min_ratio(moves.iter().map(|(d, e)| (*d, *e)));
+            let mut row = at.var_head(j, *c);
+            row.extend([sens_num(coef - used, noise(coef.abs() + used_terms)), num(coef), sens_num(increase, EPS), sens_num(decrease, EPS)]);
+            row
+        })
+        .collect();
+    Ok((
+        var_rows,
+        con_rows,
+        vec!["Cell", "Name", "Final Value", "Reduced Cost", "Objective Coefficient", "Allowable Increase", "Allowable Decrease"],
+        vec!["Cell", "Name", "Final Value", "Shadow Price", "Constraint R.H. Side", "Allowable Increase", "Allowable Decrease"],
+    ))
+}
+
+/// Whether `x` satisfies the LP (rows and bounds) and reaches its optimum `best` (minimising
+/// `c·x`), within the constraint precision.
+fn optimal_point(c: &[f64], rows: &[LpRow], (lo, hi): (&[f64], &[f64]), x: &[f64], best: f64) -> bool {
+    let near = |v: f64, size: f64| v <= super::PRECISION * (1.0 + size);
+    let dot = |a: &[f64]| a.iter().zip(x).fold((0.0, 0.0), |(s, m), (p, q)| (s + p * q, m + (p * q).abs()));
+    let rows_hold = rows.iter().all(|(a, sense, b)| {
+        let (lhs, size) = dot(a);
+        let size = size + b.abs();
+        match sense {
+            Sense::Le => near(lhs - b, size),
+            Sense::Ge => near(b - lhs, size),
+            Sense::Eq => near((lhs - b).abs(), size),
+        }
+    });
+    let bounds_hold = x.iter().enumerate().all(|(j, v)| {
+        let (l, h) = (lo.get(j).copied().unwrap_or(f64::NEG_INFINITY), hi.get(j).copied().unwrap_or(f64::INFINITY));
+        v.is_finite() && near(l - v, v.abs()) && near(v - h, v.abs())
+    });
+    let (value, size) = dot(c);
+    // Two-sided: a point just outside a bound can beat the optimum by far more than the
+    // precision when its costs are large.
+    rows_hold && bounds_hold && near((value - best).abs(), size + best.abs())
+}
+
+fn entry(tab: &FinalTableau, r: usize, col: usize) -> f64 {
+    tab.t.get(r).and_then(|row| row.get(col)).copied().unwrap_or(0.0)
+}
+
+fn real_basic(tab: &FinalTableau, r: usize) -> bool {
+    tab.basis.get(r).is_some_and(|b| *b < tab.first_art)
+}
+
+/// -1 for a row the tableau negated, else 1.
+fn flip(tab: &FinalTableau, i: usize) -> f64 {
+    if tab.flipped.get(i).copied().unwrap_or(false) { -1.0 } else { 1.0 }
+}
+
+/// Central-difference gradients, inside the bounds, of the objective and of each row's
+/// `lhs − rhs` at the solution.
+fn gradients(ev: &mut Evaluator, inp: &ReportInput) -> Result<(Vec<f64>, Vec<Vec<f64>>)> {
+    let x = &inp.solution;
+    let n = x.len();
+    let mut df = vec![0.0; n];
+    let mut dg = vec![vec![0.0; n]; ev.rows.len()];
+    for j in 0..n {
+        let xj = x.get(j).copied().unwrap_or(0.0);
+        let (l, u) = (inp.prob.lo.get(j).copied().unwrap_or(f64::NEG_INFINITY), inp.prob.hi.get(j).copied().unwrap_or(f64::INFINITY));
+        let h = 1e-6 * xj.abs().max(1.0);
+        let (a, b) = ((xj - h).max(l), (xj + h).min(u));
+        if b <= a {
+            continue;
+        }
+        let mut point = |v: f64| {
+            let mut p = x.clone();
+            if let Some(e) = p.get_mut(j) {
+                *e = v;
+            }
+            ev.eval(&p).ok_or_else(out_of_evaluations)
+        };
+        let (pa, pb) = (point(a)?, point(b)?);
+        if let Some(e) = df.get_mut(j) {
+            *e = (pb.f - pa.f) / (b - a);
+        }
+        for (row, (ga, gb)) in dg.iter_mut().zip(pa.g.iter().zip(&pb.g)) {
+            if let Some(e) = row.get_mut(j) {
+                *e = (gb - ga) / (b - a);
+            }
+        }
+    }
+    Ok((df, dg))
+}
+
+/// The rounding noise in the Lagrange multiplier of row `i`, one of the `active` rows whose
+/// multipliers `lambda` solve the normal equations of ∇f = Σ λₖ ∇gₖ over the `free` variables.
+/// Row i's equation gives λᵢ = (Σⱼ gᵢⱼ·dfⱼ − Σ_{k≠i} λₖ·Σⱼ gᵢⱼ·gₖⱼ) / Σⱼ gᵢⱼ²; the difference
+/// gradients are good to about 1E-6 of each component, and the rounding error of a sum scales
+/// with the sum of its terms' sizes, so the noise is 1E-6 of the terms' sizes over the row's
+/// own squared norm. The bound is local to the row: another row enters only through its
+/// coupling term with this one, so a row whose gradient is orthogonal to the others has a bound
+/// from its own gradient alone, and a large component of ∇f that the row's gradient doesn't
+/// touch hides nothing. Nearly parallel active gradients make the normal equations
+/// ill-conditioned and amplify every multiplier's error beyond any bound local to a row; such
+/// multipliers are ill-determined whatever the bound, and they print as computed.
+fn multiplier_noise(i: usize, active: &[usize], lambda: &[f64], free: &[usize], (df, dg): (&[f64], &[Vec<f64>])) -> f64 {
+    let grad = |i: usize, j: usize| dg.get(i).and_then(|r| r.get(j)).copied().unwrap_or(0.0);
+    let own: f64 = free.iter().map(|&j| grad(i, j) * grad(i, j)).sum();
+    let objective: f64 = free.iter().map(|&j| (grad(i, j) * df.get(j).copied().unwrap_or(0.0)).abs()).sum();
+    let coupling: f64 = active
+        .iter()
+        .zip(lambda)
+        .filter(|(k, _)| **k != i)
+        .map(|(&k, l)| free.iter().map(|&j| (grad(i, j) * grad(k, j) * l).abs()).sum::<f64>())
+        .sum();
+    1e-6 * (objective + coupling) / own.max(f64::MIN_POSITIVE)
+}
+
+/// Least squares on the free variables: ∇f = Σ λᵢ ∇gᵢ over the `active` rows. An inequality's
+/// multiplier has the sign that tightening it costs (λ ≥ 0 for `lhs ≤ rhs` when maximising); a
+/// row whose multiplier has the other sign beyond its own noise ([`multiplier_noise`], the
+/// criterion the report rounds it with) isn't binding, so it leaves and the rest are solved
+/// again; the objective's value plays no part. Returns each row's multiplier and its noise (0
+/// and 0 for a row that doesn't bind). The work is bounded by `budget` (see
+/// [`MULTIPLIER_WORK`]); past it the report is refused.
+fn multipliers(
+    rows: &[Row],
+    goal: Goal,
+    mut active: Vec<usize>,
+    free: &[usize],
+    (df, dg): (&[f64], &[Vec<f64>]),
+    budget: u64,
+) -> Result<(Vec<f64>, Vec<f64>)> {
+    let grad = |i: usize, j: usize| dg.get(i).and_then(|r| r.get(j)).copied().unwrap_or(0.0);
+    let goal_sign = if goal == Goal::Max { 1.0 } else { -1.0 };
+    let mut lambda = vec![0.0; rows.len()];
+    let mut noise = vec![0.0; rows.len()];
+    let mut work = 0u64;
+    while !active.is_empty() {
+        let (k, f) = (active.len() as u64, (free.len() as u64).max(1));
+        work = work.saturating_add(3u64.saturating_mul(k).saturating_mul(k).saturating_mul(f));
+        if work > budget {
+            return Err(EngineError::Other(
+                "The Sensitivity report of this nonlinear model has too many binding constraints to work out in reasonable time.".into(),
+            ));
+        }
+        let normal: Vec<Vec<f64>> =
+            active.iter().map(|&i| active.iter().map(|&k| free.iter().map(|&j| grad(i, j) * grad(k, j)).sum()).collect()).collect();
+        let rhs: Vec<f64> = active.iter().map(|&i| free.iter().map(|&j| grad(i, j) * df.get(j).copied().unwrap_or(0.0)).sum()).collect();
+        let solved = solve_dense(normal, rhs);
+        let bounds: Vec<f64> = active.iter().map(|&i| multiplier_noise(i, &active, &solved, free, (df, dg))).collect();
+        let wrong = active
+            .iter()
+            .zip(&solved)
+            .zip(&bounds)
+            .enumerate()
+            .filter_map(|(k, ((&i, l), bound))| {
+                let want = match rows.get(i).map_or(Sense::Eq, |r| r.sense) {
+                    Sense::Le => goal_sign,
+                    Sense::Ge => -goal_sign,
+                    Sense::Eq => return None,
+                };
+                (want * l < -bound).then_some((k, want * l))
+            })
+            .min_by(|a, b| a.1.total_cmp(&b.1));
+        match wrong {
+            Some((k, _)) => {
+                active.remove(k);
+            }
+            None => {
+                for ((i, l), bound) in active.iter().zip(solved).zip(bounds) {
+                    if let (Some(e), Some(n)) = (lambda.get_mut(*i), noise.get_mut(*i)) {
+                        *e = l;
+                        *n = bound;
+                    }
+                }
+                break;
+            }
+        }
+    }
+    Ok((lambda, noise))
+}
+
+/// Whether the solution is a Karush-Kuhn-Tucker point of the model, from the gradients the
+/// report computed: every row and bound holds, and each reduced gradient (the objective's
+/// gradient less what the binding rows' multipliers explain) is near zero relative to `terms`,
+/// the sizes of the terms it is made of, or, for a variable at a bound, points only out of that
+/// bound in the objective's direction. A binding row whose multiplier had the wrong sign was
+/// dropped by [`multipliers`], so the part of the gradient it would have explained remains in
+/// the reduced gradients and the point is refused here. Every tolerance is local to the row or
+/// variable it tests: a large right side or gradient component elsewhere in the model excuses
+/// nothing here.
+fn kkt_point(at: &At, rows: &[Row], terms: &[f64], reduced: &[f64]) -> bool {
+    let (inp, x) = (at.inp, &at.inp.solution);
+    // Feasibility: each row within ten times the constraint precision of its own right side.
+    // The nonlinear method accepts a solution within that of the model's largest constant right
+    // side, which lets a small row go by its own size; a row over by more than its own
+    // precision isn't satisfied, whatever the other rows' sizes, so the report refuses it.
+    let rows_hold = rows.iter().enumerate().all(|(i, r)| {
+        let (g, b) = (at.point.g.get(i).copied().unwrap_or(f64::NAN), at.finals.get(i).map_or(f64::NAN, |f| f.1));
+        let over = match r.sense {
+            Sense::Le => g,
+            Sense::Ge => -g,
+            Sense::Eq => g.abs(),
+        };
+        over <= 10.0 * super::PRECISION * (1.0 + b.abs())
+    });
+    if !rows_hold || !at.point.f.is_finite() {
+        return false;
+    }
+    // Stationarity: each reduced gradient relative to the terms it is made of (the objective's
+    // component and the multipliers' parts, at least 1). At optima the nonlinear method solves,
+    // a reduced gradient is below about 1E-6 of those terms; 1E-3 leaves a wide margin for the
+    // method's stopping test and for rounding in the difference gradients of badly scaled
+    // models, and still refuses a point whose objective improves visibly along a direction the
+    // constraints allow. The bound is the variable's own: with a bound from the gradient's
+    // largest component, max 1E9·x + y at (1, 0.5) passed, although raising y gains 1 per unit
+    // against a tolerance of 1E6.
+    // Moving up a variable improves a maximum by +r and a minimum by −r per unit.
+    let sign = if inp.model.goal == Goal::Max { 1.0 } else { -1.0 };
+    x.iter().enumerate().all(|(j, &v)| {
+        let (l, h) = (inp.prob.lo.get(j).copied().unwrap_or(f64::NEG_INFINITY), inp.prob.hi.get(j).copied().unwrap_or(f64::INFINITY));
+        let near = super::PRECISION * (1.0 + v.abs());
+        if !v.is_finite() || v < l - near || v > h + near {
+            return false;
+        }
+        // The same test for "at a bound" as the choice of free variables above.
+        let edge = 1e-6 * v.abs().max(1.0);
+        let (at_lo, at_hi) = (v <= l + edge, v >= h - edge);
+        let gain = reduced.get(j).map_or(f64::NAN, |r| sign * r);
+        let tol = 1e-3 * terms.get(j).copied().unwrap_or(0.0).max(1.0);
+        gain.is_finite() && (at_hi || gain <= tol) && (at_lo || gain >= -tol)
+    })
+}
+
+/// GRG Nonlinear: reduced gradients and the binding rows' Lagrange multipliers.
+fn nonlinear_body(ev: &mut Evaluator, at: &At) -> Result<Body> {
+    let (inp, m) = (at.inp, at.inp.model);
+    let x = &inp.solution;
+    let (df, dg) = gradients(ev, inp)?;
+    // Binding rows (within the nonlinear method's accuracy) and the variables strictly inside
+    // their bounds.
+    let active: Vec<usize> = ev
+        .rows
+        .iter()
+        .enumerate()
+        .filter(|(i, r)| {
+            let (g, b) = (at.point.g.get(*i).copied().unwrap_or(f64::NAN), at.finals.get(*i).map_or(0.0, |f| f.1));
+            r.sense == Sense::Eq || g.abs() <= 1e-4 * (1.0 + b.abs())
+        })
+        .map(|(i, _)| i)
+        .collect();
+    let free: Vec<usize> = (0..x.len())
+        .filter(|&j| {
+            let v = x.get(j).copied().unwrap_or(0.0);
+            let tol = 1e-6 * v.abs().max(1.0);
+            v > inp.prob.lo.get(j).copied().unwrap_or(f64::NEG_INFINITY) + tol && v < inp.prob.hi.get(j).copied().unwrap_or(f64::INFINITY) - tol
+        })
+        .collect();
+    let (lambda, lambda_noise) = multipliers(&ev.rows, m.goal, active, &free, (&df, &dg), MULTIPLIER_WORK)?;
+    let grad = |i: usize, j: usize| dg.get(i).and_then(|r| r.get(j)).copied().unwrap_or(0.0);
+    let reduced: Vec<f64> =
+        (0..x.len()).map(|j| df.get(j).copied().unwrap_or(0.0) - lambda.iter().enumerate().map(|(i, l)| l * grad(i, j)).sum::<f64>()).collect();
+    // A reduced gradient dfⱼ − Σᵢ λᵢ·gᵢⱼ is a sum of the gradient's own component and the
+    // multipliers' terms; the sizes of those terms bound its rounding error and scale the test
+    // of whether it is zero.
+    let terms: Vec<f64> = (0..x.len())
+        .map(|j| df.get(j).copied().unwrap_or(0.0).abs() + lambda.iter().enumerate().map(|(i, l)| (l * grad(i, j)).abs()).sum::<f64>())
+        .collect();
+    // The report describes the optimum, so the values it shows must be one.
+    if !kkt_point(at, &ev.rows, &terms, &reduced) {
+        return Err(not_optimal());
+    }
+    // Difference gradients are good to about 1E-6 of each component, so a reduced gradient's
+    // noise is 1E-6 of its terms' sizes, plus each multiplier's own noise
+    // ([`multiplier_noise`]) times |gᵢⱼ|. The bound is the variable's own: a small gradient
+    // keeps its value beside a large component of another variable.
+    let reduced_noise = |j: usize| {
+        let carried: f64 = lambda_noise.iter().enumerate().map(|(i, n)| if grad(i, j) == 0.0 { 0.0 } else { n * grad(i, j).abs() }).sum();
+        1e-6 * terms.get(j).copied().unwrap_or(0.0) + carried
+    };
+    let var_rows = m
+        .vars
+        .iter()
+        .enumerate()
+        .map(|(j, c)| {
+            let mut row = at.var_head(j, *c);
+            row.push(sens_num(reduced.get(j).copied().unwrap_or(f64::NAN), reduced_noise(j)));
+            row
+        })
+        .collect();
+    let con_rows = ev
+        .rows
+        .iter()
+        .enumerate()
+        .map(|(i, r)| {
+            let mut row = at.row_head(i, r);
+            row.push(sens_num(lambda.get(i).copied().unwrap_or(0.0), lambda_noise.get(i).copied().unwrap_or(0.0)));
+            row
+        })
+        .collect();
+    Ok((var_rows, con_rows, vec!["Cell", "Name", "Final Value", "Reduced Gradient"], vec!["Cell", "Name", "Final Value", "Lagrange Multiplier"]))
+}
+
+#[cfg(test)]
+mod tests {
+    use super::{Goal, MULTIPLIER_WORK, Row, Sense, multipliers};
+    use gridcraft_core::CellRef;
+
+    #[test]
+    fn multipliers_at_the_constraint_limit_stop_at_the_work_limit() {
+        // 1,000 binding `<=` rows over 200 free variables, each multiplier with the wrong sign:
+        // without a limit every pass would drop one row and solve again, about 4·10⁸ operations
+        // a pass for up to 1,000 passes. The limit (lowered in tests) refuses it at once.
+        let (k, f) = (1000usize, 200usize);
+        let rows: Vec<Row> = (0..k).map(|_| Row { lhs: CellRef::new(0, 0), sense: Sense::Le, rhs: None, rhs_const: 0.0 }).collect();
+        let dg: Vec<Vec<f64>> = (0..k).map(|i| (0..f).map(|j| if j == i % f { 1.0 } else { 0.0 }).collect()).collect();
+        let df = vec![-1.0; f];
+        let free: Vec<usize> = (0..f).collect();
+        let got = multipliers(&rows, Goal::Max, (0..k).collect(), &free, (&df, &dg), MULTIPLIER_WORK);
+        assert!(got.is_err());
+        // A small model within the limit still gets its multipliers: ∇f = (1, 2) on two
+        // binding rows with unit gradients.
+        let dg2 = vec![vec![1.0, 0.0], vec![0.0, 1.0]];
+        let (got, noise) = multipliers(&rows[..2], Goal::Max, vec![0, 1], &free[..2], (&[1.0, 2.0], &dg2), MULTIPLIER_WORK).unwrap();
+        assert_eq!(got, vec![1.0, 2.0]);
+        // Each multiplier's noise is 1E-6 of its own row's term (1 · 1 and 1 · 2 over a unit
+        // norm): the rows are orthogonal, so neither enters the other's bound.
+        assert_eq!(noise, vec![1e-6, 2e-6]);
+    }
+}

crates/engine/src/cmd/solver/tests.rs

@@ -283,6 +283,8 @@ fn refused_solver_reports_preserve_the_document_and_history() {
         (true, "simplex", json!(["A1 <= 4"]), json!(["answer"]), "workbook is protected"),
         (false, "simplex", json!(["A1 <= 4", "A1 int"]), json!(["answer", "limits"]), "integer or binary"),
         (false, "evolutionary", json!(["A1 <= 4"]), json!(["answer", "limits"]), "Simplex LP and GRG"),
+        (false, "simplex", json!(["A1 <= 4", "A1 int"]), json!(["answer", "sensitivity"]), "integer or binary"),
+        (false, "evolutionary", json!(["A1 <= 4"]), json!(["sensitivity"]), "Simplex LP and GRG"),
     ] {
         for (keep, save) in [(true, true), (false, true), (true, false), (false, false)] {
             let mut s = s();
@@ -400,6 +402,823 @@ fn solver_report_insertion_keeps_variable_edit_guards() {
     }
 }
 
+/// The numbers from column D on in the row below the one holding `text` on report `name`.
+fn report_row(s: &Session, name: &str, text: &str, cols: &[&str]) -> Vec<f64> {
+    let r = row_of(s, name, text) + 1;
+    cols.iter().map(|c| at(s, name, &format!("{c}{r}")).as_f64().unwrap()).collect()
+}
+
+fn close(got: &[f64], want: &[f64], tol: f64) -> bool {
+    got.len() == want.len() && got.iter().zip(want).all(|(g, w)| (g - w).abs() <= tol * (1.0 + w.abs()))
+}
+
+#[test]
+fn sensitivity_report_gives_the_wyndor_glass_shadow_prices_and_ranges() {
+    // Hillier & Lieberman's Wyndor Glass sensitivity analysis: shadow prices 0, 1.5 and 1
+    // for the three plants; plant 1's hours may fall by 2 and rise without limit, plants 2's
+    // and 3's by 6 either way; the doors' profit 3 may range over 0..7.5, the windows' 5 over
+    // 2 and up.
+    let mut s = s();
+    s.execute(
+        "range.setValues",
+        json!({"range": "A1", "values": [[null, "Doors", "Windows", "Used", "Available"], ["Batches", 0, 0, null, null], ["Plant 1", 1, 0, null, 4], ["Plant 2", 0, 2, null, 12], ["Plant 3", 3, 2, null, 18], ["Profit", 3, 5, null, null]]}),
+    )
+    .unwrap();
+    for r in 3..=6 {
+        set(&mut s, &[(&format!("D{r}"), &format!("=SUMPRODUCT(B{r}:C{r},$B$2:$C$2)"))]);
+    }
+    let solved = s
+        .execute(
+            "data.solver",
+            json!({"objective": "D6", "variables": "B2:C2", "constraints": ["D3:D5 <= E3:E5"], "method": "simplex", "reports": ["sensitivity"]}),
+        )
+        .unwrap();
+    assert_eq!(solved["reports"], json!(["Sensitivity Report 1"]), "{solved}");
+    let rep = "Sensitivity Report 1";
+    assert_eq!(at(&s, rep, "A1"), Value::text("GridCraft Sensitivity Report"));
+    let r = row_of(&s, rep, "Reduced Cost");
+    assert_eq!(at(&s, rep, &format!("F{}", r + 1)), Value::text("Objective Coefficient"));
+    let cols = ["D", "E", "F", "G", "H"];
+    // Final value, reduced cost, objective coefficient, allowable increase and decrease.
+    assert!(close(&report_row(&s, rep, "$B$2", &cols)[1..], &[0.0, 3.0, 4.5, 3.0], 1e-9));
+    assert!(close(&report_row(&s, rep, "$C$2", &cols), &[6.0, 0.0, 5.0, 1e30, 3.0], 1e-9));
+    // Final value, shadow price, right side, allowable increase and decrease.
+    assert!(close(&report_row(&s, rep, "$D$3", &cols), &[2.0, 0.0, 4.0, 1e30, 2.0], 1e-9));
+    assert!(close(&report_row(&s, rep, "$D$4", &cols), &[12.0, 1.5, 12.0, 6.0, 6.0], 1e-9));
+    assert!(close(&report_row(&s, rep, "$D$5", &cols), &[18.0, 1.0, 18.0, 6.0, 6.0], 1e-9));
+    // A minimisation with `>=` rows: min 2x + 3y, x + y ≥ 4, x + 3y ≥ 6 → (3, 1), cost 9.
+    // The duals 1.5 and 0.5 solve u + v = 2, u + 3v = 3; the basis holds for 2 ≤ b₁ ≤ 6,
+    // 4 ≤ b₂ ≤ 12, 1 ≤ c₁ ≤ 3 and 2 ≤ c₂ ≤ 6.
+    set(&mut s, &[("H1", "0"), ("H2", "0"), ("I1", "=2*H1+3*H2"), ("I2", "=H1+H2"), ("I3", "=H1+3*H2")]);
+    let r = s
+        .execute(
+            "data.solver",
+            json!({"objective": "I1", "goal": "min", "variables": "H1:H2", "constraints": ["I2 >= 4", "I3 >= 6"], "method": "simplex", "reports": ["sensitivity"]}),
+        )
+        .unwrap();
+    let rep = r["reports"][0].as_str().unwrap().to_string();
+    assert!(close(&report_row(&s, &rep, "$H$1", &cols), &[3.0, 0.0, 2.0, 1.0, 1.0], 1e-9), "{r}");
+    assert!(close(&report_row(&s, &rep, "$H$2", &cols), &[1.0, 0.0, 3.0, 3.0, 1.0], 1e-9));
+    assert!(close(&report_row(&s, &rep, "$I$2", &cols), &[4.0, 1.5, 4.0, 2.0, 2.0], 1e-9));
+    assert!(close(&report_row(&s, &rep, "$I$3", &cols), &[6.0, 0.5, 6.0, 6.0, 2.0], 1e-9));
+    // A variable held at zero has a reduced cost: min 2x + 3y with x + y ≥ 4 alone uses
+    // only x; y would add 3 − 2 = 1 per unit.
+    let r = s
+        .execute(
+            "data.solver",
+            json!({"objective": "I1", "goal": "min", "variables": "H1:H2", "constraints": ["I2 >= 4"], "method": "simplex", "reports": ["sensitivity"]}),
+        )
+        .unwrap();
+    let rep = r["reports"][0].as_str().unwrap().to_string();
+    assert!(close(&report_row(&s, &rep, "$H$2", &cols), &[0.0, 1.0, 3.0, 1e30, 1.0], 1e-9), "{r}");
+    // Integer models, the Evolutionary method and Value Of goals have no Sensitivity report.
+    s.execute("data.solver", json!({"objective": "D6", "variables": "B2:C2", "constraints": ["D3:D5 <= E3:E5"], "method": "simplex"})).unwrap();
+    for p in [
+        json!({"reports": ["sensitivity"], "constraints": ["D3:D5 <= E3:E5", "B2:C2 int"]}),
+        json!({"reports": ["sensitivity"], "method": "evolutionary", "constraints": ["D3:D5 <= E3:E5", "B2:C2 <= 10"]}),
+        json!({"reports": ["sensitivity"], "goal": "value", "value": 30}),
+    ] {
+        assert!(s.execute("data.solverReport", p.clone()).is_err(), "{p}");
+    }
+    // data.solver refuses a Value Of goal's Sensitivity report before it changes anything.
+    set(&mut s, &[("B2", "1"), ("C2", "1")]);
+    let undo = s.doc().unwrap().undo.len();
+    let p = json!({"objective": "D6", "goal": "value", "value": 30, "variables": "B2:C2", "constraints": ["D3:D5 <= E3:E5"], "method": "simplex", "reports": ["sensitivity"]});
+    let err = s.execute("data.solver", p).unwrap_err();
+    assert!(err.to_string().contains("maximise or minimise"), "{err}");
+    assert_eq!((at(&s, "Sheet1", "B2"), at(&s, "Sheet1", "C2")), (Value::number(1.0), Value::number(1.0)));
+    assert_eq!(s.doc().unwrap().undo.len(), undo);
+}
+
+#[test]
+fn simplex_sensitivity_ranges_of_a_redundant_equality_are_zero() {
+    // max x with x = 2 and 2x = 4: the second row repeats the first, so any change to either
+    // right side alone has no solution, and both ranges of each row are 0. With x = 2 alone the
+    // right side may rise without limit and fall by 2 (to x = 0).
+    let mut s = s();
+    set(&mut s, &[("A1", "0"), ("B1", "=A1"), ("C1", "=A1"), ("C2", "=2*A1")]);
+    let cols = ["D", "F", "G", "H"];
+    let r = s
+        .execute(
+            "data.solver",
+            json!({"objective": "B1", "variables": "A1", "constraints": ["C1 = 2", "C2 = 4"], "method": "simplex", "reports": ["sensitivity"]}),
+        )
+        .unwrap();
+    let rep = r["reports"][0].as_str().unwrap().to_string();
+    assert!(close(&report_row(&s, &rep, "$C$1", &cols), &[2.0, 2.0, 0.0, 0.0], 1e-9), "{:?}", report_row(&s, &rep, "$C$1", &cols));
+    assert!(close(&report_row(&s, &rep, "$C$2", &cols), &[4.0, 4.0, 0.0, 0.0], 1e-9), "{:?}", report_row(&s, &rep, "$C$2", &cols));
+    // The redundant row doesn't restrict the objective coefficient: x = 2 stays the only
+    // solution whatever its profit.
+    assert!(close(&report_row(&s, &rep, "$A$1", &["D", "F", "G", "H"]), &[2.0, 1.0, 1e30, 1e30], 1e-9));
+    let r = s
+        .execute(
+            "data.solver",
+            json!({"objective": "B1", "variables": "A1", "constraints": ["C1 = 2"], "method": "simplex", "reports": ["sensitivity"]}),
+        )
+        .unwrap();
+    let rep = r["reports"][0].as_str().unwrap().to_string();
+    assert!(close(&report_row(&s, &rep, "$C$1", &cols), &[2.0, 2.0, 1e30, 2.0], 1e-9), "{:?}", report_row(&s, &rep, "$C$1", &cols));
+}
+
+#[test]
+fn simplex_sensitivity_report_needs_an_optimal_point() {
+    // Wyndor Glass: (2, 6) is the optimum (36). (0, 0) is feasible but not optimal; (4, 4.8)
+    // reaches 36 but breaks plant 3's hours (21.6 > 18).
+    let mut s = s();
+    s.execute(
+        "range.setValues",
+        json!({"range": "A1", "values": [[null, "Doors", "Windows", "Used", "Available"], ["Batches", 0, 0, null, null], ["Plant 1", 1, 0, null, 4], ["Plant 2", 0, 2, null, 12], ["Plant 3", 3, 2, null, 18], ["Profit", 3, 5, null, null]]}),
+    )
+    .unwrap();
+    for r in 3..=6 {
+        set(&mut s, &[(&format!("D{r}"), &format!("=SUMPRODUCT(B{r}:C{r},$B$2:$C$2)"))]);
+    }
+    let model = json!({"objective": "D6", "variables": "B2:C2", "constraints": ["D3:D5 <= E3:E5"], "method": "simplex"});
+    let report = |b: f64, c: f64| json!({"reports": ["sensitivity"], "model": model, "values": {"B2": b, "C2": c}});
+    for (b, c) in [(0.0, 0.0), (4.0, 4.8)] {
+        let before = s.doc().unwrap().clone();
+        let err = s.execute("data.solverReport", report(b, c)).unwrap_err();
+        assert!(err.to_string().contains("not optimal"), "({b}, {c}): {err}");
+        assert!(std::sync::Arc::ptr_eq(&s.doc().unwrap().wb, &before.wb));
+    }
+    let r = s.execute("data.solverReport", report(2.0, 6.0)).unwrap();
+    assert_eq!(r["sheets"], json!(["Sensitivity Report 1"]), "{r}");
+}
+
+#[test]
+fn simplex_sensitivity_report_describes_the_vertex_the_solver_returns() {
+    // max x + y with x + y <= 4, x <= 3, y <= 3: every point from (3, 1) to (1, 3) is optimal,
+    // but the ranges belong to one vertex. The vertex Simplex LP returns is reported (from
+    // several starting values, as Solver Results passes the solve's result on); the other vertex
+    // and the midpoint are refused.
+    let model = json!({"objective": "D1", "variables": "A1:A2", "constraints": ["C1 <= 4", "A1 <= 3", "A2 <= 3"], "method": "simplex"});
+    let mut found = Vec::new();
+    for (x, y) in [("0", "0"), ("0.1", "0.7"), ("2.9", "0.3"), ("1", "1"), ("-5", "12.25")] {
+        let mut s = s();
+        set(&mut s, &[("A1", x), ("A2", y), ("C1", "=A1+A2"), ("D1", "=A1+A2")]);
+        let mut p = model.clone();
+        p["keep"] = json!(false);
+        let r = s.execute("data.solver", p).unwrap();
+        assert_eq!(r["status"], "optimal", "{r}");
+        let rep = s.execute("data.solverReport", json!({"reports": ["sensitivity"], "model": r["model"], "values": r["values"]})).unwrap();
+        assert_eq!(rep["sheets"], json!(["Sensitivity Report 1"]), "from ({x}, {y}): {rep}");
+        found.push((r["values"]["A1"].as_f64().unwrap(), r["values"]["A2"].as_f64().unwrap()));
+    }
+    let (vx, vy) = found[0];
+    assert!(found.iter().all(|v| *v == (vx, vy)), "{found:?}");
+    let other = if vx == 3.0 { (1.0, 3.0) } else { (3.0, 1.0) };
+    assert!([(3.0, 1.0), (1.0, 3.0)].contains(&(vx, vy)), "{found:?}");
+    let mut s = s();
+    set(&mut s, &[("A1", "0"), ("A2", "0"), ("C1", "=A1+A2"), ("D1", "=A1+A2")]);
+    for (x, y) in [other, (2.0, 2.0)] {
+        let before = s.doc().unwrap().clone();
+        let err = s.execute("data.solverReport", json!({"reports": ["sensitivity"], "model": model, "values": {"A1": x, "A2": y}})).unwrap_err();
+        assert!(err.to_string().contains("another optimal solution"), "({x}, {y}): {err}");
+        assert!(std::sync::Arc::ptr_eq(&s.doc().unwrap().wb, &before.wb));
+        // data.solver's Skip mode (as Solver Results passes) skips only the Sensitivity report.
+        let r = s
+            .execute(
+                "data.solverReport",
+                json!({"reports": ["answer", "sensitivity"], "model": model, "values": {"A1": x, "A2": y}, "skipNotOptimal": true}),
+            )
+            .unwrap();
+        assert_eq!(r["skippedReports"][0]["report"], "sensitivity", "{r}");
+        assert_eq!(r["sheets"].as_array().map(Vec::len), Some(1), "{r}");
+    }
+    let r = s.execute("data.solverReport", json!({"reports": ["sensitivity"], "model": model, "values": {"A1": vx, "A2": vy}})).unwrap();
+    assert!(r["sheets"][0].as_str().unwrap().starts_with("Sensitivity Report"), "{r}");
+}
+
+#[test]
+fn simplex_sensitivity_report_repeats_the_solve_from_its_start_values() {
+    // max 0.1x + 0.1y with 0.1x + 0.1y <= 0.4, 0.7x + 0.7y <= 2.8, x <= 3, y <= 3: the optimum is
+    // not unique, and the difference coefficients of these formulas differ in their last bits
+    // with the point they are taken at. From these starts the solve returns (3, 1) or (1, 3);
+    // linearised at that solution instead, the report's solve broke the tie the other way and
+    // refused the solve's own values. Built from the same start values, the program is the
+    // same to the bit, so each solve's own report is written.
+    let model =
+        json!({"objective": "D1", "variables": "A1:A2", "constraints": ["C1 <= 0.4", "E1 <= 2.8", "A1 <= 3", "A2 <= 3"], "method": "simplex"});
+    let formulas = [("C1", "=A1*0.1+A2*0.1"), ("D1", "=A1*0.1+A2*0.1"), ("E1", "=A1*0.7+A2*0.7")];
+    let program = |s: &Session, x: &[f64]| {
+        let m = super::parse_model(s, &model).unwrap();
+        let prob = super::build_problem(&m);
+        let mut ev = super::Evaluator {
+            trial: super::Trial::new(s.doc().unwrap()),
+            sheet: m.sheet,
+            vars: m.vars.clone(),
+            objective: m.objective,
+            rows: prob.rows.clone(),
+            evals: 0,
+            budget: 100,
+            target: None,
+        };
+        let p = super::lp::simplex_program(&mut ev, m.goal, x).unwrap();
+        let bits = |v: &[f64]| v.iter().map(|x| x.to_bits()).collect::<Vec<u64>>();
+        let rows: Vec<(Vec<u64>, u64)> = p.rows.iter().map(|(a, _, b)| (bits(a), b.to_bits())).collect();
+        (bits(&p.cf), bits(&p.c), rows, p.work)
+    };
+    let mut vertices = Vec::new();
+    for (x, y) in [("-3", "-5"), ("-3", "2.42"), ("-3", "6.13"), ("-3", "9.84"), ("-3", "13.55"), ("0.37", "0.3"), ("2.55", "-1.2")] {
+        let mut s = s();
+        set(&mut s, &[("A1", x), ("A2", y)]);
+        set(&mut s, &formulas);
+        let start: Vec<f64> = [x, y].iter().map(|v| v.parse().unwrap()).collect();
+        let before = program(&s, &start);
+        let r = s.execute("data.solver", model.clone()).unwrap();
+        assert_eq!(r["status"], "optimal", "{r}");
+        // The program the report builds after the solve (cells now at the solution) from the
+        // solve's `original`, as Solver Results passes it, equals the solve's.
+        let original: Vec<f64> = ["A1", "A2"].iter().map(|c| r["original"][c].as_f64().unwrap()).collect();
+        assert_eq!(original, start);
+        assert_eq!(program(&s, &original), before, "from ({x}, {y})");
+        let rep = s
+            .execute("data.solverReport", json!({"reports": ["sensitivity"], "model": r["model"], "values": r["values"], "original": r["original"]}))
+            .unwrap();
+        assert_eq!(rep["sheets"], json!(["Sensitivity Report 1"]), "from ({x}, {y}): {rep}");
+        let v = (r["values"]["A1"].as_f64().unwrap(), r["values"]["A2"].as_f64().unwrap());
+        vertices.push(if v.0 > v.1 { "(3, 1)" } else { "(1, 3)" });
+    }
+    // Both vertices occur, so the starts exercise the tie.
+    assert!(vertices.contains(&"(3, 1)") && vertices.contains(&"(1, 3)"), "{vertices:?}");
+}
+
+#[test]
+fn simplex_sensitivity_report_refuses_a_point_that_beats_the_optimum() {
+    // min 1E+9·A1 with A1 ≥ 0: the optimum is 0. A1 = -1E-7 is outside the bound by less than
+    // the precision, yet its objective is -100, far from the optimum.
+    let mut s = s();
+    set(&mut s, &[("A1", "0"), ("B1", "=1E+9*A1")]);
+    let model = json!({"objective": "B1", "goal": "min", "variables": "A1", "constraints": ["A1 <= 5"], "method": "simplex"});
+    let before = s.doc().unwrap().clone();
+    let err = s.execute("data.solverReport", json!({"reports": ["sensitivity"], "model": model, "values": {"A1": -1e-7}})).unwrap_err();
+    assert!(err.to_string().contains("not optimal"), "{err}");
+    assert!(std::sync::Arc::ptr_eq(&s.doc().unwrap().wb, &before.wb));
+    let r = s.execute("data.solverReport", json!({"reports": ["sensitivity"], "model": model, "values": {"A1": 0}})).unwrap();
+    assert_eq!(r["sheets"], json!(["Sensitivity Report 1"]), "{r}");
+}
+
+#[test]
+fn sensitivity_report_of_a_nonlinear_model_gives_lagrange_multipliers() {
+    // The largest rectangle with perimeter at most 20: the area (P/4)² grows by P/8 = 2.5 per
+    // unit of perimeter, and both reduced gradients are 0.
+    let mut s = s();
+    set(&mut s, &[("A1", "1"), ("A2", "1"), ("B1", "=A1*A2"), ("B2", "=2*A1+2*A2")]);
+    let r =
+        s.execute("data.solver", json!({"objective": "B1", "variables": "A1:A2", "constraints": ["B2 <= 20"], "reports": ["sensitivity"]})).unwrap();
+    let rep = "Sensitivity Report 1";
+    assert_eq!(r["reports"], json!([rep]), "{r}");
+    assert_eq!(at(&s, rep, &format!("E{}", row_of(&s, rep, "$A$1"))), Value::text("Reduced Gradient"));
+    assert!(close(&report_row(&s, rep, "$A$1", &["D", "E"]), &[5.0, 0.0], 1e-4));
+    assert!(close(&report_row(&s, rep, "$B$2", &["D", "E"]), &[20.0, 2.5], 1e-4));
+    // min (x + 1)² + (y − 2)² with x + y ≤ 1 and x, y ≥ 0: (0, 1). The constraint's
+    // multiplier is −2; x sits at its bound, where raising it costs 2(x + 1) + 2 = 4.
+    set(&mut s, &[("C1", "0.5"), ("C2", "0.5"), ("D1", "=(C1+1)^2+(C2-2)^2"), ("D2", "=C1+C2")]);
+    let r = s
+        .execute(
+            "data.solver",
+            json!({"objective": "D1", "goal": "min", "variables": "C1:C2", "constraints": ["D2 <= 1"], "reports": ["sensitivity"]}),
+        )
+        .unwrap();
+    let rep = r["reports"][0].as_str().unwrap().to_string();
+    assert!(close(&report_row(&s, &rep, "$C$1", &["D", "E"]), &[0.0, 4.0], 1e-4), "{r}");
+    assert!(close(&report_row(&s, &rep, "$C$2", &["D", "E"]), &[1.0, 0.0], 1e-4));
+    assert!(close(&report_row(&s, &rep, "$D$2", &["D", "E"]), &[1.0, -2.0], 1e-4));
+}
+
+#[test]
+fn nonlinear_sensitivity_report_needs_an_optimal_point() {
+    let mut s = s();
+    set(&mut s, &[("A1", "1"), ("B1", "=A1^2")]);
+    set(&mut s, &[("C1", "0"), ("D1", "=(C1-2)^2")]);
+    set(&mut s, &[("E1", "1"), ("E2", "1"), ("F1", "=E1*E2"), ("F2", "=2*E1+2*E2")]);
+    set(&mut s, &[("G1", "0"), ("G2", "1"), ("H1", "=(G1+1)^2+(G2-2)^2"), ("H2", "=G1+G2")]);
+    let square = json!({"objective": "B1", "goal": "min", "variables": "A1", "constraints": ["A1 <= 5"]});
+    let shifted = json!({"objective": "D1", "goal": "min", "variables": "C1", "constraints": ["C1 <= 5"]});
+    let rectangle = json!({"objective": "F1", "variables": "E1:E2", "constraints": ["F2 <= 20"]});
+    let bounded = json!({"objective": "H1", "goal": "min", "variables": "G1:G2", "constraints": ["H2 <= 1"]});
+    let report = |model: &serde_json::Value, values: serde_json::Value| json!({"reports": ["sensitivity"], "model": model, "values": values});
+    for (model, values, why) in [
+        // min A1² at A1 = 1: the free variable's reduced gradient is 2, not 0.
+        (&square, json!({"A1": 1}), "free gradient"),
+        // The rectangle at 6 × 6 has a larger area than the optimum but a perimeter of 24 > 20.
+        (&rectangle, json!({"E1": 6, "E2": 6}), "infeasible"),
+        // At 4 × 6 the perimeter binds, but its multiplier can't explain the gradient (6, 4).
+        (&rectangle, json!({"E1": 4, "E2": 6}), "unexplained gradient"),
+        // min (C1 − 2)² at the bound C1 = 0: raising C1 lowers the objective.
+        (&shifted, json!({"C1": 0}), "wrong-way bound"),
+        // (1, 0) in the bounded model: the binding row's multiplier would be +4, the wrong sign
+        // for a `<=` row when minimising, so the row is dropped and the gradient is unexplained.
+        (&bounded, json!({"G1": 1, "G2": 0}), "wrong-sign multiplier"),
+    ] {
+        let before = s.doc().unwrap().clone();
+        let err = s.execute("data.solverReport", report(model, values)).unwrap_err();
+        assert!(err.to_string().contains("not optimal"), "{why}: {err}");
+        let d = s.doc().unwrap();
+        assert!(std::sync::Arc::ptr_eq(&d.wb, &before.wb), "{why}");
+        assert_eq!((d.undo.len(), d.redo.len()), (before.undo.len(), before.redo.len()), "{why}");
+    }
+    // The optima are still reported.
+    for (k, (model, values)) in
+        [(&square, json!({"A1": 0})), (&shifted, json!({"C1": 2})), (&rectangle, json!({"E1": 5, "E2": 5})), (&bounded, json!({"G1": 0, "G2": 1}))]
+            .into_iter()
+            .enumerate()
+    {
+        let r = s.execute("data.solverReport", report(model, values)).unwrap();
+        assert_eq!(r["sheets"], json!([format!("Sensitivity Report {}", k + 1)]), "{r}");
+    }
+}
+
+#[test]
+fn solver_writes_the_other_reports_when_the_run_stops_before_an_optimum() {
+    // Rosenbrock's function from (-1.2, 1) with 60 evaluations: the run stops at the limit with
+    // a solution that is not the optimum (1, 1). The solve applies and the Answer report is
+    // written; the Sensitivity report is left out with its reason.
+    let mut s = s();
+    set(&mut s, &[("A1", "-1.2"), ("A2", "1"), ("B1", "=(1-A1)^2+100*(A2-A1^2)^2")]);
+    let undo = s.doc().unwrap().undo.len();
+    let r = s
+        .execute(
+            "data.solver",
+            json!({"objective": "B1", "goal": "min", "variables": "A1:A2", "nonNegative": false, "maxEvaluations": 60, "reports": ["answer", "sensitivity"]}),
+        )
+        .unwrap();
+    assert_eq!(r["status"], "iterationLimit", "{r}");
+    assert_eq!(r["changed"], true);
+    assert_eq!(n(&s, "A1"), r["values"]["A1"].as_f64().unwrap());
+    assert_eq!(r["reports"], json!(["Answer Report 1"]), "{r}");
+    assert_eq!(r["skippedReports"][0]["report"], "sensitivity", "{r}");
+    assert!(r["skippedReports"][0]["reason"].as_str().unwrap().contains("optimal"), "{r}");
+    assert_eq!(s.doc().unwrap().undo.len(), undo + 1);
+    // A solve that reaches the optimum still writes the Sensitivity report and skips nothing.
+    let r = s
+        .execute("data.solver", json!({"objective": "B1", "goal": "min", "variables": "A1:A2", "nonNegative": false, "reports": ["sensitivity"]}))
+        .unwrap();
+    assert_eq!(r["status"], "converged", "{r}");
+    assert_eq!(r["reports"], json!(["Sensitivity Report 1"]), "{r}");
+    assert!(r.get("skippedReports").is_none(), "{r}");
+}
+
+#[test]
+fn solver_report_can_skip_a_sensitivity_report_the_values_do_not_support() {
+    // min A1² at A1 = 1 is not an optimum. Called directly, data.solverReport refuses it; with
+    // skipNotOptimal (as Solver Results passes) the Answer report is still written.
+    let mut s = s();
+    set(&mut s, &[("A1", "1"), ("B1", "=A1^2")]);
+    let model = json!({"objective": "B1", "goal": "min", "variables": "A1", "constraints": ["A1 <= 5"]});
+    let p = |reports: serde_json::Value, skip: bool| json!({"reports": reports, "model": model, "values": {"A1": 1}, "skipNotOptimal": skip});
+    let err = s.execute("data.solverReport", p(json!(["answer", "sensitivity"]), false)).unwrap_err();
+    assert!(err.to_string().contains("not optimal"), "{err}");
+    let undo = s.doc().unwrap().undo.len();
+    let r = s.execute("data.solverReport", p(json!(["answer", "sensitivity"]), true)).unwrap();
+    assert_eq!(r["sheets"], json!(["Answer Report 1"]), "{r}");
+    assert!(r["skippedReports"][0]["reason"].as_str().unwrap().contains("not optimal"), "{r}");
+    assert_eq!(s.doc().unwrap().undo.len(), undo + 1);
+    // With nothing left to write, nothing changes and no Undo step is added.
+    let before = s.doc().unwrap().clone();
+    let r = s.execute("data.solverReport", p(json!(["sensitivity"]), true)).unwrap();
+    assert_eq!(r["sheets"], json!([]), "{r}");
+    let d = s.doc().unwrap();
+    assert!(std::sync::Arc::ptr_eq(&d.wb, &before.wb));
+    assert_eq!((d.undo.len(), d.redo.len()), (before.undo.len(), before.redo.len()));
+    // A model with no Sensitivity report is still refused, skipNotOptimal or not.
+    let err = s.execute("data.solverReport", json!({"reports": ["answer", "sensitivity"], "model": {"objective": "B1", "variables": "A1", "constraints": ["A1 <= 5", "A1 int"]}, "skipNotOptimal": true})).unwrap_err();
+    assert!(err.to_string().contains("integer or binary"), "{err}");
+    // In a workbook whose structure is protected, a request whose only report is skipped
+    // writes no sheet, so it succeeds; one that would write a sheet is still refused.
+    let d = s.doc_mut().unwrap();
+    std::sync::Arc::make_mut(&mut d.wb).protected_structure = true;
+    let before = s.doc().unwrap().clone();
+    let r = s.execute("data.solverReport", p(json!(["sensitivity"]), true)).unwrap();
+    assert_eq!(r["sheets"], json!([]), "{r}");
+    assert!(r["skippedReports"][0]["reason"].as_str().unwrap().contains("not optimal"), "{r}");
+    let err = s.execute("data.solverReport", p(json!(["answer", "sensitivity"]), true)).unwrap_err();
+    assert!(err.to_string().contains("protected"), "{err}");
+    // A request that must write a sheet is refused before its reports are built, so an
+    // integer model's Limits report gives the protection message, not the report's own refusal.
+    let int_model = json!({"objective": "B1", "variables": "A1", "constraints": ["A1 <= 5", "A1 int"]});
+    let err = s.execute("data.solverReport", json!({"reports": ["limits"], "model": int_model, "values": {"A1": 1}})).unwrap_err();
+    assert!(err.to_string().contains("protected"), "{err}");
+    let d = s.doc().unwrap();
+    assert!(std::sync::Arc::ptr_eq(&d.wb, &before.wb));
+    assert_eq!((d.undo.len(), d.redo.len()), (before.undo.len(), before.redo.len()));
+}
+
+#[test]
+fn a_sensitivity_refusal_of_the_model_comes_before_any_report_is_built() {
+    // The refusals that apply to the model (integer or binary constraints, the Evolutionary
+    // method, a Value Of goal) are checked before any report is built, in both modes and from
+    // both commands: a Limits report listed first would otherwise run its two optimisations per
+    // variable, and data.solver its whole solve, for an error. An integer model asking for
+    // limits and sensitivity gets the Sensitivity refusal, not the Limits report's, which the
+    // report loop reached first; a Value Of model, whose Limits report is not refused, gets the
+    // Sensitivity refusal with no sheet written and the workbook as it was.
+    let mut s = s();
+    set(&mut s, &[("A1", "1"), ("B1", "=2*A1")]);
+    let before = s.doc().unwrap().clone();
+    let int_model = json!({"objective": "B1", "variables": "A1", "constraints": ["A1 <= 5", "A1 int"]});
+    let value_model = json!({"objective": "B1", "goal": "value", "value": 4, "variables": "A1", "constraints": ["A1 <= 5"]});
+    for (model, values, why) in
+        [(&int_model, json!({"A1": 1}), "The Sensitivity report doesn't apply"), (&value_model, json!({"A1": 2}), "maximise or minimise")]
+    {
+        for skip in [false, true] {
+            let p = json!({"reports": ["limits", "sensitivity"], "model": model, "values": values, "skipNotOptimal": skip});
+            let err = s.execute("data.solverReport", p).unwrap_err();
+            assert!(err.to_string().contains(why), "{skip}: {err}");
+        }
+    }
+    let p = json!({"objective": "B1", "variables": "A1", "constraints": ["A1 <= 5", "A1 int"], "method": "simplex", "reports": ["limits", "sensitivity"]});
+    let err = s.execute("data.solver", p).unwrap_err();
+    assert!(err.to_string().contains("The Sensitivity report doesn't apply"), "{err}");
+    let d = s.doc().unwrap();
+    assert!(std::sync::Arc::ptr_eq(&d.wb, &before.wb));
+    assert_eq!((d.undo.len(), d.redo.len()), (before.undo.len(), before.redo.len()));
+}
+
+#[test]
+fn nonlinear_sensitivity_with_too_many_binding_constraints_is_refused_without_change() {
+    // 80 binding rows on two free variables pass the work limit (lowered in tests); 40 don't.
+    let mut s = s();
+    set(&mut s, &[("A1", "5"), ("A2", "5"), ("B1", "=A1*A2")]);
+    for r in 1..=80 {
+        set(&mut s, &[(&format!("D{r}"), "=A1+A2")]);
+    }
+    let report = |rows: u32| json!({"reports": ["sensitivity"], "model": {"objective": "B1", "variables": "A1:A2", "constraints": [format!("D1:D{rows} <= 10")]}});
+    let before = s.doc().unwrap().clone();
+    let err = s.execute("data.solverReport", report(80)).unwrap_err();
+    assert!(err.to_string().contains("too many binding constraints"), "{err}");
+    let d = s.doc().unwrap();
+    assert!(std::sync::Arc::ptr_eq(&d.wb, &before.wb));
+    assert_eq!((d.undo.len(), d.redo.len()), (before.undo.len(), before.redo.len()));
+    let r = s.execute("data.solverReport", report(40)).unwrap();
+    assert_eq!(r["sheets"], json!(["Sensitivity Report 1"]), "{r}");
+}
+
+#[test]
+fn sensitivity_applies_agrees_with_the_report() {
+    let mut s = s();
+    set(&mut s, &[("A1", "1"), ("B1", "=3*A1")]);
+    for (extra, ok) in [
+        (json!({}), true),
+        (json!({"goal": "min"}), true),
+        (json!({"method": "grg"}), true),
+        (json!({"goal": "value", "value": 3}), false),
+        (json!({"constraints": ["A1 <= 4", "A1 int"]}), false),
+        (json!({"constraints": ["A1 <= 4", {"lhs": "A1", "op": "bin"}]}), false),
+        (json!({"method": "evolutionary"}), false),
+    ] {
+        let mut p = json!({"objective": "B1", "variables": "A1", "constraints": ["A1 <= 4"], "method": "simplex", "keep": false});
+        for (k, v) in extra.as_object().unwrap() {
+            p[k] = v.clone();
+        }
+        assert_eq!(super::sensitivity_applies(&p), ok, "{p}");
+        let r = s.execute("data.solver", p.clone()).unwrap();
+        assert_eq!(super::sensitivity_applies(&r["model"]), ok, "{p}");
+        let rep = s.execute("data.solverReport", json!({"reports": ["sensitivity"], "model": r["model"], "values": r["values"]}));
+        assert_eq!(rep.is_ok(), ok, "{p}: {rep:?}");
+    }
+}
+
+/// How many cells of sheet `name` hold `text`.
+fn count_of(s: &Session, name: &str, text: &str) -> usize {
+    let d = s.doc().unwrap();
+    let sh = d.wb.sheet(d.wb.sheet_index(name).unwrap()).unwrap();
+    (0..60).flat_map(|r| (0..8).map(move |c| (r, c))).filter(|(r, c)| sh.value(CellRef::new(*r, *c)) == Value::text(text)).count()
+}
+
+#[test]
+fn sensitivity_report_is_the_same_whatever_constant_the_objective_adds() {
+    // A constant in the objective changes no price: max 1E9 + x with x ≤ 10 (a row) has the
+    // same report as max x with x ≤ 10. At x = 10 the row's shadow price is 1 (one more unit
+    // of right side adds one to the objective); x is basic, with reduced cost 0 and
+    // coefficient 1, which may rise without limit and fall by 1 (to 0, where the slack ties);
+    // the right side 10 may rise without limit and fall by 10 (to x = 0). A tolerance taken
+    // from the objective's value (1E-9 · 1E9 = 1) wrote the shadow price as 0.
+    let mut s = s();
+    set(&mut s, &[("A1", "0"), ("B1", "=1000000000+A1"), ("C1", "=A1"), ("E1", "0"), ("F1", "=E1"), ("G1", "=E1")]);
+    let cols = ["D", "E", "F", "G", "H"];
+    let mut rows = Vec::new();
+    for (obj, var, con, v, c) in [("B1", "A1", "C1", "$A$1", "$C$1"), ("F1", "E1", "G1", "$E$1", "$G$1")] {
+        let p = json!({"objective": obj, "variables": var, "constraints": [format!("{con} <= 10")], "method": "simplex", "reports": ["sensitivity"]});
+        let r = s.execute("data.solver", p).unwrap();
+        assert_eq!(r["status"], "optimal", "{r}");
+        let rep = r["reports"][0].as_str().unwrap().to_string();
+        rows.push((report_row(&s, &rep, v, &cols), report_row(&s, &rep, c, &cols)));
+    }
+    assert!(close(&rows[0].0, &[10.0, 0.0, 1.0, 1e30, 1.0], 1e-9), "{rows:?}");
+    assert!(close(&rows[0].1, &[10.0, 1.0, 10.0, 1e30, 10.0], 1e-9), "{rows:?}");
+    assert_eq!(rows[0], rows[1]);
+}
+
+#[test]
+fn nonlinear_sensitivity_report_is_the_same_whatever_constant_the_objective_adds() {
+    // min 1E6 + (x − 3)² with x ≤ 2.75 has the same report as min (x − 3)² with x ≤ 2.75. At
+    // x = 2.75 the objective's gradient is 2(2.75 − 3) = −0.5. As a bound on the variable cell
+    // (no constraint row) it is x's reduced gradient, −0.5; as a row y = x ≤ 2.75 it is the
+    // row's Lagrange multiplier, −0.5 (∇f = λ∇y with ∇y = 1), and x's reduced gradient is 0.
+    // A tolerance taken from the objective's value (1E-6 · 1E6 = 1) wrote both as 0. The
+    // difference gradients of the model with the constant lose about four digits to rounding,
+    // so the figures agree to 1E-3.
+    let mut s = s();
+    set(&mut s, &[("A1", "0"), ("B1", "=1000000+(A1-3)^2"), ("C1", "=A1"), ("E1", "0"), ("F1", "=(E1-3)^2"), ("G1", "=E1")]);
+    let mut figures = Vec::new();
+    for (obj, var, con, v, c) in [("B1", "A1", "C1", "$A$1", "$C$1"), ("F1", "E1", "G1", "$E$1", "$G$1")] {
+        let solve = |s: &mut Session, bound: String| {
+            let p = json!({"objective": obj, "goal": "min", "variables": var, "constraints": [bound], "reports": ["sensitivity"]});
+            let r = s.execute("data.solver", p).unwrap();
+            assert_eq!(r["status"], "converged", "{r}");
+            r["reports"][0].as_str().unwrap().to_string()
+        };
+        let rep = solve(&mut s, format!("{var} <= 2.75"));
+        let at_bound = report_row(&s, &rep, v, &["D", "E"]);
+        assert_eq!(count_of(&s, &rep, c), 0, "a bound has no constraint row");
+        let rep = solve(&mut s, format!("{con} <= 2.75"));
+        figures.push((at_bound, report_row(&s, &rep, v, &["D", "E"]), report_row(&s, &rep, c, &["D", "E"])));
+    }
+    let (at_bound, free, row) = &figures[0];
+    assert!(close(at_bound, &[2.75, -0.5], 1e-3), "{figures:?}");
+    assert!(close(free, &[2.75, 0.0], 1e-3), "{figures:?}");
+    assert!(close(row, &[2.75, -0.5], 1e-3), "{figures:?}");
+    assert!(close(&figures[1].0, at_bound, 1e-3) && close(&figures[1].1, free, 1e-3) && close(&figures[1].2, row, 1e-3), "{figures:?}");
+}
+
+#[test]
+fn simplex_sensitivity_keeps_a_small_price_beside_a_large_coefficient() {
+    // max 1E9·x + y with x ≤ 1 and y ≤ 1 as rows and x, y ≥ 0: (1, 1). Both rows bind, the
+    // basis is {x, y} with B = I, so the shadow prices are the coefficients, 1E9 and 1, and both
+    // reduced costs are 0. x's coefficient may rise without limit and fall by 1E9 (to 0, where
+    // its row's slack ties), y's without limit and by 1. Each right side may rise without limit
+    // (the other row doesn't involve the variable) and fall by 1 (to 0). The second row's price
+    // is the single term 1 · 1, so the first row's coefficient 1E9 is no rounding error in it; a
+    // noise level taken from the largest coefficient (1E-9 · 1E9 = 1) wrote it as 0.
+    let mut s = s();
+    set(&mut s, &[("A1", "0"), ("A2", "0"), ("B1", "=1000000000*A1+A2"), ("C1", "=A1"), ("C2", "=A2")]);
+    let p = json!({"objective": "B1", "variables": "A1:A2", "constraints": ["C1 <= 1", "C2 <= 1"], "method": "simplex", "reports": ["sensitivity"]});
+    let r = s.execute("data.solver", p).unwrap();
+    assert_eq!(r["status"], "optimal", "{r}");
+    let rep = r["reports"][0].as_str().unwrap().to_string();
+    let cols = ["D", "E", "F", "G", "H"];
+    let row = |c: &str| report_row(&s, &rep, c, &cols);
+    assert!(close(&row("$A$1"), &[1.0, 0.0, 1e9, 1e30, 1e9], 1e-9), "{:?}", row("$A$1"));
+    assert!(close(&row("$A$2"), &[1.0, 0.0, 1.0, 1e30, 1.0], 1e-9), "{:?}", row("$A$2"));
+    assert!(close(&row("$C$1"), &[1.0, 1e9, 1.0, 1e30, 1.0], 1e-9), "{:?}", row("$C$1"));
+    assert!(close(&row("$C$2"), &[1.0, 1.0, 1.0, 1e30, 1.0], 1e-9), "{:?}", row("$C$2"));
+}
+
+#[test]
+fn simplex_sensitivity_keeps_a_small_reduced_cost_beside_a_large_coefficient() {
+    // max 1E9·x + y − z with x ≤ 1 and y ≤ 1 as rows and x, y, z ≥ 0 (z's bound is the only
+    // constraint on it): (1, 1, 0). z is nonbasic at zero with reduced cost −1 − 0 = −1 (no
+    // row involves it): its coefficient must rise by 1 before z is made, so the Allowable
+    // Increase is 1 and the Allowable Decrease 1E+30, as for a maximum's variable at zero in
+    // `simplex_sensitivity_gives_the_dakota_furniture_report`. The other rows are those of
+    // `simplex_sensitivity_keeps_a_small_price_beside_a_large_coefficient`. The reduced cost
+    // is the coefficient alone, so x's coefficient 1E9 is no rounding error in it; a noise
+    // level taken from the largest coefficient wrote it as 0.
+    let mut s = s();
+    set(&mut s, &[("A1", "0"), ("A2", "0"), ("A3", "0"), ("B1", "=1000000000*A1+A2-A3"), ("C1", "=A1"), ("C2", "=A2")]);
+    let p = json!({"objective": "B1", "variables": "A1:A3", "constraints": ["C1 <= 1", "C2 <= 1"], "method": "simplex", "reports": ["sensitivity"]});
+    let r = s.execute("data.solver", p).unwrap();
+    assert_eq!(r["status"], "optimal", "{r}");
+    let rep = r["reports"][0].as_str().unwrap().to_string();
+    let cols = ["D", "E", "F", "G", "H"];
+    let row = |c: &str| report_row(&s, &rep, c, &cols);
+    assert!(close(&row("$A$1"), &[1.0, 0.0, 1e9, 1e30, 1e9], 1e-9), "{:?}", row("$A$1"));
+    assert!(close(&row("$A$2"), &[1.0, 0.0, 1.0, 1e30, 1.0], 1e-9), "{:?}", row("$A$2"));
+    assert!(close(&row("$A$3"), &[0.0, -1.0, -1.0, 1.0, 1e30], 1e-9), "{:?}", row("$A$3"));
+    assert!(close(&row("$C$1"), &[1.0, 1e9, 1.0, 1e30, 1.0], 1e-9), "{:?}", row("$C$1"));
+    assert!(close(&row("$C$2"), &[1.0, 1.0, 1.0, 1e30, 1.0], 1e-9), "{:?}", row("$C$2"));
+}
+
+#[test]
+fn nonlinear_sensitivity_keeps_a_small_multiplier_beside_a_large_gradient() {
+    // The model of `simplex_sensitivity_keeps_a_small_price_beside_a_large_coefficient` under
+    // GRG Nonlinear: at (1, 1) both rows bind; with ∇f = (1E9, 1) and unit row gradients the
+    // Lagrange multipliers are 1E9 and 1, and both reduced gradients are 0. The second row's
+    // multiplier is the single term 1 · 1, so the gradient's large first component is no
+    // rounding error in it; a noise level taken from the gradient's size (1E-6 · 1E9 = 1000)
+    // wrote it as 0. The difference gradient of y is the difference of two values near 1E9,
+    // where the last place is 2⁻²² ≈ 2.4E-7: the step 1E-6 rounds to 4 places, 9.54E-7, so
+    // df_y ≈ 0.954 and the multiplier is checked to 0.1. The nonlinear method stops on this
+    // model from (0, 0) with "unbounded", so the report is asked for at the optimum (1, 1).
+    let mut s = s();
+    set(&mut s, &[("A1", "0"), ("A2", "0"), ("B1", "=1000000000*A1+A2"), ("C1", "=A1"), ("C2", "=A2")]);
+    let model = json!({"objective": "B1", "variables": "A1:A2", "constraints": ["C1 <= 1", "C2 <= 1"]});
+    let r = s.execute("data.solverReport", json!({"reports": ["sensitivity"], "model": model, "values": {"A1": 1, "A2": 1}})).unwrap();
+    let rep = r["sheets"][0].as_str().unwrap().to_string();
+    let row = |c: &str| report_row(&s, &rep, c, &["D", "E"]);
+    assert!(close(&row("$A$1"), &[1.0, 0.0], 1e-3), "{:?}", row("$A$1"));
+    assert!(close(&row("$A$2"), &[1.0, 0.0], 1e-3), "{:?}", row("$A$2"));
+    assert!(close(&row("$C$1"), &[1.0, 1e9], 1e-3), "{:?}", row("$C$1"));
+    assert!(close(&row("$C$2"), &[1.0, 1.0], 0.1), "{:?}", row("$C$2"));
+}
+
+#[test]
+fn nonlinear_sensitivity_keeps_a_small_reduced_gradient_beside_a_large_gradient() {
+    // The model of `simplex_sensitivity_keeps_a_small_reduced_cost_beside_a_large_coefficient`
+    // under GRG Nonlinear: z sits at its bound 0, where its reduced gradient is its own
+    // gradient −1 (no row involves it); a noise level taken from the gradient's size (1000)
+    // wrote it as 0. Its difference gradient is one-sided at the bound, over the step 1E-6
+    // from a value near 1E9 that rounds as in
+    // `nonlinear_sensitivity_keeps_a_small_multiplier_beside_a_large_gradient`: about −0.954,
+    // checked to 0.1. As there, the report is asked for at the optimum (1, 1, 0).
+    let mut s = s();
+    set(&mut s, &[("A1", "0"), ("A2", "0"), ("A3", "0"), ("B1", "=1000000000*A1+A2-A3"), ("C1", "=A1"), ("C2", "=A2")]);
+    let model = json!({"objective": "B1", "variables": "A1:A3", "constraints": ["C1 <= 1", "C2 <= 1"]});
+    let r = s.execute("data.solverReport", json!({"reports": ["sensitivity"], "model": model, "values": {"A1": 1, "A2": 1, "A3": 0}})).unwrap();
+    let rep = r["sheets"][0].as_str().unwrap().to_string();
+    let row = |c: &str| report_row(&s, &rep, c, &["D", "E"]);
+    assert!(close(&row("$A$1"), &[1.0, 0.0], 1e-3), "{:?}", row("$A$1"));
+    assert!(close(&row("$A$2"), &[1.0, 0.0], 1e-3), "{:?}", row("$A$2"));
+    assert!(close(&row("$A$3"), &[0.0, -1.0], 0.1), "{:?}", row("$A$3"));
+    assert!(close(&row("$C$1"), &[1.0, 1e9], 1e-3), "{:?}", row("$C$1"));
+    assert!(close(&row("$C$2"), &[1.0, 1.0], 0.1), "{:?}", row("$C$2"));
+}
+
+#[test]
+fn nonlinear_sensitivity_refuses_a_point_that_is_not_optimal_in_a_small_variable() {
+    // max 1E9·x + y with x ≤ 1 and y ≤ 1 as rows, at (1, 0.5): the second row doesn't bind, y
+    // is free and its reduced gradient is its own gradient, about 1 (0.954 with the rounding
+    // in `nonlinear_sensitivity_keeps_a_small_multiplier_beside_a_large_gradient`), so raising
+    // y to 1 gains 0.5 and the point is not optimal: the report is refused and the workbook
+    // left as it was. A stationarity tolerance taken from the largest gradient component
+    // (1E-3 · 1E9 = 1E6) accepted the point and wrote the report. The optimum (1, 1) is
+    // reported.
+    let mut s = s();
+    set(&mut s, &[("A1", "0"), ("A2", "0"), ("B1", "=1000000000*A1+A2"), ("C1", "=A1"), ("C2", "=A2")]);
+    let model = json!({"objective": "B1", "variables": "A1:A2", "constraints": ["C1 <= 1", "C2 <= 1"]});
+    let report = |values: serde_json::Value| json!({"reports": ["sensitivity"], "model": model, "values": values});
+    let before = s.doc().unwrap().clone();
+    let err = s.execute("data.solverReport", report(json!({"A1": 1, "A2": 0.5}))).unwrap_err();
+    assert!(err.to_string().contains("not optimal"), "{err}");
+    let d = s.doc().unwrap();
+    assert!(std::sync::Arc::ptr_eq(&d.wb, &before.wb));
+    assert_eq!((d.undo.len(), d.redo.len()), (before.undo.len(), before.redo.len()));
+    let r = s.execute("data.solverReport", report(json!({"A1": 1, "A2": 1}))).unwrap();
+    assert_eq!(r["sheets"], json!(["Sensitivity Report 1"]), "{r}");
+}
+
+#[test]
+fn nonlinear_sensitivity_refuses_a_point_that_breaks_a_small_row_beside_a_large_one() {
+    // max y − (x − 2)² with x ≤ 1E9 and y ≤ 1 as rows: (2, 1). At (2, 1.0001) the second row
+    // is over by 1E-4, more than ten times the constraint precision of its own right side
+    // (1E-5 · (1 + 1) = 2E-5), so the point is infeasible and the report is refused; it is
+    // stationary (x's gradient is 0 and the row, within the binding tolerance, explains y's),
+    // so feasibility alone decides. A feasibility tolerance relative to the largest constant
+    // right side (1E-5 · 1E9 = 1E4) accepted the point. The optimum (2, 1) is reported.
+    let mut s = s();
+    set(&mut s, &[("A1", "0"), ("A2", "0"), ("B1", "=A2-(A1-2)^2"), ("C1", "=A1"), ("C2", "=A2")]);
+    let model = json!({"objective": "B1", "variables": "A1:A2", "constraints": ["C1 <= 1000000000", "C2 <= 1"]});
+    let report = |values: serde_json::Value| json!({"reports": ["sensitivity"], "model": model, "values": values});
+    let before = s.doc().unwrap().clone();
+    let err = s.execute("data.solverReport", report(json!({"A1": 2, "A2": 1.0001}))).unwrap_err();
+    assert!(err.to_string().contains("not optimal"), "{err}");
+    let d = s.doc().unwrap();
+    assert!(std::sync::Arc::ptr_eq(&d.wb, &before.wb));
+    assert_eq!((d.undo.len(), d.redo.len()), (before.undo.len(), before.redo.len()));
+    let r = s.execute("data.solverReport", report(json!({"A1": 2, "A2": 1}))).unwrap();
+    assert_eq!(r["sheets"], json!(["Sensitivity Report 1"]), "{r}");
+    assert!(close(&report_row(&s, "Sensitivity Report 1", "$C$2", &["D", "E"]), &[1.0, 1.0], 1e-3), "{r}");
+}
+
+#[test]
+fn simplex_sensitivity_of_a_variable_at_its_upper_bound() {
+    // max 3x + 2y with x + y ≤ 4 (a row), x ≤ 3 and y ≤ 5 (bounds on the variable cells, which
+    // have no row in the Constraints table): x = 3, y = 1, objective 11. The row's shadow price
+    // is y's coefficient 2 (y is basic in it). x, held at its upper bound, has reduced cost
+    // 3 − 1·2 = 1: one more unit of its bound adds 1. This is the convention of Frontline's
+    // "Excel Solver - Interpreting the Sensitivity Report" (solver.com): a variable's dual
+    // value is nonzero only at a bound, measures the objective's change per unit increase in
+    // the variable, and a bound on a variable cell has no row in the Constraints section, its
+    // reduced cost standing for the bound's shadow price. Its coefficient may rise without limit
+    // and fall by 1 (to 2, a tie with y): for a maximum's variable at its upper bound the
+    // Allowable Increase is 1E+30 and the Allowable Decrease the reduced cost, the mirror
+    // image of a variable at zero (`simplex_sensitivity_gives_the_dakota_furniture_report`).
+    // y's coefficient may rise by 1 (to 3, a tie with x) and fall by 2 (to 0, where the row's
+    // slack ties). The right side 4 may rise by 4 (y reaches its bound 5 at 8) and fall by 1
+    // (y reaches 0 at 3).
+    let mut s = s();
+    set(&mut s, &[("A1", "0"), ("A2", "0"), ("B1", "=3*A1+2*A2"), ("C1", "=A1+A2")]);
+    let p = json!({"objective": "B1", "variables": "A1:A2", "constraints": ["C1 <= 4", "A1 <= 3", "A2 <= 5"], "method": "simplex", "reports": ["sensitivity"]});
+    let r = s.execute("data.solver", p).unwrap();
+    assert_eq!(r["status"], "optimal", "{r}");
+    let rep = r["reports"][0].as_str().unwrap().to_string();
+    let cols = ["D", "E", "F", "G", "H"];
+    let row = |c: &str| report_row(&s, &rep, c, &cols);
+    assert!(close(&row("$A$1"), &[3.0, 1.0, 3.0, 1e30, 1.0], 1e-9), "{:?}", row("$A$1"));
+    assert!(close(&row("$A$2"), &[1.0, 0.0, 2.0, 1.0, 2.0], 1e-9), "{:?}", row("$A$2"));
+    assert!(close(&row("$C$1"), &[4.0, 2.0, 4.0, 4.0, 1.0], 1e-9), "{:?}", row("$C$1"));
+    assert_eq!((count_of(&s, &rep, "$A$1"), count_of(&s, &rep, "$A$2"), count_of(&s, &rep, "$C$1")), (1, 1, 1));
+}
+
+#[test]
+fn simplex_sensitivity_of_a_minimum_with_a_greater_row_and_a_bound() {
+    // min 2x + 3y with x + y ≥ 5 (a row) and x ≤ 2 (a bound): x = 2, y = 3, cost 13. The row's
+    // shadow price is 3: one more unit of right side costs one more y. x, at its upper bound,
+    // has reduced cost 2 − 1·3 = −1 (one more unit of its bound saves 1). Its coefficient may
+    // rise by 1 (to 3, a tie with y) and fall without limit: for a minimum's variable at its
+    // upper bound the Allowable Increase is the reduced cost's size and the Allowable Decrease
+    // 1E+30, the mirror image of a minimum's variable at zero (the reduced cost 1 with
+    // Allowable Increase 1E+30 and Allowable Decrease 1 in
+    // `sensitivity_report_gives_the_wyndor_glass_shadow_prices_and_ranges`). y's coefficient
+    // may rise without limit (x can take no more) and fall by 1 (to 2, a tie with x). The right
+    // side 5 may rise without limit and fall by 3 (y reaches 0 at 2).
+    let mut s = s();
+    set(&mut s, &[("A1", "0"), ("A2", "0"), ("B1", "=2*A1+3*A2"), ("C1", "=A1+A2")]);
+    let p = json!({"objective": "B1", "goal": "min", "variables": "A1:A2", "constraints": ["C1 >= 5", "A1 <= 2"], "method": "simplex", "reports": ["sensitivity"]});
+    let r = s.execute("data.solver", p).unwrap();
+    assert_eq!(r["status"], "optimal", "{r}");
+    let rep = r["reports"][0].as_str().unwrap().to_string();
+    let cols = ["D", "E", "F", "G", "H"];
+    let row = |c: &str| report_row(&s, &rep, c, &cols);
+    assert!(close(&row("$A$1"), &[2.0, -1.0, 2.0, 1.0, 1e30], 1e-9), "{:?}", row("$A$1"));
+    assert!(close(&row("$A$2"), &[3.0, 0.0, 3.0, 1e30, 1.0], 1e-9), "{:?}", row("$A$2"));
+    assert!(close(&row("$C$1"), &[5.0, 3.0, 5.0, 1e30, 3.0], 1e-9), "{:?}", row("$C$1"));
+}
+
+#[test]
+fn simplex_sensitivity_gives_the_dakota_furniture_report() {
+    // Winston, Operations Research: Applications and Algorithms, the Dakota Furniture problem:
+    // max 60x₁ + 30x₂ + 20x₃ (desks, tables, chairs) with 8x₁ + 6x₂ + x₃ ≤ 48 (lumber),
+    // 4x₁ + 2x₂ + 1.5x₃ ≤ 20 (finishing) and 2x₁ + 1.5x₂ + 0.5x₃ ≤ 8 (carpentry): x = (2, 0, 8),
+    // revenue 280. For the basis {s₁, x₃, x₁} the rows of B⁻¹ are (1, 2, −8), (0, 2, −4) and
+    // (0, −0.5, 1.5), so the basic values are B⁻¹b = (24, 8, 2) and the shadow prices
+    // c_B·B⁻¹ = (0, 20, 60)·B⁻¹ = (0, 10, 10). Tables, at zero, have reduced cost
+    // 30 − (6·0 + 2·10 + 1.5·10) = −5: their price must rise by 5 before they are made, so the
+    // Allowable Increase is 5 and the Allowable Decrease 1E+30, as Excel's report shows a
+    // maximum's variable at zero. A change Δ in a basic variable's coefficient moves each
+    // nonbasic reduced cost by Δ times that variable's row of the final tableau; the tables'
+    // column there is B⁻¹(6, 2, 1.5) = (−2, −2, 1.25) and the slacks' columns are B⁻¹'s. Desks
+    // (row 3): −5 − 1.25Δ ≤ 0, −10 + 0.5Δ ≤ 0 and −10 − 1.5Δ ≤ 0 hold for −4 ≤ Δ ≤ 20. Chairs
+    // (row 2): −5 + 2Δ ≤ 0, −10 − 2Δ ≤ 0 and −10 + 4Δ ≤ 0 hold for −5 ≤ Δ ≤ 2.5. A right side
+    // moves the basic values along its column of B⁻¹: lumber (24 + Δ, 8, 2) ≥ 0 for Δ ≥ −24;
+    // finishing (24 + 2Δ, 8 + 2Δ, 2 − 0.5Δ) ≥ 0 for −4 ≤ Δ ≤ 4; carpentry
+    // (24 − 8Δ, 8 − 4Δ, 2 + 1.5Δ) ≥ 0 for −4/3 ≤ Δ ≤ 2.
+    let mut s = s();
+    set(
+        &mut s,
+        &[
+            ("A1", "0"),
+            ("A2", "0"),
+            ("A3", "0"),
+            ("B1", "=60*A1+30*A2+20*A3"),
+            ("C1", "=8*A1+6*A2+A3"),
+            ("C2", "=4*A1+2*A2+1.5*A3"),
+            ("C3", "=2*A1+1.5*A2+0.5*A3"),
+        ],
+    );
+    let p = json!({"objective": "B1", "variables": "A1:A3", "constraints": ["C1 <= 48", "C2 <= 20", "C3 <= 8"], "method": "simplex", "reports": ["sensitivity"]});
+    let r = s.execute("data.solver", p).unwrap();
+    assert_eq!(r["status"], "optimal", "{r}");
+    assert_eq!(r["objective"], 280.0, "{r}");
+    let rep = r["reports"][0].as_str().unwrap().to_string();
+    let cols = ["D", "E", "F", "G", "H"];
+    let row = |c: &str| report_row(&s, &rep, c, &cols);
+    assert!(close(&row("$A$1"), &[2.0, 0.0, 60.0, 20.0, 4.0], 1e-9), "{:?}", row("$A$1"));
+    assert!(close(&row("$A$2"), &[0.0, -5.0, 30.0, 5.0, 1e30], 1e-9), "{:?}", row("$A$2"));
+    assert!(close(&row("$A$3"), &[8.0, 0.0, 20.0, 2.5, 5.0], 1e-9), "{:?}", row("$A$3"));
+    assert!(close(&row("$C$1"), &[24.0, 0.0, 48.0, 1e30, 24.0], 1e-9), "{:?}", row("$C$1"));
+    assert!(close(&row("$C$2"), &[20.0, 10.0, 20.0, 4.0, 4.0], 1e-9), "{:?}", row("$C$2"));
+    assert!(close(&row("$C$3"), &[8.0, 10.0, 8.0, 2.0, 4.0 / 3.0], 1e-9), "{:?}", row("$C$3"));
+}
+
+#[test]
+fn nonlinear_sensitivity_signs_follow_the_goal() {
+    // A Lagrange multiplier is the derivative of the optimal objective with respect to the
+    // row's right side, and a reduced gradient that with respect to the variable's bound,
+    // whichever way the goal points. Frontline's "Excel Solver - Interpreting the Sensitivity
+    // Report" (solver.com) names a binding constraint's dual value its Lagrange multiplier and a
+    // variable's its reduced gradient, each the objective's change per unit increase in the
+    // constraint's bound or the variable, valid at the optimal point alone for a nonlinear
+    // model.
+    // min (x − 3)² + (y − 1)² with x + y ≤ b: the point of the line nearest (3, 1) is
+    // (2.5, 0.5) at b = 3, f*(b) = (4 − b)²/2 and df*/db = −(4 − b) = −1, the multiplier;
+    // maximising the negated objective reaches the same point with multiplier +1. Both
+    // variables are interior, so their reduced gradients are 0.
+    let mut s = s();
+    set(&mut s, &[("A1", "0"), ("A2", "0"), ("B1", "=(A1-3)^2+(A2-1)^2"), ("B2", "=0-(A1-3)^2-(A2-1)^2"), ("C1", "=A1+A2")]);
+    for (obj, goal, lambda) in [("B1", "min", -1.0), ("B2", "max", 1.0)] {
+        set(&mut s, &[("A1", "0"), ("A2", "0")]);
+        let p = json!({"objective": obj, "goal": goal, "variables": "A1:A2", "constraints": ["C1 <= 3"], "reports": ["sensitivity"]});
+        let r = s.execute("data.solver", p).unwrap();
+        let rep = r["reports"][0].as_str().unwrap().to_string();
+        assert!(close(&report_row(&s, &rep, "$A$1", &["D", "E"]), &[2.5, 0.0], 1e-3), "{goal}: {r}");
+        assert!(close(&report_row(&s, &rep, "$A$2", &["D", "E"]), &[0.5, 0.0], 1e-3), "{goal}: {r}");
+        assert!(close(&report_row(&s, &rep, "$C$1", &["D", "E"]), &[3.0, lambda], 1e-3), "{goal}: {r}");
+    }
+    // max −(x + 1)² − (y − 2)² with x + y ≤ 1 and x, y ≥ 0 mirrors the minimum in
+    // `sensitivity_report_of_a_nonlinear_model_gives_lagrange_multipliers`: at (0, 1) the
+    // multiplier is +2 (−2 there) and x's reduced gradient −4 (4 there), since raising x from
+    // its bound lowers the maximum by 2(x + 1) + 2 = 4 per unit.
+    set(&mut s, &[("D1", "0.5"), ("D2", "0.5"), ("E1", "=0-(D1+1)^2-(D2-2)^2"), ("E2", "=D1+D2")]);
+    let p = json!({"objective": "E1", "goal": "max", "variables": "D1:D2", "constraints": ["E2 <= 1"], "reports": ["sensitivity"]});
+    let r = s.execute("data.solver", p).unwrap();
+    let rep = r["reports"][0].as_str().unwrap().to_string();
+    assert!(close(&report_row(&s, &rep, "$D$1", &["D", "E"]), &[0.0, -4.0], 1e-3), "{r}");
+    assert!(close(&report_row(&s, &rep, "$D$2", &["D", "E"]), &[1.0, 0.0], 1e-3), "{r}");
+    assert!(close(&report_row(&s, &rep, "$E$2", &["D", "E"]), &[1.0, 2.0], 1e-3), "{r}");
+}
+
 fn formula(s: &Session, a1: &str) -> String {
     let d = s.doc().unwrap();
     d.wb.active().unwrap().cell(CellRef::parse(a1).unwrap()).and_then(|c| c.formula.as_ref()).map(|f| f.text.clone()).unwrap_or_default()

crates/ui-egui/src/dialogs/analysis.rs

@@ -632,12 +632,25 @@ pub(super) fn solver_results(app: &mut SheetApp, ui: &mut egui::Ui, d: &mut Dial
     }
     // Reports describe the solution whether it is kept or not.
     let solved = d.values.get("values").and_then(Json::as_object).is_some_and(|v| !v.is_empty());
+    // Integer models, the Evolutionary method and Value Of goals have no Sensitivity report, and
+    // it describes an optimum, so a run that stopped before one (at the evaluation limit) has
+    // none either.
+    let model_allows = d.values.get("model").is_none_or(gridcraft_engine::cmd::solver::sensitivity_applies);
+    let optimum = d.values.get("status").and_then(Json::as_str).is_some_and(gridcraft_engine::cmd::solver::sensitivity_status);
+    let sensitivity = model_allows && optimum;
+    let why_not = if model_allows {
+        "Needs an optimal or converged solution: Solver stopped before reaching one."
+    } else {
+        "Needs Simplex LP or GRG Nonlinear, a Max or Min objective and no integer constraints."
+    };
+    let offered = |key: &str| key != "sensitivityReport" || sensitivity;
     if solved {
         ui.add_space(4.0);
         ui.label("Reports:");
-        for (key, label) in [("answerReport", "Answer"), ("limitsReport", "Limits")] {
+        for (key, label) in [("answerReport", "Answer"), ("sensitivityReport", "Sensitivity"), ("limitsReport", "Limits")] {
-            let mut on = d.values.get(key).and_then(Json::as_bool).unwrap_or(false);
+            let mut on = offered(key) && d.values.get(key).and_then(Json::as_bool).unwrap_or(false);
-            if ui.checkbox(&mut on, label).changed() {
+            let check = ui.add_enabled(offered(key), egui::Checkbox::new(&mut on, label)).on_disabled_hover_text(why_not);
+            if check.changed() {
                 d.values.insert(key.into(), json!(on));
             }
         }
@@ -649,15 +662,32 @@ pub(super) fn solver_results(app: &mut SheetApp, ui: &mut egui::Ui, d: &mut Dial
         if changed && restore {
             app.run_or_alert("edit.undo", json!({}));
         }
-        let reports: Vec<&str> = [("answerReport", "answer"), ("limitsReport", "limits")]
+        let reports: Vec<&str> = [("answerReport", "answer"), ("sensitivityReport", "sensitivity"), ("limitsReport", "limits")]
             .into_iter()
-            .filter(|(key, _)| d.values.get(*key).and_then(Json::as_bool) == Some(true))
+            .filter(|(key, _)| offered(key) && d.values.get(*key).and_then(Json::as_bool) == Some(true))
             .map(|(_, kind)| kind)
             .collect();
         if solved && !reports.is_empty() {
             let field = |k: &str| d.values.get(k).cloned().unwrap_or(Json::Null);
-            let params = json!({"reports": reports, "model": field("model"), "values": field("values"), "original": field("original"), "message": field("message")});
+            // A Sensitivity report the values turn out not to support is left out, and the other
-            app.run_or_alert("data.solverReport", params);
+            // reports are still written.
+            let params = json!({"reports": reports, "model": field("model"), "values": field("values"), "original": field("original"),
+                "message": field("message"), "skipNotOptimal": true});
+            match app.run("data.solverReport", params) {
+                Ok(r) => {
+                    let reasons: Vec<&str> = r
+                        .get("skippedReports")
+                        .and_then(Json::as_array)
+                        .into_iter()
+                        .flatten()
+                        .filter_map(|k| k.get("reason").and_then(Json::as_str))
+                        .collect();
+                    if !reasons.is_empty() {
+                        app.message = Some(("Solver Results".into(), reasons.join("\n")));
+                    }
+                }
+                Err(e) => app.message = Some(("Solver Results".into(), crate::clean_error(&e))),
+            }
         }
         *open = false;
     }

crates/ui-egui/tests/it/interaction.rs

@@ -1,7 +1,7 @@
 //! UI interaction tests: real egui input through the whole app (no GPU needed).
 
 use egui::{Event, Key, Modifiers, PointerButton};
-use egui_kittest::kittest::Queryable as _;
+use egui_kittest::kittest::{NodeT as _, Queryable as _};
 use gridcraft_engine::Session;
 use gridcraft_engine::core::{CellRef, Value};
 use gridcraft_ui_egui::i18n::Language;
@@ -905,6 +905,63 @@ fn solver_dialog_loads_the_saved_model_and_results_can_restore() {
     assert_eq!(value(&h, "B1"), Value::Number(0.0));
 }
 
+#[test]
+fn solver_results_skip_the_sensitivity_report_for_integer_models() {
+    let mut s = blank();
+    for (c, v) in [("A1", "0"), ("B1", "=3*A1")] {
+        s.execute("cell.set", json!({"cell": c, "input": v})).unwrap();
+    }
+    let r =
+        s.execute("data.solver", json!({"objective": "B1", "variables": "A1", "constraints": ["A1 <= 4.5", "A1 int"], "method": "simplex"})).unwrap();
+    let mut h = harness(s);
+    // Sensitivity stays off for an integer model, so the Answer report is still written.
+    h.state_mut().open_dialog("solverResults", r);
+    for key in ["answerReport", "sensitivityReport"] {
+        h.state_mut().dialog.as_mut().unwrap().values.insert(key.into(), json!(true));
+    }
+    h.run_steps(2);
+    key(&mut h, Key::Enter, Modifiers::NONE);
+    assert!(h.state().dialog.is_none());
+    assert!(h.state().message.is_none(), "{:?}", h.state().message);
+    let names: Vec<String> = h.state().session.active().unwrap().wb.sheets.iter().map(|s| s.name.clone()).collect();
+    assert_eq!(names, ["Answer Report 1", "Sheet1"]);
+}
+
+#[test]
+fn solver_results_offer_sensitivity_only_after_an_optimum() {
+    // Rosenbrock's function with 60 evaluations stops at the limit with a solution that is not
+    // an optimum: Sensitivity is disabled, and Answer with Sensitivity ticked writes the Answer
+    // report.
+    let mut s = blank();
+    for (c, v) in [("A1", "-1.2"), ("A2", "1"), ("B1", "=(1-A1)^2+100*(A2-A1^2)^2")] {
+        s.execute("cell.set", json!({"cell": c, "input": v})).unwrap();
+    }
+    let r = s
+        .execute("data.solver", json!({"objective": "B1", "goal": "min", "variables": "A1:A2", "nonNegative": false, "maxEvaluations": 60}))
+        .unwrap();
+    assert_eq!(r["status"], "iterationLimit", "{r}");
+    let mut h = harness(s);
+    h.state_mut().open_dialog("solverResults", r.clone());
+    h.run_steps(2);
+    assert!(h.get_by_label("Sensitivity").accesskit_node().is_disabled());
+    assert!(!h.get_by_label("Answer").accesskit_node().is_disabled());
+    for key in ["answerReport", "sensitivityReport"] {
+        h.state_mut().dialog.as_mut().unwrap().values.insert(key.into(), json!(true));
+    }
+    h.run_steps(2);
+    key(&mut h, Key::Enter, Modifiers::NONE);
+    assert!(h.state().dialog.is_none());
+    assert!(h.state().message.is_none(), "{:?}", h.state().message);
+    let names: Vec<String> = h.state().session.active().unwrap().wb.sheets.iter().map(|s| s.name.clone()).collect();
+    assert_eq!(names, ["Answer Report 1", "Sheet1"]);
+    // The same model after an optimum offers it.
+    let mut ok = r;
+    ok["status"] = json!("converged");
+    h.state_mut().open_dialog("solverResults", ok);
+    h.run_steps(2);
+    assert!(!h.get_by_label("Sensitivity").accesskit_node().is_disabled());
+}
+
 #[test]
 fn solver_dialog_saves_loads_and_resets_the_model() {
     let mut s = blank();
@@ -1201,9 +1258,9 @@ fn solver_results_write_the_chosen_reports() {
         )
         .unwrap();
     let mut h = harness(s);
-    // Restore Original Values with the Answer and Limits reports: the cells go back, the reports stay.
+    // Restore Original Values with the three reports: the cells go back, the reports stay.
     h.state_mut().open_dialog("solverResults", r);
-    for key in ["restore", "answerReport", "limitsReport"] {
+    for key in ["restore", "answerReport", "sensitivityReport", "limitsReport"] {
         h.state_mut().dialog.as_mut().unwrap().values.insert(key.into(), json!(true));
     }
     h.run_steps(2);
@@ -1211,7 +1268,7 @@ fn solver_results_write_the_chosen_reports() {
     assert!(h.state().dialog.is_none());
     assert!(h.state().message.is_none(), "{:?}", h.state().message);
     let names: Vec<String> = h.state().session.active().unwrap().wb.sheets.iter().map(|s| s.name.clone()).collect();
-    assert_eq!(names, ["Answer Report 1", "Limits Report 1", "Sheet1"]);
+    assert_eq!(names, ["Answer Report 1", "Sensitivity Report 1", "Limits Report 1", "Sheet1"]);
     assert_eq!(value(&h, "B1"), Value::Number(0.0));
 }
 

docs/gaps.md

@@ -1,6 +1,6 @@
 # Gaps
 
-> **Last reviewed:** 2026-10-11 · **Last updated:** 2026-11-02 · **Change:** minor (§4: chart resizing is no longer a gap; earlier: §13: Analysis ToolPak has all 19 tools; earlier: §4: chart data tables implemented with XLSX round-trip and print/PDF; earlier: §13: Solver Options; earlier: §13: Solver Load/Save Model and Reset All) · **Target:** Microsoft Excel (Microsoft 365)
+> **Last reviewed:** 2026-10-11 · **Last updated:** 2026-11-02 · **Change:** minor (§13: Solver Sensitivity report; earlier: §4: chart resizing is no longer a gap; earlier: §13: Analysis ToolPak has all 19 tools; earlier: §4: chart data tables implemented with XLSX round-trip and print/PDF; earlier: §13: Solver Options; earlier: §13: Solver Load/Save Model and Reset All) · **Target:** Microsoft Excel (Microsoft 365)
 
 Every known shortfall against Excel, one entry each, ranked by user impact (what blocks someone
 using GridCraft for real work first). This is the work list: pick from the top, check the open
@@ -121,12 +121,12 @@ No alpha blockers: every core workflow passes the [alpha gate](roadmap.md#alpha-
 - Forecast Sheet (`data.forecastSheet`) and the FORECAST.ETS family are in, including Forecast
   Start (hindcast, #785).
 - Solver (`data.solver`, #327) has Simplex LP with integer/binary branch and bound, a local
-  nonlinear method (augmented Lagrangian, not GRG), the Evolutionary method, and Answer and
+  nonlinear method (augmented Lagrangian, not GRG), the Evolutionary method, and Answer,
-  Limits reports. Report refusal leaves the solution, saved model and history unchanged; a
+  Sensitivity and Limits reports. Report refusal leaves the solution, saved model and history unchanged; a
   successful solve and its reports share one Undo step (#1050). Load/Save Model writes the
   model on cells of its sheet and reads it back; Reset All removes it. Solver Options sets the
   evaluation limit, Integer Optimality and the Evolutionary options, which are saved with the
-  model. AllDifferent constraints and the Sensitivity and Population reports remain.
+  model. AllDifferent constraints and the Population report remain.
 - The Analysis ToolPak (`data.analysisToolPak`) has all 19 of Excel's tools; their chart output
   remains. **Estimate:** 12–20 h.
 
@@ -185,6 +185,7 @@ No alpha blockers: every core workflow passes the [alpha gate](roadmap.md#alpha-
 
 | Date | Change | Summary |
 |---|---|---|
+| 2026-11-02 | minor | §13: Solver Sensitivity report: shadow prices, reduced costs and allowable ranges for Simplex LP; Lagrange multipliers and reduced gradients for GRG Nonlinear |
 | 2026-11-02 | minor | §4: charts, like every floating object, move and resize from eight handles (#190), so resizing is no longer listed as missing; the evidence is one paragraph again, without the stale "PRs #158, #161 open" |
 | 2026-11-01 | minor | §13: the Analysis ToolPak has all 19 of Excel's tools (two-factor Anova with and without replication, z-Test, Fourier Analysis, Random Number Generation, Sampling); its chart output remains. The stale "`data.solver`, `data.analysisToolPak` missing" and "no Forecast Start" sentences removed |
 | 2026-11-01 | minor | §4: chart data tables implemented for column, bar, line, area and compatible combo charts, with legend keys and border settings, live source values, XLSX round-trip and print/PDF; other chart gaps remain. |

docs/target-app-parity.md

@@ -1,6 +1,6 @@
 # Target-app parity: GridCraft vs Microsoft Excel
 
-> **Last reviewed:** 2026-10-11 · **Last updated:** 2026-11-02 · **Change:** minor (chart resizing no longer listed as missing; no score change; earlier: what-if and analysis: Analysis ToolPak has all 19 tools; earlier: page layout and print: Page Setup's page order; no score change; earlier: charts: data tables with XLSX round-trip and print/PDF; no score change; earlier: what-if and analysis: Solver Options; the Insert Function dialog is named as in Excel; no score change; earlier: Solver Load/Save Model and Reset All; automation: macro interface; catalog count 276 / 292, as `cargo xtask parity` measures; automation: typed macros in the engine; several windows: a native window per workbook window, New Window on the same workbook, Arrange All; no score change) · **Target:** Microsoft Excel (Microsoft 365)
+> **Last reviewed:** 2026-10-11 · **Last updated:** 2026-11-02 · **Change:** minor (what-if and analysis: Solver Sensitivity report; no score change; earlier: chart resizing no longer listed as missing; no score change; earlier: what-if and analysis: Analysis ToolPak has all 19 tools; earlier: page layout and print: Page Setup's page order; no score change; earlier: charts: data tables with XLSX round-trip and print/PDF; no score change; earlier: what-if and analysis: Solver Options; the Insert Function dialog is named as in Excel; no score change; earlier: Solver Load/Save Model and Reset All; automation: macro interface; catalog count 276 / 292, as `cargo xtask parity` measures; automation: typed macros in the engine; several windows: a native window per workbook window, New Window on the same workbook, Arrange All; no score change) · **Target:** Microsoft Excel (Microsoft 365)
 
 This is the authoritative parity assessment. [`ROADMAP.md`](../ROADMAP.md) summarizes it,
 [`gaps.md`](gaps.md) lists every shortfall, and the area checklists hold the detail:
@@ -89,7 +89,7 @@ performance 8, stability 8, localization 4, hardware 2, platforms 2, ecosystem 2
 | Page layout, print, PDF | 5 | ~80% | ~55% | 12–20 | Page setup, print titles, breaks, page order, Page Break Preview, headers/footers, PDF export and printing. Paper sizes match Japanese Excel: B4/B5/B6 are the JIS sizes, ISO B4/B5/B6 stay reachable by name, and Japanese Postcard (はがき), Japanese Double Postcard (往復はがき), the chou/kaku envelopes and A6 are available. No true Page Layout view; the PDF writer uses only the 14 standard fonts without embedding (no CJK or other non-Latin text in PDF); web print only exports (#141, #143) |
 | Charts and sparklines | 8 | ~50% | ~30% | 40–60 | 25 chart kinds, all 2-D. The chart model holds title, legend position, a data-labels flag, gridlines flag, axis titles and chart data tables (legend keys, borders, XLSX round-trip and print/PDF); no axis scaling, error bars, per-point formatting, chart sheets, 3-D, surface, map, Pareto. Kinds survive XLSX round trip (#158, #161). [chart-parity.md](chart-parity.md) |
 | Page layout, print, PDF | 5 | ~80% | ~55% | 12–20 | Page setup, print titles, breaks, page order, Page Break Preview, headers/footers, PDF export and printing. No true Page Layout view; the PDF writer uses only the 14 standard fonts without embedding (no CJK or other non-Latin text in PDF); web print only exports (#141, #143) |
-| What-if and analysis | 4 | ~50% | ~40% | 22–35 | Goal Seek, Scenario Manager, data tables, Analysis ToolPak (all 19 tools, no charts), Solver (Simplex LP, a local nonlinear method and Evolutionary; Answer and Limits reports; Load/Save Model, Reset All and Solver Options; no Sensitivity report), Forecast Sheet. Missing: Analyze Data |
+| What-if and analysis | 4 | ~50% | ~40% | 22–35 | Goal Seek, Scenario Manager, data tables, Analysis ToolPak (all 19 tools, no charts), Solver (Simplex LP, a local nonlinear method and Evolutionary; Answer, Sensitivity and Limits reports; Load/Save Model, Reset All and Solver Options; no Population report), Forecast Sheet. Missing: Analyze Data |
 | Get & Transform, external data, data model | 5 | ~15% | ~10% | 60–100 | Import of CSV/TSV, JSON, HTML. No Power Query editor, no From Web, database or folder sources, no connections or refresh, no relationships or data model, no external workbook links (dropped with a warning) |
 | Insert objects | 3 | ~55% | ~45% | 25–40 | Pictures, 8 shape kinds (Excel has ~170), text boxes, own icon set, ink with Ink to Shape, hyperlinks, checkboxes. Missing: SmartArt, WordArt, equations, signature line, objects, 3-D models, screenshot, pictures in cells (PR #105) |
 | Review, comments, protection, collaboration | 4 | ~55% | ~40% | 35–55 | Threaded comments, notes, spelling (system dictionary), thesaurus, sheet and workbook protection, Allow Edit Ranges with Excel-compatible password hashes (#307). Missing: comment navigation, Show Changes, notes-to-comments conversion, accessibility checker, co-authoring, opening password-encrypted files (refused with a message since #17) |
@@ -334,6 +334,7 @@ The inventory carried over from the 2026-10-07 ROADMAP.md, updated for what land
 
 | Date | Change | Summary |
 |---|---|---|
+| 2026-11-02 | minor | What-if and analysis: Solver Sensitivity report (no score change) |
 | 2026-11-02 | minor | Charts move and resize from eight handles like every floating object (#190): the Essentials chart work is polish only (5–8 h), and the chart rows and the interaction discount no longer list resizing as missing; the Mainstream chart row no longer lists trendlines as missing (no score change) |
 | 2026-11-01 | minor | What-if and analysis: Analysis ToolPak has all 19 tools (was 13 of 19; no charts); no score change |
 | 2026-11-01 | minor | Charts: data tables for column, bar, line, area and compatible combo charts, with legend keys and borders; live source values, native XLSX round-trip and shared print/PDF rendering. Readiness estimates unchanged. |