TIME-X Test Hackathon: Finish test for `SwitchEstimator`

[lisawim]

During the TIME-X test hackathon at TU Darmstadt, I started with writing a test for SwitchEstimator. A previous version of the test is available here. The file provides a combination of unit-tests and integration tests (?).

1. Here, it is checked when the event does occur exactly at the boundary (i.e., on a collocation node at the boundary) whether a restart is executed or not. When the event does occur on the boundary no restart should be executed since here the resolution of the numerical solution for a specific accuracy is already given.

2. Another test checks whether the adapt_interpolation_info does what it is supposed to do. This method is intended to adapt the interpolation x-axis t_interp with corresponding y-values state_function. Using quad_type='RADAU-RIGHT' leads to len(t_interp) = self.coll.num_nodes while len(state_function) = self.coll.num_nodes + 1. The reason is that the initial time L.time is not included here because it is no collocation node but can be used for the interpolation anyway. Hence, it should be added to t_interp in adapt_interpolation_info. For quad_type='LOBATTO', we have len(t_interp) = self.coll.num_nodes + 1 but len(state_function) = self.coll.num_nodes. Here, the initial time of the subinterval is also a collocation nodes, hence t_interp[0] = t_interp[1] and the first entry has to be removed.

   • [ ] quad_type='RADAU-RIGHT': After adding L.time to t_interp test assert t_interp[0] == L.time.

3. In a third test SwitchEstimator is applied to a problem for different tolerances tol passed to the switch estimator. This tolerance is used to decide whether the event is already close to the boundary or if the value of the state function at the event does already satisfy the tolerance. Intuitively, the smaller tol is the smaller the event error could be expected.

For all these tests, a specific dummy problem ExactDiscontinuousTestODE is implemented. This class inherits from DiscontinuousTestODE is does nothing than returning the exact dynamics of the parent class:

• eval_f returns the derivative of the exact solution u_exact,
• solve_system returns the exact solution.

Since this class returns the exact dynamics, for switch estimation the event should be estimated in a very accurately way, i.e., the event to be found can be assumed to be very close to the exact event time t_switch_exact. Discussion with @brownbaerchen led to this idea. The question is: Does this dummy problem make sense to test the correct framework of SwitchEstimator?

For the test, also piecewise linear interpolation for SwitchEstimator is implemented, and first observations does indicate a very small rate of restarts as compared to Lagrange interpolation. This could be observed for number of collocation nodes $M=3,4,5$, but seems to be not true for $M=2$.

Final question(s) to this issue:

• Are tests sensitive enough to throw errors if the switch estimator does not work correctly anymore after an update?
• Do the tests provide an answer to the question of whether the switch estimator works well in a specifically intended sense?

[brownbaerchen]

The test I was suggesting is slightly different, actually. One source of error in the switch estimator is the polynomial interpolation. But you can remove this source of error when the interpolation is exact. That is the idea behind the polynomial_testequation. Maybe you can modify this class such that you can force it to have a root within some interval. Then, you can setup a controller, but instead of doing the whole time-stepping, you can just execute the specific parts of the switch estimator that you want to test. I do a similar thing here. Now, the switch estimator has to find the "switch event" i.e. the root of the state function exactly. Does this make sense? Of course, a root of the state function is not the same as a root of the right hand side. So maybe you need to write a completely different test problem. But to my understanding, the switch estimator would be well tested if you verify that for exact interpolation, the correct step size is chosen to very high accuracy and if you switch to Thibaut's interpolation, which is already tested. Ofc, if you do some other kind of interpolation, you would have to test that separately.

[lisawim]

Note for myself: The choice of the finite difference (FD) approximation to approximate the derivative of the state function inherently influences the performance of the switch estimator. Battery models do not like centered FDs but backward ones, whereas for the discontinuous test problems centered FDs seem to be a good choice. It seems like that it does make a difference if the state function changes its sign from positive to negative or vice versa. The choice of the step size inside the FD is important as well.

UPDATE: In general, backward FDs works better for the numerical solutions, whereas for the classes with exact dynamics in tests, centered FDs provide good results.

[brownbaerchen]

If I understand correctly, you are using the derivative of the state function in a Newton solver to find roots. But the state function is not smooth, is it? So if you are using a backward stencil you are maybe not including the non-smooth part? That would be my guess why centered stencils sometimes perform poorly. If you want to be a bit fancy, you can use the generic FD framework inside pySDC to generate a suitable stencil for you. You can use this function and experiment with high order centered stencils, lopsided stencils or whatever.

[lisawim]

So would only be a note for myself, so thanks a lot for your comment! Yess that's true, the state function is not smooth so your comment might make sense for me! This inherently improves the performance of the switch estimator.