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.analysis.function; 19 20 import org.apache.commons.math3.analysis.DifferentiableUnivariateFunction; 21 import org.apache.commons.math3.analysis.FunctionUtils; 22 import org.apache.commons.math3.analysis.ParametricUnivariateFunction; 23 import org.apache.commons.math3.analysis.UnivariateFunction; 24 import org.apache.commons.math3.analysis.differentiation.DerivativeStructure; 25 import org.apache.commons.math3.analysis.differentiation.UnivariateDifferentiableFunction; 26 import org.apache.commons.math3.exception.DimensionMismatchException; 27 import org.apache.commons.math3.exception.NullArgumentException; 28 import org.apache.commons.math3.util.FastMath; 29 30 /** 31 * <a href="http://en.wikipedia.org/wiki/Harmonic_oscillator"> 32 * simple harmonic oscillator</a> function. 33 * 34 * @since 3.0 35 */ 36 public class HarmonicOscillator implements UnivariateDifferentiableFunction, DifferentiableUnivariateFunction { 37 /** Amplitude. */ 38 private final double amplitude; 39 /** Angular frequency. */ 40 private final double omega; 41 /** Phase. */ 42 private final double phase; 43 44 /** 45 * Harmonic oscillator function. 46 * 47 * @param amplitude Amplitude. 48 * @param omega Angular frequency. 49 * @param phase Phase. 50 */ HarmonicOscillator(double amplitude, double omega, double phase)51 public HarmonicOscillator(double amplitude, 52 double omega, 53 double phase) { 54 this.amplitude = amplitude; 55 this.omega = omega; 56 this.phase = phase; 57 } 58 59 /** {@inheritDoc} */ value(double x)60 public double value(double x) { 61 return value(omega * x + phase, amplitude); 62 } 63 64 /** {@inheritDoc} 65 * @deprecated as of 3.1, replaced by {@link #value(DerivativeStructure)} 66 */ 67 @Deprecated derivative()68 public UnivariateFunction derivative() { 69 return FunctionUtils.toDifferentiableUnivariateFunction(this).derivative(); 70 } 71 72 /** 73 * Parametric function where the input array contains the parameters of 74 * the harmonic oscillator function, ordered as follows: 75 * <ul> 76 * <li>Amplitude</li> 77 * <li>Angular frequency</li> 78 * <li>Phase</li> 79 * </ul> 80 */ 81 public static class Parametric implements ParametricUnivariateFunction { 82 /** 83 * Computes the value of the harmonic oscillator at {@code x}. 84 * 85 * @param x Value for which the function must be computed. 86 * @param param Values of norm, mean and standard deviation. 87 * @return the value of the function. 88 * @throws NullArgumentException if {@code param} is {@code null}. 89 * @throws DimensionMismatchException if the size of {@code param} is 90 * not 3. 91 */ value(double x, double ... param)92 public double value(double x, double ... param) 93 throws NullArgumentException, 94 DimensionMismatchException { 95 validateParameters(param); 96 return HarmonicOscillator.value(x * param[1] + param[2], param[0]); 97 } 98 99 /** 100 * Computes the value of the gradient at {@code x}. 101 * The components of the gradient vector are the partial 102 * derivatives of the function with respect to each of the 103 * <em>parameters</em> (amplitude, angular frequency and phase). 104 * 105 * @param x Value at which the gradient must be computed. 106 * @param param Values of amplitude, angular frequency and phase. 107 * @return the gradient vector at {@code x}. 108 * @throws NullArgumentException if {@code param} is {@code null}. 109 * @throws DimensionMismatchException if the size of {@code param} is 110 * not 3. 111 */ gradient(double x, double ... param)112 public double[] gradient(double x, double ... param) 113 throws NullArgumentException, 114 DimensionMismatchException { 115 validateParameters(param); 116 117 final double amplitude = param[0]; 118 final double omega = param[1]; 119 final double phase = param[2]; 120 121 final double xTimesOmegaPlusPhase = omega * x + phase; 122 final double a = HarmonicOscillator.value(xTimesOmegaPlusPhase, 1); 123 final double p = -amplitude * FastMath.sin(xTimesOmegaPlusPhase); 124 final double w = p * x; 125 126 return new double[] { a, w, p }; 127 } 128 129 /** 130 * Validates parameters to ensure they are appropriate for the evaluation of 131 * the {@link #value(double,double[])} and {@link #gradient(double,double[])} 132 * methods. 133 * 134 * @param param Values of norm, mean and standard deviation. 135 * @throws NullArgumentException if {@code param} is {@code null}. 136 * @throws DimensionMismatchException if the size of {@code param} is 137 * not 3. 138 */ validateParameters(double[] param)139 private void validateParameters(double[] param) 140 throws NullArgumentException, 141 DimensionMismatchException { 142 if (param == null) { 143 throw new NullArgumentException(); 144 } 145 if (param.length != 3) { 146 throw new DimensionMismatchException(param.length, 3); 147 } 148 } 149 } 150 151 /** 152 * @param xTimesOmegaPlusPhase {@code omega * x + phase}. 153 * @param amplitude Amplitude. 154 * @return the value of the harmonic oscillator function at {@code x}. 155 */ value(double xTimesOmegaPlusPhase, double amplitude)156 private static double value(double xTimesOmegaPlusPhase, 157 double amplitude) { 158 return amplitude * FastMath.cos(xTimesOmegaPlusPhase); 159 } 160 161 /** {@inheritDoc} 162 * @since 3.1 163 */ value(final DerivativeStructure t)164 public DerivativeStructure value(final DerivativeStructure t) 165 throws DimensionMismatchException { 166 final double x = t.getValue(); 167 double[] f = new double[t.getOrder() + 1]; 168 169 final double alpha = omega * x + phase; 170 f[0] = amplitude * FastMath.cos(alpha); 171 if (f.length > 1) { 172 f[1] = -amplitude * omega * FastMath.sin(alpha); 173 final double mo2 = - omega * omega; 174 for (int i = 2; i < f.length; ++i) { 175 f[i] = mo2 * f[i - 2]; 176 } 177 } 178 179 return t.compose(f); 180 181 } 182 183 } 184