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