Numerical issues close to steady-state in ODE mode with default nonlinear solver

[Francesco Casella]

Please check the PowerGrids.IPSLValidation.Examples.SmallCase.OpenLoopLoadStepPowerGrids model. One can monitor the status of the simulation by looking at gen_pwGenerator4W__GEN______SM.CePu which is the electrical torque in per unit. The model starts in steady-state, then a step change is applied to the load at time = 1. What should ensue is a damped oscillation which eventually dies out to a new steady state.

After #5049 was fixed, the model run successfully in ODE mode, so the -daeMode=new custom annotation was removed from it.

However, if one simulates OpenLoopLoadStepPowerGrids with the standard solver, StopTime = 400 and -noEquidistantTimeGrid, in order to see the actual solver steps in the solution, it is apparent how the step adaptation mechanism is not working correctly. The state variables of the system completely stop oscillating around time = 40, and reach a steady state around time = 120. One would expect a high-order implicit solver such as DASSL to increase the step size significantly between time = 40 and time = 120, and then to take increasingly long step sizes as the solution has reached the steady-state.

This in fact does not happen at all, as shown in the attached step-sizes-standard-solver.png figure. In fact, the solver starts increasing the step size around time = 80, but then goes back to short time steps around 0.01 seconds. Later on, it manages to increase the step size to 10-20 s between time = 150 and time = 220, but then goes back to shorter time steps.

The culprit is the default nonlinear solver, which causes numerous convergence failures (1480 for a 400 s simulation). If kinsol is selected instead, there are 0 convergence failures and the step sizes increases as expected once the solution gets close to steady-state, see step-sizes-kinsol.png.

By comparing the two step-size diagrams, it is apparent how the standard solver has issues as soon as the step size goes above 0.1 (which happens around time = 20); when this happens, sooner or later there are solver failures which are only recovered when the ODE solver step size is reduced below 0.01 s.

Although the test case under consideration contains non linear equations that need to be solved to causalize the system, all the involved variables have very small variations after time = 20, so that the system is essentially working in a linear transient regime. Thus, there should be no convergence issues at all (which is in fact what happens when kinsol is used).

The fact that failures are likely to occur when the step-size is increased is quite difficult to explain, as the solution is close to steady-state, so there is no dependency between the step length and the solution of the system required to get the right-hand side of the ODEs f(x,t).

One possible connection is the extrapolation routine that computes the initial guess for the solver; however, @wbraun notices that we are using the same routine for all solvers (including kinsol), so there is probably some other issue with the default nonlinear solver, unles the issue is in how the interpolation routine is called. This may possibly involve scaling, as the failing system of equations has unknowns spanning 5 orders of magnitude.

[Vitalij Ruge]

It would be interesting whether the tolerance has some impact or the stop criterion inside the nonlinear solver

[Francesco Casella]

From what I understand, it doesn't.

[Francesco Casella]

Rescheduled to 1.14.0 after 1.13.0 releasee

[Francesco Casella]

Releasing 1.14.0 which is stable and has many improvements w.r.t. 1.13.2. This issue is rescheduled to 1.16.0

[Francesco Casella]

Retargeted to 1.17.0 after 1.16.0 release

[Francesco Casella]

Retargeted to 1.18.0 because of 1.17.0 timed release.

[Francesco Casella]

Ticket retargeted after milestone closed