View Javadoc
1   /*
2    * Licensed to the Hipparchus project under one or more
3    * contributor license agreements.  See the NOTICE file distributed with
4    * this work for additional information regarding copyright ownership.
5    * The ASF licenses this file to You under the Apache License, Version 2.0
6    * (the "License"); you may not use this file except in compliance with
7    * the License.  You may obtain a copy of the License at
8    *
9    *      https://www.apache.org/licenses/LICENSE-2.0
10   *
11   * Unless required by applicable law or agreed to in writing, software
12   * distributed under the License is distributed on an "AS IS" BASIS,
13   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14   * See the License for the specific language governing permissions and
15   * limitations under the License.
16   */
17  package org.hipparchus.ode.events;
18  
19  import org.hipparchus.Field;
20  import org.hipparchus.analysis.solvers.BracketedRealFieldUnivariateSolver;
21  import org.hipparchus.analysis.solvers.FieldBracketingNthOrderBrentSolver;
22  import org.hipparchus.ode.FieldExpandableODE;
23  import org.hipparchus.ode.FieldODEIntegrator;
24  import org.hipparchus.ode.FieldODEState;
25  import org.hipparchus.ode.FieldODEStateAndDerivative;
26  import org.hipparchus.ode.FieldOrdinaryDifferentialEquation;
27  import org.hipparchus.ode.nonstiff.DormandPrince853FieldIntegrator;
28  import org.hipparchus.ode.nonstiff.LutherFieldIntegrator;
29  import org.hipparchus.ode.sampling.FieldODEStateInterpolator;
30  import org.hipparchus.ode.sampling.FieldODEStepHandler;
31  import org.hipparchus.util.Binary64;
32  import org.hipparchus.util.Binary64Field;
33  import org.hipparchus.util.FastMath;
34  import org.junit.jupiter.api.Test;
35  
36  import java.util.ArrayList;
37  import java.util.Arrays;
38  import java.util.List;
39  
40  import static org.junit.jupiter.api.Assertions.assertEquals;
41  import static org.junit.jupiter.api.Assertions.assertFalse;
42  import static org.junit.jupiter.api.Assertions.assertSame;
43  import static org.junit.jupiter.api.Assertions.assertTrue;
44  import static org.junit.jupiter.api.Assertions.fail;
45  
46  /**
47   * Check events are detected correctly when the event times are close.
48   *
49   * @author Evan Ward
50   */
51  class FieldCloseEventsTest {
52  
53      /** type of field. */
54      private static final Field<Binary64> field = Binary64Field.getInstance();
55      private static final Binary64 zero = field.getZero();
56      private static final Binary64 one = field.getOne();
57      private static final FieldODEState<Binary64> initialState =
58              new FieldODEState<>(zero, new Binary64[]{zero, zero});
59  
60      @Test
61      void testEventAtFinalTime() {
62          final DormandPrince853FieldIntegrator<Binary64> integrator = new DormandPrince853FieldIntegrator<>(field, 10, 100.0, 1e-7, 1e-7);
63          final IntervalDetector detector = new IntervalDetector(720, 1e-10, 100, Action.RESET_DERIVATIVES, 16296.238, 17016.238);
64          integrator.addEventDetector(detector);
65          integrator.integrate(new Equation(), initialState, zero.add(17016.238));
66          assertEquals(2, detector.getEvents().size());
67      }
68  
69      @Test
70      void testCloseEventsFirstOneIsReset() {
71          // setup
72          // a fairly rare state to reproduce this bug. Two dates, d1 < d2, that
73          // are very close. Event triggers on d1 will reset state to break out of
74          // event handling loop in AbstractIntegrator.acceptStep(). At this point
75          // detector2 has g0Positive == true but the event time is set to just
76          // before the event so g(t0) is negative. Now on processing the
77          // next step the root solver checks the sign of the start, midpoint,
78          // and end of the interval so we need another event less than half a max
79          // check interval after d2 so that the g function will be negative at
80          // all three times. Then we get a non bracketing exception.
81          FieldODEIntegrator<Binary64> integrator =
82                  new DormandPrince853FieldIntegrator<>(field, 10, 100.0, 1e-7, 1e-7);
83  
84          TimeDetector detector1 = new TimeDetector(10, 1e-9, 100, Action.RESET_DERIVATIVES, 9);
85          integrator.addEventDetector(detector1);
86          TimeDetector detector2 = new TimeDetector(11, 1e-9, 100, 9 + 1e-15, 9 + 4.9);
87          integrator.addEventDetector(detector2);
88  
89          // action
90          integrator.integrate(new Equation(), initialState, zero.add(20));
91  
92          // verify
93          List<Event> events1 = detector1.getEvents();
94          assertEquals(1, events1.size());
95          assertEquals(9, events1.get(0).getT(), 0.0);
96          List<Event> events2 = detector2.getEvents();
97          assertEquals(0, events2.size());
98      }
99  
100     @Test
101     void testCloseEvents() {
102         // setup
103         double e = 1e-15;
104         FieldODEIntegrator<Binary64> integrator =
105                 new DormandPrince853FieldIntegrator<>(field, e, 100.0, 1e-7, 1e-7);
106 
107         TimeDetector detector1 = new TimeDetector(10, 1, 100, 5);
108         integrator.addEventDetector(detector1);
109         TimeDetector detector2 = new TimeDetector(10, 1, 100, 5.5);
110         integrator.addEventDetector(detector2);
111 
112         // action
113         integrator.integrate(new Equation(), initialState, one.add(20));
114 
115         // verify
116         List<Event> events1 = detector1.getEvents();
117         assertEquals(1, events1.size());
118         assertEquals(5, events1.get(0).getT(), 0.0);
119         List<Event> events2 = detector2.getEvents();
120         assertEquals(1, events2.size());
121         assertEquals(5.5, events2.get(0).getT(), 0.0);
122     }
123 
124     @Test
125     void testSimultaneousEvents() {
126         // setup
127         FieldODEIntegrator<Binary64> integrator =
128                 new DormandPrince853FieldIntegrator<>(field, 10, 100.0, 1e-7, 1e-7);
129 
130         TimeDetector detector1 = new TimeDetector(10, 1, 100, 5);
131         integrator.addEventDetector(detector1);
132         TimeDetector detector2 = new TimeDetector(10, 1, 100, 5);
133         integrator.addEventDetector(detector2);
134 
135         // action
136         integrator.integrate(new Equation(), initialState, zero.add(20));
137 
138         // verify
139         List<Event> events1 = detector1.getEvents();
140         assertEquals(1, events1.size());
141         assertEquals(5, events1.get(0).getT(), 0.0);
142         List<Event> events2 = detector2.getEvents();
143         assertEquals(1, events2.size());
144         assertEquals(5, events2.get(0).getT(), 0.0);
145     }
146 
147     /**
148      * Previously there were some branches when tryAdvance() returned false but did not
149      * set {@code t0 = t}. This allowed the order of events to not be chronological and to
150      * detect events that should not have occurred, both of which are problems.
151      */
152     @Test
153     void testSimultaneousEventsResetReverse() {
154         // setup
155         double tol = 1e-10;
156         FieldODEIntegrator<Binary64> integrator =
157                 new DormandPrince853FieldIntegrator<>(field, 10, 100.0, 1e-7, 1e-7);
158         boolean[] firstEventOccurred = {false};
159         List<Event> events = new ArrayList<>();
160 
161         TimeDetector detector1 = new TimeDetector(10, tol, 100, events, -5) {
162             @Override
163             public FieldODEEventHandler<Binary64> getHandler() {
164                 return (state, detector, increasing) -> {
165                 firstEventOccurred[0] = true;
166                 super.getHandler().eventOccurred(state, detector, increasing);
167                 return Action.RESET_STATE;
168                 };
169             }
170         };
171         integrator.addEventDetector(detector1);
172         // this detector changes it's g function definition when detector1 fires
173         TimeDetector detector2 = new TimeDetector(1, tol, 100, events, -1, -3, -5) {
174             @Override
175             public Binary64 g(final FieldODEStateAndDerivative<Binary64> state) {
176                 if (firstEventOccurred[0]) {
177                     return super.g(state);
178                 }
179                 return new TimeDetector(10, tol, 100, -5).g(state);
180             }
181         };
182         integrator.addEventDetector(detector2);
183 
184         // action
185         integrator.integrate(new Equation(), initialState, zero.add(-20));
186 
187         // verify
188         // order is important to make sure the test checks what it is supposed to
189         assertEquals(-5, events.get(0).getT(), 0.0);
190         assertTrue(events.get(0).isIncreasing());
191         assertEquals(detector1, events.get(0).getDetector());
192         assertEquals(-5, events.get(1).getT(), 0.0);
193         assertTrue(events.get(1).isIncreasing());
194         assertEquals(detector2, events.get(1).getDetector());
195         assertEquals(2, events.size());
196     }
197 
198     /**
199      * When two event detectors have a discontinuous event caused by a {@link
200      * Action#RESET_STATE} or {@link Action#RESET_DERIVATIVES}. The two event detectors
201      * would each say they had an event that had to be handled before the other one, but
202      * neither would actually back up at all. For Hipparchus GitHub #91.
203      */
204     @Test
205     void testSimultaneousDiscontinuousEventsAfterReset() {
206         // setup
207         double t = FastMath.PI;
208         double tol = 1e-10;
209         FieldODEIntegrator<Binary64> integrator =
210                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
211         List<Event> events = new ArrayList<>();
212 
213         TimeDetector resetDetector =
214                 new ResetDetector(10, tol, 100, events, new FieldODEState<>(zero.add(t), new Binary64[]{zero.add(1e100), zero}), t);
215         integrator.addEventDetector(resetDetector);
216         List<BaseDetector> detectors = new ArrayList<>();
217         for (int i = 0; i < 2; i++) {
218             BaseDetector detector1 = new StateDetector(10, tol, 100, events, 0.0);
219             integrator.addEventDetector(detector1);
220             detectors.add(detector1);
221         }
222 
223         // action
224         integrator.integrate(new Equation(), new FieldODEState<>(zero, new Binary64[]{zero.add(-1e100), zero}), zero.add(10));
225 
226         // verify
227         assertEquals(t, events.get(0).getT(), tol);
228         assertTrue(events.get(0).isIncreasing());
229         assertEquals(resetDetector, events.get(0).getDetector());
230         // next two events can occur in either order
231         assertEquals(t, events.get(1).getT(), tol);
232         assertTrue(events.get(1).isIncreasing());
233         assertEquals(detectors.get(0), events.get(1).getDetector());
234         assertEquals(t, events.get(2).getT(), tol);
235         assertTrue(events.get(2).isIncreasing());
236         assertEquals(detectors.get(1), events.get(2).getDetector());
237         assertEquals(3, events.size());
238     }
239 
240     /**
241      * test the g function switching with a period shorter than the tolerance. We don't
242      * need to find any of the events, but we do need to not crash. And we need to
243      * preserve the alternating increasing / decreasing sequence.
244      */
245     @Test
246     void testFastSwitching() {
247         // setup
248         // step size of 10 to land in between two events we would otherwise miss
249         FieldODEIntegrator<Binary64> integrator =
250                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
251 
252         TimeDetector detector1 = new TimeDetector(10, 0.2, 100, 9.9, 10.1, 12);
253         integrator.addEventDetector(detector1);
254 
255         // action
256         integrator.integrate(new Equation(), initialState, zero.add(20));
257 
258         //verify
259         // finds one or three events. Not 2.
260         List<Event> events1 = detector1.getEvents();
261         assertEquals(1, events1.size());
262         assertEquals(9.9, events1.get(0).getT(), 0.1);
263         assertTrue(events1.get(0).isIncreasing());
264     }
265 
266     /** "A Tricky Problem" from bug #239. */
267     @Test
268     void testTrickyCaseLower() {
269         // setup
270         double maxCheck = 10;
271         double tolerance = 1e-6;
272         double t1 = 1.0, t2 = 15, t3 = 16, t4 = 17, t5 = 18;
273         // shared event list so we know the order in which they occurred
274         List<Event> events = new ArrayList<>();
275         TimeDetector detectorA = new TimeDetector(maxCheck, tolerance, 100, events, t3);
276         TimeDetector detectorB = new TimeDetector(maxCheck, tolerance, 100, events, -10, t1, t2, t5);
277         TimeDetector detectorC = new TimeDetector(maxCheck, tolerance, 100, Action.RESET_DERIVATIVES, events, t4);
278 
279         FieldODEIntegrator<Binary64> integrator =
280                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
281         integrator.addEventDetector(detectorA);
282         integrator.addEventDetector(detectorB);
283         integrator.addEventDetector(detectorC);
284 
285         // action
286         integrator.integrate(new Equation(), initialState, zero.add(30.0));
287 
288         //verify
289         // really we only care that the Rules of Event Handling are not violated,
290         // but I only know one way to do that in this case.
291         assertEquals(5, events.size());
292         assertEquals(t1, events.get(0).getT(), tolerance);
293         assertFalse(events.get(0).isIncreasing());
294         assertEquals(t2, events.get(1).getT(), tolerance);
295         assertTrue(events.get(1).isIncreasing());
296         assertEquals(t3, events.get(2).getT(), tolerance);
297         assertTrue(events.get(2).isIncreasing());
298         assertEquals(t4, events.get(3).getT(), tolerance);
299         assertTrue(events.get(3).isIncreasing());
300         assertEquals(t5, events.get(4).getT(), tolerance);
301         assertFalse(events.get(4).isIncreasing());
302     }
303 
304     /**
305      * Test case for two event detectors. DetectorA has event at t2, DetectorB at t3, but
306      * due to the root finding tolerance DetectorB's event occurs at t1. With t1 < t2 <
307      * t3.
308      */
309     @Test
310     void testRootFindingTolerance() {
311         //setup
312         double maxCheck = 10;
313         double t2 = 11, t3 = t2 + 1e-5;
314         List<Event> events = new ArrayList<>();
315         TimeDetector detectorA = new TimeDetector(maxCheck, 1e-6, 100, events, t2);
316         TimeDetector detectorB = new FlatDetector(maxCheck, 0.5, 100, events, t3);
317         FieldODEIntegrator<Binary64> integrator =
318                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
319         integrator.addEventDetector(detectorA);
320         integrator.addEventDetector(detectorB);
321 
322         // action
323         integrator.integrate(new Equation(), initialState, zero.add(30.0));
324 
325         // verify
326         // if these fail the event finding did its job,
327         // but this test isn't testing what it is supposed to be
328         assertSame(detectorB, events.get(0).getDetector());
329         assertSame(detectorA, events.get(1).getDetector());
330         assertTrue(events.get(0).getT() < events.get(1).getT());
331 
332         // check event detection worked
333         assertEquals(2, events.size());
334         assertEquals(t3, events.get(0).getT(), 0.5);
335         assertTrue(events.get(0).isIncreasing());
336         assertEquals(t2, events.get(1).getT(), 1e-6);
337         assertTrue(events.get(1).isIncreasing());
338     }
339 
340     /** check when g(t < root) < 0,  g(root + convergence) < 0. */
341     @Test
342     void testRootPlusToleranceHasWrongSign() {
343         // setup
344         double maxCheck = 10;
345         double tolerance = 1e-6;
346         final double toleranceB = 0.3;
347         double t1 = 11, t2 = 11.1, t3 = 11.2;
348         // shared event list so we know the order in which they occurred
349         List<Event> events = new ArrayList<>();
350         TimeDetector detectorA = new TimeDetector(maxCheck, 1e-6, 100, events, t2);
351         TimeDetector detectorB = new TimeDetector(maxCheck, toleranceB, 100, events, t1, t3);
352         FieldODEIntegrator<Binary64> integrator =
353                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
354         integrator.addEventDetector(detectorA);
355         integrator.addEventDetector(detectorB);
356 
357         // action
358         integrator.integrate(new Equation(), initialState, zero.add(30.0));
359 
360         // verify
361         // we only care that the rules are satisfied, there are other solutions
362         assertEquals(3, events.size());
363         assertEquals(t1, events.get(0).getT(), toleranceB);
364         assertTrue(events.get(0).isIncreasing());
365         assertSame(detectorB, events.get(0).getDetector());
366         assertEquals(t2, events.get(1).getT(), tolerance);
367         assertTrue(events.get(1).isIncreasing());
368         assertSame(detectorA, events.get(1).getDetector());
369         assertEquals(t3, events.get(2).getT(), toleranceB);
370         assertFalse(events.get(2).isIncreasing());
371         assertSame(detectorB, events.get(2).getDetector());
372         // chronological
373         for (int i = 1; i < events.size(); i++) {
374             assertTrue(events.get(i).getT() >= events.get(i - 1).getT());
375         }
376     }
377 
378     /** check when g(t < root) < 0,  g(root + convergence) < 0. */
379     @Test
380     void testRootPlusToleranceHasWrongSignAndLessThanTb() {
381         // setup
382         // test is fragile w.r.t. implementation and these parameters
383         double maxCheck = 10;
384         double tolerance = 0.5;
385         double t1 = 11, t2 = 11.4, t3 = 12.0;
386         // shared event list so we know the order in which they occurred
387         List<Event> events = new ArrayList<>();
388         TimeDetector detectorB = new FlatDetector(maxCheck, tolerance, 100, events, t1, t2, t3);
389         FieldODEIntegrator<Binary64> integrator =
390                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
391         integrator.addEventDetector(detectorB);
392 
393         // action
394         integrator.integrate(new Equation(), initialState, zero.add(30.0));
395 
396         // verify
397         // allowed to find t1 or t3.
398         assertEquals(1, events.size());
399         assertEquals(t1, events.get(0).getT(), tolerance);
400         assertTrue(events.get(0).isIncreasing());
401         assertSame(detectorB, events.get(0).getDetector());
402     }
403 
404     /**
405      * Check when g(t) has a multiple root. e.g. g(t < root) < 0, g(root) = 0, g(t > root)
406      * < 0.
407      */
408     @Test
409     void testDoubleRoot() {
410         // setup
411         double maxCheck = 10;
412         double tolerance = 1e-6;
413         double t1 = 11;
414         // shared event list so we know the order in which they occurred
415         List<Event> events = new ArrayList<>();
416         TimeDetector detectorA = new TimeDetector(maxCheck, tolerance, 100, events, t1);
417         TimeDetector detectorB = new TimeDetector(maxCheck, tolerance, 100, events, t1, t1);
418         FieldODEIntegrator<Binary64> integrator =
419                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
420         integrator.addEventDetector(detectorA);
421         integrator.addEventDetector(detectorB);
422 
423         // action
424         integrator.integrate(new Equation(), initialState, zero.add(30.0));
425 
426         // verify
427         assertEquals(1, events.size());
428         assertEquals(t1, events.get(0).getT(), 0.0);
429         assertTrue(events.get(0).isIncreasing());
430         assertSame(detectorA, events.get(0).getDetector());
431         // detector worked correctly
432         assertEquals(0.0, detectorB.g(state(t1)).getReal());
433         assertTrue(detectorB.g(state(t1 - 1e-6)).getReal() < 0);
434         assertTrue(detectorB.g(state(t1 + 1e-6)).getReal() < 0);
435     }
436 
437     /**
438      * Check when g(t) has a multiple root. e.g. g(t < root) > 0, g(root) = 0, g(t > root)
439      * > 0.
440      */
441     @Test
442     void testDoubleRootOppositeSign() {
443         // setup
444         double maxCheck = 10;
445         double tolerance = 1e-6;
446         double t1 = 11;
447         // shared event list so we know the order in which they occurred
448         List<Event> events = new ArrayList<>();
449         TimeDetector detectorA = new TimeDetector(maxCheck, tolerance, 100, events, t1);
450         TimeDetector detectorB = new ContinuousDetector(maxCheck, tolerance, 100, events, -20, t1, t1);
451         detectorB.g(state(t1));
452         FieldODEIntegrator<Binary64> integrator =
453                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
454         integrator.addEventDetector(detectorA);
455         integrator.addEventDetector(detectorB);
456 
457         // action
458         integrator.integrate(new Equation(), initialState, zero.add(30.0));
459 
460         // verify
461         assertEquals(1, events.size());
462         assertEquals(t1, events.get(0).getT(), 0.0);
463         assertTrue(events.get(0).isIncreasing());
464         assertSame(detectorA, events.get(0).getDetector());
465         // detector worked correctly
466         assertEquals(0.0, detectorB.g(state(t1)).getReal(), 0.0);
467         assertTrue(detectorB.g(state(t1 - 1e-6)).getReal() > 0);
468         assertTrue(detectorB.g(state(t1 + 1e-6)).getReal() > 0);
469     }
470 
471     /** check root finding when zero at both ends. */
472     @Test
473     void testZeroAtBeginningAndEndOfInterval() {
474         // setup
475         double maxCheck = 10;
476         double tolerance = 1e-6;
477         double t1 = 10, t2 = 20;
478         // shared event list so we know the order in which they occurred
479         List<Event> events = new ArrayList<>();
480         TimeDetector detectorA = new ContinuousDetector(maxCheck, tolerance, 100, events, t1, t2);
481         FieldODEIntegrator<Binary64> integrator =
482                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
483         integrator.addEventDetector(detectorA);
484 
485         // action
486         integrator.integrate(new Equation(), initialState, zero.add(30.0));
487 
488         // verify
489         assertEquals(2, events.size());
490         assertEquals(t1, events.get(0).getT(), 0.0);
491         assertTrue(events.get(0).isIncreasing());
492         assertSame(detectorA, events.get(0).getDetector());
493         assertEquals(t2, events.get(1).getT(), 0.0);
494         assertFalse(events.get(1).isIncreasing());
495         assertSame(detectorA, events.get(1).getDetector());
496     }
497 
498     /** check root finding when zero at both ends. */
499     @Test
500     void testZeroAtBeginningAndEndOfIntervalOppositeSign() {
501         // setup
502         double maxCheck = 10;
503         double tolerance = 1e-6;
504         double t1 = 10, t2 = 20;
505         // shared event list so we know the order in which they occurred
506         List<Event> events = new ArrayList<>();
507         TimeDetector detectorA = new ContinuousDetector(maxCheck, tolerance, 100, events, -10, t1, t2);
508         FieldODEIntegrator<Binary64> integrator =
509                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
510         integrator.addEventDetector(detectorA);
511 
512         // action
513         integrator.integrate(new Equation(), initialState, zero.add(30.0));
514 
515         // verify
516         assertEquals(2, events.size());
517         assertEquals(t1, events.get(0).getT(), 0.0);
518         assertFalse(events.get(0).isIncreasing());
519         assertSame(detectorA, events.get(0).getDetector());
520         assertEquals(t2, events.get(1).getT(), 0.0);
521         assertTrue(events.get(1).isIncreasing());
522         assertSame(detectorA, events.get(1).getDetector());
523     }
524 
525     /** Test where an event detector has to back up multiple times. */
526     @Test
527     void testMultipleBackups() {
528         // setup
529         double maxCheck = 10;
530         double tolerance = 1e-6;
531         double t1 = 11.0, t2 = 12, t3 = 13, t4 = 14, t5 = 15, t6 = 16, t7 = 17;
532         // shared event list so we know the order in which they occurred
533         List<Event> events = new ArrayList<>();
534         TimeDetector detectorA = new ContinuousDetector(maxCheck, tolerance, 100, events, t6);
535         TimeDetector detectorB = new ContinuousDetector(maxCheck, tolerance, 100, events, t1, t3, t4, t7);
536         TimeDetector detectorC = new ContinuousDetector(maxCheck, tolerance, 100, events, t2, t5);
537 
538         FieldODEIntegrator<Binary64> integrator =
539                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
540         integrator.addEventDetector(detectorA);
541         integrator.addEventDetector(detectorB);
542         integrator.addEventDetector(detectorC);
543 
544         // action
545         integrator.integrate(new Equation(), initialState, zero.add(30.0));
546 
547         //verify
548         // really we only care that the Rules of Event Handling are not violated,
549         assertEquals(5, events.size());
550         assertEquals(t1, events.get(0).getT(), tolerance);
551         assertTrue(events.get(0).isIncreasing());
552         assertEquals(detectorB, events.get(0).getDetector());
553         assertEquals(t2, events.get(1).getT(), tolerance);
554         assertTrue(events.get(1).isIncreasing());
555         assertEquals(detectorC, events.get(1).getDetector());
556         // reporting t3 and t4 is optional, seeing them is not.
557         // we know a root was found at t3 because events are reported at t2 and t5.
558         /*
559         Assertions.assertEquals(t3, events.get(2).getT(), tolerance);
560         Assertions.assertEquals(false, events.get(2).isIncreasing());
561         Assertions.assertEquals(detectorB, events.get(2).getDetector());
562         Assertions.assertEquals(t4, events.get(3).getT(), tolerance);
563         Assertions.assertEquals(true, events.get(3).isIncreasing());
564         Assertions.assertEquals(detectorB, events.get(3).getDetector());
565         */
566         assertEquals(t5, events.get(2).getT(), tolerance);
567         assertFalse(events.get(2).isIncreasing());
568         assertEquals(detectorC, events.get(2).getDetector());
569         assertEquals(t6, events.get(3).getT(), tolerance);
570         assertTrue(events.get(3).isIncreasing());
571         assertEquals(detectorA, events.get(3).getDetector());
572         assertEquals(t7, events.get(4).getT(), tolerance);
573         assertFalse(events.get(4).isIncreasing());
574         assertEquals(detectorB, events.get(4).getDetector());
575     }
576 
577     /** Test a reset event triggering another event at the same time. */
578     @Test
579     void testEventCausedByStateReset() {
580         // setup
581         double maxCheck = 10;
582         double tolerance = 1e-6;
583         double t1 = 15.0;
584         final FieldODEState<Binary64> newState = new FieldODEState<>(
585                 zero.add(t1), new Binary64[]{zero.add(-20), zero});
586         // shared event list so we know the order in which they occurred
587         List<Event> events = new ArrayList<>();
588         TimeDetector detectorA = new ResetDetector(maxCheck, tolerance, 100, events, newState, t1);
589         BaseDetector detectorB = new StateDetector(maxCheck, tolerance, 100, events, -1);
590 
591         FieldODEIntegrator<Binary64> integrator =
592                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
593         integrator.addEventDetector(detectorA);
594         integrator.addEventDetector(detectorB);
595 
596         // action
597         integrator.integrate(new Equation(), initialState, zero.add(40.0));
598 
599         //verify
600         // really we only care that the Rules of Event Handling are not violated,
601         assertEquals(3, events.size());
602         assertEquals(t1, events.get(0).getT(), tolerance);
603         assertTrue(events.get(0).isIncreasing());
604         assertEquals(detectorA, events.get(0).getDetector());
605         assertEquals(t1, events.get(1).getT(), tolerance);
606         assertFalse(events.get(1).isIncreasing());
607         assertEquals(detectorB, events.get(1).getDetector());
608         assertEquals(t1 + 19, events.get(2).getT(), tolerance);
609         assertTrue(events.get(2).isIncreasing());
610         assertEquals(detectorB, events.get(2).getDetector());
611     }
612 
613     /** check when t + tolerance == t. */
614     @Test
615     void testConvergenceTooTight() {
616         // setup
617         double maxCheck = 10;
618         double tolerance = 1e-18;
619         double t1 = 15;
620         // shared event list so we know the order in which they occurred
621         List<Event> events = new ArrayList<>();
622         TimeDetector detectorA = new ContinuousDetector(maxCheck, tolerance, 100, events, t1);
623         FieldODEIntegrator<Binary64> integrator =
624                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
625         integrator.addEventDetector(detectorA);
626 
627         // action
628         integrator.integrate(new Equation(), initialState, zero.add(30.0));
629 
630         // verify
631         assertEquals(1, events.size());
632         assertEquals(t1, events.get(0).getT(), 0.0);
633         assertTrue(events.get(0).isIncreasing());
634         assertSame(detectorA, events.get(0).getDetector());
635     }
636 
637     /** check when root finding tolerance > event finding tolerance. */
638     @Test
639     void testToleranceMismatch() {
640         // setup
641         double maxCheck = 10;
642         double tolerance = 1e-18;
643         double t1 = 15.1;
644         // shared event list so we know the order in which they occurred
645         List<Event> events = new ArrayList<>();
646         TimeDetector detectorA = new FlatDetector(maxCheck, tolerance, 100, events, t1);
647         FieldODEIntegrator<Binary64> integrator =
648                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
649         integrator.addEventDetector(detectorA);
650 
651         // action
652         integrator.integrate(new Equation(), initialState, zero.add(30.0));
653 
654         // verify
655         assertEquals(1, events.size());
656         // use root finder tolerance instead of event finder tolerance.
657         assertEquals(t1, events.get(0).getT(), 1e-3);
658         assertTrue(events.get(0).isIncreasing());
659         assertSame(detectorA, events.get(0).getDetector());
660     }
661 
662     /**
663      * test when one event detector changes the definition of another's g function before
664      * the end of the step as a result of a continue action. Not sure if this should be
665      * officially supported, but it is used in Orekit's DateDetector, it's useful, and not
666      * too hard to implement.
667      */
668     @Test
669     void testEventChangesGFunctionDefinition() {
670         // setup
671         double maxCheck = 5;
672         double tolerance = 1e-6;
673         double t1 = 11, t2 = 19;
674         // shared event list so we know the order in which they occurred
675         List<Event> events = new ArrayList<>();
676         // mutable boolean
677         boolean[] swap = new boolean[1];
678         final TimeDetector detectorA = new ContinuousDetector(maxCheck, tolerance, 100, events, t1) {
679             @Override
680             public FieldODEEventHandler<Binary64> getHandler() {
681                 return (state, detector, increasing) -> {
682                     swap[0] = true;
683                     return super.getHandler().eventOccurred(state, detector, increasing);
684                 };
685             }
686         };
687         final TimeDetector detectorB = new ContinuousDetector(maxCheck, tolerance, 100, events, t2);
688         BaseDetector detectorC = new BaseDetector(maxCheck, tolerance, 100, Action.CONTINUE, events) {
689 
690             @Override
691             public Binary64 g(FieldODEStateAndDerivative<Binary64> state) {
692                 if (swap[0]) {
693                     return detectorB.g(state);
694                 } else {
695                     return zero.add(-1);
696                 }
697             }
698 
699         };
700         FieldODEIntegrator<Binary64> integrator =
701                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
702         integrator.addEventDetector(detectorA);
703         integrator.addEventDetector(detectorC);
704 
705         // action
706         integrator.integrate(new Equation(), initialState, zero.add(30.0));
707 
708         // verify
709         assertEquals(2, events.size());
710         assertEquals(t1, events.get(0).getT(), tolerance);
711         assertTrue(events.get(0).isIncreasing());
712         assertSame(detectorA, events.get(0).getDetector());
713         assertEquals(t2, events.get(1).getT(), tolerance);
714         assertTrue(events.get(1).isIncreasing());
715         assertSame(detectorC, events.get(1).getDetector());
716     }
717 
718 
719     /**
720      * test when one event detector changes the definition of another's g function before
721      * the end of the step as a result of an event occurring. In this case the change
722      * cancels the occurrence of the event.
723      */
724     @Test
725     void testEventChangesGFunctionDefinitionCancel() {
726         // setup
727         double maxCheck = 5;
728         double tolerance = 1e-6;
729         double t1 = 11, t2 = 11.1;
730         // shared event list so we know the order in which they occurred
731         List<Event> events = new ArrayList<>();
732         // mutable boolean
733         boolean[] swap = new boolean[1];
734         final TimeDetector detectorA =
735                         new ContinuousDetector(maxCheck, tolerance, 100, Action.RESET_EVENTS, events, t1) {
736             @Override
737             public FieldODEEventHandler<Binary64> getHandler() {
738                 return (state, detector, increasing) -> {
739                     swap[0] = true;
740                     return super.getHandler().eventOccurred(state, detector, increasing);
741                 };
742             }
743         };
744         final TimeDetector detectorB = new ContinuousDetector(maxCheck, tolerance, 100, events, t2);
745         BaseDetector detectorC = new BaseDetector(maxCheck, tolerance, 100, Action.CONTINUE, events) {
746 
747             @Override
748             public Binary64 g(FieldODEStateAndDerivative<Binary64> state) {
749                 if (!swap[0]) {
750                     return detectorB.g(state);
751                 } else {
752                     return zero.add(-1);
753                 }
754             }
755 
756         };
757         FieldODEIntegrator<Binary64> integrator =
758                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
759         integrator.addEventDetector(detectorA);
760         integrator.addEventDetector(detectorC);
761 
762         // action
763         integrator.integrate(new Equation(), initialState, zero.add(30.0));
764 
765         // verify
766         assertEquals(1, events.size());
767         assertEquals(t1, events.get(0).getT(), tolerance);
768         assertTrue(events.get(0).isIncreasing());
769         assertSame(detectorA, events.get(0).getDetector());
770     }
771 
772     /**
773      * test when one event detector changes the definition of another's g function before
774      * the end of the step as a result of an event occurring. In this case the change
775      * delays the occurrence of the event.
776      */
777     @Test
778     void testEventChangesGFunctionDefinitionDelay() {
779         // setup
780         double maxCheck = 5;
781         double tolerance = 1e-6;
782         double t1 = 11, t2 = 11.1, t3 = 11.2;
783         // shared event list so we know the order in which they occurred
784         List<Event> events = new ArrayList<>();
785         // mutable boolean
786         boolean[] swap = new boolean[1];
787         final TimeDetector detectorA =
788                         new ContinuousDetector(maxCheck, tolerance, 100, Action.RESET_EVENTS, events, t1) {
789             @Override
790             public FieldODEEventHandler<Binary64> getHandler() {
791                 return (state, detector, increasing) -> {
792                     swap[0] = true;
793                     return super.getHandler().eventOccurred(state, detector, increasing);
794                 };
795             }
796         };
797         final TimeDetector detectorB = new ContinuousDetector(maxCheck, tolerance, 100, events, t2);
798         final TimeDetector detectorD = new ContinuousDetector(maxCheck, tolerance, 100, events, t3);
799         BaseDetector detectorC = new BaseDetector(maxCheck, tolerance, 100, Action.CONTINUE, events) {
800 
801             @Override
802             public Binary64 g(FieldODEStateAndDerivative<Binary64> state) {
803                 if (!swap[0]) {
804                     return detectorB.g(state);
805                 } else {
806                     return detectorD.g(state);
807                 }
808             }
809 
810         };
811         FieldODEIntegrator<Binary64> integrator =
812                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
813         integrator.addEventDetector(detectorA);
814         integrator.addEventDetector(detectorC);
815 
816         // action
817         integrator.integrate(new Equation(), initialState, zero.add(30.0));
818 
819         // verify
820         assertEquals(2, events.size());
821         assertEquals(t1, events.get(0).getT(), tolerance);
822         assertTrue(events.get(0).isIncreasing());
823         assertSame(detectorA, events.get(0).getDetector());
824         assertEquals(t3, events.get(1).getT(), tolerance);
825         assertTrue(events.get(1).isIncreasing());
826         assertSame(detectorC, events.get(1).getDetector());
827     }
828 
829     /**
830      * test when one event detector changes the definition of another's g function before
831      * the end of the step as a result of an event occurring. In this case the change
832      * causes the event to happen sooner than originally expected.
833      */
834     @Test
835     void testEventChangesGFunctionDefinitionAccelerate() {
836         // setup
837         double maxCheck = 5;
838         double tolerance = 1e-6;
839         double t1 = 11, t2 = 11.1, t3 = 11.2;
840         // shared event list so we know the order in which they occurred
841         List<Event> events = new ArrayList<>();
842         // mutable boolean
843         boolean[] swap = new boolean[1];
844         final TimeDetector detectorA =
845                         new ContinuousDetector(maxCheck, tolerance, 100, Action.RESET_EVENTS, events, t1) {
846             @Override
847             public FieldODEEventHandler<Binary64> getHandler() {
848                 return (state, detector, increasing) -> {
849                     swap[0] = true;
850                     return super.getHandler().eventOccurred(state, detector, increasing);
851                 };
852             }
853         };
854         final TimeDetector detectorB = new ContinuousDetector(maxCheck, tolerance, 100, events, t2);
855         final TimeDetector detectorD = new ContinuousDetector(maxCheck, tolerance, 100, events, t3);
856         BaseDetector detectorC = new BaseDetector(maxCheck, tolerance, 100, Action.CONTINUE, events) {
857 
858             @Override
859             public Binary64 g(FieldODEStateAndDerivative<Binary64> state) {
860                 if (swap[0]) {
861                     return detectorB.g(state);
862                 } else {
863                     return detectorD.g(state);
864                 }
865             }
866 
867         };
868         FieldODEIntegrator<Binary64> integrator =
869                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
870         integrator.addEventDetector(detectorA);
871         integrator.addEventDetector(detectorC);
872 
873         // action
874         integrator.integrate(new Equation(), initialState, zero.add(30.0));
875 
876         // verify
877         assertEquals(2, events.size());
878         assertEquals(t1, events.get(0).getT(), tolerance);
879         assertTrue(events.get(0).isIncreasing());
880         assertSame(detectorA, events.get(0).getDetector());
881         assertEquals(t2, events.get(1).getT(), tolerance);
882         assertTrue(events.get(1).isIncreasing());
883         assertSame(detectorC, events.get(1).getDetector());
884     }
885 
886     /** check when root finding tolerance > event finding tolerance. */
887     @Test
888     void testToleranceStop() {
889         // setup
890         double maxCheck = 10;
891         double tolerance = 1e-18; // less than 1 ulp
892         double t1 = 15.1;
893         // shared event list so we know the order in which they occurred
894         List<Event> events = new ArrayList<>();
895         TimeDetector detectorA = new FlatDetector(maxCheck, tolerance, 100, Action.STOP, events, t1);
896         FieldODEIntegrator<Binary64> integrator =
897                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
898         integrator.addEventDetector(detectorA);
899 
900         // action
901         FieldODEStateAndDerivative<Binary64> finalState =
902                 integrator.integrate(new Equation(), initialState, zero.add(30.0));
903 
904         // verify
905         assertEquals(1, events.size());
906         // use root finder tolerance instead of event finder tolerance.
907         assertEquals(t1, events.get(0).getT(), tolerance);
908         assertTrue(events.get(0).isIncreasing());
909         assertSame(detectorA, events.get(0).getDetector());
910         assertEquals(t1, finalState.getTime().getReal(), tolerance);
911 
912         // try to resume propagation
913         finalState = integrator.integrate(new Equation(), finalState, zero.add(30.0));
914 
915         // verify it got to the end
916         assertEquals(30.0, finalState.getTime().getReal(), 0.0);
917     }
918 
919     /**
920      * Test when g function is initially zero for longer than the tolerance. Can occur
921      * when restarting after a stop and cancellation occurs in the g function.
922      */
923     @Test
924     void testLongInitialZero() {
925         // setup
926         double maxCheck = 10;
927         double tolerance = 1;
928         // shared event list so we know the order in which they occurred
929         List<Event> events = new ArrayList<>();
930         TimeDetector detectorA = new TimeDetector(maxCheck, tolerance, 100, Action.STOP, events, -50) {
931             @Override
932             public Binary64 g(FieldODEStateAndDerivative<Binary64> state) {
933                 if (state.getTime().getReal() < 2) {
934                     return zero;
935                 } else {
936                     return super.g(state);
937                 }
938             }
939         };
940         FieldODEIntegrator<Binary64> integrator =
941                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
942         integrator.addEventDetector(detectorA);
943 
944         // action
945         integrator.integrate(new Equation(), initialState, zero.add(30.0));
946 
947         // verify
948         assertEquals(0, events.size());
949     }
950 
951     /**
952      * The root finder requires the start point to be in the interval (a, b) which is hard
953      * when there aren't many numbers between a and b. This test uses a second event
954      * detector to force a very small window for the first event detector.
955      */
956     @Test
957     void testShortBracketingInterval() {
958         // setup
959         double maxCheck = 10;
960         double tolerance = 1e-6;
961         final double t1 = FastMath.nextUp(10.0), t2 = 10.5;
962         // shared event list so we know the order in which they occurred
963         List<Event> events = new ArrayList<>();
964         // never zero so there is no easy way out
965         TimeDetector detectorA = new TimeDetector(maxCheck, tolerance, 100, events) {
966             @Override
967             public Binary64 g(FieldODEStateAndDerivative<Binary64> state) {
968                 final Binary64 t = state.getTime();
969                 if (t.getReal() < t1) {
970                     return one.negate();
971                 } else if (t.getReal() < t2) {
972                     return one;
973                 } else {
974                     return one.negate();
975                 }
976             }
977         };
978         TimeDetector detectorB = new TimeDetector(maxCheck, tolerance, 100, events, t1);
979         FieldODEIntegrator<Binary64> integrator =
980                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
981         integrator.addEventDetector(detectorA);
982         integrator.addEventDetector(detectorB);
983 
984         // action
985         integrator.integrate(new Equation(), initialState, zero.add(30.0));
986 
987         // verify
988         assertEquals(3, events.size());
989         assertEquals(t1, events.get(0).getT(), tolerance);
990         assertTrue(events.get(0).isIncreasing());
991         assertSame(detectorA, events.get(0).getDetector());
992         assertEquals(t1, events.get(1).getT(), tolerance);
993         assertTrue(events.get(1).isIncreasing());
994         assertSame(detectorB, events.get(1).getDetector());
995         assertEquals(t2, events.get(2).getT(), tolerance);
996         assertFalse(events.get(2).isIncreasing());
997         assertSame(detectorA, events.get(2).getDetector());
998     }
999 
1000     /** Check that steps are restricted correctly with a continue event. */
1001     @Test
1002     void testEventStepHandler() {
1003         // setup
1004         double tolerance = 1e-18;
1005         FieldODEIntegrator<Binary64> integrator =
1006                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1007         integrator.addEventDetector(new TimeDetector(100.0, tolerance, 100, 5));
1008         StepHandler stepHandler = new StepHandler();
1009         integrator.addStepHandler(stepHandler);
1010 
1011         // action
1012         FieldODEStateAndDerivative<Binary64> finalState = integrator
1013                 .integrate(new Equation(), initialState, zero.add(10));
1014 
1015         // verify
1016         assertEquals(10.0, finalState.getTime().getReal(), tolerance);
1017         assertEquals(0.0,
1018                 stepHandler.initialState.getTime().getReal(), tolerance);
1019         assertEquals(10.0, stepHandler.finalTime.getReal(), tolerance);
1020         assertEquals(10.0,
1021                 stepHandler.finalState.getTime().getReal(), tolerance);
1022         FieldODEStateInterpolator<Binary64> interpolator = stepHandler.interpolators.get(0);
1023         assertEquals(0.0,
1024                 interpolator.getPreviousState().getTime().getReal(), tolerance);
1025         assertEquals(5.0,
1026                 interpolator.getCurrentState().getTime().getReal(), tolerance);
1027         interpolator = stepHandler.interpolators.get(1);
1028         assertEquals(5.0,
1029                 interpolator.getPreviousState().getTime().getReal(), tolerance);
1030         assertEquals(10.0,
1031                 interpolator.getCurrentState().getTime().getReal(), tolerance);
1032         assertEquals(2, stepHandler.interpolators.size());
1033     }
1034 
1035     /** Test resetState(...) returns {@code null}. */
1036     @Test
1037     void testEventCausedByDerivativesReset() {
1038         // setup
1039         TimeDetector detectorA = new TimeDetector(10, 1e-6, 100, Action.RESET_STATE, 15.0) {
1040             @Override
1041             public FieldODEEventHandler<Binary64> getHandler() {
1042                 return new FieldODEEventHandler<Binary64>() {
1043                     @Override
1044                     public Action eventOccurred(FieldODEStateAndDerivative<Binary64> state,
1045                                                 FieldODEEventDetector<Binary64> detector,
1046                                                 boolean increasing) {
1047                         return Action.RESET_STATE;
1048                     }
1049                     @Override
1050                     public FieldODEState<Binary64> resetState(FieldODEEventDetector<Binary64> detector,
1051                                                                FieldODEStateAndDerivative<Binary64> state) {
1052                         return null;
1053                     }
1054                 };
1055             }
1056         };
1057         FieldODEIntegrator<Binary64> integrator =
1058                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1059         integrator.addEventDetector(detectorA);
1060 
1061         try {
1062             // action
1063             integrator.integrate(new Equation(), initialState, zero.add(20.0));
1064             fail("Expected Exception");
1065         } catch (NullPointerException e) {
1066             // expected
1067         }
1068     }
1069 
1070     @Test
1071     void testResetChangesSign() {
1072         FieldOrdinaryDifferentialEquation<Binary64> equation = new FieldOrdinaryDifferentialEquation<Binary64>() {
1073             public int getDimension() { return 1; }
1074             public Binary64[] computeDerivatives(Binary64 t, Binary64[] y) { return new Binary64[] { new Binary64(1.0) }; }
1075         };
1076 
1077         LutherFieldIntegrator<Binary64> integrator = new LutherFieldIntegrator<>(Binary64Field.getInstance(), new Binary64(20.0));
1078         final double small = 1.0e-10;
1079         ResetChangesSignGenerator eventsGenerator = new ResetChangesSignGenerator(6.0, 9.0, -0.5 * small, 8.0, small, 1000);
1080         integrator.addEventDetector(eventsGenerator);
1081         final FieldODEStateAndDerivative<Binary64> end = integrator.integrate(new FieldExpandableODE<>(equation),
1082                 new FieldODEState<>(new Binary64(0.0),
1083                         new Binary64[] { new Binary64(0.0) }),
1084                 new Binary64(100.0));
1085         assertEquals(2,                 eventsGenerator.getCount());
1086         assertEquals(9.0,               end.getCompleteState()[0].getReal(), 1.0e-12);
1087         assertEquals(9.0 + 0.5 * small, end.getTime().getReal(),             1.0e-12);
1088     }
1089 
1090     /* The following tests are copies of the above tests, except that they propagate in
1091      * the reverse direction and all the signs on the time values are negated.
1092      */
1093 
1094 
1095     @Test
1096     void testCloseEventsFirstOneIsResetReverse() {
1097         // setup
1098         // a fairly rare state to reproduce this bug. Two dates, d1 < d2, that
1099         // are very close. Event triggers on d1 will reset state to break out of
1100         // event handling loop in AbstractIntegrator.acceptStep(). At this point
1101         // detector2 has g0Positive == true but the event time is set to just
1102         // before the event so g(t0) is negative. Now on processing the
1103         // next step the root solver checks the sign of the start, midpoint,
1104         // and end of the interval so we need another event less than half a max
1105         // check interval after d2 so that the g function will be negative at
1106         // all three times. Then we get a non bracketing exception.
1107         FieldODEIntegrator<Binary64> integrator =
1108                 new DormandPrince853FieldIntegrator<>(field, 10, 100.0, 1e-7, 1e-7);
1109 
1110         // switched for 9 to 1 to be close to the start of the step
1111         double t1 = -1;
1112         TimeDetector detector1 = new TimeDetector(10, 1e-9, 100, Action.RESET_DERIVATIVES, t1);
1113         integrator.addEventDetector(detector1);
1114         TimeDetector detector2 = new TimeDetector(11, 1e-9, 100, t1 - 1e-15, t1 - 4.9);
1115         integrator.addEventDetector(detector2);
1116 
1117         // action
1118         integrator.integrate(new Equation(), initialState, zero.add(-20));
1119 
1120         // verify
1121         List<Event> events1 = detector1.getEvents();
1122         assertEquals(1, events1.size());
1123         assertEquals(t1, events1.get(0).getT(), 0.0);
1124         List<Event> events2 = detector2.getEvents();
1125         assertEquals(0, events2.size());
1126     }
1127 
1128     @Test
1129     void testCloseEventsReverse() {
1130         // setup
1131         double e = 1e-15;
1132         FieldODEIntegrator<Binary64> integrator =
1133                 new DormandPrince853FieldIntegrator<>(field, e, 100.0, 1e-7, 1e-7);
1134 
1135         TimeDetector detector1 = new TimeDetector(10, 1, 100, -5);
1136         integrator.addEventDetector(detector1);
1137         TimeDetector detector2 = new TimeDetector(10, 1, 100, -5.5);
1138         integrator.addEventDetector(detector2);
1139 
1140         // action
1141         integrator.integrate(new Equation(), initialState, zero.add(-20));
1142 
1143         // verify
1144         List<Event> events1 = detector1.getEvents();
1145         assertEquals(1, events1.size());
1146         assertEquals(-5, events1.get(0).getT(), 0.0);
1147         List<Event> events2 = detector2.getEvents();
1148         assertEquals(1, events2.size());
1149         assertEquals(-5.5, events2.get(0).getT(), 0.0);
1150     }
1151 
1152     @Test
1153     void testSimultaneousEventsReverse() {
1154         // setup
1155         FieldODEIntegrator<Binary64> integrator =
1156                 new DormandPrince853FieldIntegrator<>(field, 10, 100.0, 1e-7, 1e-7);
1157 
1158         TimeDetector detector1 = new TimeDetector(10, 1, 100, -5);
1159         integrator.addEventDetector(detector1);
1160         TimeDetector detector2 = new TimeDetector(10, 1, 100, -5);
1161         integrator.addEventDetector(detector2);
1162 
1163         // action
1164         integrator.integrate(new Equation(), initialState, zero.add(-20));
1165 
1166         // verify
1167         List<Event> events1 = detector1.getEvents();
1168         assertEquals(1, events1.size());
1169         assertEquals(-5, events1.get(0).getT(), 0.0);
1170         List<Event> events2 = detector2.getEvents();
1171         assertEquals(1, events2.size());
1172         assertEquals(-5, events2.get(0).getT(), 0.0);
1173     }
1174 
1175     /**
1176      * Previously there were some branches when tryAdvance() returned false but did not
1177      * set {@code t0 = t}. This allowed the order of events to not be chronological and to
1178      * detect events that should not have occurred, both of which are problems.
1179      */
1180     @Test
1181     void testSimultaneousEventsReset() {
1182         // setup
1183         double tol = 1e-10;
1184         FieldODEIntegrator<Binary64> integrator =
1185                 new DormandPrince853FieldIntegrator<>(field, 10, 100.0, 1e-7, 1e-7);
1186         boolean[] firstEventOccurred = {false};
1187         List<Event> events = new ArrayList<>();
1188 
1189         TimeDetector detector1 = new TimeDetector(10, tol, 100, events, 5) {
1190             @Override
1191             public FieldODEEventHandler<Binary64> getHandler() {
1192                 return (state, detector, increasing) -> {
1193                     firstEventOccurred[0] = true;
1194                     super.getHandler().eventOccurred(state, detector, increasing);
1195                     return Action.RESET_STATE;
1196                 };
1197             }
1198         };
1199         integrator.addEventDetector(detector1);
1200         // this detector changes it's g function definition when detector1 fires
1201         TimeDetector detector2 = new TimeDetector(1, tol, 100, events, 1, 3, 5) {
1202             @Override
1203             public Binary64 g(final FieldODEStateAndDerivative<Binary64> state) {
1204                 if (firstEventOccurred[0]) {
1205                     return super.g(state);
1206                 }
1207                 return new TimeDetector(1, tol, 100, 5).g(state);
1208             }
1209         };
1210         integrator.addEventDetector(detector2);
1211 
1212         // action
1213         integrator.integrate(new Equation(), initialState, zero.add(20));
1214 
1215         // verify
1216         // order is important to make sure the test checks what it is supposed to
1217         assertEquals(5, events.get(0).getT(), 0.0);
1218         assertTrue(events.get(0).isIncreasing());
1219         assertEquals(detector1, events.get(0).getDetector());
1220         assertEquals(5, events.get(1).getT(), 0.0);
1221         assertTrue(events.get(1).isIncreasing());
1222         assertEquals(detector2, events.get(1).getDetector());
1223         assertEquals(2, events.size());
1224     }
1225 
1226     /**
1227      * When two event detectors have a discontinuous event caused by a {@link
1228      * Action#RESET_STATE} or {@link Action#RESET_DERIVATIVES}. The two event detectors
1229      * would each say they had an event that had to be handled before the other one, but
1230      * neither would actually back up at all. For Hipparchus GitHub #91.
1231      */
1232     @Test
1233     void testSimultaneousDiscontinuousEventsAfterResetReverse() {
1234         // setup
1235         double t = -FastMath.PI;
1236         double tol = 1e-10;
1237         FieldODEIntegrator<Binary64> integrator =
1238                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1239         List<Event> events = new ArrayList<>();
1240 
1241         TimeDetector resetDetector =
1242                 new ResetDetector(10, tol, 100, events, new FieldODEState<>(zero.add(t), new Binary64[]{zero.add(1e100), zero}), t);
1243         integrator.addEventDetector(resetDetector);
1244         List<BaseDetector> detectors = new ArrayList<>();
1245         for (int i = 0; i < 2; i++) {
1246             BaseDetector detector1 = new StateDetector(10, tol, 100, events, 0.0);
1247             integrator.addEventDetector(detector1);
1248             detectors.add(detector1);
1249         }
1250 
1251         // action
1252         integrator.integrate(new Equation(), new FieldODEState<>(zero, new Binary64[]{zero.add(-1e100), zero}), zero.add(-10));
1253 
1254         // verify
1255         assertEquals(t, events.get(0).getT(), tol);
1256         assertTrue(events.get(0).isIncreasing());
1257         assertEquals(resetDetector, events.get(0).getDetector());
1258         // next two events can occur in either order
1259         assertEquals(t, events.get(1).getT(), tol);
1260         assertFalse(events.get(1).isIncreasing());
1261         assertEquals(detectors.get(0), events.get(1).getDetector());
1262         assertEquals(t, events.get(2).getT(), tol);
1263         assertFalse(events.get(2).isIncreasing());
1264         assertEquals(detectors.get(1), events.get(2).getDetector());
1265         assertEquals(3, events.size());
1266     }
1267 
1268     /**
1269      * test the g function switching with a period shorter than the tolerance. We don't
1270      * need to find any of the events, but we do need to not crash. And we need to
1271      * preserve the alternating increasing / decreasing sequence.
1272      */
1273     @Test
1274     void testFastSwitchingReverse() {
1275         // setup
1276         // step size of 10 to land in between two events we would otherwise miss
1277         FieldODEIntegrator<Binary64> integrator =
1278                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1279 
1280         TimeDetector detector1 = new TimeDetector(10, 0.2, 100, -9.9, -10.1, -12);
1281         integrator.addEventDetector(detector1);
1282 
1283         // action
1284         integrator.integrate(new Equation(), initialState, zero.add(-20));
1285 
1286         //verify
1287         // finds one or three events. Not 2.
1288         List<Event> events1 = detector1.getEvents();
1289         assertEquals(1, events1.size());
1290         assertEquals(-9.9, events1.get(0).getT(), 0.2);
1291         assertTrue(events1.get(0).isIncreasing());
1292     }
1293 
1294     /** "A Tricky Problem" from bug #239. */
1295     @Test
1296     void testTrickyCaseLowerReverse() {
1297         // setup
1298         double maxCheck = 10;
1299         double tolerance = 1e-6;
1300         double t1 = -1.0, t2 = -15, t3 = -16, t4 = -17, t5 = -18;
1301         // shared event list so we know the order in which they occurred
1302         List<Event> events = new ArrayList<>();
1303         TimeDetector detectorA = new TimeDetector(maxCheck, tolerance, 100, events, t3);
1304         TimeDetector detectorB = new TimeDetector(maxCheck, tolerance, 100, events, -50, t1, t2, t5);
1305         TimeDetector detectorC = new TimeDetector(maxCheck, tolerance, 100, Action.RESET_DERIVATIVES, events, t4);
1306 
1307         FieldODEIntegrator<Binary64> integrator =
1308                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1309         integrator.addEventDetector(detectorA);
1310         integrator.addEventDetector(detectorB);
1311         integrator.addEventDetector(detectorC);
1312 
1313         // action
1314         integrator.integrate(new Equation(), initialState, zero.add(-30.0));
1315 
1316         //verify
1317         // really we only care that the Rules of Event Handling are not violated,
1318         // but I only know one way to do that in this case.
1319         assertEquals(5, events.size());
1320         assertEquals(t1, events.get(0).getT(), tolerance);
1321         assertFalse(events.get(0).isIncreasing());
1322         assertEquals(t2, events.get(1).getT(), tolerance);
1323         assertTrue(events.get(1).isIncreasing());
1324         assertEquals(t3, events.get(2).getT(), tolerance);
1325         assertTrue(events.get(2).isIncreasing());
1326         assertEquals(t4, events.get(3).getT(), tolerance);
1327         assertTrue(events.get(3).isIncreasing());
1328         assertEquals(t5, events.get(4).getT(), tolerance);
1329         assertFalse(events.get(4).isIncreasing());
1330     }
1331 
1332     /**
1333      * Test case for two event detectors. DetectorA has event at t2, DetectorB at t3, but
1334      * due to the root finding tolerance DetectorB's event occurs at t1. With t1 < t2 <
1335      * t3.
1336      */
1337     @Test
1338     void testRootFindingToleranceReverse() {
1339         //setup
1340         double maxCheck = 10;
1341         double t2 = -11, t3 = t2 - 1e-5;
1342         List<Event> events = new ArrayList<>();
1343         TimeDetector detectorA = new TimeDetector(maxCheck, 1e-6, 100, events, t2);
1344         FlatDetector detectorB = new FlatDetector(maxCheck, 0.5, 100, events, t3);
1345         FieldODEIntegrator<Binary64> integrator =
1346                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1347         integrator.addEventDetector(detectorA);
1348         integrator.addEventDetector(detectorB);
1349 
1350         // action
1351         integrator.integrate(new Equation(), initialState, zero.add(-30.0));
1352 
1353         // verify
1354         // if these fail the event finding did its job,
1355         // but this test isn't testing what it is supposed to be
1356         assertSame(detectorB, events.get(0).getDetector());
1357         assertSame(detectorA, events.get(1).getDetector());
1358         assertTrue(events.get(0).getT() > events.get(1).getT());
1359 
1360         // check event detection worked
1361         assertEquals(2, events.size());
1362         assertEquals(t3, events.get(0).getT(), 0.5);
1363         assertTrue(events.get(0).isIncreasing());
1364         assertEquals(t2, events.get(1).getT(), 1e-6);
1365         assertTrue(events.get(1).isIncreasing());
1366     }
1367 
1368     /** check when g(t < root) < 0,  g(root + convergence) < 0. */
1369     @Test
1370     void testRootPlusToleranceHasWrongSignReverse() {
1371         // setup
1372         double maxCheck = 10;
1373         double tolerance = 1e-6;
1374         final double toleranceB = 0.3;
1375         double t1 = -11, t2 = -11.1, t3 = -11.2;
1376         // shared event list so we know the order in which they occurred
1377         List<Event> events = new ArrayList<>();
1378         TimeDetector detectorA = new TimeDetector(maxCheck, tolerance, 100, events, t2);
1379         TimeDetector detectorB = new TimeDetector(maxCheck, toleranceB, 100, events, -50, t1, t3);
1380         FieldODEIntegrator<Binary64> integrator =
1381                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1382         integrator.addEventDetector(detectorA);
1383         integrator.addEventDetector(detectorB);
1384 
1385         // action
1386         integrator.integrate(new Equation(), initialState, zero.add(-30.0));
1387 
1388         // verify
1389         // we only care that the rules are satisfied. There are multiple solutions.
1390         assertEquals(3, events.size());
1391         assertEquals(t1, events.get(0).getT(), toleranceB);
1392         assertTrue(events.get(0).isIncreasing());
1393         assertSame(detectorB, events.get(0).getDetector());
1394         assertEquals(t3, events.get(1).getT(), toleranceB);
1395         assertFalse(events.get(1).isIncreasing());
1396         assertSame(detectorB, events.get(1).getDetector());
1397         assertEquals(t2, events.get(2).getT(), tolerance);
1398         assertTrue(events.get(2).isIncreasing());
1399         assertSame(detectorA, events.get(2).getDetector());
1400         // ascending order
1401         assertTrue(events.get(0).getT() >= events.get(1).getT());
1402         assertTrue(events.get(1).getT() >= events.get(2).getT());
1403     }
1404 
1405     /** check when g(t < root) < 0,  g(root + convergence) < 0. */
1406     @Test
1407     void testRootPlusToleranceHasWrongSignAndLessThanTbReverse() {
1408         // setup
1409         // test is fragile w.r.t. implementation and these parameters
1410         double maxCheck = 10;
1411         double tolerance = 0.5;
1412         double t1 = -11, t2 = -11.4, t3 = -12.0;
1413         // shared event list so we know the order in which they occurred
1414         List<Event> events = new ArrayList<>();
1415         TimeDetector detectorB = new FlatDetector(maxCheck, tolerance, 100, events, t1, t2, t3);
1416         FieldODEIntegrator<Binary64> integrator =
1417                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1418         integrator.addEventDetector(detectorB);
1419 
1420         // action
1421         integrator.integrate(new Equation(), initialState, zero.add(-30.0));
1422 
1423         // verify
1424         // allowed to report t1 or t3.
1425         assertEquals(1, events.size());
1426         assertEquals(t1, events.get(0).getT(), tolerance);
1427         assertTrue(events.get(0).isIncreasing());
1428         assertSame(detectorB, events.get(0).getDetector());
1429     }
1430 
1431     /**
1432      * Check when g(t) has a multiple root. e.g. g(t < root) < 0, g(root) = 0, g(t > root)
1433      * < 0.
1434      */
1435     @Test
1436     void testDoubleRootReverse() {
1437         // setup
1438         double maxCheck = 10;
1439         double tolerance = 1e-6;
1440         double t1 = -11;
1441         // shared event list so we know the order in which they occurred
1442         List<Event> events = new ArrayList<>();
1443         TimeDetector detectorA = new TimeDetector(maxCheck, tolerance, 100, events, t1);
1444         TimeDetector detectorB = new TimeDetector(maxCheck, tolerance, 100, events, t1, t1);
1445         FieldODEIntegrator<Binary64> integrator =
1446                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1447         integrator.addEventDetector(detectorA);
1448         integrator.addEventDetector(detectorB);
1449 
1450         // action
1451         integrator.integrate(new Equation(), initialState, zero.add(-30.0));
1452 
1453         // verify
1454         assertEquals(1, events.size());
1455         assertEquals(t1, events.get(0).getT(), 0.0);
1456         assertTrue(events.get(0).isIncreasing());
1457         assertSame(detectorA, events.get(0).getDetector());
1458         // detector worked correctly
1459         assertEquals(0.0, detectorB.g(state(t1)).getReal());
1460         assertTrue(detectorB.g(state(t1 + 1e-6)).getReal() < 0);
1461         assertTrue(detectorB.g(state(t1 - 1e-6)).getReal() < 0);
1462     }
1463 
1464     /**
1465      * Check when g(t) has a multiple root. e.g. g(t < root) > 0, g(root) = 0, g(t > root)
1466      * > 0.
1467      */
1468     @Test
1469     void testDoubleRootOppositeSignReverse() {
1470         // setup
1471         double maxCheck = 10;
1472         double tolerance = 1e-6;
1473         double t1 = -11;
1474         // shared event list so we know the order in which they occurred
1475         List<Event> events = new ArrayList<>();
1476         TimeDetector detectorA = new TimeDetector(maxCheck, tolerance, 100, events, t1);
1477         TimeDetector detectorB = new ContinuousDetector(maxCheck, tolerance, 100, events, -50, t1, t1);
1478         detectorB.g(state(t1));
1479         FieldODEIntegrator<Binary64> integrator =
1480                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1481         integrator.addEventDetector(detectorA);
1482         integrator.addEventDetector(detectorB);
1483 
1484         // action
1485         integrator.integrate(new Equation(), initialState, zero.add(-30.0));
1486 
1487         // verify
1488         assertEquals(1, events.size());
1489         assertEquals(t1, events.get(0).getT(), 0.0);
1490         assertTrue(events.get(0).isIncreasing());
1491         assertSame(detectorA, events.get(0).getDetector());
1492         // detector worked correctly
1493         assertEquals(0.0, detectorB.g(state(t1)).getReal(), 0.0);
1494         assertTrue(detectorB.g(state(t1 + 1e-6)).getReal() > 0);
1495         assertTrue(detectorB.g(state(t1 - 1e-6)).getReal() > 0);
1496     }
1497 
1498     /** check root finding when zero at both ends. */
1499     @Test
1500     void testZeroAtBeginningAndEndOfIntervalReverse() {
1501         // setup
1502         double maxCheck = 10;
1503         double tolerance = 1e-6;
1504         double t1 = -10, t2 = -20;
1505         // shared event list so we know the order in which they occurred
1506         List<Event> events = new ArrayList<>();
1507         TimeDetector detectorA = new ContinuousDetector(maxCheck, tolerance, 100, events, -50, t1, t2);
1508         FieldODEIntegrator<Binary64> integrator =
1509                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1510         integrator.addEventDetector(detectorA);
1511 
1512         // action
1513         integrator.integrate(new Equation(), initialState, zero.add(-30.0));
1514 
1515         // verify
1516         assertEquals(2, events.size());
1517         assertEquals(t1, events.get(0).getT(), 0.0);
1518         assertTrue(events.get(0).isIncreasing());
1519         assertSame(detectorA, events.get(0).getDetector());
1520         assertEquals(t2, events.get(1).getT(), 0.0);
1521         assertFalse(events.get(1).isIncreasing());
1522         assertSame(detectorA, events.get(1).getDetector());
1523     }
1524 
1525     /** check root finding when zero at both ends. */
1526     @Test
1527     void testZeroAtBeginningAndEndOfIntervalOppositeSignReverse() {
1528         // setup
1529         double maxCheck = 10;
1530         double tolerance = 1e-6;
1531         double t1 = -10, t2 = -20;
1532         // shared event list so we know the order in which they occurred
1533         List<Event> events = new ArrayList<>();
1534         TimeDetector detectorA = new ContinuousDetector(maxCheck, tolerance, 100, events, t1, t2);
1535         FieldODEIntegrator<Binary64> integrator =
1536                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1537         integrator.addEventDetector(detectorA);
1538 
1539         // action
1540         integrator.integrate(new Equation(), initialState, zero.add(-30.0));
1541 
1542         // verify
1543         assertEquals(2, events.size());
1544         assertEquals(t1, events.get(0).getT(), 0.0);
1545         assertFalse(events.get(0).isIncreasing());
1546         assertSame(detectorA, events.get(0).getDetector());
1547         assertEquals(t2, events.get(1).getT(), 0.0);
1548         assertTrue(events.get(1).isIncreasing());
1549         assertSame(detectorA, events.get(1).getDetector());
1550     }
1551 
1552     /** Test where an event detector has to back up multiple times. */
1553     @Test
1554     void testMultipleBackupsReverse() {
1555         // setup
1556         double maxCheck = 10;
1557         double tolerance = 1e-6;
1558         double t1 = -11.0, t2 = -12, t3 = -13, t4 = -14, t5 = -15, t6 = -16, t7 = -17;
1559         // shared event list so we know the order in which they occurred
1560         List<Event> events = new ArrayList<>();
1561         TimeDetector detectorA = new ContinuousDetector(maxCheck, tolerance, 100, events, t6);
1562         TimeDetector detectorB = new ContinuousDetector(maxCheck, tolerance, 100, events, -50, t1, t3, t4, t7);
1563         TimeDetector detectorC = new ContinuousDetector(maxCheck, tolerance, 100, events, -50, t2, t5);
1564 
1565         FieldODEIntegrator<Binary64> integrator =
1566                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1567         integrator.addEventDetector(detectorA);
1568         integrator.addEventDetector(detectorB);
1569         integrator.addEventDetector(detectorC);
1570 
1571         // action
1572         integrator.integrate(new Equation(), initialState, zero.add(-30.0));
1573 
1574         //verify
1575         // really we only care that the Rules of Event Handling are not violated,
1576         assertEquals(5, events.size());
1577         assertEquals(t1, events.get(0).getT(), tolerance);
1578         assertTrue(events.get(0).isIncreasing());
1579         assertEquals(detectorB, events.get(0).getDetector());
1580         assertEquals(t2, events.get(1).getT(), tolerance);
1581         assertTrue(events.get(1).isIncreasing());
1582         assertEquals(detectorC, events.get(1).getDetector());
1583         // reporting t3 and t4 is optional, seeing them is not.
1584         // we know a root was found at t3 because events are reported at t2 and t5.
1585         /*
1586         Assertions.assertEquals(t3, events.get(2).getT(), tolerance);
1587         Assertions.assertEquals(false, events.get(2).isIncreasing());
1588         Assertions.assertEquals(detectorB, events.get(2).getDetector());
1589         Assertions.assertEquals(t4, events.get(3).getT(), tolerance);
1590         Assertions.assertEquals(true, events.get(3).isIncreasing());
1591         Assertions.assertEquals(detectorB, events.get(3).getDetector());
1592         */
1593         assertEquals(t5, events.get(2).getT(), tolerance);
1594         assertFalse(events.get(2).isIncreasing());
1595         assertEquals(detectorC, events.get(2).getDetector());
1596         assertEquals(t6, events.get(3).getT(), tolerance);
1597         assertTrue(events.get(3).isIncreasing());
1598         assertEquals(detectorA, events.get(3).getDetector());
1599         assertEquals(t7, events.get(4).getT(), tolerance);
1600         assertFalse(events.get(4).isIncreasing());
1601         assertEquals(detectorB, events.get(4).getDetector());
1602     }
1603 
1604     /** Test a reset event triggering another event at the same time. */
1605     @Test
1606     void testEventCausedByStateResetReverse() {
1607         // setup
1608         double maxCheck = 10;
1609         double tolerance = 1e-6;
1610         double t1 = -15.0;
1611         final FieldODEState<Binary64> newState =
1612                 new FieldODEState<>(zero.add(t1), new Binary64[]{zero.add(20), zero});
1613         // shared event list so we know the order in which they occurred
1614         List<Event> events = new ArrayList<>();
1615         TimeDetector detectorA = new ResetDetector(maxCheck, tolerance, 100, events, newState, t1);
1616         BaseDetector detectorB = new StateDetector(maxCheck, tolerance, 100, events, 1);
1617 
1618         FieldODEIntegrator<Binary64> integrator =
1619                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1620         integrator.addEventDetector(detectorA);
1621         integrator.addEventDetector(detectorB);
1622 
1623         // action
1624         integrator.integrate(new Equation(), initialState, zero.add(-40.0));
1625 
1626         //verify
1627         // really we only care that the Rules of Event Handling are not violated,
1628         assertEquals(3, events.size());
1629         assertEquals(t1, events.get(0).getT(), tolerance);
1630         assertTrue(events.get(0).isIncreasing());
1631         assertEquals(detectorA, events.get(0).getDetector());
1632         assertEquals(t1, events.get(1).getT(), tolerance);
1633         assertFalse(events.get(1).isIncreasing());
1634         assertEquals(detectorB, events.get(1).getDetector());
1635         assertEquals(t1 - 19, events.get(2).getT(), tolerance);
1636         assertTrue(events.get(2).isIncreasing());
1637         assertEquals(detectorB, events.get(2).getDetector());
1638     }
1639 
1640     /** check when t + tolerance == t. */
1641     @Test
1642     void testConvergenceTooTightReverse() {
1643         // setup
1644         double maxCheck = 10;
1645         double tolerance = 1e-18;
1646         double t1 = -15;
1647         // shared event list so we know the order in which they occurred
1648         List<Event> events = new ArrayList<>();
1649         TimeDetector detectorA = new ContinuousDetector(maxCheck, tolerance, 100, events, t1);
1650         FieldODEIntegrator<Binary64> integrator =
1651                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1652         integrator.addEventDetector(detectorA);
1653 
1654         // action
1655         integrator.integrate(new Equation(), initialState, zero.add(-30.0));
1656 
1657         // verify
1658         assertEquals(1, events.size());
1659         assertEquals(t1, events.get(0).getT(), 0.0);
1660         assertTrue(events.get(0).isIncreasing());
1661         assertSame(detectorA, events.get(0).getDetector());
1662     }
1663 
1664     /** check when root finding tolerance > event finding tolerance. */
1665     @Test
1666     void testToleranceMismatchReverse() {
1667         // setup
1668         double maxCheck = 10;
1669         double tolerance = 1e-18;
1670         double t1 = -15.1;
1671         // shared event list so we know the order in which they occurred
1672         List<Event> events = new ArrayList<>();
1673         TimeDetector detectorA = new FlatDetector(maxCheck, tolerance, 100, events, t1);
1674         FieldODEIntegrator<Binary64> integrator =
1675                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1676         integrator.addEventDetector(detectorA);
1677 
1678         // action
1679         integrator.integrate(new Equation(), initialState, zero.add(-30.0));
1680 
1681         // verify
1682         assertEquals(1, events.size());
1683         // use root finding tolerance since it is larger
1684         assertEquals(t1, events.get(0).getT(), 1e-3);
1685         assertTrue(events.get(0).isIncreasing());
1686         assertSame(detectorA, events.get(0).getDetector());
1687     }
1688 
1689     /**
1690      * test when one event detector changes the definition of another's g function before
1691      * the end of the step as a result of a continue action. Not sure if this should be
1692      * officially supported, but it is used in Orekit's DateDetector, it's useful, and not
1693      * too hard to implement.
1694      */
1695     @Test
1696     void testEventChangesGFunctionDefinitionReverse() {
1697         // setup
1698         double maxCheck = 5;
1699         double tolerance = 1e-6;
1700         double t1 = -11, t2 = -19;
1701         // shared event list so we know the order in which they occurred
1702         List<Event> events = new ArrayList<>();
1703         // mutable boolean
1704         boolean[] swap = new boolean[1];
1705         final TimeDetector detectorA = new ContinuousDetector(maxCheck, tolerance, 100, events, t1) {
1706             @Override
1707             public FieldODEEventHandler<Binary64> getHandler() {
1708                 return (state, detector, increasing) -> {
1709                     swap[0] = true;
1710                     return super.getHandler().eventOccurred(state, detector, increasing);
1711                 };
1712             }
1713         };
1714         final TimeDetector detectorB = new ContinuousDetector(maxCheck, tolerance, 100, events, t2);
1715         BaseDetector detectorC = new BaseDetector(maxCheck, tolerance, 100, Action.CONTINUE, events) {
1716 
1717             @Override
1718             public Binary64 g(FieldODEStateAndDerivative<Binary64> state) {
1719                 if (swap[0]) {
1720                     return detectorB.g(state);
1721                 } else {
1722                     return one;
1723                 }
1724             }
1725 
1726         };
1727         FieldODEIntegrator<Binary64> integrator =
1728                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1729         integrator.addEventDetector(detectorA);
1730         integrator.addEventDetector(detectorC);
1731 
1732         // action
1733         integrator.integrate(new Equation(), initialState, zero.add(-30.0));
1734 
1735         // verify
1736         assertEquals(2, events.size());
1737         assertEquals(t1, events.get(0).getT(), tolerance);
1738         assertTrue(events.get(0).isIncreasing());
1739         assertSame(detectorA, events.get(0).getDetector());
1740         assertEquals(t2, events.get(1).getT(), tolerance);
1741         assertTrue(events.get(1).isIncreasing());
1742         assertSame(detectorC, events.get(1).getDetector());
1743     }
1744 
1745 
1746     /**
1747      * test when one event detector changes the definition of another's g function before
1748      * the end of the step as a result of an event occurring. In this case the change
1749      * cancels the occurrence of the event.
1750      */
1751     @Test
1752     void testEventChangesGFunctionDefinitionCancelReverse() {
1753         // setup
1754         double maxCheck = 5;
1755         double tolerance = 1e-6;
1756         double t1 = -11, t2 = -11.1;
1757         // shared event list so we know the order in which they occurred
1758         List<Event> events = new ArrayList<>();
1759         // mutable boolean
1760         boolean[] swap = new boolean[1];
1761         final TimeDetector detectorA =
1762                         new ContinuousDetector(maxCheck, tolerance, 100, Action.RESET_EVENTS, events, t1) {
1763             @Override
1764             public FieldODEEventHandler<Binary64> getHandler() {
1765                 return (state, detector, increasing) -> {
1766                     swap[0] = true;
1767                     return super.getHandler().eventOccurred(state, detector, increasing);
1768                 };
1769             }
1770         };
1771         final TimeDetector detectorB = new ContinuousDetector(maxCheck, tolerance, 100, events, t2);
1772         BaseDetector detectorC = new BaseDetector(maxCheck, tolerance, 100, Action.CONTINUE, events) {
1773 
1774             @Override
1775             public Binary64 g(FieldODEStateAndDerivative<Binary64> state) {
1776                 if (!swap[0]) {
1777                     return detectorB.g(state);
1778                 } else {
1779                     return zero.add(1);
1780                 }
1781             }
1782 
1783         };
1784         FieldODEIntegrator<Binary64> integrator =
1785                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1786         integrator.addEventDetector(detectorA);
1787         integrator.addEventDetector(detectorC);
1788 
1789         // action
1790         integrator.integrate(new Equation(), initialState, zero.add(-30.0));
1791 
1792         // verify
1793         assertEquals(1, events.size());
1794         assertEquals(t1, events.get(0).getT(), tolerance);
1795         assertTrue(events.get(0).isIncreasing());
1796         assertSame(detectorA, events.get(0).getDetector());
1797     }
1798 
1799 
1800     /**
1801      * test when one event detector changes the definition of another's g function before
1802      * the end of the step as a result of an event occurring. In this case the change
1803      * delays the occurrence of the event.
1804      */
1805     @Test
1806     void testEventChangesGFunctionDefinitionDelayReverse() {
1807         // setup
1808         double maxCheck = 5;
1809         double tolerance = 1e-6;
1810         double t1 = -11, t2 = -11.1, t3 = -11.2;
1811         // shared event list so we know the order in which they occurred
1812         List<Event> events = new ArrayList<>();
1813         // mutable boolean
1814         boolean[] swap = new boolean[1];
1815         final TimeDetector detectorA =
1816                         new ContinuousDetector(maxCheck, tolerance, 100, Action.RESET_EVENTS, events, t1) {
1817             @Override
1818             public FieldODEEventHandler<Binary64> getHandler() {
1819                 return (state, detector, increasing) -> {
1820                     swap[0] = true;
1821                     return super.getHandler().eventOccurred(state, detector, increasing);
1822                 };
1823             }
1824         };
1825         final TimeDetector detectorB = new ContinuousDetector(maxCheck, tolerance, 100, events, t2);
1826         final TimeDetector detectorD = new ContinuousDetector(maxCheck, tolerance, 100, events, t3);
1827         BaseDetector detectorC = new BaseDetector(maxCheck, tolerance, 100, Action.CONTINUE, events) {
1828 
1829             @Override
1830             public Binary64 g(FieldODEStateAndDerivative<Binary64> state) {
1831                 if (!swap[0]) {
1832                     return detectorB.g(state);
1833                 } else {
1834                     return detectorD.g(state);
1835                 }
1836             }
1837 
1838         };
1839         FieldODEIntegrator<Binary64> integrator =
1840                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1841         integrator.addEventDetector(detectorA);
1842         integrator.addEventDetector(detectorC);
1843 
1844         // action
1845         integrator.integrate(new Equation(), initialState, zero.add(-30.0));
1846 
1847         // verify
1848         assertEquals(2, events.size());
1849         assertEquals(t1, events.get(0).getT(), tolerance);
1850         assertTrue(events.get(0).isIncreasing());
1851         assertSame(detectorA, events.get(0).getDetector());
1852         assertEquals(t3, events.get(1).getT(), tolerance);
1853         assertTrue(events.get(1).isIncreasing());
1854         assertSame(detectorC, events.get(1).getDetector());
1855     }
1856 
1857     /**
1858      * test when one event detector changes the definition of another's g function before
1859      * the end of the step as a result of an event occurring. In this case the change
1860      * causes the event to happen sooner than originally expected.
1861      */
1862     @Test
1863     void testEventChangesGFunctionDefinitionAccelerateReverse() {
1864         // setup
1865         double maxCheck = 5;
1866         double tolerance = 1e-6;
1867         double t1 = -11, t2 = -11.1, t3 = -11.2;
1868         // shared event list so we know the order in which they occurred
1869         List<Event> events = new ArrayList<>();
1870         // mutable boolean
1871         boolean[] swap = new boolean[1];
1872         final TimeDetector detectorA =
1873                         new ContinuousDetector(maxCheck, tolerance, 100, Action.RESET_EVENTS, events, t1) {
1874             @Override
1875             public FieldODEEventHandler<Binary64> getHandler() {
1876                 return (state, detector, increasing) -> {
1877                     swap[0] = true;
1878                     return super.getHandler().eventOccurred(state, detector, increasing);
1879                 };
1880             }
1881         };
1882         final TimeDetector detectorB = new ContinuousDetector(maxCheck, tolerance, 100, events, t2);
1883         final TimeDetector detectorD = new ContinuousDetector(maxCheck, tolerance, 100, events, t3);
1884         BaseDetector detectorC = new BaseDetector(maxCheck, tolerance, 100, Action.CONTINUE, events) {
1885 
1886             @Override
1887             public Binary64 g(FieldODEStateAndDerivative<Binary64> state) {
1888                 if (swap[0]) {
1889                     return detectorB.g(state);
1890                 } else {
1891                     return detectorD.g(state);
1892                 }
1893             }
1894 
1895         };
1896         FieldODEIntegrator<Binary64> integrator =
1897                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1898         integrator.addEventDetector(detectorA);
1899         integrator.addEventDetector(detectorC);
1900 
1901         // action
1902         integrator.integrate(new Equation(), initialState, zero.add(-30.0));
1903 
1904         // verify
1905         assertEquals(2, events.size());
1906         assertEquals(t1, events.get(0).getT(), tolerance);
1907         assertTrue(events.get(0).isIncreasing());
1908         assertSame(detectorA, events.get(0).getDetector());
1909         assertEquals(t2, events.get(1).getT(), tolerance);
1910         assertTrue(events.get(1).isIncreasing());
1911         assertSame(detectorC, events.get(1).getDetector());
1912     }
1913 
1914     /** check when root finding tolerance > event finding tolerance. */
1915     @Test
1916     void testToleranceStopReverse() {
1917         // setup
1918         double maxCheck = 10;
1919         double tolerance = 1e-18; // less than 1 ulp
1920         double t1 = -15.1;
1921         // shared event list so we know the order in which they occurred
1922         List<Event> events = new ArrayList<>();
1923         TimeDetector detectorA = new FlatDetector(maxCheck, tolerance, 100, Action.STOP, events, t1);
1924         FieldODEIntegrator<Binary64> integrator =
1925                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1926         integrator.addEventDetector(detectorA);
1927 
1928         // action
1929         FieldODEStateAndDerivative<Binary64> finalState =
1930                 integrator.integrate(new Equation(), initialState, zero.add(-30.0));
1931 
1932         // verify
1933         assertEquals(1, events.size());
1934         // use root finder tolerance instead of event finder tolerance.
1935         assertEquals(t1, events.get(0).getT(), tolerance);
1936         assertTrue(events.get(0).isIncreasing());
1937         assertSame(detectorA, events.get(0).getDetector());
1938         assertEquals(t1, finalState.getTime().getReal(), tolerance);
1939 
1940         // try to resume propagation
1941         finalState = integrator.integrate(new Equation(), finalState, zero.add(-30.0));
1942 
1943         // verify it got to the end
1944         assertEquals(-30.0, finalState.getTime().getReal(), 0.0);
1945     }
1946 
1947     /**
1948      * Test when g function is initially zero for longer than the tolerance. Can occur
1949      * when restarting after a stop and cancellation occurs in the g function.
1950      */
1951     @Test
1952     void testLongInitialZeroReverse() {
1953         // setup
1954         double maxCheck = 10;
1955         double tolerance = 1;
1956         // shared event list so we know the order in which they occurred
1957         List<Event> events = new ArrayList<>();
1958         TimeDetector detectorA = new TimeDetector(maxCheck, tolerance, 100, Action.STOP, events, 50) {
1959             @Override
1960             public Binary64 g(FieldODEStateAndDerivative<Binary64> state) {
1961                 if (state.getTime().getReal() > -2) {
1962                     return zero;
1963                 } else {
1964                     return super.g(state);
1965                 }
1966             }
1967         };
1968         FieldODEIntegrator<Binary64> integrator =
1969                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
1970         integrator.addEventDetector(detectorA);
1971 
1972         // action
1973         integrator.integrate(new Equation(), initialState, zero.add(-30.0));
1974 
1975         // verify
1976         assertEquals(0, events.size());
1977     }
1978 
1979     /**
1980      * The root finder requires the start point to be in the interval (a, b) which is hard
1981      * when there aren't many numbers between a and b. This test uses a second event
1982      * detector to force a very small window for the first event detector.
1983      */
1984     @Test
1985     void testShortBracketingIntervalReverse() {
1986         // setup
1987         double maxCheck = 10;
1988         double tolerance = 1e-6;
1989         final double t1 = FastMath.nextDown(-10.0), t2 = -10.5;
1990         // shared event list so we know the order in which they occurred
1991         List<Event> events = new ArrayList<>();
1992         // never zero so there is no easy way out
1993         TimeDetector detectorA = new TimeDetector(maxCheck, tolerance, 100, events) {
1994             @Override
1995             public Binary64 g(FieldODEStateAndDerivative<Binary64> state) {
1996                 final Binary64 t = state.getTime();
1997                 if (t.getReal() > t1) {
1998                     return one.negate();
1999                 } else if (t.getReal() > t2) {
2000                     return one;
2001                 } else {
2002                     return one.negate();
2003                 }
2004             }
2005         };
2006         TimeDetector detectorB = new TimeDetector(maxCheck, tolerance, 100, events, t1);
2007         FieldODEIntegrator<Binary64> integrator =
2008                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
2009         integrator.addEventDetector(detectorA);
2010         integrator.addEventDetector(detectorB);
2011 
2012         // action
2013         integrator.integrate(new Equation(), initialState, zero.add(-30.0));
2014 
2015         // verify
2016         assertEquals(3, events.size());
2017         assertEquals(t1, events.get(0).getT(), tolerance);
2018         assertFalse(events.get(0).isIncreasing());
2019         assertSame(detectorA, events.get(0).getDetector());
2020         assertEquals(t1, events.get(1).getT(), tolerance);
2021         assertTrue(events.get(1).isIncreasing());
2022         assertSame(detectorB, events.get(1).getDetector());
2023         assertEquals(t2, events.get(2).getT(), tolerance);
2024         assertTrue(events.get(2).isIncreasing());
2025         assertSame(detectorA, events.get(2).getDetector());
2026     }
2027 
2028     /** Check that steps are restricted correctly with a continue event. */
2029     @Test
2030     void testEventStepHandlerReverse() {
2031         // setup
2032         double tolerance = 1e-18;
2033         FieldODEIntegrator<Binary64> integrator =
2034                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
2035         integrator.addEventDetector(new TimeDetector(100.0, tolerance, 100, -5));
2036         StepHandler stepHandler = new StepHandler();
2037         integrator.addStepHandler(stepHandler);
2038 
2039         // action
2040         FieldODEStateAndDerivative<Binary64> finalState = integrator
2041                 .integrate(new Equation(), initialState, zero.add(-10));
2042 
2043         // verify
2044         assertEquals(-10.0, finalState.getTime().getReal(), tolerance);
2045         assertEquals(0.0,
2046                 stepHandler.initialState.getTime().getReal(), tolerance);
2047         assertEquals(-10.0, stepHandler.finalTime.getReal(), tolerance);
2048         assertEquals(-10.0,
2049                 stepHandler.finalState.getTime().getReal(), tolerance);
2050         FieldODEStateInterpolator<Binary64> interpolator = stepHandler.interpolators.get(0);
2051         assertEquals(0.0,
2052                 interpolator.getPreviousState().getTime().getReal(), tolerance);
2053         assertEquals(-5.0,
2054                 interpolator.getCurrentState().getTime().getReal(), tolerance);
2055         interpolator = stepHandler.interpolators.get(1);
2056         assertEquals(-5.0,
2057                 interpolator.getPreviousState().getTime().getReal(), tolerance);
2058         assertEquals(-10.0,
2059                 interpolator.getCurrentState().getTime().getReal(), tolerance);
2060         assertEquals(2, stepHandler.interpolators.size());
2061     }
2062 
2063     /** Test resetState(...) returns {@code null}. */
2064     @Test
2065     void testEventCausedByDerivativesResetReverse() {
2066         // setup
2067         TimeDetector detectorA = new TimeDetector(10, 1e-6, 100, Action.RESET_STATE, -15.0) {
2068             @Override
2069             public FieldODEEventHandler<Binary64> getHandler() {
2070                 return new FieldODEEventHandler<Binary64>() {
2071                     @Override
2072                     public Action eventOccurred(FieldODEStateAndDerivative<Binary64> state,
2073                                                 FieldODEEventDetector<Binary64> detector,
2074                                                 boolean increasing) {
2075                         return Action.RESET_STATE;
2076                     }
2077                     @Override
2078                     public FieldODEState<Binary64> resetState(FieldODEEventDetector<Binary64> detector,
2079                                                                FieldODEStateAndDerivative<Binary64> state) {
2080                         return null;
2081                     }
2082                 };
2083             }
2084         };
2085         FieldODEIntegrator<Binary64> integrator =
2086                 new DormandPrince853FieldIntegrator<>(field, 10, 10, 1e-7, 1e-7);
2087         integrator.addEventDetector(detectorA);
2088 
2089         try {
2090             // action
2091             integrator.integrate(new Equation(), initialState, zero.add(-20.0));
2092             fail("Expected Exception");
2093         } catch (NullPointerException e) {
2094             // expected
2095         }
2096     }
2097 
2098     @Test
2099     void testResetChangesSignReverse() {
2100         FieldOrdinaryDifferentialEquation<Binary64> equation = new FieldOrdinaryDifferentialEquation<Binary64>() {
2101             public int getDimension() { return 1; }
2102             public Binary64[] computeDerivatives(Binary64 t, Binary64[] y) { return new Binary64[] { new Binary64(1.0) }; }
2103         };
2104 
2105         LutherFieldIntegrator<Binary64> integrator = new LutherFieldIntegrator<>(Binary64Field.getInstance(), new Binary64(20.0));
2106         final double small = 1.0e-10;
2107         ResetChangesSignGenerator eventsGenerator = new ResetChangesSignGenerator(-6.0, -9.0, +0.5 * small, 8.0, small, 1000);
2108         integrator.addEventDetector(eventsGenerator);
2109         final FieldODEStateAndDerivative<Binary64> end = integrator.integrate(new FieldExpandableODE<>(equation),
2110                                                                                new FieldODEState<>(new Binary64(0.0),
2111                                                                                                    new Binary64[] { new Binary64(0.0) }),
2112                                                                                new Binary64(-100.0));
2113         assertEquals(2,                  eventsGenerator.getCount());
2114         assertEquals(-9.0,               end.getCompleteState()[0].getReal(), 1.0e-12);
2115         assertEquals(-9.0 - 0.5 * small, end.getTime().getReal(),             1.0e-12);
2116     }
2117 
2118     /* utility classes and methods */
2119 
2120     /**
2121      * Create a state at a time.
2122      *
2123      * @param t time of state.
2124      * @return new state.
2125      */
2126     private FieldODEStateAndDerivative<Binary64> state(double t) {
2127         return new FieldODEStateAndDerivative<>(
2128                 zero.add(t), new Binary64[0], new Binary64[0]);
2129     }
2130 
2131     /** Base class to record events that occurred. */
2132     private static abstract class BaseDetector implements FieldODEEventDetector<Binary64> {
2133 
2134         private final FieldAdaptableInterval<Binary64>             maxCheck;
2135         private final int                                          maxIter;
2136         private final BracketedRealFieldUnivariateSolver<Binary64> solver;
2137         protected final Action                                     action;
2138 
2139         /** times the event was actually triggered. */
2140         private final List<Event> events;
2141 
2142         public BaseDetector(final double maxCheck, final double threshold, final int maxIter,
2143                             Action action, List<Event> events) {
2144             this.maxCheck  = (s, isForward) -> maxCheck;
2145             this.maxIter   = maxIter;
2146             this.solver    = new FieldBracketingNthOrderBrentSolver<>(new Binary64(0),
2147                                                                       new Binary64(threshold),
2148                                                                       new Binary64(0),
2149                                                                       5);
2150             this.action    = action;
2151             this.events    = events;
2152         }
2153 
2154         public FieldAdaptableInterval<Binary64> getMaxCheckInterval() {
2155             return maxCheck;
2156         }
2157 
2158         public int getMaxIterationCount() {
2159             return maxIter;
2160         }
2161 
2162         public BracketedRealFieldUnivariateSolver<Binary64> getSolver() {
2163             return solver;
2164         }
2165 
2166         /**
2167          * Get all events that occurred.
2168          *
2169          * @return list of events in the order they occurred.
2170          */
2171         public List<Event> getEvents() {
2172             return events;
2173         }
2174 
2175         @Override
2176         public FieldODEEventHandler<Binary64> getHandler() {
2177             return new FieldODEEventHandler<Binary64>() {
2178                 @Override
2179                 public Action eventOccurred(FieldODEStateAndDerivative<Binary64> state,
2180                                             FieldODEEventDetector<Binary64> detector,
2181                                             boolean increasing) {
2182                     events.add(new Event(state, detector, increasing));
2183                     return action;
2184                 }
2185             };
2186         }
2187 
2188     }
2189 
2190     /** Trigger an event at a given time interval. */
2191     private static class IntervalDetector extends BaseDetector {
2192 
2193         /** Interval start. */
2194         protected final double start;
2195 
2196         /** Interval end. */
2197         protected final double end;
2198 
2199         public IntervalDetector(final double maxCheck, final double threshold, final int maxIter,
2200                                 final Action action, final double startTime, final double endTime) {
2201             super(maxCheck, threshold, maxIter, action, new ArrayList<>());
2202             this.start = startTime;
2203             this.end = endTime;
2204         }
2205 
2206         @Override
2207         public Binary64 g(final FieldODEStateAndDerivative<Binary64> state) {
2208             final Binary64 t = state.getTime();
2209             return t.subtract(start).multiply(t.negate().add(end));
2210         }
2211 
2212     }
2213 
2214     /** Trigger an event at a particular time. */
2215     private static class TimeDetector extends BaseDetector {
2216 
2217         /** time of the event to trigger. */
2218         protected final double[] eventTs;
2219 
2220         /**
2221          * Create a new detector.
2222          *
2223          * @param eventTs the time to trigger an event.
2224          */
2225         public TimeDetector(final double maxCheck, final double threshold, final int maxIter,
2226                             final double... eventTs) {
2227             this(maxCheck, threshold, maxIter, Action.CONTINUE, eventTs);
2228         }
2229 
2230         public TimeDetector(final double maxCheck, final double threshold, final int maxIter,
2231                             final List<Event> events, final double... eventTs) {
2232             this(maxCheck, threshold, maxIter, Action.CONTINUE, events, eventTs);
2233         }
2234 
2235         public TimeDetector(final double maxCheck, final double threshold, final int maxIter,
2236                             final Action action, final double... eventTs) {
2237             this(maxCheck, threshold, maxIter, action, new ArrayList<Event>(), eventTs);
2238         }
2239 
2240         public TimeDetector(final double maxCheck, final double threshold, final int maxIter,
2241                             final Action action, final List<Event> events, final double... eventTs) {
2242             super(maxCheck, threshold, maxIter, action, events);
2243             this.eventTs = eventTs.clone();
2244             Arrays.sort(this.eventTs);
2245         }
2246 
2247         @Override
2248         public Binary64 g(final FieldODEStateAndDerivative<Binary64> state) {
2249             final Binary64 t = state.getTime();
2250             int i = 0;
2251             while (i < eventTs.length && t.getReal() > eventTs[i]) {
2252                 i++;
2253             }
2254             i--;
2255             if (i < 0) {
2256                 return t.subtract(eventTs[0]);
2257             } else {
2258                 int sign = (i % 2) * 2 - 1;
2259                 return t.subtract(eventTs[i]).multiply(-sign);
2260             }
2261         }
2262 
2263     }
2264 
2265     private static class Event {
2266 
2267         private final FieldODEStateAndDerivative<Binary64> state;
2268         private final boolean increasing;
2269         private final FieldODEEventDetector<Binary64> detector;
2270 
2271         public Event(FieldODEStateAndDerivative<Binary64> state,
2272                      FieldODEEventDetector<Binary64> detector,
2273                      boolean increasing) {
2274             this.increasing = increasing;
2275             this.state = state;
2276             this.detector = detector;
2277         }
2278 
2279         public boolean isIncreasing() {
2280             return increasing;
2281         }
2282 
2283         public double getT() {
2284             return state.getTime().getReal();
2285         }
2286 
2287         public FieldODEEventDetector<Binary64> getDetector() {
2288             return detector;
2289         }
2290 
2291         @Override
2292         public String toString() {
2293             return "Event{" +
2294                     "increasing=" + increasing +
2295                     ", time=" + state.getTime() +
2296                     ", state=" + Arrays.toString(state.getCompleteState()) +
2297                     '}';
2298         }
2299     }
2300 
2301     /**
2302      * Same as {@link TimeDetector} except that it has a very flat g function which makes
2303      * root finding hard.
2304      */
2305     private static class FlatDetector extends TimeDetector {
2306 
2307         public FlatDetector(final double maxCheck, final double threshold, final int maxIter,
2308                             final List<Event> events, final double... eventTs) {
2309             super(maxCheck, threshold, maxIter, events, eventTs);
2310         }
2311 
2312         public FlatDetector(final double maxCheck, final double threshold, final int maxIter,
2313                             final Action action, final List<Event> events, final double... eventTs) {
2314             super(maxCheck, threshold, maxIter, action, events, eventTs);
2315         }
2316 
2317         @Override
2318         public Binary64 g(FieldODEStateAndDerivative<Binary64> state) {
2319             final Binary64 g = super.g(state);
2320             return g.sign();
2321         }
2322 
2323     }
2324 
2325     /** Linear on both ends, parabolic in the middle. */
2326     private static class ContinuousDetector extends TimeDetector {
2327 
2328         public ContinuousDetector(final double maxCheck, final double threshold, final int maxIter,
2329                                   final List<Event> events, final double... eventTs) {
2330             super(maxCheck, threshold, maxIter, events, eventTs);
2331         }
2332 
2333         public ContinuousDetector(final double maxCheck, final double threshold, final int maxIter,
2334                                   final Action action, final List<Event> events, final double... eventTs) {
2335             super(maxCheck, threshold, maxIter, action, events, eventTs);
2336         }
2337 
2338         @Override
2339         public Binary64 g(FieldODEStateAndDerivative<Binary64> state) {
2340             final Binary64 t = state.getTime();
2341             int i = 0;
2342             while (i < eventTs.length && t.getReal() > eventTs[i]) {
2343                 i++;
2344             }
2345             i--;
2346             if (i < 0) {
2347                 return t.subtract(eventTs[0]);
2348             } else if (i < eventTs.length - 1) {
2349                 int sign = (i % 2) * 2 - 1;
2350                 return t.subtract(eventTs[i]).multiply(t.negate().add(eventTs[i + 1])).multiply(-sign);
2351             } else {
2352                 int sign = (i % 2) * 2 - 1;
2353                 return t.subtract(eventTs[i]).multiply(-sign);
2354             }
2355         }
2356     }
2357 
2358     /** Reset the state at a particular time. */
2359     private static class ResetDetector extends TimeDetector {
2360 
2361         private final FieldODEState<Binary64> resetState;
2362 
2363         public ResetDetector(final double maxCheck, final double threshold, final int maxIter,
2364                              final List<Event> events, final FieldODEState<Binary64> state, final double eventT) {
2365             super(maxCheck, threshold, maxIter, Action.RESET_STATE, events, eventT);
2366             this.resetState = state;
2367         }
2368 
2369         @Override
2370         public FieldODEEventHandler<Binary64> getHandler() {
2371             return new FieldODEEventHandler<Binary64>() {
2372                 @Override
2373                 public Action eventOccurred(FieldODEStateAndDerivative<Binary64> state,
2374                                             FieldODEEventDetector<Binary64> detector, boolean increasing) {
2375                     return ResetDetector.super.getHandler().eventOccurred(state, detector, increasing);
2376                 }
2377                 @Override
2378                 public FieldODEState<Binary64> resetState(FieldODEEventDetector<Binary64> detector,
2379                                                            FieldODEStateAndDerivative<Binary64> state) {
2380                     assertEquals(eventTs[0], state.getTime().getReal(), 0);
2381                     return resetState;
2382                 }
2383             };
2384         }
2385 
2386     }
2387 
2388     /** Switching function based on the first primary state. */
2389     private static class StateDetector extends BaseDetector {
2390 
2391         private final double triggerState;
2392 
2393         public StateDetector(final double maxCheck, final double threshold, final int maxIter,
2394                              final List<Event> events, final double triggerState) {
2395             super(maxCheck, threshold, maxIter, Action.CONTINUE, events);
2396             this.triggerState = triggerState;
2397         }
2398 
2399         @Override
2400         public Binary64 g(FieldODEStateAndDerivative<Binary64> state) {
2401             return state.getPrimaryState()[0].subtract(this.triggerState);
2402         }
2403     }
2404 
2405     /** Some basic equations to integrate. */
2406     public static class Equation extends FieldExpandableODE<Binary64> {
2407         public Equation() {
2408             super(new EquationODE());
2409         }
2410     }
2411 
2412     /** Some basic equations to integrate. */
2413     public static class EquationODE implements FieldOrdinaryDifferentialEquation<Binary64> {
2414 
2415         public int getDimension() {
2416             return 2;
2417         }
2418 
2419         @Override
2420         public Binary64[] computeDerivatives(Binary64 t, Binary64[] y) {
2421             return new Binary64[]{one, one.multiply(2)};
2422         }
2423 
2424     }
2425 
2426     private static class StepHandler implements FieldODEStepHandler<Binary64> {
2427 
2428         private FieldODEStateAndDerivative<Binary64> initialState;
2429         private Binary64 finalTime;
2430         private List<FieldODEStateInterpolator<Binary64>> interpolators = new ArrayList<>();
2431         private FieldODEStateAndDerivative<Binary64> finalState;
2432 
2433         @Override
2434         public void init(FieldODEStateAndDerivative<Binary64> initialState,
2435                          Binary64 finalTime) {
2436             this.initialState = initialState;
2437             this.finalTime = finalTime;
2438         }
2439 
2440         @Override
2441         public void handleStep(FieldODEStateInterpolator<Binary64> interpolator) {
2442             this.interpolators.add(interpolator);
2443         }
2444 
2445         @Override
2446         public void finish(FieldODEStateAndDerivative<Binary64> finalState) {
2447             this.finalState = finalState;
2448         }
2449     }
2450 
2451     private class ResetChangesSignGenerator implements FieldODEEventDetector<Binary64> {
2452 
2453         private final FieldAdaptableInterval<Binary64>             maxCheck;
2454         private final int                                          maxIter;
2455         private final BracketedRealFieldUnivariateSolver<Binary64> solver;
2456         final double                                               y1;
2457         final double                                               y2;
2458         final double                                               change;
2459         int                                                        count;
2460 
2461         public ResetChangesSignGenerator(final double y1, final double y2, final double change,
2462                                          final double maxCheck, final double threshold, final int maxIter) {
2463             this.maxCheck  = (s, isForward) -> maxCheck;
2464             this.maxIter   = maxIter;
2465             this.solver    = new FieldBracketingNthOrderBrentSolver<>(new Binary64(0),
2466                             new Binary64(threshold),
2467                             new Binary64(0),
2468                             5);
2469             this.y1        = y1;
2470             this.y2        = y2;
2471             this.change    = change;
2472             this.count     = 0;
2473         }
2474 
2475         public FieldAdaptableInterval<Binary64> getMaxCheckInterval() {
2476             return maxCheck;
2477         }
2478 
2479         public int getMaxIterationCount() {
2480             return maxIter;
2481         }
2482 
2483         public BracketedRealFieldUnivariateSolver<Binary64> getSolver() {
2484             return solver;
2485         }
2486 
2487         public FieldODEEventHandler<Binary64> getHandler() {
2488             return new FieldODEEventHandler<Binary64>() {
2489                 public Action eventOccurred(FieldODEStateAndDerivative<Binary64> s,
2490                                             FieldODEEventDetector<Binary64> detector,
2491                                             boolean increasing) {
2492                     return ++count < 2 ? Action.RESET_STATE : Action.STOP;
2493                 }
2494 
2495                 public FieldODEState<Binary64> resetState(FieldODEEventDetector<Binary64> detector,
2496                                                            FieldODEStateAndDerivative<Binary64> s) {
2497                     return new FieldODEState<>(s.getTime(), new Binary64[] { s.getCompleteState()[0].add(change) });
2498                 }
2499             };
2500         }
2501 
2502         public Binary64 g(FieldODEStateAndDerivative<Binary64> s) {
2503             return s.getCompleteState()[0].subtract(y1).multiply(s.getCompleteState()[0].subtract(y2));
2504         }
2505 
2506         public int getCount() {
2507             return count;
2508         }
2509 
2510     }
2511 
2512 }