1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
|
/*
* Licensed to the Apache Software Foundation (ASF) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* The ASF licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.apache.commons.math3.ode.events;
import java.util.Arrays;
/** Wrapper used to detect only increasing or decreasing events.
*
* <p>General {@link EventHandler events} are defined implicitly
* by a {@link EventHandler#g(double, double[]) g function} crossing
* zero. This function needs to be continuous in the event neighborhood,
* and its sign must remain consistent between events. This implies that
* during an ODE integration, events triggered are alternately events
* for which the function increases from negative to positive values,
* and events for which the function decreases from positive to
* negative values.
* </p>
*
* <p>Sometimes, users are only interested in one type of event (say
* increasing events for example) and not in the other type. In these
* cases, looking precisely for all events location and triggering
* events that will later be ignored is a waste of computing time.</p>
*
* <p>Users can wrap a regular {@link EventHandler event handler} in
* an instance of this class and provide this wrapping instance to
* the {@link org.apache.commons.math3.ode.FirstOrderIntegrator ODE solver}
* in order to avoid wasting time looking for uninteresting events.
* The wrapper will intercept the calls to the {@link
* EventHandler#g(double, double[]) g function} and to the {@link
* EventHandler#eventOccurred(double, double[], boolean)
* eventOccurred} method in order to ignore uninteresting events. The
* wrapped regular {@link EventHandler event handler} will the see only
* the interesting events, i.e. either only {@code increasing} events or
* {@code decreasing} events. the number of calls to the {@link
* EventHandler#g(double, double[]) g function} will also be reduced.</p>
*
* @since 3.2
*/
public class EventFilter implements EventHandler {
/** Number of past transformers updates stored. */
private static final int HISTORY_SIZE = 100;
/** Wrapped event handler. */
private final EventHandler rawHandler;
/** Filter to use. */
private final FilterType filter;
/** Transformers of the g function. */
private final Transformer[] transformers;
/** Update time of the transformers. */
private final double[] updates;
/** Indicator for forward integration. */
private boolean forward;
/** Extreme time encountered so far. */
private double extremeT;
/** Wrap an {@link EventHandler event handler}.
* @param rawHandler event handler to wrap
* @param filter filter to use
*/
public EventFilter(final EventHandler rawHandler, final FilterType filter) {
this.rawHandler = rawHandler;
this.filter = filter;
this.transformers = new Transformer[HISTORY_SIZE];
this.updates = new double[HISTORY_SIZE];
}
/** {@inheritDoc} */
public void init(double t0, double[] y0, double t) {
// delegate to raw handler
rawHandler.init(t0, y0, t);
// initialize events triggering logic
forward = t >= t0;
extremeT = forward ? Double.NEGATIVE_INFINITY : Double.POSITIVE_INFINITY;
Arrays.fill(transformers, Transformer.UNINITIALIZED);
Arrays.fill(updates, extremeT);
}
/** {@inheritDoc} */
public double g(double t, double[] y) {
final double rawG = rawHandler.g(t, y);
// search which transformer should be applied to g
if (forward) {
final int last = transformers.length - 1;
if (extremeT < t) {
// we are at the forward end of the history
// check if a new rough root has been crossed
final Transformer previous = transformers[last];
final Transformer next = filter.selectTransformer(previous, rawG, forward);
if (next != previous) {
// there is a root somewhere between extremeT and t.
// the new transformer is valid for t (this is how we have just computed
// it above), but it is in fact valid on both sides of the root, so
// it was already valid before t and even up to previous time. We store
// the switch at extremeT for safety, to ensure the previous transformer
// is not applied too close of the root
System.arraycopy(updates, 1, updates, 0, last);
System.arraycopy(transformers, 1, transformers, 0, last);
updates[last] = extremeT;
transformers[last] = next;
}
extremeT = t;
// apply the transform
return next.transformed(rawG);
} else {
// we are in the middle of the history
// select the transformer
for (int i = last; i > 0; --i) {
if (updates[i] <= t) {
// apply the transform
return transformers[i].transformed(rawG);
}
}
return transformers[0].transformed(rawG);
}
} else {
if (t < extremeT) {
// we are at the backward end of the history
// check if a new rough root has been crossed
final Transformer previous = transformers[0];
final Transformer next = filter.selectTransformer(previous, rawG, forward);
if (next != previous) {
// there is a root somewhere between extremeT and t.
// the new transformer is valid for t (this is how we have just computed
// it above), but it is in fact valid on both sides of the root, so
// it was already valid before t and even up to previous time. We store
// the switch at extremeT for safety, to ensure the previous transformer
// is not applied too close of the root
System.arraycopy(updates, 0, updates, 1, updates.length - 1);
System.arraycopy(transformers, 0, transformers, 1, transformers.length - 1);
updates[0] = extremeT;
transformers[0] = next;
}
extremeT = t;
// apply the transform
return next.transformed(rawG);
} else {
// we are in the middle of the history
// select the transformer
for (int i = 0; i < updates.length - 1; ++i) {
if (t <= updates[i]) {
// apply the transform
return transformers[i].transformed(rawG);
}
}
return transformers[updates.length - 1].transformed(rawG);
}
}
}
/** {@inheritDoc} */
public Action eventOccurred(double t, double[] y, boolean increasing) {
// delegate to raw handler, fixing increasing status on the fly
return rawHandler.eventOccurred(t, y, filter.getTriggeredIncreasing());
}
/** {@inheritDoc} */
public void resetState(double t, double[] y) {
// delegate to raw handler
rawHandler.resetState(t, y);
}
}
|
