1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23 package org.hipparchus.ode;
24
25 import java.util.ArrayList;
26 import java.util.Collections;
27 import java.util.Comparator;
28 import java.util.List;
29 import java.util.PriorityQueue;
30 import java.util.Queue;
31 import java.util.stream.Collectors;
32 import java.util.stream.Stream;
33
34 import org.hipparchus.exception.MathIllegalArgumentException;
35 import org.hipparchus.exception.MathIllegalStateException;
36 import org.hipparchus.ode.events.Action;
37 import org.hipparchus.ode.events.DetectorBasedEventState;
38 import org.hipparchus.ode.events.EventOccurrence;
39 import org.hipparchus.ode.events.EventState;
40 import org.hipparchus.ode.events.ODEEventDetector;
41 import org.hipparchus.ode.events.ODEStepEndHandler;
42 import org.hipparchus.ode.events.StepEndEventState;
43 import org.hipparchus.ode.sampling.AbstractODEStateInterpolator;
44 import org.hipparchus.ode.sampling.ODEStateInterpolator;
45 import org.hipparchus.ode.sampling.ODEStepHandler;
46 import org.hipparchus.util.FastMath;
47 import org.hipparchus.util.Incrementor;
48
49
50
51
52 public abstract class AbstractIntegrator implements ODEIntegrator {
53
54
55 private List<ODEStepHandler> stepHandlers;
56
57
58 private ODEStateAndDerivative stepStart;
59
60
61 private double stepSize;
62
63
64 private boolean isLastStep;
65
66
67 private boolean resetOccurred;
68
69
70 private List<DetectorBasedEventState> detectorBasedEventsStates;
71
72
73 private List<StepEndEventState> stepEndEventsStates;
74
75
76 private boolean statesInitialized;
77
78
79 private final String name;
80
81
82 private Incrementor evaluations;
83
84
85 private ExpandableODE equations;
86
87
88
89
90 protected AbstractIntegrator(final String name) {
91 this.name = name;
92 stepHandlers = new ArrayList<>();
93 stepStart = null;
94 stepSize = Double.NaN;
95 detectorBasedEventsStates = new ArrayList<>();
96 stepEndEventsStates = new ArrayList<>();
97 statesInitialized = false;
98 evaluations = new Incrementor();
99 }
100
101
102 @Override
103 public String getName() {
104 return name;
105 }
106
107
108 @Override
109 public void addStepHandler(final ODEStepHandler handler) {
110 stepHandlers.add(handler);
111 }
112
113
114 @Override
115 public List<ODEStepHandler> getStepHandlers() {
116 return Collections.unmodifiableList(stepHandlers);
117 }
118
119
120 @Override
121 public void clearStepHandlers() {
122 stepHandlers.clear();
123 }
124
125
126 @Override
127 public void addEventDetector(final ODEEventDetector detector) {
128 detectorBasedEventsStates.add(new DetectorBasedEventState(detector));
129 }
130
131
132 @Override
133 public List<ODEEventDetector> getEventDetectors() {
134 return detectorBasedEventsStates.stream().map(DetectorBasedEventState::getEventDetector).collect(Collectors.toList());
135 }
136
137
138 @Override
139 public void clearEventDetectors() {
140 detectorBasedEventsStates.clear();
141 }
142
143
144 @Override
145 public void addStepEndHandler(ODEStepEndHandler handler) {
146 stepEndEventsStates.add(new StepEndEventState(handler));
147 }
148
149
150 @Override
151 public List<ODEStepEndHandler> getStepEndHandlers() {
152 return stepEndEventsStates.stream().map(StepEndEventState::getHandler).collect(Collectors.toList());
153 }
154
155
156 @Override
157 public void clearStepEndHandlers() {
158 stepEndEventsStates.clear();
159 }
160
161
162 @Override
163 public double getCurrentSignedStepsize() {
164 return stepSize;
165 }
166
167
168 @Override
169 public void setMaxEvaluations(int maxEvaluations) {
170 evaluations = evaluations.withMaximalCount((maxEvaluations < 0) ? Integer.MAX_VALUE : maxEvaluations);
171 }
172
173
174 @Override
175 public int getMaxEvaluations() {
176 return evaluations.getMaximalCount();
177 }
178
179
180 @Override
181 public int getEvaluations() {
182 return evaluations.getCount();
183 }
184
185
186
187
188
189
190
191
192
193 protected ODEStateAndDerivative initIntegration(final ExpandableODE eqn,
194 final ODEState s0, final double t) {
195
196 this.equations = eqn;
197 evaluations = evaluations.withCount(0);
198
199
200 eqn.init(s0, t);
201
202
203 final double t0 = s0.getTime();
204 final double[] y0 = s0.getCompleteState();
205 final double[] y0Dot = computeDerivatives(t0, y0);
206
207
208 final ODEStateAndDerivative s0WithDerivatives =
209 eqn.getMapper().mapStateAndDerivative(t0, y0, y0Dot);
210
211
212 detectorBasedEventsStates.forEach(s -> {
213 s.init(s0WithDerivatives, t);
214 s.getEventDetector().getHandler().init(s0WithDerivatives, t, s.getEventDetector());
215 });
216
217
218 stepEndEventsStates.forEach(s -> {
219 s.init(s0WithDerivatives, t);
220 s.getHandler().init(s0WithDerivatives, t);
221 });
222
223
224 for (ODEStepHandler handler : stepHandlers) {
225 handler.init(s0WithDerivatives, t);
226 }
227
228 setStateInitialized(false);
229
230 return s0WithDerivatives;
231
232 }
233
234
235
236
237 protected ExpandableODE getEquations() {
238 return equations;
239 }
240
241
242
243
244 protected Incrementor getEvaluationsCounter() {
245 return evaluations;
246 }
247
248
249
250
251
252
253
254
255
256
257 public double[] computeDerivatives(final double t, final double[] y)
258 throws MathIllegalArgumentException, MathIllegalStateException, NullPointerException {
259 evaluations.increment();
260 return equations.computeDerivatives(t, y);
261 }
262
263
264
265
266
267 protected void incrementEvaluations(final int nTimes) {
268 evaluations.increment(nTimes);
269 }
270
271
272
273
274
275
276
277 protected void setStateInitialized(final boolean stateInitialized) {
278 this.statesInitialized = stateInitialized;
279 }
280
281
282
283
284
285
286
287
288
289
290 protected ODEStateAndDerivative acceptStep(final AbstractODEStateInterpolator interpolator,
291 final double tEnd)
292 throws MathIllegalArgumentException, MathIllegalStateException {
293
294 ODEStateAndDerivative previousState = interpolator.getGlobalPreviousState();
295 final ODEStateAndDerivative currentState = interpolator.getGlobalCurrentState();
296 ODEStateInterpolator restricted = interpolator;
297
298
299
300 if (!statesInitialized) {
301
302 detectorBasedEventsStates.forEach(s -> s.reinitializeBegin(interpolator));
303 statesInitialized = true;
304 }
305
306
307 stepEndEventsStates.forEach(s -> s.setStepEnd(currentState.getTime()));
308
309
310 final int orderingSign = interpolator.isForward() ? +1 : -1;
311 final Queue<EventState> occurringEvents = new PriorityQueue<>(new Comparator<EventState>() {
312
313 @Override
314 public int compare(final EventState es0, final EventState es1) {
315 return orderingSign * Double.compare(es0.getEventTime(), es1.getEventTime());
316 }
317 });
318
319 resetOccurred = false;
320 boolean doneWithStep = false;
321 resetEvents:
322 do {
323
324
325 occurringEvents.clear();
326 final ODEStateInterpolator finalRestricted = restricted;
327 Stream.concat(detectorBasedEventsStates.stream(), stepEndEventsStates.stream()).
328 forEach(s -> { if (s.evaluateStep(finalRestricted)) {
329
330 occurringEvents.add(s);
331 }
332 });
333
334 do {
335
336 eventLoop:
337 while (!occurringEvents.isEmpty()) {
338
339
340 final EventState currentEvent = occurringEvents.poll();
341
342
343 ODEStateAndDerivative eventState = restricted.getInterpolatedState(currentEvent.getEventTime());
344
345
346 restricted = restricted.restrictStep(previousState, eventState);
347
348
349 for (final DetectorBasedEventState state : detectorBasedEventsStates) {
350 if (state != currentEvent && state.tryAdvance(eventState, interpolator)) {
351
352
353 occurringEvents.remove(state);
354
355 occurringEvents.add(state);
356
357 occurringEvents.add(currentEvent);
358 continue eventLoop;
359 }
360 }
361
362
363
364 for (final ODEStepHandler handler : stepHandlers) {
365 handler.handleStep(restricted);
366 }
367
368
369 final EventOccurrence occurrence = currentEvent.doEvent(eventState);
370 final Action action = occurrence.getAction();
371 isLastStep = action == Action.STOP;
372
373 if (isLastStep) {
374
375
376
377
378 final ODEStateAndDerivative savedState = eventState;
379 eventState = interpolator.getInterpolatedState(occurrence.getStopTime());
380 restricted = interpolator.restrictStep(savedState, eventState);
381
382
383 for (final ODEStepHandler handler : stepHandlers) {
384 handler.handleStep(restricted);
385 handler.finish(restricted.getCurrentState());
386 }
387
388 }
389
390 if (isLastStep) {
391
392 return eventState;
393 }
394
395 if (action == Action.RESET_DERIVATIVES || action == Action.RESET_STATE) {
396
397
398 final ODEState newState = occurrence.getNewState();
399 final double[] y = newState.getCompleteState();
400 final double[] yDot = computeDerivatives(newState.getTime(), y);
401 resetOccurred = true;
402 final ODEStateAndDerivative newStateAndDerivative = equations.getMapper().mapStateAndDerivative(newState.getTime(),
403 y, yDot);
404 detectorBasedEventsStates.forEach(e -> e.getEventDetector().reset(newStateAndDerivative, tEnd));
405 return newStateAndDerivative;
406 }
407
408
409
410 previousState = eventState;
411 restricted = restricted.restrictStep(eventState, currentState);
412
413 if (action == Action.RESET_EVENTS) {
414 continue resetEvents;
415 }
416
417
418
419 if (currentEvent.evaluateStep(restricted)) {
420
421 occurringEvents.add(currentEvent);
422 }
423
424 }
425
426
427
428
429
430
431
432
433 for (final DetectorBasedEventState state : detectorBasedEventsStates) {
434 if (state.tryAdvance(currentState, interpolator)) {
435 occurringEvents.add(state);
436 }
437 }
438
439 } while (!occurringEvents.isEmpty());
440
441 doneWithStep = true;
442 } while (!doneWithStep);
443
444 isLastStep = isLastStep || FastMath.abs(currentState.getTime() - tEnd) < FastMath.ulp(tEnd);
445
446
447 for (ODEStepHandler handler : stepHandlers) {
448 handler.handleStep(restricted);
449 if (isLastStep) {
450 handler.finish(restricted.getCurrentState());
451 }
452 }
453
454 return currentState;
455
456 }
457
458
459
460
461
462
463
464
465 protected void sanityChecks(final ODEState initialState, final double t)
466 throws MathIllegalArgumentException {
467
468 final double threshold = 1000 * FastMath.ulp(FastMath.max(FastMath.abs(initialState.getTime()),
469 FastMath.abs(t)));
470 final double dt = FastMath.abs(initialState.getTime() - t);
471 if (dt < threshold) {
472 throw new MathIllegalArgumentException(LocalizedODEFormats.TOO_SMALL_INTEGRATION_INTERVAL,
473 dt, threshold, false);
474 }
475
476 }
477
478
479
480
481 protected boolean resetOccurred() {
482 return resetOccurred;
483 }
484
485
486
487
488 protected void setStepSize(final double stepSize) {
489 this.stepSize = stepSize;
490 }
491
492
493
494
495 protected double getStepSize() {
496 return stepSize;
497 }
498
499
500
501 protected void setStepStart(final ODEStateAndDerivative stepStart) {
502 this.stepStart = stepStart;
503 }
504
505
506 @Override
507 public ODEStateAndDerivative getStepStart() {
508 return stepStart;
509 }
510
511
512
513
514 protected void setIsLastStep(final boolean isLastStep) {
515 this.isLastStep = isLastStep;
516 }
517
518
519
520
521 protected boolean isLastStep() {
522 return isLastStep;
523 }
524
525 }