1 /*
2  * Licensed to the Apache Software Foundation (ASF) under one or more
3  * contributor license agreements.  See the NOTICE file distributed with
4  * this work for additional information regarding copyright ownership.
5  * The ASF licenses this file to You under the Apache License, Version 2.0
6  * (the "License"); you may not use this file except in compliance with
7  * the License.  You may obtain a copy of the License at
8  *
9  *      http://www.apache.org/licenses/LICENSE-2.0
10  *
11  * Unless required by applicable law or agreed to in writing, software
12  * distributed under the License is distributed on an "AS IS" BASIS,
13  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14  * See the License for the specific language governing permissions and
15  * limitations under the License.
16  */
17 
18 package org.apache.commons.math3.ode.nonstiff;
19 
20 import org.apache.commons.math3.Field;
21 import org.apache.commons.math3.RealFieldElement;
22 import org.apache.commons.math3.ode.FieldEquationsMapper;
23 import org.apache.commons.math3.ode.FieldODEStateAndDerivative;
24 import org.apache.commons.math3.ode.sampling.AbstractFieldStepInterpolator;
25 import org.apache.commons.math3.util.MathArrays;
26 
27 /** This class represents an interpolator over the last step during an
28  * ODE integration for Runge-Kutta and embedded Runge-Kutta integrators.
29  *
30  * @see RungeKuttaFieldIntegrator
31  * @see EmbeddedRungeKuttaFieldIntegrator
32  *
33  * @param <T> the type of the field elements
34  * @since 3.6
35  */
36 
37 abstract class RungeKuttaFieldStepInterpolator<T extends RealFieldElement<T>>
38     extends AbstractFieldStepInterpolator<T> {
39 
40     /** Field to which the time and state vector elements belong. */
41     private final Field<T> field;
42 
43     /** Slopes at the intermediate points. */
44     private final T[][] yDotK;
45 
46     /** Simple constructor.
47      * @param field field to which the time and state vector elements belong
48      * @param forward integration direction indicator
49      * @param yDotK slopes at the intermediate points
50      * @param globalPreviousState start of the global step
51      * @param globalCurrentState end of the global step
52      * @param softPreviousState start of the restricted step
53      * @param softCurrentState end of the restricted step
54      * @param mapper equations mapper for the all equations
55      */
RungeKuttaFieldStepInterpolator(final Field<T> field, final boolean forward, final T[][] yDotK, final FieldODEStateAndDerivative<T> globalPreviousState, final FieldODEStateAndDerivative<T> globalCurrentState, final FieldODEStateAndDerivative<T> softPreviousState, final FieldODEStateAndDerivative<T> softCurrentState, final FieldEquationsMapper<T> mapper)56     protected RungeKuttaFieldStepInterpolator(final Field<T> field, final boolean forward,
57                                               final T[][] yDotK,
58                                               final FieldODEStateAndDerivative<T> globalPreviousState,
59                                               final FieldODEStateAndDerivative<T> globalCurrentState,
60                                               final FieldODEStateAndDerivative<T> softPreviousState,
61                                               final FieldODEStateAndDerivative<T> softCurrentState,
62                                               final FieldEquationsMapper<T> mapper) {
63         super(forward, globalPreviousState, globalCurrentState, softPreviousState, softCurrentState, mapper);
64         this.field = field;
65         this.yDotK = MathArrays.buildArray(field, yDotK.length, -1);
66         for (int i = 0; i < yDotK.length; ++i) {
67             this.yDotK[i] = yDotK[i].clone();
68         }
69     }
70 
71     /** {@inheritDoc} */
72     @Override
create(boolean newForward, FieldODEStateAndDerivative<T> newGlobalPreviousState, FieldODEStateAndDerivative<T> newGlobalCurrentState, FieldODEStateAndDerivative<T> newSoftPreviousState, FieldODEStateAndDerivative<T> newSoftCurrentState, FieldEquationsMapper<T> newMapper)73     protected RungeKuttaFieldStepInterpolator<T> create(boolean newForward,
74                                                         FieldODEStateAndDerivative<T> newGlobalPreviousState,
75                                                         FieldODEStateAndDerivative<T> newGlobalCurrentState,
76                                                         FieldODEStateAndDerivative<T> newSoftPreviousState,
77                                                         FieldODEStateAndDerivative<T> newSoftCurrentState,
78                                                         FieldEquationsMapper<T> newMapper) {
79         return create(field, newForward, yDotK,
80                       newGlobalPreviousState, newGlobalCurrentState,
81                       newSoftPreviousState, newSoftCurrentState,
82                       newMapper);
83     }
84 
85     /** Create a new instance.
86      * @param newField field to which the time and state vector elements belong
87      * @param newForward integration direction indicator
88      * @param newYDotK slopes at the intermediate points
89      * @param newGlobalPreviousState start of the global step
90      * @param newGlobalCurrentState end of the global step
91      * @param newSoftPreviousState start of the restricted step
92      * @param newSoftCurrentState end of the restricted step
93      * @param newMapper equations mapper for the all equations
94      * @return a new instance
95      */
create(Field<T> newField, boolean newForward, T[][] newYDotK, FieldODEStateAndDerivative<T> newGlobalPreviousState, FieldODEStateAndDerivative<T> newGlobalCurrentState, FieldODEStateAndDerivative<T> newSoftPreviousState, FieldODEStateAndDerivative<T> newSoftCurrentState, FieldEquationsMapper<T> newMapper)96     protected abstract RungeKuttaFieldStepInterpolator<T> create(Field<T> newField, boolean newForward, T[][] newYDotK,
97                                                                  FieldODEStateAndDerivative<T> newGlobalPreviousState,
98                                                                  FieldODEStateAndDerivative<T> newGlobalCurrentState,
99                                                                  FieldODEStateAndDerivative<T> newSoftPreviousState,
100                                                                  FieldODEStateAndDerivative<T> newSoftCurrentState,
101                                                                  FieldEquationsMapper<T> newMapper);
102 
103     /** Compute a state by linear combination added to previous state.
104      * @param coefficients coefficients to apply to the method staged derivatives
105      * @return combined state
106      */
previousStateLinearCombination(final T ... coefficients)107     protected final T[] previousStateLinearCombination(final T ... coefficients) {
108         return combine(getPreviousState().getState(),
109                        coefficients);
110     }
111 
112     /** Compute a state by linear combination added to current state.
113      * @param coefficients coefficients to apply to the method staged derivatives
114      * @return combined state
115      */
currentStateLinearCombination(final T ... coefficients)116     protected T[] currentStateLinearCombination(final T ... coefficients) {
117         return combine(getCurrentState().getState(),
118                        coefficients);
119     }
120 
121     /** Compute a state derivative by linear combination.
122      * @param coefficients coefficients to apply to the method staged derivatives
123      * @return combined state
124      */
derivativeLinearCombination(final T ... coefficients)125     protected T[] derivativeLinearCombination(final T ... coefficients) {
126         return combine(MathArrays.buildArray(field, yDotK[0].length), coefficients);
127     }
128 
129     /** Linearly combine arrays.
130      * @param a array to add to
131      * @param coefficients coefficients to apply to the method staged derivatives
132      * @return a itself, as a convenience for fluent API
133      */
combine(final T[] a, final T ... coefficients)134     private T[] combine(final T[] a, final T ... coefficients) {
135         for (int i = 0; i < a.length; ++i) {
136             for (int k = 0; k < coefficients.length; ++k) {
137                 a[i] = a[i].add(coefficients[k].multiply(yDotK[k][i]));
138             }
139         }
140         return a;
141     }
142 
143 }
144