You are viewing a plain text version of this content. The canonical link for it is here.
Posted to dev@commons.apache.org by ps...@apache.org on 2006/11/19 22:36:23 UTC
svn commit: r476930 [13/22] - in
/jakarta/commons/proper/math/trunk/src/mantissa: ./ src/ src/org/
src/org/spaceroots/ src/org/spaceroots/mantissa/
src/org/spaceroots/mantissa/algebra/
src/org/spaceroots/mantissa/estimation/ src/org/spaceroots/mantissa...
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/CostFunction.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/CostFunction.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/CostFunction.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/CostFunction.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,35 @@
+// 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.spaceroots.mantissa.optimization;
+
+/** This interface represents a cost function to be minimized.
+ * @author Luc Maisonobe
+ * @version $Id: CostFunction.java 1705 2006-09-17 19:57:39Z luc $
+ */
+public interface CostFunction {
+
+
+ /** Compute the cost associated to the given parameters array.
+ * @param x parameters array
+ * @return cost associated to the parameters array
+ * @exception CostException if no cost can be computed for the parameters
+ * @see PointCostPair
+ */
+ public double cost(double[] x) throws CostException;
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/CostFunction.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/DirectSearchOptimizer.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/DirectSearchOptimizer.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/DirectSearchOptimizer.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/DirectSearchOptimizer.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,583 @@
+// 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.spaceroots.mantissa.optimization;
+
+import org.spaceroots.mantissa.random.RandomVectorGenerator;
+import org.spaceroots.mantissa.random.UncorrelatedRandomVectorGenerator;
+import org.spaceroots.mantissa.random.CorrelatedRandomVectorGenerator;
+import org.spaceroots.mantissa.random.UniformRandomGenerator;
+import org.spaceroots.mantissa.random.VectorialSampleStatistics;
+import org.spaceroots.mantissa.random.NotPositiveDefiniteMatrixException;
+
+import java.util.Arrays;
+import java.util.Comparator;
+
+/** This class implements simplex-based direct search optimization
+ * algorithms.
+
+ * <p>Direct search method only use cost function values, they don't
+ * need derivatives and don't either try to compute approximation of
+ * the derivatives. According to a 1996 paper by Margaret H. Wright
+ * (<a href="http://cm.bell-labs.com/cm/cs/doc/96/4-02.ps.gz">Direct
+ * Search Methods: Once Scorned, Now Respectable</a>), they are used
+ * when either the computation of the derivative is impossible (noisy
+ * functions, unpredictable dicontinuities) or difficult (complexity,
+ * computation cost). In the first cases, rather than an optimum, a
+ * <em>not too bad</em> point is desired. In the latter cases, an
+ * optimum is desired but cannot be reasonably found. In all cases
+ * direct search methods can be useful.</p>
+
+ * <p>Simplex-based direct search methods are based on comparison of
+ * the cost function values at the vertices of a simplex (which is a
+ * set of n+1 points in dimension n) that is updated by the algorithms
+ * steps.</p>
+
+ * <p>The instances can be built either in single-start or in
+ * multi-start mode. Multi-start is a traditional way to try to avoid
+ * beeing trapped in a local minimum and miss the global minimum of a
+ * function. It can also be used to verify the convergence of an
+ * algorithm. In multi-start mode, the {@link #minimizes(CostFunction,
+ * int, ConvergenceChecker, double[], double[]) minimizes}
+ * method returns the best minimum found after all starts, and the
+ * {@link #getMinima getMinima} method can be used to retrieve all
+ * minima from all starts (including the one already provided by the
+ * {@link #minimizes(CostFunction, int, ConvergenceChecker, double[],
+ * double[]) minimizes} method).</p>
+
+ * <p>This class is the base class performing the boilerplate simplex
+ * initialization and handling. The simplex update by itself is
+ * performed by the derived classes according to the implemented
+ * algorithms.</p>
+
+ * @author Luc Maisonobe
+ * @version $Id: DirectSearchOptimizer.java 1705 2006-09-17 19:57:39Z luc $
+ * @see CostFunction
+ * @see NelderMead
+ * @see MultiDirectional
+ */
+public abstract class DirectSearchOptimizer {
+
+ /** Simple constructor.
+ */
+ protected DirectSearchOptimizer() {
+ }
+
+ /** Minimizes a cost function.
+ * <p>The initial simplex is built from two vertices that are
+ * considered to represent two opposite vertices of a box parallel
+ * to the canonical axes of the space. The simplex is the subset of
+ * vertices encountered while going from vertexA to vertexB
+ * travelling along the box edges only. This can be seen as a scaled
+ * regular simplex using the projected separation between the given
+ * points as the scaling factor along each coordinate axis.</p>
+ * <p>The optimization is performed in single-start mode.</p>
+ * @param f cost function
+ * @param maxEvaluations maximal number of function calls for each
+ * start (note that the number will be checked <em>after</em>
+ * complete simplices have been evaluated, this means that in some
+ * cases this number will be exceeded by a few units, depending on
+ * the dimension of the problem)
+ * @param checker object to use to check for convergence
+ * @param vertexA first vertex
+ * @param vertexB last vertex
+ * @return the point/cost pairs giving the minimal cost
+ * @exception CostException if the cost function throws one during
+ * the search
+ * @exception NoConvergenceException if none of the starts did
+ * converge (it is not thrown if at least one start did converge)
+ */
+ public PointCostPair minimizes(CostFunction f, int maxEvaluations,
+ ConvergenceChecker checker,
+ double[] vertexA, double[] vertexB)
+ throws CostException, NoConvergenceException {
+
+ // set up optimizer
+ buildSimplex(vertexA, vertexB);
+ setSingleStart();
+
+ // compute minimum
+ return minimizes(f, maxEvaluations, checker);
+
+ }
+
+ /** Minimizes a cost function.
+ * <p>The initial simplex is built from two vertices that are
+ * considered to represent two opposite vertices of a box parallel
+ * to the canonical axes of the space. The simplex is the subset of
+ * vertices encountered while going from vertexA to vertexB
+ * travelling along the box edges only. This can be seen as a scaled
+ * regular simplex using the projected separation between the given
+ * points as the scaling factor along each coordinate axis.</p>
+ * <p>The optimization is performed in multi-start mode.</p>
+ * @param f cost function
+ * @param maxEvaluations maximal number of function calls for each
+ * start (note that the number will be checked <em>after</em>
+ * complete simplices have been evaluated, this means that in some
+ * cases this number will be exceeded by a few units, depending on
+ * the dimension of the problem)
+ * @param checker object to use to check for convergence
+ * @param vertexA first vertex
+ * @param vertexB last vertex
+ * @param starts number of starts to perform (including the
+ * first one), multi-start is disabled if value is less than or
+ * equal to 1
+ * @param seed seed for the random vector generator
+ * @return the point/cost pairs giving the minimal cost
+ * @exception CostException if the cost function throws one during
+ * the search
+ * @exception NoConvergenceException if none of the starts did
+ * converge (it is not thrown if at least one start did converge)
+ */
+ public PointCostPair minimizes(CostFunction f, int maxEvaluations,
+ ConvergenceChecker checker,
+ double[] vertexA, double[] vertexB,
+ int starts, long seed)
+ throws CostException, NoConvergenceException {
+
+ // set up the simplex travelling around the box
+ buildSimplex(vertexA, vertexB);
+
+ // we consider the simplex could have been produced by a generator
+ // having its mean value at the center of the box, the standard
+ // deviation along each axe beeing the corresponding half size
+ double[] mean = new double[vertexA.length];
+ double[] standardDeviation = new double[vertexA.length];
+ for (int i = 0; i < vertexA.length; ++i) {
+ mean[i] = 0.5 * (vertexA[i] + vertexB[i]);
+ standardDeviation[i] = 0.5 * Math.abs(vertexA[i] - vertexB[i]);
+ }
+
+ RandomVectorGenerator rvg =
+ new UncorrelatedRandomVectorGenerator(mean, standardDeviation,
+ new UniformRandomGenerator(seed));
+ setMultiStart(starts, rvg);
+
+ // compute minimum
+ return minimizes(f, maxEvaluations, checker);
+
+ }
+
+ /** Minimizes a cost function.
+ * <p>The simplex is built from all its vertices.</p>
+ * <p>The optimization is performed in single-start mode.</p>
+ * @param f cost function
+ * @param maxEvaluations maximal number of function calls for each
+ * start (note that the number will be checked <em>after</em>
+ * complete simplices have been evaluated, this means that in some
+ * cases this number will be exceeded by a few units, depending on
+ * the dimension of the problem)
+ * @param checker object to use to check for convergence
+ * @param vertices array containing all vertices of the simplex
+ * @return the point/cost pairs giving the minimal cost
+ * @exception CostException if the cost function throws one during
+ * the search
+ * @exception NoConvergenceException if none of the starts did
+ * converge (it is not thrown if at least one start did converge)
+ */
+ public PointCostPair minimizes(CostFunction f, int maxEvaluations,
+ ConvergenceChecker checker,
+ double[][] vertices)
+ throws CostException, NoConvergenceException {
+
+ // set up optimizer
+ buildSimplex(vertices);
+ setSingleStart();
+
+ // compute minimum
+ return minimizes(f, maxEvaluations, checker);
+
+ }
+
+ /** Minimizes a cost function.
+ * <p>The simplex is built from all its vertices.</p>
+ * <p>The optimization is performed in multi-start mode.</p>
+ * @param f cost function
+ * @param maxEvaluations maximal number of function calls for each
+ * start (note that the number will be checked <em>after</em>
+ * complete simplices have been evaluated, this means that in some
+ * cases this number will be exceeded by a few units, depending on
+ * the dimension of the problem)
+ * @param checker object to use to check for convergence
+ * @param vertices array containing all vertices of the simplex
+ * @param starts number of starts to perform (including the
+ * first one), multi-start is disabled if value is less than or
+ * equal to 1
+ * @param seed seed for the random vector generator
+ * @return the point/cost pairs giving the minimal cost
+ * @exception NotPositiveDefiniteMatrixException if the vertices
+ * array is degenerated
+ * @exception CostException if the cost function throws one during
+ * the search
+ * @exception NoConvergenceException if none of the starts did
+ * converge (it is not thrown if at least one start did converge)
+ */
+ public PointCostPair minimizes(CostFunction f, int maxEvaluations,
+ ConvergenceChecker checker,
+ double[][] vertices,
+ int starts, long seed)
+ throws NotPositiveDefiniteMatrixException,
+ CostException, NoConvergenceException {
+
+ // store the points into the simplex
+ buildSimplex(vertices);
+
+ // compute the statistical properties of the simplex points
+ VectorialSampleStatistics statistics = new VectorialSampleStatistics();
+ for (int i = 0; i < vertices.length; ++i) {
+ statistics.add(vertices[i]);
+ }
+
+ RandomVectorGenerator rvg =
+ new CorrelatedRandomVectorGenerator(statistics.getMean(null),
+ statistics.getCovarianceMatrix(null),
+ new UniformRandomGenerator(seed));
+ setMultiStart(starts, rvg);
+
+ // compute minimum
+ return minimizes(f, maxEvaluations, checker);
+
+ }
+
+ /** Minimizes a cost function.
+ * <p>The simplex is built randomly.</p>
+ * <p>The optimization is performed in single-start mode.</p>
+ * @param f cost function
+ * @param maxEvaluations maximal number of function calls for each
+ * start (note that the number will be checked <em>after</em>
+ * complete simplices have been evaluated, this means that in some
+ * cases this number will be exceeded by a few units, depending on
+ * the dimension of the problem)
+ * @param checker object to use to check for convergence
+ * @param generator random vector generator
+ * @return the point/cost pairs giving the minimal cost
+ * @exception CostException if the cost function throws one during
+ * the search
+ * @exception NoConvergenceException if none of the starts did
+ * converge (it is not thrown if at least one start did converge)
+ */
+ public PointCostPair minimizes(CostFunction f, int maxEvaluations,
+ ConvergenceChecker checker,
+ RandomVectorGenerator generator)
+ throws CostException, NoConvergenceException {
+
+ // set up optimizer
+ buildSimplex(generator);
+ setSingleStart();
+
+ // compute minimum
+ return minimizes(f, maxEvaluations, checker);
+
+ }
+
+ /** Minimizes a cost function.
+ * <p>The simplex is built randomly.</p>
+ * <p>The optimization is performed in multi-start mode.</p>
+ * @param f cost function
+ * @param maxEvaluations maximal number of function calls for each
+ * start (note that the number will be checked <em>after</em>
+ * complete simplices have been evaluated, this means that in some
+ * cases this number will be exceeded by a few units, depending on
+ * the dimension of the problem)
+ * @param checker object to use to check for convergence
+ * @param generator random vector generator
+ * @param starts number of starts to perform (including the
+ * first one), multi-start is disabled if value is less than or
+ * equal to 1
+ * @return the point/cost pairs giving the minimal cost
+ * @exception CostException if the cost function throws one during
+ * the search
+ * @exception NoConvergenceException if none of the starts did
+ * converge (it is not thrown if at least one start did converge)
+ */
+ public PointCostPair minimizes(CostFunction f, int maxEvaluations,
+ ConvergenceChecker checker,
+ RandomVectorGenerator generator,
+ int starts)
+ throws CostException, NoConvergenceException {
+
+ // set up optimizer
+ buildSimplex(generator);
+ setMultiStart(starts, generator);
+
+ // compute minimum
+ return minimizes(f, maxEvaluations, checker);
+
+ }
+
+ /** Build a simplex from two extreme vertices.
+ * <p>The two vertices are considered to represent two opposite
+ * vertices of a box parallel to the canonical axes of the
+ * space. The simplex is the subset of vertices encountered while
+ * going from vertexA to vertexB travelling along the box edges
+ * only. This can be seen as a scaled regular simplex using the
+ * projected separation between the given points as the scaling
+ * factor along each coordinate axis.</p>
+ * @param vertexA first vertex
+ * @param vertexB last vertex
+ */
+ private void buildSimplex(double[] vertexA, double[] vertexB) {
+
+ int n = vertexA.length;
+ simplex = new PointCostPair[n + 1];
+
+ // set up the simplex travelling around the box
+ for (int i = 0; i <= n; ++i) {
+ double[] vertex = new double[n];
+ if (i > 0) {
+ System.arraycopy(vertexB, 0, vertex, 0, i);
+ }
+ if (i < n) {
+ System.arraycopy(vertexA, i, vertex, i, n - i);
+ }
+ simplex[i] = new PointCostPair(vertex);
+ }
+
+ }
+
+ /** Build a simplex from all its points.
+ * @param vertices array containing all vertices of the simplex
+ */
+ private void buildSimplex(double[][] vertices) {
+ int n = vertices.length - 1;
+ simplex = new PointCostPair[n + 1];
+ for (int i = 0; i <= n; ++i) {
+ simplex[i] = new PointCostPair(vertices[i]);
+ }
+ }
+
+ /** Build a simplex randomly.
+ * @param generator random vector generator
+ */
+ private void buildSimplex(RandomVectorGenerator generator) {
+
+ // use first vector size to compute the number of points
+ double[] vertex = generator.nextVector();
+ int n = vertex.length;
+ simplex = new PointCostPair[n + 1];
+ simplex[0] = new PointCostPair(vertex);
+
+ // fill up the vertex
+ for (int i = 1; i <= n; ++i) {
+ simplex[i] = new PointCostPair(generator.nextVector());
+ }
+
+ }
+
+ /** Set up single-start mode.
+ */
+ private void setSingleStart() {
+ starts = 1;
+ generator = null;
+ minima = null;
+ }
+
+ /** Set up multi-start mode.
+ * @param starts number of starts to perform (including the
+ * first one), multi-start is disabled if value is less than or
+ * equal to 1
+ * @param generator random vector generator to use for restarts
+ */
+ public void setMultiStart(int starts, RandomVectorGenerator generator) {
+ if (starts < 2) {
+ this.starts = 1;
+ this.generator = null;
+ minima = null;
+ } else {
+ this.starts = starts;
+ this.generator = generator;
+ minima = null;
+ }
+ }
+
+ /** Get all the minima found during the last call to {@link
+ * #minimizes(CostFunction, int, ConvergenceChecker, double[], double[])
+ * minimizes}.
+ * <p>The optimizer stores all the minima found during a set of
+ * restarts when multi-start mode is enabled. The {@link
+ * #minimizes(CostFunction, int, ConvergenceChecker, double[], double[])
+ * minimizes} method returns the best point only. This method
+ * returns all the points found at the end of each starts, including
+ * the best one already returned by the {@link #minimizes(CostFunction,
+ * int, ConvergenceChecker, double[], double[]) minimizes} method.
+ * The array as one element for each start as specified in the constructor
+ * (it has one element only if optimizer has been set up for single-start).</p>
+ * <p>The array containing the minima is ordered with the results
+ * from the runs that did converge first, sorted from lowest to
+ * highest minimum cost, and null elements corresponding to the runs
+ * that did not converge (all elements will be null if the {@link
+ * #minimizes(CostFunction, int, ConvergenceChecker, double[], double[])
+ * minimizes} method throwed a {@link NoConvergenceException
+ * NoConvergenceException}).</p>
+ * @return array containing the minima, or null if {@link
+ * #minimizes(CostFunction, int, ConvergenceChecker, double[], double[])
+ * minimizes} has not been called
+ */
+ public PointCostPair[] getMinima() {
+ return minima;
+ }
+
+ /** Minimizes a cost function.
+ * @param f cost function
+ * @param maxEvaluations maximal number of function calls for each
+ * start (note that the number will be checked <em>after</em>
+ * complete simplices have been evaluated, this means that in some
+ * cases this number will be exceeded by a few units, depending on
+ * the dimension of the problem)
+ * @param checker object to use to check for convergence
+ * @return the point/cost pairs giving the minimal cost
+ * @exception CostException if the cost function throws one during
+ * the search
+ * @exception NoConvergenceException if none of the starts did
+ * converge (it is not thrown if at least one start did converge)
+ */
+ private PointCostPair minimizes(CostFunction f, int maxEvaluations,
+ ConvergenceChecker checker)
+ throws CostException, NoConvergenceException {
+
+ this.f = f;
+ minima = new PointCostPair[starts];
+
+ // multi-start loop
+ for (int i = 0; i < starts; ++i) {
+
+ evaluations = 0;
+ evaluateSimplex();
+
+ for (boolean loop = true; loop;) {
+ if (checker.converged(simplex)) {
+ // we have found a minimum
+ minima[i] = simplex[0];
+ loop = false;
+ } else if (evaluations >= maxEvaluations) {
+ // this start did not converge, try a new one
+ minima[i] = null;
+ loop = false;
+ } else {
+ iterateSimplex();
+ }
+ }
+
+ if (i < (starts - 1)) {
+ // restart
+ buildSimplex(generator);
+ }
+
+ }
+
+ // sort the minima from lowest cost to highest cost, followed by
+ // null elements
+ Arrays.sort(minima, pointCostPairComparator);
+
+ // return the found point given the lowest cost
+ if (minima[0] == null) {
+ throw new NoConvergenceException("none of the {0} start points"
+ + " lead to convergence",
+ new String[] {
+ Integer.toString(starts)
+ });
+ }
+ return minima[0];
+
+ }
+
+ /** Compute the next simplex of the algorithm.
+ */
+ protected abstract void iterateSimplex()
+ throws CostException;
+
+ /** Evaluate the cost on one point.
+ * <p>A side effect of this method is to count the number of
+ * function evaluations</p>
+ * @param x point on which the cost function should be evaluated
+ * @return cost at the given point
+ * @exception CostException if no cost can be computed for the parameters
+ */
+ protected double evaluateCost(double[] x)
+ throws CostException {
+ evaluations++;
+ return f.cost(x);
+ }
+
+ /** Evaluate all the non-evaluated points of the simplex.
+ * @exception CostException if no cost can be computed for the parameters
+ */
+ protected void evaluateSimplex()
+ throws CostException {
+
+ // evaluate the cost at all non-evaluated simplex points
+ for (int i = 0; i < simplex.length; ++i) {
+ PointCostPair pair = simplex[i];
+ if (! pair.isEvaluated()) {
+ pair.setCost(evaluateCost(pair.getPoint()));
+ }
+ }
+
+ // sort the simplex from lowest cost to highest cost
+ Arrays.sort(simplex, pointCostPairComparator);
+
+ }
+
+ /** Replace the worst point of the simplex by a new point.
+ * @param pointCostPair point to insert
+ */
+ protected void replaceWorstPoint(PointCostPair pointCostPair) {
+ int n = simplex.length - 1;
+ for (int i = 0; i < n; ++i) {
+ if (simplex[i].getCost() > pointCostPair.getCost()) {
+ PointCostPair tmp = simplex[i];
+ simplex[i] = pointCostPair;
+ pointCostPair = tmp;
+ }
+ }
+ simplex[n] = pointCostPair;
+ }
+
+ /** Comparator for {@link PointCostPair PointCostPair} objects. */
+ private static Comparator pointCostPairComparator = new Comparator() {
+ public int compare(Object o1, Object o2) {
+ if (o1 == null) {
+ return (o2 == null) ? 0 : +1;
+ } else if (o2 == null) {
+ return -1;
+ } else {
+ double cost1 = ((PointCostPair) o1).getCost();
+ double cost2 = ((PointCostPair) o2).getCost();
+ return (cost1 < cost2) ? -1 : ((o1 == o2) ? 0 : +1);
+ }
+ }
+ };
+
+ /** Simplex. */
+ protected PointCostPair[] simplex;
+
+ /** Cost function. */
+ private CostFunction f;
+
+ /** Number of evaluations already performed. */
+ private int evaluations;
+
+ /** Number of starts to go. */
+ private int starts;
+
+ /** Random generator for multi-start. */
+ private RandomVectorGenerator generator;
+
+ /** Found minima. */
+ private PointCostPair[] minima;
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/DirectSearchOptimizer.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/MultiDirectional.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/MultiDirectional.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/MultiDirectional.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/MultiDirectional.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,123 @@
+// 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.spaceroots.mantissa.optimization;
+
+/** This class implements the multi-directional direct search method.
+
+ * @author Luc Maisonobe
+ * @version $Id: MultiDirectional.java 1705 2006-09-17 19:57:39Z luc $
+ * @see NelderMead
+ */
+public class MultiDirectional
+ extends DirectSearchOptimizer {
+
+ /** Build a multi-directional optimizer with default coefficients.
+ * <p>The default values are 2.0 for khi and 0.5 for gamma.</p>
+ */
+ public MultiDirectional() {
+ super();
+ this.khi = 2.0;
+ this.gamma = 0.5;
+ }
+
+ /** Build a multi-directional optimizer with specified coefficients.
+ * @param khi expansion coefficient
+ * @param gamma contraction coefficient
+ */
+ public MultiDirectional(double khi, double gamma) {
+ super();
+ this.khi = khi;
+ this.gamma = gamma;
+ }
+
+ /** Compute the next simplex of the algorithm.
+ */
+ protected void iterateSimplex()
+ throws CostException {
+
+ while (true) {
+
+ // save the original vertex
+ PointCostPair[] original = simplex;
+ double originalCost = original[0].getCost();
+
+ // perform a reflection step
+ double reflectedCost = evaluateNewSimplex(original, 1.0);
+ if (reflectedCost < originalCost) {
+
+ // compute the expanded simplex
+ PointCostPair[] reflected = simplex;
+ double expandedCost = evaluateNewSimplex(original, khi);
+ if (reflectedCost <= expandedCost) {
+ // accept the reflected simplex
+ simplex = reflected;
+ }
+
+ return;
+
+ }
+
+ // compute the contracted simplex
+ double contractedCost = evaluateNewSimplex(original, gamma);
+ if (contractedCost < originalCost) {
+ // accept the contracted simplex
+ return;
+ }
+
+ }
+
+ }
+
+ /** Compute and evaluate a new simplex.
+ * @param original original simplex (to be preserved)
+ * @param coeff linear coefficient
+ * @return smallest cost in the transformed simplex
+ * @exception CostException if the function cannot be evaluated at
+ * some point
+ */
+ private double evaluateNewSimplex(PointCostPair[] original, double coeff)
+ throws CostException {
+
+ double[] xSmallest = original[0].getPoint();
+ int n = xSmallest.length;
+
+ // create the linearly transformed simplex
+ simplex = new PointCostPair[n + 1];
+ simplex[0] = original[0];
+ for (int i = 1; i <= n; ++i) {
+ double[] xOriginal = original[i].getPoint();
+ double[] xTransformed = new double[n];
+ for (int j = 0; j < n; ++j) {
+ xTransformed[j] = xSmallest[j] + coeff * (xSmallest[j] - xOriginal[j]);
+ }
+ simplex[i] = new PointCostPair(xTransformed);
+ }
+
+ // evaluate it
+ evaluateSimplex();
+ return simplex[0].getCost();
+
+ }
+
+ /** Expansion coefficient. */
+ private double khi;
+
+ /** Contraction coefficient. */
+ private double gamma;
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/MultiDirectional.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/NelderMead.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/NelderMead.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/NelderMead.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/NelderMead.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,184 @@
+// 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.spaceroots.mantissa.optimization;
+
+/** This class implements the Nelder-Mead direct search method.
+
+ * @author Luc Maisonobe
+ * @version $Id: NelderMead.java 1705 2006-09-17 19:57:39Z luc $
+ * @see MultiDirectional
+ */
+public class NelderMead
+ extends DirectSearchOptimizer {
+
+ /** Build a Nelder-Mead optimizer with default coefficients.
+ * <p>The default coefficients are 1.0 for rho, 2.0 for khi and 0.5
+ * for both gamma and sigma.</p>
+ */
+ public NelderMead() {
+ super();
+ this.rho = 1.0;
+ this.khi = 2.0;
+ this.gamma = 0.5;
+ this.sigma = 0.5;
+ }
+
+ /** Build a Nelder-Mead optimizer with specified coefficients.
+ * @param rho reflection coefficient
+ * @param khi expansion coefficient
+ * @param gamma contraction coefficient
+ * @param sigma shrinkage coefficient
+ */
+ public NelderMead(double rho, double khi, double gamma, double sigma) {
+ super();
+ this.rho = rho;
+ this.khi = khi;
+ this.gamma = gamma;
+ this.sigma = sigma;
+ }
+
+ /** Compute the next simplex of the algorithm.
+ */
+ protected void iterateSimplex()
+ throws CostException {
+
+ // the simplex has n+1 point if dimension is n
+ int n = simplex.length - 1;
+
+ // interesting costs
+ double smallest = simplex[0].getCost();
+ double secondLargest = simplex[n-1].getCost();
+ double largest = simplex[n].getCost();
+ double[] xLargest = simplex[n].getPoint();
+
+ // compute the centroid of the best vertices
+ // (dismissing the worst point at index n)
+ double[] centroid = new double[n];
+ for (int i = 0; i < n; ++i) {
+ double[] x = simplex[i].getPoint();
+ for (int j = 0; j < n; ++j) {
+ centroid[j] += x[j];
+ }
+ }
+ double scaling = 1.0 / n;
+ for (int j = 0; j < n; ++j) {
+ centroid[j] *= scaling;
+ }
+
+ // compute the reflection point
+ double[] xR = new double[n];
+ for (int j = 0; j < n; ++j) {
+ xR[j] = centroid[j] + rho * (centroid[j] - xLargest[j]);
+ }
+ double costR = evaluateCost(xR);
+
+ if ((smallest <= costR) && (costR < secondLargest)) {
+
+ // accept the reflected point
+ PointCostPair r = new PointCostPair(xR);
+ r.setCost(costR);
+ replaceWorstPoint(r);
+
+ } else if (costR < smallest) {
+
+ // compute the expansion point
+ double[] xE = new double[n];
+ for (int j = 0; j < n; ++j) {
+ xE[j] = centroid[j] + khi * (xR[j] - centroid[j]);
+ }
+ double costE = evaluateCost(xE);
+
+ if (costE < costR) {
+ // accept the expansion point
+ PointCostPair e = new PointCostPair(xE);
+ e.setCost(costE);
+ replaceWorstPoint(e);
+ } else {
+ // accept the reflected point
+ PointCostPair r = new PointCostPair(xR);
+ r.setCost(costR);
+ replaceWorstPoint(r);
+ }
+
+ } else {
+
+ if (costR < largest) {
+
+ // perform an outside contraction
+ double[] xC = new double[n];
+ for (int j = 0; j < n; ++j) {
+ xC[j] = centroid[j] + gamma * (xR[j] - centroid[j]);
+ }
+ double costC = evaluateCost(xC);
+
+ if (costC <= costR) {
+ // accept the contraction point
+ PointCostPair c = new PointCostPair(xC);
+ c.setCost(costC);
+ replaceWorstPoint(c);
+ return;
+ }
+
+ } else {
+
+ // perform an inside contraction
+ double[] xC = new double[n];
+ for (int j = 0; j < n; ++j) {
+ xC[j] = centroid[j] - gamma * (centroid[j] - xLargest[j]);
+ }
+ double costC = evaluateCost(xC);
+
+ if (costC < largest) {
+ // accept the contraction point
+ PointCostPair c = new PointCostPair(xC);
+ c.setCost(costC);
+ replaceWorstPoint(c);
+ return;
+ }
+
+ }
+
+ // perform a shrink
+ double[] xSmallest = simplex[0].getPoint();
+ for (int i = 1; i < simplex.length; ++i) {
+ PointCostPair pair = simplex[i];
+ double[] x = pair.getPoint();
+ for (int j = 0; j < n; ++j) {
+ x[j] = xSmallest[j] + sigma * (x[j] - xSmallest[j]);
+ }
+ pair.setCost(Double.NaN);
+ }
+ evaluateSimplex();
+
+ }
+
+ }
+
+ /** Reflection coefficient. */
+ private double rho;
+
+ /** Expansion coefficient. */
+ private double khi;
+
+ /** Contraction coefficient. */
+ private double gamma;
+
+ /** Shrinkage coefficient. */
+ private double sigma;
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/NelderMead.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/NoConvergenceException.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/NoConvergenceException.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/NoConvergenceException.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/NoConvergenceException.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,43 @@
+// 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.spaceroots.mantissa.optimization;
+
+import org.spaceroots.mantissa.MantissaException;
+
+/** This class represents exceptions thrown by optimization algorithms.
+
+ * @version $Id: NoConvergenceException.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class NoConvergenceException
+ extends MantissaException {
+
+ /** Simple constructor.
+ * Build an exception by translating and formating a message
+ * @param specifier format specifier (to be translated)
+ * @param parts to insert in the format (no translation)
+ */
+ public NoConvergenceException(String specifier, String[] parts) {
+ super(specifier, parts);
+ }
+
+ private static final long serialVersionUID = 4854864422540042859L;
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/NoConvergenceException.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/PointCostPair.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/PointCostPair.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/PointCostPair.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/PointCostPair.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,89 @@
+// 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.spaceroots.mantissa.optimization;
+
+/** This class holds a point and its associated cost.
+ * <p>A cost/point pair is not evaluated at build time. Its associated
+ * cost set to <code>Double.NaN</code> until it is evaluated.</p>
+ * @author Luc Maisonobe
+ * @version $Id: PointCostPair.java 1705 2006-09-17 19:57:39Z luc $
+ * @see CostFunction
+ */
+public class PointCostPair {
+
+ /** Build a point/cost pair with non-evaluated cost.
+ * @param point point coordinates
+ */
+ public PointCostPair(double[] point) {
+ this.point = point;
+ cost = Double.NaN;
+ }
+
+ /** Reset the point coordinates.
+ * <p>Resetting the points coordinates automatically reset the cost
+ * to non-evaluated</p>
+ * @param point new point coordinates
+ * @return old point coordinates (this can be re-used to put the
+ * coordinates of another point without re-allocating an array)
+ */
+ public double[] setPoint(double[] point) {
+ double[] oldPoint = this.point;
+ this.point = point;
+ cost = Double.NaN;
+ return oldPoint;
+ }
+
+ /** Get the point coordinates.
+ * @return point coordinates
+ */
+ public double[] getPoint() {
+ return point;
+ }
+
+ /** Set the cost.
+ * @param cost cost to store in the instance (can be
+ * <code>Double.NaN</code> to reset the instance to non-evaluated)
+ */
+ public void setCost(double cost) {
+ this.cost = cost;
+ }
+
+ /** Get the cost.
+ * @return cost associated to the point (or <code>Double.NaN</code>
+ * if the instance is not evaluated)
+ */
+ public double getCost() {
+ return cost;
+ }
+
+ /** Check if the cost has been evaluated.
+ * <p>The cost is considered to be non-evaluated if it is
+ * <code>Double.isNaN(pair.getCost())</code> would return true</p>
+ * @return true if the cost has been evaluated
+ */
+ public boolean isEvaluated() {
+ return ! Double.isNaN(cost);
+ }
+
+ /** Point coordinates. */
+ private double[] point;
+
+ /** Cost associated to the point. */
+ private double cost;
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/optimization/PointCostPair.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/ComputableFunctionIntegrator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/ComputableFunctionIntegrator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/ComputableFunctionIntegrator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/ComputableFunctionIntegrator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,49 @@
+// 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.spaceroots.mantissa.quadrature.scalar;
+
+import org.spaceroots.mantissa.functions.scalar.ComputableFunction;
+import org.spaceroots.mantissa.functions.FunctionException;
+
+/** This interface represents an integrator for scalar functions.
+
+ * <p>The classes which are devoted to integrate scalar functions
+ * should implement this interface. The functions which can be handled
+ * should implement the {@link
+ * org.spaceroots.mantissa.functions.scalar.ComputableFunction
+ * ComputableFunction} interface.</p>
+
+ * @see org.spaceroots.mantissa.functions.scalar.ComputableFunction
+
+ * @version $Id: ComputableFunctionIntegrator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public interface ComputableFunctionIntegrator {
+ /** Integrate a function over a defined range.
+ * @param f function to integrate
+ * @param a first bound of the range (can be lesser or greater than b)
+ * @param b second bound of the range (can be lesser or greater than a)
+ * @return value of the integral over the range
+ * @exception FunctionException if the underlying function throws one
+ */
+ public double integrate(ComputableFunction f, double a, double b)
+ throws FunctionException;
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/ComputableFunctionIntegrator.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/EnhancedSimpsonIntegrator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/EnhancedSimpsonIntegrator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/EnhancedSimpsonIntegrator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/EnhancedSimpsonIntegrator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,57 @@
+// 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.spaceroots.mantissa.quadrature.scalar;
+
+import org.spaceroots.mantissa.functions.FunctionException;
+import org.spaceroots.mantissa.functions.ExhaustedSampleException;
+import org.spaceroots.mantissa.functions.scalar.SampledFunctionIterator;
+
+/** This class implements an enhanced Simpson-like integrator.
+
+ * <p>A traditional Simpson integrator is based on a quadratic
+ * approximation of the function on three equally spaced points. This
+ * integrator does the same thing but can handle non-equally spaced
+ * points. If it is used on a regular sample, it behaves exactly as a
+ * traditional Simpson integrator.</p>
+
+ * @version $Id: EnhancedSimpsonIntegrator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class EnhancedSimpsonIntegrator
+ implements SampledFunctionIntegrator {
+ public double integrate(SampledFunctionIterator iter)
+ throws ExhaustedSampleException, FunctionException {
+
+ EnhancedSimpsonIntegratorSampler sampler
+ = new EnhancedSimpsonIntegratorSampler(iter);
+ double sum = 0.0;
+
+ try {
+ while (true) {
+ sum = sampler.nextSamplePoint().getY();
+ }
+ } catch(ExhaustedSampleException e) {
+ }
+
+ return sum;
+
+ }
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/EnhancedSimpsonIntegrator.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/EnhancedSimpsonIntegratorSampler.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/EnhancedSimpsonIntegratorSampler.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/EnhancedSimpsonIntegratorSampler.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/EnhancedSimpsonIntegratorSampler.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,100 @@
+// 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.spaceroots.mantissa.quadrature.scalar;
+
+import org.spaceroots.mantissa.functions.scalar.*;
+import org.spaceroots.mantissa.functions.FunctionException;
+import org.spaceroots.mantissa.functions.ExhaustedSampleException;
+
+/** This class implements an enhanced Simpson integrator as a sample.
+
+ * <p>A traditional Simpson integrator is based on a quadratic
+ * approximation of the function on three equally spaced points. This
+ * integrator does the same thing but can handle non-equally spaced
+ * points. If it is used on a regular sample, it behaves exactly as a
+ * traditional Simpson integrator.</p>
+
+ * @see EnhancedSimpsonIntegrator
+
+ * @version $Id: EnhancedSimpsonIntegratorSampler.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class EnhancedSimpsonIntegratorSampler
+ implements SampledFunctionIterator {
+
+ /** Underlying sampled function iterator. */
+ private SampledFunctionIterator iter;
+
+ /** Next point. */
+ private ScalarValuedPair next;
+
+ /** Current running sum. */
+ private double sum;
+
+ /** Constructor.
+ * Build an integrator from an underlying sample iterator.
+ * @param iter iterator over the base function
+ */
+ public EnhancedSimpsonIntegratorSampler(SampledFunctionIterator iter)
+ throws ExhaustedSampleException, FunctionException {
+
+ this.iter = iter;
+
+ // get the first point
+ next = iter.nextSamplePoint();
+
+ // initialize the sum
+ sum = 0.0;
+
+ }
+
+ public boolean hasNext() {
+ return iter.hasNext();
+ }
+
+ public ScalarValuedPair nextSamplePoint()
+ throws ExhaustedSampleException, FunctionException {
+ // performs one step of an enhanced Simpson scheme
+ ScalarValuedPair previous = next;
+ ScalarValuedPair current = iter.nextSamplePoint();
+
+ try {
+ next = iter.nextSamplePoint();
+
+ double h1 = current.getX() - previous.getX();
+ double h2 = next.getX() - current.getX();
+ double cP = (h1 + h2) * (2 * h1 - h2) / (6 * h1);
+ double cC = (h1 + h2) * (h1 + h2) * (h1 + h2) / (6 * h1 * h2);
+ double cN = (h1 + h2) * (2 * h2 - h1) / (6 * h2);
+
+ sum += cP * previous.getY() + cC * current.getY() + cN * next.getY();
+
+ } catch(ExhaustedSampleException e) {
+ // we have an incomplete step at the end of the sample
+ // we use a trapezoid scheme for this last step
+ sum += 0.5 * (current.getX() - previous.getX()) * (previous.getY() + current.getY());
+ return new ScalarValuedPair(current.getX(), sum);
+ }
+
+ return new ScalarValuedPair(next.getX(), sum);
+
+ }
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/EnhancedSimpsonIntegratorSampler.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/GaussLegendreIntegrator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/GaussLegendreIntegrator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/GaussLegendreIntegrator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/GaussLegendreIntegrator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,151 @@
+// 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.spaceroots.mantissa.quadrature.scalar;
+
+import org.spaceroots.mantissa.functions.scalar.ComputableFunction;
+import org.spaceroots.mantissa.functions.FunctionException;
+
+/** This class implements a Gauss-Legendre integrator.
+
+ * <p>Gauss-Legendre integrators are efficient integrators that can
+ * accurately integrate functions with few functions evaluations. A
+ * Gauss-Legendre integrator using an n-points quadrature formula can
+ * integrate exactly 2n-1 degree polynoms.</p>
+
+ * <p>These integrators evaluate the function on n carefully chosen
+ * points in each step interval. These points are not evenly
+ * spaced. The function is <emph>never</emph> evaluated at the
+ * boundary points, which means it can be undefined at these
+ * points.</p>
+
+ * @version $Id: GaussLegendreIntegrator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class GaussLegendreIntegrator
+ implements ComputableFunctionIntegrator {
+ /** Build a Gauss-Legendre integrator.
+
+ * <p>A Gauss-Legendre integrator is a formula like:
+ * <pre>
+ * int (f) from -1 to +1 = Sum (ai * f(xi))
+ * </pre>
+ * </p>
+ *
+ * <p>The coefficients of the formula are computed as follow:
+ * <pre>
+ * let n be the desired number of points
+ * the xi are the roots of the degree n Legendre polynomial
+ * the ai are the integrals int (Li^2) from -1 to +1
+ * where Li (x) = Prod (x-xk)/(xi-xk) for k != i
+ * </pre>
+ * </p>
+ *
+ * <p>A formula in n points can integrate exactly polynoms of degree
+ * up to 2n-1.</p>
+ *
+ * @param minPoints minimal number of points desired
+ * @param rawStep raw integration step (the precise step will be
+ * adjusted in order to have an integer number of steps in the
+ * integration range).
+ * */
+ public GaussLegendreIntegrator(int minPoints, double rawStep) {
+ if (minPoints <= 2) {
+ weightedRoots = new double[][] {
+ { 1.0, -1.0 / Math.sqrt(3.0) },
+ { 1.0, 1.0 / Math.sqrt(3.0) }
+ };
+ } else if (minPoints <= 3) {
+ weightedRoots = new double[][] {
+ { 5.0 / 9.0, -Math.sqrt(0.6) },
+ { 8.0 / 9.0, 0.0 },
+ { 5.0 / 9.0, Math.sqrt(0.6) }
+ };
+ } else if (minPoints <= 4) {
+ weightedRoots = new double[][] {
+ { (90.0 - 5.0 * Math.sqrt(30.0)) / 180.0,
+ -Math.sqrt((15.0 + 2.0 * Math.sqrt(30.0)) / 35.0) },
+ { (90.0 + 5.0 * Math.sqrt(30.0)) / 180.0,
+ -Math.sqrt((15.0 - 2.0 * Math.sqrt(30.0)) / 35.0) },
+ { (90.0 + 5.0 * Math.sqrt(30.0)) / 180.0,
+ Math.sqrt((15.0 - 2.0 * Math.sqrt(30.0)) / 35.0) },
+ { (90.0 - 5.0 * Math.sqrt(30.0)) / 180.0,
+ Math.sqrt((15.0 + 2.0 * Math.sqrt(30.0)) / 35.0) }
+ };
+ } else {
+ weightedRoots = new double[][] {
+ { (322.0 - 13.0 * Math.sqrt(70.0)) / 900.0,
+ -Math.sqrt((35.0 + 2.0 * Math.sqrt(70.0)) / 63.0) },
+ { (322.0 + 13.0 * Math.sqrt(70.0)) / 900.0,
+ -Math.sqrt((35.0 - 2.0 * Math.sqrt(70.0)) / 63.0) },
+ { 128.0 / 225.0,
+ 0.0 },
+ { (322.0 + 13.0 * Math.sqrt(70.0)) / 900.0,
+ Math.sqrt((35.0 - 2.0 * Math.sqrt(70.0)) / 63.0) },
+ { (322.0 - 13.0 * Math.sqrt(70.0)) / 900.0,
+ Math.sqrt((35.0 + 2.0 * Math.sqrt(70.0)) / 63.0) }
+ };
+ }
+
+ this.rawStep = rawStep;
+
+ }
+
+ /** Get the number of functions evaluation per step.
+ * @return number of function evaluation per step
+ */
+ public int getEvaluationsPerStep() {
+ return weightedRoots.length;
+ }
+
+ public double integrate(ComputableFunction f, double a, double b)
+ throws FunctionException {
+
+ // swap the bounds if they are not in ascending order
+ if (b < a) {
+ double tmp = b;
+ b = a;
+ a = tmp;
+ }
+
+ // adjust the step according to the bounds
+ long n = Math.round(0.5 + (b - a) / rawStep);
+ double step = (b - a) / n;
+
+ // integrate over all elementary steps
+ double halfStep = step / 2.0;
+ double midPoint = a + halfStep;
+ double sum = 0.0;
+ for (long i = 0; i < n; ++i) {
+ for (int j = 0; j < weightedRoots.length; ++j) {
+ sum += weightedRoots[j][0]
+ * f.valueAt(midPoint + halfStep * weightedRoots[j][1]);
+ }
+ midPoint += step;
+ }
+
+ return halfStep * sum;
+
+ }
+
+ double[][] weightedRoots;
+
+ double rawStep;
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/GaussLegendreIntegrator.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/RiemannIntegrator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/RiemannIntegrator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/RiemannIntegrator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/RiemannIntegrator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,62 @@
+// 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.spaceroots.mantissa.quadrature.scalar;
+
+import org.spaceroots.mantissa.functions.FunctionException;
+import org.spaceroots.mantissa.functions.ExhaustedSampleException;
+import org.spaceroots.mantissa.functions.scalar.SampledFunctionIterator;
+
+/** This class implements a Riemann integrator.
+
+ * <p>A Riemann integrator is a very simple one that assumes the
+ * function is constant over the integration step. Since it is very
+ * simple, this algorithm needs very small steps to achieve high
+ * accuracy, and small steps lead to numerical errors and
+ * instabilities.</p>
+
+ * <p>This algorithm is almost never used and has been included in
+ * this package only as a simple template for more useful
+ * integrators.</p>
+
+ * @see TrapezoidIntegrator
+
+ * @version $Id: RiemannIntegrator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class RiemannIntegrator
+ implements SampledFunctionIntegrator {
+ public double integrate(SampledFunctionIterator iter)
+ throws ExhaustedSampleException, FunctionException {
+
+ RiemannIntegratorSampler sampler = new RiemannIntegratorSampler(iter);
+ double sum = 0.0;
+
+ try {
+ while (true) {
+ sum = sampler.nextSamplePoint().getY();
+ }
+ } catch(ExhaustedSampleException e) {
+ }
+
+ return sum;
+
+ }
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/RiemannIntegrator.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/RiemannIntegratorSampler.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/RiemannIntegratorSampler.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/RiemannIntegratorSampler.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/RiemannIntegratorSampler.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,87 @@
+// 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.spaceroots.mantissa.quadrature.scalar;
+
+import org.spaceroots.mantissa.functions.scalar.*;
+import org.spaceroots.mantissa.functions.FunctionException;
+import org.spaceroots.mantissa.functions.ExhaustedSampleException;
+
+/** This class implements a Riemann integrator as a sample.
+
+ * <p>A Riemann integrator is a very simple one that assumes the
+ * function is constant over the integration step. Since it is very
+ * simple, this algorithm needs very small steps to achieve high
+ * accuracy, and small steps lead to numerical errors and
+ * instabilities.</p>
+
+ * <p>This algorithm is almost never used and has been included in
+ * this package only as a simple template for more useful
+ * integrators.</p>
+
+ * @see RiemannIntegrator
+
+ * @version $Id: RiemannIntegratorSampler.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class RiemannIntegratorSampler
+ implements SampledFunctionIterator {
+
+ /** Underlying sample iterator. */
+ private SampledFunctionIterator iter;
+
+ /** Current point. */
+ private ScalarValuedPair current;
+
+ /** Current running sum. */
+ private double sum;
+
+ /** Constructor.
+ * Build an integrator from an underlying sample iterator.
+ * @param iter iterator over the base function
+ */
+ public RiemannIntegratorSampler(SampledFunctionIterator iter)
+ throws ExhaustedSampleException, FunctionException {
+
+ this.iter = iter;
+
+ // get the first point
+ current = iter.nextSamplePoint();
+
+ // initialize the sum
+ sum = 0.0;
+
+ }
+
+ public boolean hasNext() {
+ return iter.hasNext();
+ }
+
+ public ScalarValuedPair nextSamplePoint()
+ throws ExhaustedSampleException, FunctionException {
+ // performs one step of a Riemann scheme
+ ScalarValuedPair previous = current;
+ current = iter.nextSamplePoint();
+ sum += (current.getX() - previous.getX()) * previous.getY();
+
+ return new ScalarValuedPair(current.getX(), sum);
+
+ }
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/RiemannIntegratorSampler.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/SampledFunctionIntegrator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/SampledFunctionIntegrator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/SampledFunctionIntegrator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/SampledFunctionIntegrator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,48 @@
+// 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.spaceroots.mantissa.quadrature.scalar;
+
+import org.spaceroots.mantissa.functions.scalar.SampledFunctionIterator;
+import org.spaceroots.mantissa.functions.ExhaustedSampleException;
+import org.spaceroots.mantissa.functions.FunctionException;
+
+/** This interface represents an integrator for scalar samples.
+
+ * <p>The classes which are devoted to integrate scalar samples
+ * should implement this interface.</p>
+
+ * @see org.spaceroots.mantissa.functions.scalar.SampledFunctionIterator
+ * @see ComputableFunctionIntegrator
+
+ * @version $Id: SampledFunctionIntegrator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public interface SampledFunctionIntegrator {
+ /** Integrate a sample over its overall range
+ * @param iter iterator over the sample to integrate
+ * @return value of the integral over the sample range
+ * @exception ExhaustedSampleException if the sample does not have
+ * enough points for the integration scheme
+ * @exception FunctionException if the underlying sampled function throws one
+ */
+ public double integrate(SampledFunctionIterator iter)
+ throws ExhaustedSampleException, FunctionException;
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/SampledFunctionIntegrator.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/TrapezoidIntegrator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/TrapezoidIntegrator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/TrapezoidIntegrator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/TrapezoidIntegrator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,53 @@
+// 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.spaceroots.mantissa.quadrature.scalar;
+
+import org.spaceroots.mantissa.functions.FunctionException;
+import org.spaceroots.mantissa.functions.ExhaustedSampleException;
+import org.spaceroots.mantissa.functions.scalar.SampledFunctionIterator;
+
+/** This class implements a trapezoid integrator.
+
+ * <p>A trapezoid integrator is a very simple one that assumes the
+ * function is linear over the integration step.</p>
+
+ * @version $Id: TrapezoidIntegrator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class TrapezoidIntegrator
+ implements SampledFunctionIntegrator {
+ public double integrate(SampledFunctionIterator iter)
+ throws ExhaustedSampleException, FunctionException {
+
+ TrapezoidIntegratorSampler sampler = new TrapezoidIntegratorSampler(iter);
+ double sum = 0.0;
+
+ try {
+ while (true) {
+ sum = sampler.nextSamplePoint().getY();
+ }
+ } catch (ExhaustedSampleException e) {
+ }
+
+ return sum;
+
+ }
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/TrapezoidIntegrator.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/TrapezoidIntegratorSampler.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/TrapezoidIntegratorSampler.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/TrapezoidIntegratorSampler.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/TrapezoidIntegratorSampler.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,89 @@
+// 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.spaceroots.mantissa.quadrature.scalar;
+
+import org.spaceroots.mantissa.functions.scalar.*;
+import org.spaceroots.mantissa.functions.FunctionException;
+import org.spaceroots.mantissa.functions.ExhaustedSampleException;
+
+/** This class implements a trapezoid integrator as a sample.
+
+ * <p>A trapezoid integrator is a very simple one that assumes the
+ * function is constant over the integration step. Since it is very
+ * simple, this algorithm needs very small steps to achieve high
+ * accuracy, and small steps lead to numerical errors and
+ * instabilities.</p>
+
+ * <p>This algorithm is almost never used and has been included in
+ * this package only as a simple template for more useful
+ * integrators.</p>
+
+ * @see TrapezoidIntegrator
+
+ * @version $Id: TrapezoidIntegratorSampler.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class TrapezoidIntegratorSampler
+ implements SampledFunctionIterator {
+
+ /** Underlying sample iterator. */
+ private SampledFunctionIterator iter;
+
+ /** Current point. */
+ private ScalarValuedPair current;
+
+ /** Current running sum. */
+ private double sum;
+
+ /** Constructor.
+ * Build an integrator from an underlying sample iterator.
+ * @param iter iterator over the base function
+ */
+ public TrapezoidIntegratorSampler(SampledFunctionIterator iter)
+ throws ExhaustedSampleException, FunctionException {
+
+ this.iter = iter;
+
+ // get the first point
+ current = iter.nextSamplePoint();
+
+ // initialize the sum
+ sum = 0.0;
+
+ }
+
+ public boolean hasNext() {
+ return iter.hasNext ();
+ }
+
+ public ScalarValuedPair nextSamplePoint()
+ throws ExhaustedSampleException, FunctionException {
+ // performs one step of a trapezoid scheme
+ ScalarValuedPair previous = current;
+ current = iter.nextSamplePoint();
+ sum += 0.5
+ * (current.getX() - previous.getX())
+ * (previous.getY() + current.getY());
+
+ return new ScalarValuedPair(current.getX(), sum);
+
+ }
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/scalar/TrapezoidIntegratorSampler.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/ComputableFunctionIntegrator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/ComputableFunctionIntegrator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/ComputableFunctionIntegrator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/ComputableFunctionIntegrator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,49 @@
+// 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.spaceroots.mantissa.quadrature.vectorial;
+
+import org.spaceroots.mantissa.functions.vectorial.ComputableFunction;
+import org.spaceroots.mantissa.functions.FunctionException;
+
+/** This interface represents an integrator for vectorial functions.
+
+ * <p>The classes which are devoted to integrate vectorial functions
+ * should implement this interface. The functions which can be handled
+ * should implement the {@link
+ * org.spaceroots.mantissa.functions.vectorial.ComputableFunction
+ * ComputableFunction} interface.</p>
+
+ * @see org.spaceroots.mantissa.functions.vectorial.ComputableFunction
+
+ * @version $Id: ComputableFunctionIntegrator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public interface ComputableFunctionIntegrator {
+ /** Integrate a function over a defined range.
+ * @param f function to integrate
+ * @param a first bound of the range (can be lesser or greater than b)
+ * @param b second bound of the range (can be lesser or greater than a)
+ * @return value of the integral over the range
+ * @exception FunctionException if the underlying function throws one
+ */
+ public double[] integrate(ComputableFunction f, double a, double b)
+ throws FunctionException;
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/ComputableFunctionIntegrator.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/EnhancedSimpsonIntegrator.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/EnhancedSimpsonIntegrator.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/EnhancedSimpsonIntegrator.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/EnhancedSimpsonIntegrator.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,57 @@
+// 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.spaceroots.mantissa.quadrature.vectorial;
+
+import org.spaceroots.mantissa.functions.FunctionException;
+import org.spaceroots.mantissa.functions.ExhaustedSampleException;
+import org.spaceroots.mantissa.functions.vectorial.SampledFunctionIterator;
+
+/** This class implements an enhanced Simpson-like integrator.
+
+ * <p>A traditional Simpson integrator is based on a quadratic
+ * approximation of the function on three equally spaced points. This
+ * integrator does the same thing but can handle non-equally spaced
+ * points. If it is used on a regular sample, it behaves exactly as a
+ * traditional Simpson integrator.</p>
+
+ * @version $Id: EnhancedSimpsonIntegrator.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class EnhancedSimpsonIntegrator
+ implements SampledFunctionIntegrator {
+ public double[] integrate(SampledFunctionIterator iter)
+ throws ExhaustedSampleException, FunctionException {
+
+ EnhancedSimpsonIntegratorSampler sampler =
+ new EnhancedSimpsonIntegratorSampler(iter);
+ double[] sum = null;
+
+ try {
+ while (true) {
+ sum = sampler.nextSamplePoint().getY();
+ }
+ } catch(ExhaustedSampleException e) {
+ }
+
+ return sum;
+
+ }
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/EnhancedSimpsonIntegrator.java
------------------------------------------------------------------------------
svn:eol-style = native
Added: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/EnhancedSimpsonIntegratorSampler.java
URL: http://svn.apache.org/viewvc/jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/EnhancedSimpsonIntegratorSampler.java?view=auto&rev=476930
==============================================================================
--- jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/EnhancedSimpsonIntegratorSampler.java (added)
+++ jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/EnhancedSimpsonIntegratorSampler.java Sun Nov 19 13:36:15 2006
@@ -0,0 +1,119 @@
+// 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.spaceroots.mantissa.quadrature.vectorial;
+
+import org.spaceroots.mantissa.functions.vectorial.*;
+import org.spaceroots.mantissa.functions.FunctionException;
+import org.spaceroots.mantissa.functions.ExhaustedSampleException;
+
+/** This class implements an enhanced Simpson integrator as a sample.
+
+ * <p>A traditional Simpson integrator is based on a quadratic
+ * approximation of the function on three equally spaced points. This
+ * integrator does the same thing but can handle non-equally spaced
+ * points. If it is used on a regular sample, it behaves exactly as a
+ * traditional Simpson integrator.</p>
+
+ * @see EnhancedSimpsonIntegrator
+
+ * @version $Id: EnhancedSimpsonIntegratorSampler.java 1705 2006-09-17 19:57:39Z luc $
+ * @author L. Maisonobe
+
+ */
+
+public class EnhancedSimpsonIntegratorSampler
+ implements SampledFunctionIterator {
+
+ /** Underlying sample iterator. */
+ private SampledFunctionIterator iter;
+
+ /** Next point. */
+ private VectorialValuedPair next;
+
+ /** Current running sum. */
+ private double[] sum;
+
+ /** Constructor.
+ * Build an integrator from an underlying sample iterator.
+ * @param iter iterator over the base function
+ */
+ public EnhancedSimpsonIntegratorSampler(SampledFunctionIterator iter)
+ throws ExhaustedSampleException, FunctionException {
+
+ this.iter = iter;
+
+ // get the first point
+ next = iter.nextSamplePoint();
+
+ // initialize the sum
+ sum = new double[iter.getDimension()];
+ for (int i = 0; i < sum.length; ++i) {
+ sum[i] = 0.0;
+ }
+
+ }
+
+ public boolean hasNext() {
+ return iter.hasNext();
+ }
+
+ public int getDimension() {
+ return iter.getDimension();
+ }
+
+ public VectorialValuedPair nextSamplePoint()
+ throws ExhaustedSampleException, FunctionException {
+ // performs one step of an enhanced Simpson scheme
+ VectorialValuedPair previous = next;
+ VectorialValuedPair current = iter.nextSamplePoint();
+
+ try {
+ next = iter.nextSamplePoint();
+
+ double h1 = current.getX() - previous.getX();
+ double h2 = next.getX() - current.getX();
+ double cP = (h1 + h2) * (2 * h1 - h2) / (6 * h1);
+ double cC = (h1 + h2) * (h1 + h2) * (h1 + h2) / (6 * h1 * h2);
+ double cN = (h1 + h2) * (2 * h2 - h1) / (6 * h2);
+
+ double[] pY = previous.getY();
+ double[] cY = current.getY();
+ double[] nY = next.getY();
+ for (int i = 0; i < sum.length; ++i) {
+ sum [i] += cP * pY[i] + cC * cY[i] + cN * nY[i];
+ }
+
+ } catch(ExhaustedSampleException e) {
+ // we have an incomplete step at the end of the sample
+ // we use a trapezoid scheme for this last step
+ double halfDx = 0.5 * (current.getX() - previous.getX());
+ double[] pY = previous.getY();
+ double[] cY = current.getY();
+ for (int i = 0; i < sum.length; ++i) {
+ sum [i] += halfDx * (pY[i] + cY[i]);
+ }
+ return new VectorialValuedPair(current.getX(), sum);
+ }
+
+ double[] values = new double[sum.length];
+ System.arraycopy(sum, 0, values, 0, sum.length);
+ return new VectorialValuedPair(next.getX(), values);
+
+ }
+
+}
Propchange: jakarta/commons/proper/math/trunk/src/mantissa/src/org/spaceroots/mantissa/quadrature/vectorial/EnhancedSimpsonIntegratorSampler.java
------------------------------------------------------------------------------
svn:eol-style = native
---------------------------------------------------------------------
To unsubscribe, e-mail: commons-dev-unsubscribe@jakarta.apache.org
For additional commands, e-mail: commons-dev-help@jakarta.apache.org