Thank you, Barry. I got a working code following your suggestion. See below. -Ling
From: Barry Smith <[email protected]> Date: Wednesday, October 1, 2025 at 11:36 AM To: Zou, Ling <[email protected]> Cc: PETSc <[email protected]> Subject: Re: [petsc-users] Proper way for exception handling No. Your callback function should not be returning a PETSC_ERR_NOT_CONVERGED since your function doesn't know anything about the status of the nonlinear solver and is not in the business of deciding if SNES is converging or not. What are you ZjQcmQRYFpfptBannerStart This Message Is From an External Sender This message came from outside your organization. ZjQcmQRYFpfptBannerEnd No. Your callback function should not be returning a PETSC_ERR_NOT_CONVERGED since your function doesn't know anything about the status of the nonlinear solver and is not in the business of deciding if SNES is converging or not. • You are right. I implemented the code as proposed in the email. It did not go the way I expected. What are you trying to convey back to PETSc if your function is not "all good"? There are generally two possibilities, 1) just give up and stop the entire program immediately, then return PETSC_ERR_USER <- Not what I am looking for. 2) indicate that SNES solve is asking your function to be evaluated at a point that is not in the domain of your function. For example, say your function is the square root and x is -1. In this situation, depending on the exact nonlinear solver being used SNES may be able to continue to try to solve the nonlinear system by changing its x value. To indicate this call SNESSetFunctionDomainError(snes); and then do a usual return PETSC_SUCCESS; If SNES can continue it will, if it cannot then IT will generate a negative SNESConvergedReason indicating that SNES did not converge. <- Yes. This is what I am looking for. To provide some context, it is a flow problem solver that solves for pressure (p) and temperature (T) nonlinear variables. Density and other properties are evaluated based on p and T. In some case, p and T go out of the physical domain, e.g., negative value, equation of state package will throw an exception, instead of error out, because the code wants to try, for example, a smaller time step size. Following your suggestion, the following code works now, thank you! ================================================================ PetscErrorCode SNESFormFunction(SNES, Vec, Vec, void*); double my_function(); ================================================================ double my_function() { if (all_good) // e.g., pressure > 0 return 1; else // e.g., pressure < 0 { throw 199; return 0; } } PetscErrorCode SNESFormFunction(SNES snes, Vec u, Vec r, void* ctx) { Try { my_value = my_function(); } Catch (int err) { SNESSetFunctionDomainError(snes); return PETSC_SUCCESS; } // compute residuals return PETSC_SUCCESS; } int main(int argc, char **argv) { Initialize_PETSc(); double dt = 1, dt_min = 0.001; while (dt > dt_min) { SNESSolve(AppCtx.snes, NULL, AppCtx.u); SNESGetConvergedReason(AppCtx.snes, &(AppCtx.snes_converged_reason)); if (not_converged) dt *= 0.5; } } PetscFinalize(); } ================================================================ Barry On Oct 1, 2025, at 11:46 AM, Zou, Ling via petsc-users <[email protected]> wrote: Although I haven’t tried yet. Does it make sense and should it work to change the code this way? ================================================================ PetscErrorCode SNESFormFunction(SNES, Vec, Vec, void*); double my_function(); ================================================================ double my_function() { if (all_good) return 1; else { throw 199; return 0; } } PetscErrorCode SNESFormFunction(SNES snes, Vec u, Vec r, void* ctx) { Try { my_value = my_function(); } Catch (int err) { return PETSC_ERR_NOT_CONVERGED; } // compute residuals return PETSC_SUCCESS; } int main(int argc, char **argv) { Initialize_PETSc(); double dt = 1, dt_min = 0.001; while (dt > dt_min) { SNESSolve(AppCtx.snes, NULL, AppCtx.u); SNESGetConvergedReason(AppCtx.snes, &(AppCtx.snes_converged_reason)); if (not_converged) dt *= 0.5; } } PetscFinalize(); } ================================================================ From: petsc-users <[email protected]<mailto:[email protected]>> on behalf of Zou, Ling via petsc-users <[email protected]<mailto:[email protected]>> Date: Wednesday, October 1, 2025 at 9:33 AM To: PETSc <[email protected]<mailto:[email protected]>> Subject: [petsc-users] Proper way for exception handling Hi, After updating to PETSc 3.23 (from a quite old version, ~3.8), I found that my old way of exception handling not working any more, on my Mac. I would like to learn the proper way to handle exceptions in PETSc code. Here is my pseudo code: ================================================================ PetscErrorCode SNESFormFunction(SNES, Vec, Vec, void*); double my_function(); ================================================================ double my_function() { if (all_good) return 1; else { throw 199; return 0; } } PetscErrorCode SNESFormFunction(SNES snes, Vec u, Vec r, void* ctx) { double my_value = my_function(); // compute residuals } int main(int argc, char **argv) { Initialize_PETSc(); double dt = 1, dt_min = 0.001; while (dt > dt_min) { try { SNESSolve(AppCtx.snes, NULL, AppCtx.u); } catch (int err) { dt *= 0.5; } } PetscFinalize(); } ================================================================ This piece of logic used to work well, but now giving me the following error: Current time (the starting time of this time step) = 0.. NL step = 0, SNES Function norm = 8.60984E+03 libc++abi: terminating due to uncaught exception of type int Abort trap: 6 Q1: why this exception catch logic not working any more? Q2: is there any good example of PETSc exception handling I can follow? Best, -Ling
