From 8c351e45204e5ec1a1e5760d60a6e749707ee8f5 Mon Sep 17 00:00:00 2001 From: Ben Dudson Date: Wed, 11 Feb 2026 10:49:18 -0800 Subject: [PATCH 01/13] snes: Change target_its to BoutReal Allows the target to be set to push timestep rather than stall once reaching a set number of iterations. --- src/solver/impls/snes/snes.hxx | 7 +++++-- 1 file changed, 5 insertions(+), 2 deletions(-) diff --git a/src/solver/impls/snes/snes.hxx b/src/solver/impls/snes/snes.hxx index 5d007f8ad7..15f97c42a7 100644 --- a/src/solver/impls/snes/snes.hxx +++ b/src/solver/impls/snes/snes.hxx @@ -162,8 +162,11 @@ private: ///< PID controller parameters bool pid_controller; ///< Use PID controller? - int target_its; ///< Target number of nonlinear iterations for the PID controller. - ///< Use with caution! Not tested values. + + /// Target number of nonlinear iterations for the PID controller. + /// This can be non-integer to push timestep more aggressively + BoutReal target_its; + BoutReal kP; ///< (0.6 - 0.8) Proportional parameter (main response to current step) BoutReal kI; ///< (0.2 - 0.4) Integral parameter (smooths history of changes) BoutReal kD; ///< (0.1 - 0.3) Derivative (dampens oscillation - optional) From 20613a4ad6419218da21e785addf01a88159aa8b Mon Sep 17 00:00:00 2001 From: Ben Dudson Date: Wed, 11 Feb 2026 10:53:55 -0800 Subject: [PATCH 02/13] snes: last_failure_weight input for tuning recent_failure_rate Calculates a moving average of recent SNES failures. Enables timestep controllers to be more cautious when SNES is repeatedly failing. --- src/solver/impls/snes/snes.cxx | 7 +++++-- src/solver/impls/snes/snes.hxx | 2 +- 2 files changed, 6 insertions(+), 3 deletions(-) diff --git a/src/solver/impls/snes/snes.cxx b/src/solver/impls/snes/snes.cxx index da66092e84..7684c79a86 100644 --- a/src/solver/impls/snes/snes.cxx +++ b/src/solver/impls/snes/snes.cxx @@ -545,6 +545,9 @@ SNESSolver::SNESSolver(Options* opts) (*options)["pid_consider_failures"] .doc("Reduce timestep increases if recent solves have failed") .withDefault(false)), + last_failure_weight((*options)["last_failure_weight"] + .doc("Weighting of last timestep in recent failure rate") + .withDefault(0.1)), diagnose( (*options)["diagnose"].doc("Print additional diagnostics").withDefault(false)), diagnose_failures((*options)["diagnose_failures"] @@ -979,12 +982,12 @@ int SNESSolver::run() { SNESGetLinearSolveIterations(snes, &lin_its); // Rolling average of recent failures - recent_failure_rate *= 1. - inv_failure_window; + recent_failure_rate *= 1. - last_failure_weight; if ((ierr != PETSC_SUCCESS) or (reason < 0)) { // Diverged or SNES failed - recent_failure_rate += inv_failure_window; + recent_failure_rate += last_failure_weight; ++snes_failures; steps_since_snes_failure = 0; diff --git a/src/solver/impls/snes/snes.hxx b/src/solver/impls/snes/snes.hxx index 15f97c42a7..0cf9eef4b0 100644 --- a/src/solver/impls/snes/snes.hxx +++ b/src/solver/impls/snes/snes.hxx @@ -172,7 +172,7 @@ private: BoutReal kD; ///< (0.1 - 0.3) Derivative (dampens oscillation - optional) bool pid_consider_failures; ///< Reduce timestep increases if recent solves have failed BoutReal recent_failure_rate; ///< Rolling average of recent failure rate - const BoutReal inv_failure_window = 0.1; ///< 1 / number of recent solves + BoutReal last_failure_weight; ///< 1 / number of recent solves int nl_its_prev; int nl_its_prev2; From a4c48fb4db25d7957f580db68b2ea76ac3510bf9 Mon Sep 17 00:00:00 2001 From: Ben Dudson Date: Wed, 11 Feb 2026 10:57:09 -0800 Subject: [PATCH 03/13] snes: updateGlobalTimestep for global timestep control Separate out the global timestep control into a function --- src/solver/impls/snes/snes.cxx | 89 +++++++++++++--------------------- src/solver/impls/snes/snes.hxx | 4 ++ 2 files changed, 39 insertions(+), 54 deletions(-) diff --git a/src/solver/impls/snes/snes.cxx b/src/solver/impls/snes/snes.cxx index 7684c79a86..32a9902543 100644 --- a/src/solver/impls/snes/snes.cxx +++ b/src/solver/impls/snes/snes.cxx @@ -851,7 +851,6 @@ int SNESSolver::run() { bool looping = true; int snes_failures = 0; // Count SNES convergence failures - int steps_since_snes_failure = 0; int saved_jacobian_lag = 0; int loop_count = 0; recent_failure_rate = 0.0; @@ -990,7 +989,6 @@ int SNESSolver::run() { recent_failure_rate += last_failure_weight; ++snes_failures; - steps_since_snes_failure = 0; if (diagnose_failures or (snes_failures == max_snes_failures)) { // Print diagnostics to help identify source of the problem @@ -1123,7 +1121,6 @@ int SNESSolver::run() { } simtime += dt; - ++steps_since_snes_failure; if (diagnose) { // Gather and print diagnostic information @@ -1145,67 +1142,25 @@ int SNESSolver::run() { } if (equation_form == BoutSnesEquationForm::pseudo_transient) { - if (pid_controller) { - // Adjust pseudo_alpha based on nonlinear iterations - pseudo_alpha = pid(pseudo_alpha, nl_its, max_timestep * atol * 100); - } + // Adjust pseudo_alpha to globally scale timesteps + pseudo_alpha = updateGlobalTimestep(pseudo_alpha, nl_its, recent_failure_rate, max_timestep * atol * 100); // Adjust local timesteps PetscCall(updatePseudoTimestepping(snes_x)); - } else if (pid_controller) { - // Changing the timestep using a PID controller. - // Note: The preconditioner depends on the timestep, - // so we recalculate the jacobian and the preconditioner - // every time the timestep changes - - timestep = pid(timestep, nl_its, max_timestep); + } else { + // Adjust timestep + timestep = updateGlobalTimestep(timestep, nl_its, recent_failure_rate, max_timestep); + } - // NOTE(malamast): Do we really need this? - // Recompute Jacobian (for now) + if (static_cast(lin_its) / nl_its > 0.5 * maxl) { + // Recompute Jacobian if number of linear iterations is too high if (saved_jacobian_lag == 0) { SNESGetLagJacobian(snes, &saved_jacobian_lag); SNESSetLagJacobian(snes, 1); } - - } else { - // Consider changing the timestep. - // Note: The preconditioner depends on the timestep, - // so if it is not recalculated the it will be less - // effective. - if ((nl_its <= lower_its) && (timestep < max_timestep) - && (steps_since_snes_failure > 2)) { - // Increase timestep slightly - timestep *= timestep_factor_on_lower_its; - - timestep = std::min(timestep, max_timestep); - - // Note: Setting the SNESJacobianFn to NULL retains - // previously set evaluation function. - // - // The SNES Jacobian is a combination of the RHS Jacobian - // and a factor involving the timestep. - // Depends on equation_form - // -> Probably call SNESSetJacobian(snes, Jfd, Jfd, NULL, fdcoloring); - if (static_cast(lin_its) / nl_its > 4) { - // Recompute Jacobian (for now) - if (saved_jacobian_lag == 0) { - SNESGetLagJacobian(snes, &saved_jacobian_lag); - SNESSetLagJacobian(snes, 1); - } - } - - } else if (nl_its >= upper_its) { - // Reduce timestep slightly - timestep *= timestep_factor_on_upper_its; - - // Recompute Jacobian - if (saved_jacobian_lag == 0) { - SNESGetLagJacobian(snes, &saved_jacobian_lag); - SNESSetLagJacobian(snes, 1); - } - } } + snes_failures = 0; } while (looping); @@ -1274,6 +1229,32 @@ int SNESSolver::run() { return 0; } +BoutReal SNESSolver::updateGlobalTimestep(BoutReal timestep, int nl_its, + BoutReal recent_failure_rate, BoutReal max_dt) { + // Note: The preconditioner depends on the timestep, + // so if it is not recalculated the it will be less + // effective. + + if (pid_controller) { + // Changing the timestep using a PID controller. + return pid(timestep, nl_its, max_dt); + } + + // Consider changing the timestep. + if ((nl_its <= lower_its) && (timestep < max_timestep) && (recent_failure_rate < 0.5)) { + // Increase timestep slightly + timestep *= timestep_factor_on_lower_its; + + return std::min(timestep, max_timestep); + } + + if (nl_its >= upper_its) { + // Reduce timestep slightly + return timestep * timestep_factor_on_upper_its; + } + return timestep; // No change +} + PetscErrorCode SNESSolver::initPseudoTimestepping() { // Storage for per-variable timestep PetscCall(VecDuplicate(snes_x, &dt_vec)); diff --git a/src/solver/impls/snes/snes.hxx b/src/solver/impls/snes/snes.hxx index 0cf9eef4b0..5b8ed14f7c 100644 --- a/src/solver/impls/snes/snes.hxx +++ b/src/solver/impls/snes/snes.hxx @@ -144,6 +144,10 @@ private: BoutReal pseudo_max_ratio; ///< Maximum timestep ratio between neighboring cells Vec dt_vec; ///< Each quantity can have its own timestep + /// Adjust the global timestep + BoutReal updateGlobalTimestep(BoutReal timestep, int nl_its, + BoutReal recent_failure_rate, BoutReal max_dt); + /// Initialize the Pseudo-Transient Continuation method PetscErrorCode initPseudoTimestepping(); /// Update dt_vec based on new solution x From 39b259594ebbe92c719716a79be6d4ff7a04f0ba Mon Sep 17 00:00:00 2001 From: Ben Dudson Date: Wed, 11 Feb 2026 12:39:18 -0800 Subject: [PATCH 04/13] snes: Split residual and pseudo-timestep calculation Now always calculates local and global residuals. If pseudo timestepping is used then another loop over variables is used to update per-cell timesteps and set `dt_vec`. When `diagnose=true`, the local residual will be saved as `snes_local_residual`, and the global residual will be printed after each internal timestep. --- src/solver/impls/snes/snes.cxx | 168 +++++++++++++++++++++++---------- src/solver/impls/snes/snes.hxx | 17 +++- 2 files changed, 133 insertions(+), 52 deletions(-) diff --git a/src/solver/impls/snes/snes.cxx b/src/solver/impls/snes/snes.cxx index 32a9902543..4513cb5aeb 100644 --- a/src/solver/impls/snes/snes.cxx +++ b/src/solver/impls/snes/snes.cxx @@ -670,6 +670,10 @@ int SNESSolver::init() { if (equation_form == BoutSnesEquationForm::pseudo_transient) { PetscCall(initPseudoTimestepping()); } + // Per-cell residuals + local_residual = 0.0; + local_residual_2d = 0.0; + global_residual = 0.0; // Nonlinear solver interface (SNES) output_info.write("Create SNES\n"); @@ -845,6 +849,15 @@ int SNESSolver::run() { PetscCall(VecRestoreArray(snes_x, &xdata)); } + // Initialise residuals + local_residual = 0.0; + local_residual_2d = 0.0; + global_residual = 0.0; + PetscCall(updateResiduals(snes_x)); + if (diagnose) { + output.write("\n Residual: {}\n", global_residual); + } + BoutReal target = simtime; for (int s = 0; s < getNumberOutputSteps(); s++) { target += getOutputTimestep(); @@ -1122,6 +1135,9 @@ int SNESSolver::run() { simtime += dt; + // Update local and global residuals + PetscCall(updateResiduals(snes_x)); + if (diagnose) { // Gather and print diagnostic information @@ -1132,6 +1148,8 @@ int SNESSolver::run() { output.write(", SNES failures: {}", snes_failures); } output.write("\n"); + // Print residual on separate line, so post-processing isn't affected + output.write(" Residual: {}\n", global_residual); } // MatFilter and MatEliminateZeros(Mat, bool) require PETSc >= 3.20 @@ -1146,7 +1164,7 @@ int SNESSolver::run() { pseudo_alpha = updateGlobalTimestep(pseudo_alpha, nl_its, recent_failure_rate, max_timestep * atol * 100); // Adjust local timesteps - PetscCall(updatePseudoTimestepping(snes_x)); + PetscCall(updatePseudoTimestepping()); } else { // Adjust timestep @@ -1255,6 +1273,95 @@ BoutReal SNESSolver::updateGlobalTimestep(BoutReal timestep, int nl_its, return timestep; // No change } +PetscErrorCode SNESSolver::updateResiduals(Vec x) { + // Push residuals to previous residuals + // Copying so that data is not shared + local_residual_prev = copy(local_residual); + local_residual_2d_prev = copy(local_residual_2d); + global_residual_prev = global_residual; + + // Call RHS function to get time derivatives + PetscCall(rhs_function(x, snes_f, false)); + + // Reading the residual vectors + const BoutReal* current_residual = nullptr; + PetscCall(VecGetArrayRead(snes_f, ¤t_residual)); + + // Note: The ordering of quantities in the PETSc vectors + // depends on the Solver::loop_vars function + Mesh* mesh = bout::globals::mesh; + int idx = 0; // Index into PETSc Vecs + + // Boundary cells + for (const auto& i2d : mesh->getRegion2D("RGN_BNDRY")) { + // Field2D quantities evolved together + BoutReal residual = 0.0; + int count = 0; + + for (const auto& f : f2d) { + if (!f.evolve_bndry) { + continue; // Not evolving boundary => Skip + } + residual += SQ(current_residual[idx++]); + ++count; + } + if (count > 0) { + local_residual_2d[i2d] = sqrt(residual / count); + } + + // Field3D quantities evolved together within a cell + for (int jz = mesh->zstart; jz <= mesh->zend; jz++) { + count = 0; + residual = 0.0; + for (const auto& f : f3d) { + if (!f.evolve_bndry) { + continue; // Not evolving boundary => Skip + } + residual += SQ(current_residual[idx++]); + ++count; + } + if (count > 0) { + auto i3d = mesh->ind2Dto3D(i2d, jz); + local_residual[i3d] = sqrt(residual / count); + } + } + } + + // Bulk of domain. + // These loops don't check the boundary flags + for (const auto& i2d : mesh->getRegion2D("RGN_NOBNDRY")) { + // Field2D quantities evolved together + if (!f2d.empty()) { + BoutReal residual = 0.0; + for (std::size_t i = 0; i != f2d.size(); ++i) { + residual += SQ(current_residual[idx++]); + } + local_residual_2d[i2d] = sqrt(residual / static_cast(f2d.size())); + } + + // Field3D quantities evolved together within a cell + if (!f3d.empty()) { + for (int jz = mesh->zstart; jz <= mesh->zend; jz++) { + auto i3d = mesh->ind2Dto3D(i2d, jz); + + BoutReal residual = 0.0; + for (std::size_t i = 0; i != f3d.size(); ++i) { + residual += SQ(current_residual[idx++]); + } + local_residual[i3d] = sqrt(residual / static_cast(f3d.size())); + } + } + } + + // Restore Vec data arrays + PetscCall(VecRestoreArrayRead(snes_f, ¤t_residual)); + + // Global residual metric (RMS) + global_residual = std::sqrt(mean(SQ(local_residual), true)); + + return PETSC_SUCCESS; +} + PetscErrorCode SNESSolver::initPseudoTimestepping() { // Storage for per-variable timestep PetscCall(VecDuplicate(snes_x, &dt_vec)); @@ -1263,22 +1370,14 @@ PetscErrorCode SNESSolver::initPseudoTimestepping() { // Diagnostic outputs pseudo_timestep = timestep; - pseudo_residual = 0.0; - pseudo_residual_2d = 0.0; return PETSC_SUCCESS; } -PetscErrorCode SNESSolver::updatePseudoTimestepping(Vec x) { - // Call RHS function to get time derivatives - PetscCall(rhs_function(x, snes_f, false)); - +PetscErrorCode SNESSolver::updatePseudoTimestepping() { // Use a per-cell timestep so that e.g density and pressure // evolve in a way consistent with the equation of state. - // Reading the residual vectors - const BoutReal* current_residual = nullptr; - PetscCall(VecGetArrayRead(snes_f, ¤t_residual)); // Modifying the dt_vec values BoutReal* dt_data = nullptr; PetscCall(VecGetArray(dt_vec, &dt_data)); @@ -1292,47 +1391,34 @@ PetscErrorCode SNESSolver::updatePseudoTimestepping(Vec x) { for (const auto& i2d : mesh->getRegion2D("RGN_BNDRY")) { // Field2D quantities evolved together int count = 0; - BoutReal residual = 0.0; - for (const auto& f : f2d) { if (!f.evolve_bndry) { continue; // Not evolving boundary => Skip } - residual += SQ(current_residual[idx + count]); ++count; } if (count > 0) { - residual = sqrt(residual / count); - // Adjust timestep for these quantities - BoutReal new_timestep = - updatePseudoTimestep(dt_data[idx], pseudo_residual_2d[i2d], residual); + BoutReal new_timestep = updatePseudoTimestep(dt_data[idx], local_residual_2d_prev[i2d], + local_residual_2d[i2d]); for (int i = 0; i != count; ++i) { dt_data[idx++] = new_timestep; } - pseudo_residual_2d[i2d] = residual; } // Field3D quantities evolved together within a cell for (int jz = mesh->zstart; jz <= mesh->zend; jz++) { count = 0; - residual = 0.0; for (const auto& f : f3d) { if (!f.evolve_bndry) { continue; // Not evolving boundary => Skip } - residual += SQ(current_residual[idx + count]); ++count; } if (count > 0) { - residual = sqrt(residual / count); - auto i3d = mesh->ind2Dto3D(i2d, jz); - const BoutReal new_timestep = - updatePseudoTimestep(dt_data[idx], pseudo_residual[i3d], residual); - - pseudo_residual[i3d] = residual; - pseudo_timestep[i3d] = new_timestep; + const BoutReal new_timestep = updatePseudoTimestep( + dt_data[idx], local_residual_prev[i3d], local_residual[i3d]); for (int i = 0; i != count; ++i) { dt_data[idx++] = new_timestep; @@ -1346,19 +1432,12 @@ PetscErrorCode SNESSolver::updatePseudoTimestepping(Vec x) { for (const auto& i2d : mesh->getRegion2D("RGN_NOBNDRY")) { // Field2D quantities evolved together if (!f2d.empty()) { - BoutReal residual = 0.0; - for (std::size_t i = 0; i != f2d.size(); ++i) { - residual += SQ(current_residual[idx + i]); - } - residual = sqrt(residual / f2d.size()); - // Adjust timestep for these quantities - const BoutReal new_timestep = - updatePseudoTimestep(dt_data[idx], pseudo_residual_2d[i2d], residual); + const BoutReal new_timestep = updatePseudoTimestep( + dt_data[idx], local_residual_2d_prev[i2d], local_residual_2d[i2d]); for (std::size_t i = 0; i != f2d.size(); ++i) { dt_data[idx++] = new_timestep; } - pseudo_residual_2d[i2d] = residual; } // Field3D quantities evolved together within a cell @@ -1366,15 +1445,8 @@ PetscErrorCode SNESSolver::updatePseudoTimestepping(Vec x) { for (int jz = mesh->zstart; jz <= mesh->zend; jz++) { auto i3d = mesh->ind2Dto3D(i2d, jz); - BoutReal residual = 0.0; - for (std::size_t i = 0; i != f3d.size(); ++i) { - residual += SQ(current_residual[idx + i]); - } - residual = sqrt(residual / f3d.size()); - BoutReal new_timestep = - updatePseudoTimestep(dt_data[idx], pseudo_residual[i3d], residual); - - pseudo_residual[i3d] = residual; + BoutReal new_timestep = updatePseudoTimestep( + dt_data[idx], local_residual_prev[i3d], local_residual[i3d]); // Compare to neighbors BoutReal min_neighboring_dt = max_timestep; @@ -1404,7 +1476,6 @@ PetscErrorCode SNESSolver::updatePseudoTimestepping(Vec x) { } // Restore Vec data arrays - PetscCall(VecRestoreArrayRead(snes_f, ¤t_residual)); PetscCall(VecRestoreArray(dt_vec, &dt_data)); // Need timesteps on neighboring processors @@ -1801,11 +1872,12 @@ void SNESSolver::outputVars(Options& output_options, bool save_repeat) { return; // Don't save diagnostics to restart files } + output_options["snes_local_residual"].assignRepeat(local_residual, "t", save_repeat, + "SNESSolver"); + if (equation_form == BoutSnesEquationForm::pseudo_transient) { output_options["snes_pseudo_alpha"].assignRepeat(pseudo_alpha, "t", save_repeat, "SNESSolver"); - output_options["snes_pseudo_residual"].assignRepeat(pseudo_residual, "t", save_repeat, - "SNESSolver"); output_options["snes_pseudo_timestep"].assignRepeat(pseudo_timestep, "t", save_repeat, "SNESSolver"); } diff --git a/src/solver/impls/snes/snes.hxx b/src/solver/impls/snes/snes.hxx index 5b8ed14f7c..b36212ade2 100644 --- a/src/solver/impls/snes/snes.hxx +++ b/src/solver/impls/snes/snes.hxx @@ -148,10 +148,20 @@ private: BoutReal updateGlobalTimestep(BoutReal timestep, int nl_its, BoutReal recent_failure_rate, BoutReal max_dt); + /// Calculate per-cell and global residuals + /// given an input system state `x` + PetscErrorCode updateResiduals(Vec x); + Field3D local_residual; ///< Residual of Field3D quantities in each cell + Field2D local_residual_2d; ///< Residual of Field2D quantities in each cell + BoutReal global_residual; ///< Global residual measure + Field3D local_residual_prev; ///< Previous Field3D local residuals + Field2D local_residual_2d_prev; ///< Previous Field2D local residuals + BoutReal global_residual_prev; ///< Previous global residual + /// Initialize the Pseudo-Transient Continuation method PetscErrorCode initPseudoTimestepping(); - /// Update dt_vec based on new solution x - PetscErrorCode updatePseudoTimestepping(Vec x); + /// Update dt_vec based on residuals + PetscErrorCode updatePseudoTimestepping(); /// Decide the next pseudo-timestep. Called by updatePseudoTimestepping BoutReal updatePseudoTimestep(BoutReal previous_timestep, BoutReal previous_residual, BoutReal current_residual); @@ -160,8 +170,7 @@ private: BoutReal updatePseudoTimestep_history_based(BoutReal previous_timestep, BoutReal previous_residual, BoutReal current_residual); - Field3D pseudo_residual; ///< Diagnostic output - Field2D pseudo_residual_2d; + Field3D pseudo_timestep; ///< PID controller parameters From b2102d058a9a9cbe89c69abc05b30f4744c71715 Mon Sep 17 00:00:00 2001 From: Ben Dudson Date: Wed, 11 Feb 2026 15:04:12 -0800 Subject: [PATCH 05/13] snes: output_trigger option Selects the trigger for writing an output. The default is `fixed_time_interval` that outputs at regular time intervals. Other choice is `residual_ratio` that outputs when the global residual falls by a factor `output_residual_ratio`. --- src/solver/impls/snes/snes.cxx | 93 +++++++++++++++++++++------------- src/solver/impls/snes/snes.hxx | 27 ++++++++-- 2 files changed, 80 insertions(+), 40 deletions(-) diff --git a/src/solver/impls/snes/snes.cxx b/src/solver/impls/snes/snes.cxx index 4513cb5aeb..1766c94186 100644 --- a/src/solver/impls/snes/snes.cxx +++ b/src/solver/impls/snes/snes.cxx @@ -474,6 +474,14 @@ PetscErrorCode SNESSolver::FDJrestoreFromPruning() { SNESSolver::SNESSolver(Options* opts) : Solver(opts), + output_trigger( + (*options)["output_trigger"] + .doc("Decides when to save outputs. fixed_time_interval, residual_ratio") + .withDefault(BoutSnesOutput::fixed_time_interval)), + output_residual_ratio( + (*options)["output_residual_ratio"] + .doc("Trigger an output when residual falls by this ratio") + .withDefault(0.5)), timestep( (*options)["timestep"].doc("Initial backward Euler timestep").withDefault(1.0)), dt_min_reset((*options)["dt_min_reset"] @@ -535,8 +543,13 @@ SNESSolver::SNESSolver(Options* opts) pseudo_max_ratio((*options)["pseudo_max_ratio"] .doc("PTC maximum timestep ratio between neighbors") .withDefault(2.)), - pid_controller( - (*options)["pid_controller"].doc("Use PID controller?").withDefault(false)), + timestep_control((*options)["timestep_control"] + .doc("Timestep control method") + .withDefault(BoutSnesTimestep::pid_nonlinear_its)), + timestep_factor((*options)["timestep_factor"] + .doc("When timestep_control=residual_ratio, multiply timestep " + "by this factor each step") + .withDefault(1.1)), target_its((*options)["target_its"].doc("Target snes iterations").withDefault(7)), kP((*options)["kP"].doc("Proportional PID parameter").withDefault(0.7)), kI((*options)["kI"].doc("Integral PID parameter").withDefault(0.3)), @@ -733,16 +746,10 @@ int SNESSolver::init() { // Note: If the 'Amat' Jacobian is matrix free, SNESComputeJacobian // always updates its reference 'u' vector every nonlinear iteration SNESSetLagJacobian(snes, lag_jacobian); - if (pid_controller) { - nl_its_prev = target_its; - nl_its_prev2 = target_its; - SNESSetLagJacobianPersists(snes, PETSC_FALSE); - SNESSetLagPreconditionerPersists(snes, PETSC_FALSE); - } else { - // Set Jacobian and preconditioner to persist across time steps - SNESSetLagJacobianPersists(snes, PETSC_TRUE); - SNESSetLagPreconditionerPersists(snes, PETSC_TRUE); - } + nl_its_prev = target_its; + nl_its_prev2 = target_its; + SNESSetLagJacobianPersists(snes, PETSC_FALSE); + SNESSetLagPreconditionerPersists(snes, PETSC_FALSE); SNESSetLagPreconditioner(snes, 1); // Rebuild when Jacobian is rebuilt } @@ -859,6 +866,7 @@ int SNESSolver::run() { } BoutReal target = simtime; + recent_failure_rate = 0.0; for (int s = 0; s < getNumberOutputSteps(); s++) { target += getOutputTimestep(); @@ -866,10 +874,14 @@ int SNESSolver::run() { int snes_failures = 0; // Count SNES convergence failures int saved_jacobian_lag = 0; int loop_count = 0; - recent_failure_rate = 0.0; + + BoutReal start_global_residual = global_residual; do { - if (simtime >= target) + if ((output_trigger == BoutSnesOutput::fixed_time_interval && (simtime >= target)) + || (output_trigger == BoutSnesOutput::residual_ratio + && (global_residual <= start_global_residual * output_residual_ratio))) break; // Could happen if step over multiple outputs + if (scale_vars) { // Individual variable scaling // Note: If variables are rescaled then the Jacobian columns @@ -930,8 +942,7 @@ int SNESSolver::run() { PetscCall(VecMin(dt_vec, nullptr, ×tep)); dt = timestep; - looping = true; - if (simtime + timestep >= target) { + if (output_trigger == BoutSnesOutput::fixed_time_interval && simtime + timestep >= target) { looping = false; } } else { @@ -954,8 +965,7 @@ int SNESSolver::run() { // Set the timestep dt = timestep; - looping = true; - if (simtime + dt >= target) { + if (output_trigger == BoutSnesOutput::fixed_time_interval && simtime + dt >= target) { // Note: When the timestep is changed the preconditioner needs to be updated // => Step over the output time and interpolate if not matrix free @@ -975,9 +985,7 @@ int SNESSolver::run() { VecAXPBY(snes_x, -beta, (1. + beta), x1); } - if (pid_controller) { - SNESSetLagJacobian(snes, lag_jacobian); - } + SNESSetLagJacobian(snes, lag_jacobian); } // Run the solver @@ -1182,7 +1190,8 @@ int SNESSolver::run() { snes_failures = 0; } while (looping); - if (!matrix_free) { + BoutReal output_time = simtime; + if (output_trigger == BoutSnesOutput::fixed_time_interval && !matrix_free) { ASSERT2(simtime >= target); ASSERT2(simtime - dt <= target); // Stepped over output timestep => Interpolate @@ -1197,6 +1206,7 @@ int SNESSolver::run() { // output_x <- alpha * x0 + (1 - alpha) * output_x VecAXPBY(output_x, alpha, 1. - alpha, x0); + output_time = target; } else { // Timestep was adjusted to hit target output time @@ -1231,7 +1241,7 @@ int SNESSolver::run() { PetscCall(VecRestoreArrayRead(scaled_x, &xdata)); try { - run_rhs(target); // Run RHS to calculate auxilliary variables + run_rhs(output_time); // Run RHS to calculate auxilliary variables } catch (BoutException& e) { output_error.write("ERROR: BoutException thrown: {}\n", e.what()); // Abort simulation. There is no way to recover unless @@ -1239,7 +1249,7 @@ int SNESSolver::run() { BoutComm::abort(1); } - if (call_monitors(target, s, getNumberOutputSteps()) != 0) { + if (call_monitors(output_time, s, getNumberOutputSteps()) != 0) { break; // User signalled to quit } } @@ -1253,22 +1263,35 @@ BoutReal SNESSolver::updateGlobalTimestep(BoutReal timestep, int nl_its, // so if it is not recalculated the it will be less // effective. - if (pid_controller) { + switch (timestep_control) { + case BoutSnesTimestep::pid_nonlinear_its: // Changing the timestep using a PID controller. return pid(timestep, nl_its, max_dt); - } - // Consider changing the timestep. - if ((nl_its <= lower_its) && (timestep < max_timestep) && (recent_failure_rate < 0.5)) { - // Increase timestep slightly - timestep *= timestep_factor_on_lower_its; + case BoutSnesTimestep::threshold_nonlinear_its: + // Consider changing the timestep, based on thresholds in NL iterations + if ((nl_its <= lower_its) && (timestep < max_timestep) && (recent_failure_rate < 0.5)) { + // Increase timestep slightly + timestep *= timestep_factor_on_lower_its; + return std::min(timestep, max_timestep); + } - return std::min(timestep, max_timestep); - } + if (nl_its >= upper_its) { + // Reduce timestep slightly + return timestep * timestep_factor_on_upper_its; + } + return timestep; // No change - if (nl_its >= upper_its) { - // Reduce timestep slightly - return timestep * timestep_factor_on_upper_its; + case BoutSnesTimestep::residual_ratio: + // Use ratio of previous and current global residual + // Intended to be the same as https://petsc.org/release/manualpages/TS/TSPSEUDO/ + // (Note the PETSc manual has the expression for dt_n upside down) + + return std::min({timestep_factor * timestep * global_residual_prev / global_residual, + max_timestep}); + + case BoutSnesTimestep::fixed: + break; } return timestep; // No change } diff --git a/src/solver/impls/snes/snes.hxx b/src/solver/impls/snes/snes.hxx index b36212ade2..768ba33520 100644 --- a/src/solver/impls/snes/snes.hxx +++ b/src/solver/impls/snes/snes.hxx @@ -4,7 +4,7 @@ * using PETSc for the SNES interface * ************************************************************************** - * Copyright 2015-2025 BOUT++ contributors + * Copyright 2015-2026 BOUT++ contributors * * Contact: Ben Dudson, dudson2@llnl.gov * @@ -59,6 +59,16 @@ BOUT_ENUM_CLASS(BoutPTCStrategy, history_based, ///< Grow/shrink dt based on residual decrease/increase hybrid); ///< Combine inverse_residual and history_based strategies +BOUT_ENUM_CLASS(BoutSnesTimestep, + pid_nonlinear_its, ///< PID controller on nonlinear iterations + threshold_nonlinear_its, ///< Use thresholds on nonlinear iterations + residual_ratio, ///< Use ratio of previous and current residual + fixed); ///< Fixed timestep (no adaptation) + +BOUT_ENUM_CLASS(BoutSnesOutput, + fixed_time_interval, ///< Output at fixed time intervals + residual_ratio); ///< When the residual is reduced by a given ratio + /// Uses PETSc's SNES interface to find a steady state solution to a /// nonlinear ODE by integrating in time with Backward Euler class SNESSolver : public Solver { @@ -112,6 +122,10 @@ private: /// @param[in] linear Specifies that the SNES solver is in a linear (KSP) inner loop PetscErrorCode rhs_function(Vec x, Vec f, bool linear); + BoutSnesOutput output_trigger; ///< Sets when outputs are written + + BoutReal output_residual_ratio; ///< Trigger an output when residual falls by this ratio + BoutReal timestep; ///< Internal timestep BoutReal dt; ///< Current timestep used in snes_function. BoutReal dt_min_reset; ///< If dt falls below this, reset solve @@ -173,8 +187,11 @@ private: Field3D pseudo_timestep; - ///< PID controller parameters - bool pid_controller; ///< Use PID controller? + /// Timestep controller method + BoutSnesTimestep timestep_control; + + /// When using BoutSnesTimestep::residual_ratio + BoutReal timestep_factor; ///< Multiply timestep by this each step /// Target number of nonlinear iterations for the PID controller. /// This can be non-integer to push timestep more aggressively @@ -187,8 +204,8 @@ private: BoutReal recent_failure_rate; ///< Rolling average of recent failure rate BoutReal last_failure_weight; ///< 1 / number of recent solves - int nl_its_prev; - int nl_its_prev2; + BoutReal nl_its_prev; + BoutReal nl_its_prev2; BoutReal pid(BoutReal timestep, int nl_its, BoutReal max_dt); ///< Updates the timestep From 4e134e9ab0080cf859802f31f893a7cb51ded0e2 Mon Sep 17 00:00:00 2001 From: Ben Dudson Date: Wed, 11 Feb 2026 15:33:23 -0800 Subject: [PATCH 06/13] snes: Update manual Describe new settings `timestep_control` and `output_trigger` that select different strategies to adjust global timesteps and output simulation state. --- manual/sphinx/user_docs/time_integration.rst | 70 ++++++++++++++++---- 1 file changed, 57 insertions(+), 13 deletions(-) diff --git a/manual/sphinx/user_docs/time_integration.rst b/manual/sphinx/user_docs/time_integration.rst index 7658602bc8..70b9c3df92 100644 --- a/manual/sphinx/user_docs/time_integration.rst +++ b/manual/sphinx/user_docs/time_integration.rst @@ -444,29 +444,73 @@ on nonlinear iteration count. dt_min_reset = 1e-6 # Reset the solver when timestep < this # Timestep adaptation + timestep_control = pid_nonlinear_its + target_its = 7 # Target number of nonlinear iterations + kP = 0.7 # Proportional gain + kI = 0.3 # Integral gain + kD = 0.2 # Derivative gain + +This uses a PID controller that adjusts the timestep to maintain approximately ``target_its`` +nonlinear iterations per solve. + +Residual Ratio +^^^^^^^^^^^^^^ + +This adjusts the timestep using the ratio of global residuals and a timestep factor: + +.. math:: + + dt_n = r dt_{n-1} \frac{||F(X_{n-1})||}{||F(X_{n})|| + +so that as the residual falls the timestep :math:`dt` is increased. +The :math:`r` parameter is input option ``timestep_factor`` that has +default value 1.1. + +.. code-block:: ini + + [solver] + timestep_control = residual_ratio # Use global residual + timestep_factor = 1.1 # Constant timestep factor + +Threshold Controller +^^^^^^^^^^^^^^^^^^^^ + +An alternative adaptive strategy uses thresholds in nonlinear iterations +to adjust the timestep: + +.. code-block:: ini + + [solver] + timestep_control = threshold_nonlinear_its lower_its = 3 # Increase dt if iterations < this upper_its = 10 # Decrease dt if iterations > this timestep_factor_on_lower_its = 1.4 # Growth factor timestep_factor_on_upper_its = 0.9 # Reduction factor timestep_factor_on_failure = 0.5 # Reduction on convergence failure -PID Controller -^^^^^^^^^^^^^^ +The adjustments are less smooth than the default PID method, but the +timestep is changed less frequently. This may enable the Jacobian and +preconditioner to be used for more iterations. -An alternative adaptive strategy using a PID controller: +Output trigger +~~~~~~~~~~~~~~ + +The default behavior is to save outputs at a regular time interval, as +BOUT++ solvers do. This is desirable when performing time-dependent +simulations, but for simulations that are trying to get to steady +state a better measure of progress is reduction of the global residual +(norm of the time-derivatives of the system). .. code-block:: ini [solver] - pid_controller = true - target_its = 7 # Target number of nonlinear iterations - kP = 0.7 # Proportional gain - kI = 0.3 # Integral gain - kD = 0.2 # Derivative gain + output_trigger = residual_ratio # Trigger an output based on the ratio of residuals + output_residual_ratio = 0.5 # Output when global residual is multiplied by this -The PID controller adjusts the timestep to maintain approximately ``target_its`` -nonlinear iterations per solve, providing smoother adaptation than threshold-based -methods. +With this choice, each output has a global residual that is less than +or equal to `output_residual_ratio` times the last output global +residual. This provides a way of measuring progress to steady state +that is independent of time integration accuracy. Pseudo-Transient Continuation and Switched Evolution Relaxation ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ @@ -518,7 +562,7 @@ state. timestep = 1.0 # Initial timestep # SER parameters - pid_controller = true # Scale timesteps based on iterations + timestep_control = pid_controller # Scale timesteps based on iterations pseudo_max_ratio = 2.0 # Limit neighbor timestep ratio # Tolerances @@ -588,7 +632,7 @@ adjust ``pseudo_alpha`` depending on the nonlinearity of the system: .. code-block:: ini [solver] - pid_controller = true + timestep_control = pid_controller # Scale global timestep using PID controller target_its = 7 # Target number of nonlinear iterations kP = 0.7 # Proportional gain kI = 0.3 # Integral gain From e5bb430240546531d2180e6bbad01543bd0809da Mon Sep 17 00:00:00 2001 From: Ben Dudson Date: Wed, 11 Feb 2026 15:39:06 -0800 Subject: [PATCH 07/13] snes: Add snes_global_residual scalar to dmp files Save the global residual diagnostic --- src/solver/impls/snes/snes.cxx | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/solver/impls/snes/snes.cxx b/src/solver/impls/snes/snes.cxx index 1766c94186..75798bfca5 100644 --- a/src/solver/impls/snes/snes.cxx +++ b/src/solver/impls/snes/snes.cxx @@ -1897,6 +1897,8 @@ void SNESSolver::outputVars(Options& output_options, bool save_repeat) { output_options["snes_local_residual"].assignRepeat(local_residual, "t", save_repeat, "SNESSolver"); + output_options["snes_global_residual"].assignRepeat(global_residual, "t", save_repeat, + "SNESSolver"); if (equation_form == BoutSnesEquationForm::pseudo_transient) { output_options["snes_pseudo_alpha"].assignRepeat(pseudo_alpha, "t", save_repeat, From b9b038b0243d3a54e8f89448fb0eb6664f7821d9 Mon Sep 17 00:00:00 2001 From: Ben Dudson Date: Wed, 11 Feb 2026 16:11:01 -0800 Subject: [PATCH 08/13] snes: Fix manual Uses setting `pid_controller` rather than `pid_nonlinear_its`. --- manual/sphinx/user_docs/time_integration.rst | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/manual/sphinx/user_docs/time_integration.rst b/manual/sphinx/user_docs/time_integration.rst index 70b9c3df92..37bf9a6abf 100644 --- a/manual/sphinx/user_docs/time_integration.rst +++ b/manual/sphinx/user_docs/time_integration.rst @@ -562,7 +562,7 @@ state. timestep = 1.0 # Initial timestep # SER parameters - timestep_control = pid_controller # Scale timesteps based on iterations + timestep_control = pid_nonlinear_its # Scale timesteps based on iterations pseudo_max_ratio = 2.0 # Limit neighbor timestep ratio # Tolerances @@ -632,7 +632,7 @@ adjust ``pseudo_alpha`` depending on the nonlinearity of the system: .. code-block:: ini [solver] - timestep_control = pid_controller # Scale global timestep using PID controller + timestep_control = pid_nonlinear_its # Scale global timestep using PID controller target_its = 7 # Target number of nonlinear iterations kP = 0.7 # Proportional gain kI = 0.3 # Integral gain From 2c2402ed59c91164eaff46dee0f79386b7347e9e Mon Sep 17 00:00:00 2001 From: bendudson <219233+bendudson@users.noreply.github.com> Date: Thu, 12 Feb 2026 00:14:36 +0000 Subject: [PATCH 09/13] [bot] Apply format changes --- src/solver/impls/snes/snes.cxx | 19 ++++++++++++------- 1 file changed, 12 insertions(+), 7 deletions(-) diff --git a/src/solver/impls/snes/snes.cxx b/src/solver/impls/snes/snes.cxx index 75798bfca5..420adaf046 100644 --- a/src/solver/impls/snes/snes.cxx +++ b/src/solver/impls/snes/snes.cxx @@ -942,7 +942,8 @@ int SNESSolver::run() { PetscCall(VecMin(dt_vec, nullptr, ×tep)); dt = timestep; - if (output_trigger == BoutSnesOutput::fixed_time_interval && simtime + timestep >= target) { + if (output_trigger == BoutSnesOutput::fixed_time_interval + && simtime + timestep >= target) { looping = false; } } else { @@ -965,7 +966,8 @@ int SNESSolver::run() { // Set the timestep dt = timestep; - if (output_trigger == BoutSnesOutput::fixed_time_interval && simtime + dt >= target) { + if (output_trigger == BoutSnesOutput::fixed_time_interval + && simtime + dt >= target) { // Note: When the timestep is changed the preconditioner needs to be updated // => Step over the output time and interpolate if not matrix free @@ -1169,14 +1171,16 @@ int SNESSolver::run() { if (equation_form == BoutSnesEquationForm::pseudo_transient) { // Adjust pseudo_alpha to globally scale timesteps - pseudo_alpha = updateGlobalTimestep(pseudo_alpha, nl_its, recent_failure_rate, max_timestep * atol * 100); + pseudo_alpha = updateGlobalTimestep(pseudo_alpha, nl_its, recent_failure_rate, + max_timestep * atol * 100); // Adjust local timesteps PetscCall(updatePseudoTimestepping()); } else { // Adjust timestep - timestep = updateGlobalTimestep(timestep, nl_its, recent_failure_rate, max_timestep); + timestep = + updateGlobalTimestep(timestep, nl_its, recent_failure_rate, max_timestep); } if (static_cast(lin_its) / nl_its > 0.5 * maxl) { @@ -1270,7 +1274,8 @@ BoutReal SNESSolver::updateGlobalTimestep(BoutReal timestep, int nl_its, case BoutSnesTimestep::threshold_nonlinear_its: // Consider changing the timestep, based on thresholds in NL iterations - if ((nl_its <= lower_its) && (timestep < max_timestep) && (recent_failure_rate < 0.5)) { + if ((nl_its <= lower_its) && (timestep < max_timestep) + && (recent_failure_rate < 0.5)) { // Increase timestep slightly timestep *= timestep_factor_on_lower_its; return std::min(timestep, max_timestep); @@ -1422,8 +1427,8 @@ PetscErrorCode SNESSolver::updatePseudoTimestepping() { } if (count > 0) { // Adjust timestep for these quantities - BoutReal new_timestep = updatePseudoTimestep(dt_data[idx], local_residual_2d_prev[i2d], - local_residual_2d[i2d]); + BoutReal new_timestep = updatePseudoTimestep( + dt_data[idx], local_residual_2d_prev[i2d], local_residual_2d[i2d]); for (int i = 0; i != count; ++i) { dt_data[idx++] = new_timestep; } From d5ad38cb0506762b72ae513ea26341afcc51a61b Mon Sep 17 00:00:00 2001 From: bendudson <219233+bendudson@users.noreply.github.com> Date: Thu, 12 Feb 2026 00:14:37 +0000 Subject: [PATCH 10/13] [bot] Add last format changes commit to ignore file --- .git-blame-ignore-revs | 1 + 1 file changed, 1 insertion(+) diff --git a/.git-blame-ignore-revs b/.git-blame-ignore-revs index a8bd688f08..3c401e4195 100644 --- a/.git-blame-ignore-revs +++ b/.git-blame-ignore-revs @@ -4,3 +4,4 @@ ed2117e6d6826a98b6988e2f18c0c34e408563b6 # Added by the bot 4b010b7634aee1045743be80c268d4644522cd29 a71cad2dd6ace5741a754e2ca7daacd4bb094e0e +2c2402ed59c91164eaff46dee0f79386b7347e9e From 74ce9aad802cf355ce7fdabab55e7c90cd0f5c9b Mon Sep 17 00:00:00 2001 From: Ben Dudson Date: Thu, 12 Feb 2026 10:59:13 -0800 Subject: [PATCH 11/13] snes: jacobian_persists option False by default. If true, re-uses the Jacobian and preconditioner across SNES calls. --- src/solver/impls/snes/snes.cxx | 31 ++++++++++++++++++++----------- src/solver/impls/snes/snes.hxx | 1 + 2 files changed, 21 insertions(+), 11 deletions(-) diff --git a/src/solver/impls/snes/snes.cxx b/src/solver/impls/snes/snes.cxx index 420adaf046..a05a19b02c 100644 --- a/src/solver/impls/snes/snes.cxx +++ b/src/solver/impls/snes/snes.cxx @@ -256,8 +256,10 @@ PetscErrorCode SNESSolver::FDJinitialise() { for (int i = 0; i < nlocal; ++i) { // Assume all elements in the z direction are potentially coupled - d_nnz.emplace_back(d_nnz_map3d[i].size() * mesh->LocalNz + d_nnz_map2d[i].size()); - o_nnz.emplace_back(o_nnz_map3d[i].size() * mesh->LocalNz + o_nnz_map2d[i].size()); + d_nnz.emplace_back((d_nnz_map3d[i].size() * mesh->LocalNz) + + d_nnz_map2d[i].size()); + o_nnz.emplace_back((o_nnz_map3d[i].size() * mesh->LocalNz) + + o_nnz_map2d[i].size()); } } @@ -369,7 +371,7 @@ PetscErrorCode SNESSolver::FDJinitialise() { // Values are 0 or 1 so tolerance (1e-5) shouldn't matter PetscBool symmetric; PetscCall(MatIsSymmetric(Jfd, 1e-5, &symmetric)); - if (!symmetric) { + if (!static_cast(symmetric)) { output_warn.write("Jacobian pattern is not symmetric\n"); } } @@ -597,6 +599,10 @@ SNESSolver::SNESSolver(Options* opts) lag_jacobian((*options)["lag_jacobian"] .doc("Re-use the Jacobian this number of SNES iterations") .withDefault(50)), + jacobian_persists( + (*options)["jacobian_persists"] + .doc("Re-use Jacobian and preconditioner across nonlinear solves") + .withDefault(false)), use_coloring((*options)["use_coloring"] .doc("Use matrix coloring to calculate Jacobian?") .withDefault(true)), @@ -630,7 +636,7 @@ int SNESSolver::init() { // Get total problem size int ntmp; if (bout::globals::mpi->MPI_Allreduce(&nlocal, &ntmp, 1, MPI_INT, MPI_SUM, - BoutComm::get())) { + BoutComm::get()) != 0) { throw BoutException("MPI_Allreduce failed!"); } neq = ntmp; @@ -748,9 +754,11 @@ int SNESSolver::init() { SNESSetLagJacobian(snes, lag_jacobian); nl_its_prev = target_its; nl_its_prev2 = target_its; - SNESSetLagJacobianPersists(snes, PETSC_FALSE); - SNESSetLagPreconditionerPersists(snes, PETSC_FALSE); - SNESSetLagPreconditioner(snes, 1); // Rebuild when Jacobian is rebuilt + PetscCall( + SNESSetLagJacobianPersists(snes, static_cast(jacobian_persists))); + PetscCall(SNESSetLagPreconditionerPersists( + snes, static_cast(jacobian_persists))); + PetscCall(SNESSetLagPreconditioner(snes, 1)); // Rebuild when Jacobian is rebuilt } // Set tolerances @@ -826,10 +834,10 @@ int SNESSolver::init() { SNESType snestype; SNESGetType(snes, &snestype); output_info.write("SNES Type : {}\n", snestype); - if (ksptype) { + if (ksptype != nullptr) { output_info.write("KSP Type : {}\n", ksptype); } - if (pctype) { + if (pctype != nullptr) { output_info.write("PC Type : {}\n", pctype); } } @@ -879,8 +887,9 @@ int SNESSolver::run() { do { if ((output_trigger == BoutSnesOutput::fixed_time_interval && (simtime >= target)) || (output_trigger == BoutSnesOutput::residual_ratio - && (global_residual <= start_global_residual * output_residual_ratio))) + && (global_residual <= start_global_residual * output_residual_ratio))) { break; // Could happen if step over multiple outputs + } if (scale_vars) { // Individual variable scaling @@ -1033,7 +1042,7 @@ int SNESSolver::run() { } if (snes_failures == max_snes_failures) { - output.write("Too many SNES failures ({}). Aborting."); + output.write("Too many SNES failures ({}). Aborting.", snes_failures); return 1; } diff --git a/src/solver/impls/snes/snes.hxx b/src/solver/impls/snes/snes.hxx index 768ba33520..ca0ca50d36 100644 --- a/src/solver/impls/snes/snes.hxx +++ b/src/solver/impls/snes/snes.hxx @@ -243,6 +243,7 @@ private: bool matrix_free; ///< Use matrix free Jacobian bool matrix_free_operator; ///< Use matrix free Jacobian in the operator? int lag_jacobian; ///< Re-use Jacobian + bool jacobian_persists; ///< Re-use Jacobian and preconditioner across nonlinear solves bool use_coloring; ///< Use matrix coloring bool jacobian_recalculated; ///< Flag set when Jacobian is recalculated From 05b7c571544c3bcb153fce67d12b9ac48947fc2d Mon Sep 17 00:00:00 2001 From: bendudson <219233+bendudson@users.noreply.github.com> Date: Thu, 12 Feb 2026 19:02:18 +0000 Subject: [PATCH 12/13] [bot] Apply format changes --- src/solver/impls/snes/snes.cxx | 3 ++- src/solver/impls/snes/snes.hxx | 2 +- 2 files changed, 3 insertions(+), 2 deletions(-) diff --git a/src/solver/impls/snes/snes.cxx b/src/solver/impls/snes/snes.cxx index a05a19b02c..2e4744eb44 100644 --- a/src/solver/impls/snes/snes.cxx +++ b/src/solver/impls/snes/snes.cxx @@ -636,7 +636,8 @@ int SNESSolver::init() { // Get total problem size int ntmp; if (bout::globals::mpi->MPI_Allreduce(&nlocal, &ntmp, 1, MPI_INT, MPI_SUM, - BoutComm::get()) != 0) { + BoutComm::get()) + != 0) { throw BoutException("MPI_Allreduce failed!"); } neq = ntmp; diff --git a/src/solver/impls/snes/snes.hxx b/src/solver/impls/snes/snes.hxx index ca0ca50d36..40412f83b7 100644 --- a/src/solver/impls/snes/snes.hxx +++ b/src/solver/impls/snes/snes.hxx @@ -243,7 +243,7 @@ private: bool matrix_free; ///< Use matrix free Jacobian bool matrix_free_operator; ///< Use matrix free Jacobian in the operator? int lag_jacobian; ///< Re-use Jacobian - bool jacobian_persists; ///< Re-use Jacobian and preconditioner across nonlinear solves + bool jacobian_persists; ///< Re-use Jacobian and preconditioner across nonlinear solves bool use_coloring; ///< Use matrix coloring bool jacobian_recalculated; ///< Flag set when Jacobian is recalculated From 2d51378948980d1bdce739d9dd7112a77f1c137f Mon Sep 17 00:00:00 2001 From: bendudson <219233+bendudson@users.noreply.github.com> Date: Thu, 12 Feb 2026 19:02:20 +0000 Subject: [PATCH 13/13] [bot] Add last format changes commit to ignore file --- .git-blame-ignore-revs | 1 + 1 file changed, 1 insertion(+) diff --git a/.git-blame-ignore-revs b/.git-blame-ignore-revs index 3c401e4195..93c768a148 100644 --- a/.git-blame-ignore-revs +++ b/.git-blame-ignore-revs @@ -5,3 +5,4 @@ ed2117e6d6826a98b6988e2f18c0c34e408563b6 4b010b7634aee1045743be80c268d4644522cd29 a71cad2dd6ace5741a754e2ca7daacd4bb094e0e 2c2402ed59c91164eaff46dee0f79386b7347e9e +05b7c571544c3bcb153fce67d12b9ac48947fc2d