001/* 002 * Licensed to the Apache Software Foundation (ASF) under one or more 003 * contributor license agreements. See the NOTICE file distributed with 004 * this work for additional information regarding copyright ownership. 005 * The ASF licenses this file to You under the Apache License, Version 2.0 006 * (the "License"); you may not use this file except in compliance with 007 * the License. You may obtain a copy of the License at 008 * 009 * http://www.apache.org/licenses/LICENSE-2.0 010 * 011 * Unless required by applicable law or agreed to in writing, software 012 * distributed under the License is distributed on an "AS IS" BASIS, 013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 014 * See the License for the specific language governing permissions and 015 * limitations under the License. 016 */ 017 018package org.apache.commons.math4.legacy.ode.events; 019 020import java.util.Arrays; 021 022/** Wrapper used to detect only increasing or decreasing events. 023 * 024 * <p>General {@link EventHandler events} are defined implicitly 025 * by a {@link EventHandler#g(double, double[]) g function} crossing 026 * zero. This function needs to be continuous in the event neighborhood, 027 * and its sign must remain consistent between events. This implies that 028 * during an ODE integration, events triggered are alternately events 029 * for which the function increases from negative to positive values, 030 * and events for which the function decreases from positive to 031 * negative values. 032 * </p> 033 * 034 * <p>Sometimes, users are only interested in one type of event (say 035 * increasing events for example) and not in the other type. In these 036 * cases, looking precisely for all events location and triggering 037 * events that will later be ignored is a waste of computing time.</p> 038 * 039 * <p>Users can wrap a regular {@link EventHandler event handler} in 040 * an instance of this class and provide this wrapping instance to 041 * the {@link org.apache.commons.math4.legacy.ode.FirstOrderIntegrator ODE solver} 042 * in order to avoid wasting time looking for uninteresting events. 043 * The wrapper will intercept the calls to the {@link 044 * EventHandler#g(double, double[]) g function} and to the {@link 045 * EventHandler#eventOccurred(double, double[], boolean) 046 * eventOccurred} method in order to ignore uninteresting events. The 047 * wrapped regular {@link EventHandler event handler} will the see only 048 * the interesting events, i.e. either only {@code increasing} events or 049 * {@code decreasing} events. the number of calls to the {@link 050 * EventHandler#g(double, double[]) g function} will also be reduced.</p> 051 * 052 * @since 3.2 053 */ 054 055public class EventFilter implements EventHandler { 056 057 /** Number of past transformers updates stored. */ 058 private static final int HISTORY_SIZE = 100; 059 060 /** Wrapped event handler. */ 061 private final EventHandler rawHandler; 062 063 /** Filter to use. */ 064 private final FilterType filter; 065 066 /** Transformers of the g function. */ 067 private final Transformer[] transformers; 068 069 /** Update time of the transformers. */ 070 private final double[] updates; 071 072 /** Indicator for forward integration. */ 073 private boolean forward; 074 075 /** Extreme time encountered so far. */ 076 private double extremeT; 077 078 /** Wrap an {@link EventHandler event handler}. 079 * @param rawHandler event handler to wrap 080 * @param filter filter to use 081 */ 082 public EventFilter(final EventHandler rawHandler, final FilterType filter) { 083 this.rawHandler = rawHandler; 084 this.filter = filter; 085 this.transformers = new Transformer[HISTORY_SIZE]; 086 this.updates = new double[HISTORY_SIZE]; 087 } 088 089 /** {@inheritDoc} */ 090 @Override 091 public void init(double t0, double[] y0, double t) { 092 093 // delegate to raw handler 094 rawHandler.init(t0, y0, t); 095 096 // initialize events triggering logic 097 forward = t >= t0; 098 extremeT = forward ? Double.NEGATIVE_INFINITY : Double.POSITIVE_INFINITY; 099 Arrays.fill(transformers, Transformer.UNINITIALIZED); 100 Arrays.fill(updates, extremeT); 101 } 102 103 /** {@inheritDoc} */ 104 @Override 105 public double g(double t, double[] y) { 106 107 final double rawG = rawHandler.g(t, y); 108 109 // search which transformer should be applied to g 110 if (forward) { 111 final int last = transformers.length - 1; 112 if (extremeT < t) { 113 // we are at the forward end of the history 114 115 // check if a new rough root has been crossed 116 final Transformer previous = transformers[last]; 117 final Transformer next = filter.selectTransformer(previous, rawG, forward); 118 if (next != previous) { 119 // there is a root somewhere between extremeT and t. 120 // the new transformer is valid for t (this is how we have just computed 121 // it above), but it is in fact valid on both sides of the root, so 122 // it was already valid before t and even up to previous time. We store 123 // the switch at extremeT for safety, to ensure the previous transformer 124 // is not applied too close of the root 125 System.arraycopy(updates, 1, updates, 0, last); 126 System.arraycopy(transformers, 1, transformers, 0, last); 127 updates[last] = extremeT; 128 transformers[last] = next; 129 } 130 131 extremeT = t; 132 133 // apply the transform 134 return next.transformed(rawG); 135 } else { 136 // we are in the middle of the history 137 138 // select the transformer 139 for (int i = last; i > 0; --i) { 140 if (updates[i] <= t) { 141 // apply the transform 142 return transformers[i].transformed(rawG); 143 } 144 } 145 146 return transformers[0].transformed(rawG); 147 } 148 } else { 149 if (t < extremeT) { 150 // we are at the backward end of the history 151 152 // check if a new rough root has been crossed 153 final Transformer previous = transformers[0]; 154 final Transformer next = filter.selectTransformer(previous, rawG, forward); 155 if (next != previous) { 156 // there is a root somewhere between extremeT and t. 157 // the new transformer is valid for t (this is how we have just computed 158 // it above), but it is in fact valid on both sides of the root, so 159 // it was already valid before t and even up to previous time. We store 160 // the switch at extremeT for safety, to ensure the previous transformer 161 // is not applied too close of the root 162 System.arraycopy(updates, 0, updates, 1, updates.length - 1); 163 System.arraycopy(transformers, 0, transformers, 1, transformers.length - 1); 164 updates[0] = extremeT; 165 transformers[0] = next; 166 } 167 168 extremeT = t; 169 170 // apply the transform 171 return next.transformed(rawG); 172 } else { 173 // we are in the middle of the history 174 175 // select the transformer 176 for (int i = 0; i < updates.length - 1; ++i) { 177 if (t <= updates[i]) { 178 // apply the transform 179 return transformers[i].transformed(rawG); 180 } 181 } 182 183 return transformers[updates.length - 1].transformed(rawG); 184 } 185 } 186 } 187 188 /** {@inheritDoc} */ 189 @Override 190 public Action eventOccurred(double t, double[] y, boolean increasing) { 191 // delegate to raw handler, fixing increasing status on the fly 192 return rawHandler.eventOccurred(t, y, filter.getTriggeredIncreasing()); 193 } 194 195 /** {@inheritDoc} */ 196 @Override 197 public void resetState(double t, double[] y) { 198 // delegate to raw handler 199 rawHandler.resetState(t, y); 200 } 201}