Use the new differentiation framework in optimization package.

As a temporary hack for backward compatibility, some optimizers do
provide "optimize" methods which accept the new functions signatures
(i.e. MultivariateDifferentiableXxxFunction), but they do *not*
advertise it in their "implements" clause. The reason for that is that
Java forbid implementing a parameterized interface with two different
parameters. These optimizers already implement the interface with the
older functions signatures, and this cannot be removed as of 3.1.

When 4.0 will be started, the old signatures will be removed, and the
"implements" declarations will be updated.

git-svn-id: https://svn.apache.org/repos/asf/commons/proper/math/trunk@1384907 13f79535-47bb-0310-9956-ffa450edef68
This commit is contained in:
Luc Maisonobe 2012-09-14 20:17:00 +00:00
parent 8bac3361f2
commit c8b8e61243
36 changed files with 1140 additions and 1130 deletions

View File

@ -22,6 +22,7 @@ package org.apache.commons.math3.analysis;
* multivariate real function.
* @version $Id$
* @since 2.0
* @deprecated as of 3.1 replaced by {@link org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableFunction}
*/
public interface DifferentiableMultivariateFunction extends MultivariateFunction {

View File

@ -17,12 +17,12 @@
package org.apache.commons.math3.analysis;
/**
* Extension of {@link MultivariateVectorFunction} representing a differentiable
* multivariate vectorial function.
* @version $Id$
* @since 2.0
* @deprecated as of 3.1 replaced by {@link org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableVectorFunction}
*/
public interface DifferentiableMultivariateVectorFunction
extends MultivariateVectorFunction {

View File

@ -25,7 +25,7 @@ import org.apache.commons.math3.analysis.MultivariateFunction;
* the following interfaces:
* <ul>
* <li>{@link org.apache.commons.math3.optimization.MultivariateOptimizer}</li>
* <li>{@link org.apache.commons.math3.optimization.DifferentiableMultivariateOptimizer}</li>
* <li>{@link org.apache.commons.math3.optimization.MultivariateDifferentiableOptimizer}</li>
* </ul>
*
* @param <FUNC> Type of the objective function to be optimized.

View File

@ -25,7 +25,7 @@ import org.apache.commons.math3.analysis.MultivariateFunction;
* the following interfaces:
* <ul>
* <li>{@link org.apache.commons.math3.optimization.MultivariateOptimizer}</li>
* <li>{@link org.apache.commons.math3.optimization.DifferentiableMultivariateOptimizer}</li>
* <li>{@link org.apache.commons.math3.optimization.MultivariateDifferentiableOptimizer}</li>
* </ul>
*
* @param <FUNC> Type of the objective function to be optimized.

View File

@ -23,8 +23,8 @@ package org.apache.commons.math3.optimization;
* the following interfaces:
* <ul>
* <li>{@link org.apache.commons.math3.optimization.MultivariateOptimizer}</li>
* <li>{@link org.apache.commons.math3.optimization.DifferentiableMultivariateOptimizer}</li>
* <li>{@link org.apache.commons.math3.optimization.DifferentiableMultivariateVectorOptimizer}</li>
* <li>{@link org.apache.commons.math3.optimization.MultivariateDifferentiableOptimizer}</li>
* <li>{@link org.apache.commons.math3.optimization.MultivariateDifferentiableVectorOptimizer}</li>
* <li>{@link org.apache.commons.math3.optimization.univariate.UnivariateOptimizer}</li>
* </ul>
*

View File

@ -30,6 +30,7 @@ import org.apache.commons.math3.random.RandomVectorGenerator;
*
* @version $Id$
* @since 2.0
* @deprecated as of 3.1 replaced by {@link MultivariateDifferentiableMultiStartOptimizer}
*/
public class DifferentiableMultivariateMultiStartOptimizer
extends BaseMultivariateMultiStartOptimizer<DifferentiableMultivariateFunction>
@ -44,8 +45,8 @@ public class DifferentiableMultivariateMultiStartOptimizer
* @param generator Random vector generator to use for restarts.
*/
public DifferentiableMultivariateMultiStartOptimizer(final DifferentiableMultivariateOptimizer optimizer,
final int starts,
final RandomVectorGenerator generator) {
final int starts,
final RandomVectorGenerator generator) {
super(optimizer, starts, generator);
}
}

View File

@ -31,6 +31,8 @@ import org.apache.commons.math3.analysis.DifferentiableMultivariateFunction;
*
* @version $Id$
* @since 2.0
* @deprecated as of 3.1 replaced by {@link MultivariateDifferentiableOptimizer}
*/
@Deprecated
public interface DifferentiableMultivariateOptimizer
extends BaseMultivariateOptimizer<DifferentiableMultivariateFunction> {}

View File

@ -30,7 +30,9 @@ import org.apache.commons.math3.random.RandomVectorGenerator;
*
* @version $Id$
* @since 2.0
* @deprecated as of 3.1 replaced by {@link MultivariateDifferentiableVectorMultiStartOptimizer}
*/
@Deprecated
public class DifferentiableMultivariateVectorMultiStartOptimizer
extends BaseMultivariateVectorMultiStartOptimizer<DifferentiableMultivariateVectorFunction>
implements DifferentiableMultivariateVectorOptimizer {

View File

@ -26,6 +26,8 @@ import org.apache.commons.math3.analysis.DifferentiableMultivariateVectorFunctio
*
* @version $Id$
* @since 3.0
* @deprecated as of 3.1 replaced by {@link MultivariateDifferentiableVectorOptimizer}
*/
@Deprecated
public interface DifferentiableMultivariateVectorOptimizer
extends BaseMultivariateVectorOptimizer<DifferentiableMultivariateVectorFunction> {}

View File

@ -0,0 +1,51 @@
/*
* Licensed to the Apache Software Foundation (ASF) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* The ASF licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.apache.commons.math3.optimization;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableFunction;
import org.apache.commons.math3.random.RandomVectorGenerator;
/**
* Special implementation of the {@link MultivariateDifferentiableOptimizer}
* interface adding multi-start features to an existing optimizer.
*
* This class wraps a classical optimizer to use it several times in
* turn with different starting points in order to avoid being trapped
* into a local extremum when looking for a global one.
*
* @version $Id$
* @since 3.1
*/
public class MultivariateDifferentiableMultiStartOptimizer
extends BaseMultivariateMultiStartOptimizer<MultivariateDifferentiableFunction>
implements MultivariateDifferentiableOptimizer {
/**
* Create a multi-start optimizer from a single-start optimizer.
*
* @param optimizer Single-start optimizer to wrap.
* @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 MultivariateDifferentiableMultiStartOptimizer(final MultivariateDifferentiableOptimizer optimizer,
final int starts,
final RandomVectorGenerator generator) {
super(optimizer, starts, generator);
}
}

View File

@ -0,0 +1,36 @@
/*
* Licensed to the Apache Software Foundation (ASF) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* The ASF licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.apache.commons.math3.optimization;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableFunction;
/**
* This interface represents an optimization algorithm for
* {@link MultivariateDifferentiableFunction scalar differentiable objective
* functions}.
* Optimization algorithms find the input point set that either {@link GoalType
* maximize or minimize} an objective function.
*
* @see MultivariateOptimizer
* @see MultivariateDifferentiableVectorOptimizer
*
* @version $Id$
* @since 3.1
*/
public interface MultivariateDifferentiableOptimizer
extends BaseMultivariateOptimizer<MultivariateDifferentiableFunction> {}

View File

@ -0,0 +1,52 @@
/*
* Licensed to the Apache Software Foundation (ASF) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* The ASF licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.apache.commons.math3.optimization;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableVectorFunction;
import org.apache.commons.math3.random.RandomVectorGenerator;
/**
* Special implementation of the {@link MultivariateDifferentiableVectorOptimizer}
* interface adding multi-start features to an existing optimizer.
*
* This class wraps a classical optimizer to use it several times in
* turn with different starting points in order to avoid being trapped
* into a local extremum when looking for a global one.
*
* @version $Id$
* @since 3.1
*/
public class MultivariateDifferentiableVectorMultiStartOptimizer
extends BaseMultivariateVectorMultiStartOptimizer<MultivariateDifferentiableVectorFunction>
implements MultivariateDifferentiableVectorOptimizer {
/**
* Create a multi-start optimizer from a single-start optimizer.
*
* @param optimizer Single-start optimizer to wrap.
* @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 MultivariateDifferentiableVectorMultiStartOptimizer(
final MultivariateDifferentiableVectorOptimizer optimizer,
final int starts,
final RandomVectorGenerator generator) {
super(optimizer, starts, generator);
}
}

View File

@ -0,0 +1,31 @@
/*
* Licensed to the Apache Software Foundation (ASF) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* The ASF licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.apache.commons.math3.optimization;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableVectorFunction;
/**
* This interface represents an optimization algorithm for
* {@link MultivariateDifferentiableVectorFunction differentiable vectorial
* objective functions}.
*
* @version $Id$
* @since 3.1
*/
public interface MultivariateDifferentiableVectorOptimizer
extends BaseMultivariateVectorOptimizer<MultivariateDifferentiableVectorFunction> {}

View File

@ -25,8 +25,8 @@ import org.apache.commons.math3.analysis.MultivariateFunction;
* <p>Optimization algorithms find the input point set that either {@link GoalType
* maximize or minimize} an objective function.</p>
*
* @see DifferentiableMultivariateOptimizer
* @see DifferentiableMultivariateVectorOptimizer
* @see MultivariateDifferentiableOptimizer
* @see MultivariateDifferentiableVectorOptimizer
* @version $Id$
* @since 2.0
*/

View File

@ -102,6 +102,28 @@ public abstract class BaseAbstractMultivariateOptimizer<FUNC extends Multivariat
/** {@inheritDoc} */
public PointValuePair optimize(int maxEval, FUNC f, GoalType goalType,
double[] startPoint) {
return optimizeInternal(maxEval, f, goalType, startPoint);
}
/**
* Optimize an objective function.
*
* @param f Objective function.
* @param goalType Type of optimization goal: either
* {@link GoalType#MAXIMIZE} or {@link GoalType#MINIMIZE}.
* @param startPoint Start point for optimization.
* @param maxEval Maximum number of function evaluations.
* @return the point/value pair giving the optimal value for objective
* function.
* @throws org.apache.commons.math3.exception.DimensionMismatchException
* if the start point dimension is wrong.
* @throws org.apache.commons.math3.exception.TooManyEvaluationsException
* if the maximal number of evaluations is exceeded.
* @throws org.apache.commons.math3.exception.NullArgumentException if
* any argument is {@code null}.
*/
protected PointValuePair optimizeInternal(int maxEval, MultivariateFunction f, GoalType goalType,
double[] startPoint) {
// Checks.
if (f == null) {
throw new NullArgumentException();

View File

@ -103,7 +103,33 @@ public abstract class BaseAbstractMultivariateVectorOptimizer<FUNC extends Multi
/** {@inheritDoc} */
public PointVectorValuePair optimize(int maxEval, FUNC f, double[] t, double[] w,
double[] startPoint) {
double[] startPoint) {
return optimizeInternal(maxEval, f, t, w, startPoint);
}
/**
* Optimize an objective function.
* Optimization is considered to be a weighted least-squares minimization.
* The cost function to be minimized is
* <code>&sum;weight<sub>i</sub>(objective<sub>i</sub> - target<sub>i</sub>)<sup>2</sup></code>
*
* @param f Objective function.
* @param target Target value for the objective functions at optimum.
* @param weight Weights for the least squares cost computation.
* @param startPoint Start point for optimization.
* @return the point/value pair giving the optimal value for objective
* function.
* @param maxEval Maximum number of function evaluations.
* @throws org.apache.commons.math3.exception.DimensionMismatchException
* if the start point dimension is wrong.
* @throws org.apache.commons.math3.exception.TooManyEvaluationsException
* if the maximal number of evaluations is exceeded.
* @throws org.apache.commons.math3.exception.NullArgumentException if
* any argument is {@code null}.
*/
protected PointVectorValuePair optimizeInternal(final int maxEval, final MultivariateVectorFunction f,
final double[] t, final double[] w,
final double[] startPoint) {
// Checks.
if (f == null) {
throw new NullArgumentException();
@ -133,6 +159,7 @@ public abstract class BaseAbstractMultivariateVectorOptimizer<FUNC extends Multi
// Perform computation.
return doOptimize();
}
/**

View File

@ -21,9 +21,12 @@ import java.util.ArrayList;
import java.util.List;
import org.apache.commons.math3.analysis.DifferentiableMultivariateVectorFunction;
import org.apache.commons.math3.analysis.ParametricUnivariateFunction;
import org.apache.commons.math3.analysis.MultivariateMatrixFunction;
import org.apache.commons.math3.analysis.ParametricUnivariateFunction;
import org.apache.commons.math3.analysis.differentiation.DerivativeStructure;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableVectorFunction;
import org.apache.commons.math3.optimization.DifferentiableMultivariateVectorOptimizer;
import org.apache.commons.math3.optimization.MultivariateDifferentiableVectorOptimizer;
import org.apache.commons.math3.optimization.PointVectorValuePair;
/** Fitter for parametric univariate real functions y = f(x).
@ -44,17 +47,37 @@ import org.apache.commons.math3.optimization.PointVectorValuePair;
* @since 2.0
*/
public class CurveFitter<T extends ParametricUnivariateFunction> {
/** Optimizer to use for the fitting.
* @deprecated as of 3.1 replaced by {@link #optimizer}
*/
@Deprecated
private final DifferentiableMultivariateVectorOptimizer oldOptimizer;
/** Optimizer to use for the fitting. */
private final DifferentiableMultivariateVectorOptimizer optimizer;
private final MultivariateDifferentiableVectorOptimizer optimizer;
/** Observed points. */
private final List<WeightedObservedPoint> observations;
/** Simple constructor.
* @param optimizer optimizer to use for the fitting
* @deprecated as of 3.1 replaced by {@link #CurveFitter(MultivariateDifferentiableVectorOptimizer)}
*/
public CurveFitter(final DifferentiableMultivariateVectorOptimizer optimizer) {
this.optimizer = optimizer;
observations = new ArrayList<WeightedObservedPoint>();
this.oldOptimizer = optimizer;
this.optimizer = null;
observations = new ArrayList<WeightedObservedPoint>();
}
/** Simple constructor.
* @param optimizer optimizer to use for the fitting
* @since 3.1
*/
public CurveFitter(final MultivariateDifferentiableVectorOptimizer optimizer) {
this.oldOptimizer = null;
this.optimizer = optimizer;
observations = new ArrayList<WeightedObservedPoint>();
}
/** Add an observed (x,y) point to the sample with unit weight.
@ -158,16 +181,23 @@ public class CurveFitter<T extends ParametricUnivariateFunction> {
}
// perform the fit
PointVectorValuePair optimum =
optimizer.optimize(maxEval, new TheoreticalValuesFunction(f),
target, weights, initialGuess);
final PointVectorValuePair optimum;
if (optimizer == null) {
// to be removed in 4.0
optimum = oldOptimizer.optimize(maxEval, new OldTheoreticalValuesFunction(f),
target, weights, initialGuess);
} else {
optimum = optimizer.optimize(maxEval, new TheoreticalValuesFunction(f),
target, weights, initialGuess);
}
// extract the coefficients
return optimum.getPointRef();
}
/** Vectorial function computing function theoretical values. */
private class TheoreticalValuesFunction
@Deprecated
private class OldTheoreticalValuesFunction
implements DifferentiableMultivariateVectorFunction {
/** Function to fit. */
private final ParametricUnivariateFunction f;
@ -175,7 +205,7 @@ public class CurveFitter<T extends ParametricUnivariateFunction> {
/** Simple constructor.
* @param f function to fit.
*/
public TheoreticalValuesFunction(final ParametricUnivariateFunction f) {
public OldTheoreticalValuesFunction(final ParametricUnivariateFunction f) {
this.f = f;
}
@ -207,4 +237,60 @@ public class CurveFitter<T extends ParametricUnivariateFunction> {
return values;
}
}
/** Vectorial function computing function theoretical values. */
private class TheoreticalValuesFunction implements MultivariateDifferentiableVectorFunction {
/** Function to fit. */
private final ParametricUnivariateFunction f;
/** Simple constructor.
* @param f function to fit.
*/
public TheoreticalValuesFunction(final ParametricUnivariateFunction f) {
this.f = f;
}
/** {@inheritDoc} */
public double[] value(double[] point) {
// compute the residuals
final double[] values = new double[observations.size()];
int i = 0;
for (WeightedObservedPoint observed : observations) {
values[i++] = f.value(observed.getX(), point);
}
return values;
}
/** {@inheritDoc} */
public DerivativeStructure[] value(DerivativeStructure[] point) {
// extract parameters
final double[] parameters = new double[point.length];
for (int k = 0; k < point.length; ++k) {
parameters[k] = point[k].getValue();
}
// compute the residuals
final DerivativeStructure[] values = new DerivativeStructure[observations.size()];
int i = 0;
for (WeightedObservedPoint observed : observations) {
// build the DerivativeStructure by adding first the value as a constant
// and then adding derivatives
DerivativeStructure vi = new DerivativeStructure(point.length, 1, f.value(observed.getX(), parameters));
for (int k = 0; k < point.length; ++k) {
vi = vi.add(new DerivativeStructure(point.length, 1, k, 0.0));
}
values[i++] = vi;
}
return values;
}
}
}

View File

@ -0,0 +1,76 @@
/*
* Licensed to the Apache Software Foundation (ASF) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* The ASF licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.apache.commons.math3.optimization.general;
import org.apache.commons.math3.analysis.MultivariateVectorFunction;
import org.apache.commons.math3.analysis.differentiation.GradientFunction;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableFunction;
import org.apache.commons.math3.optimization.ConvergenceChecker;
import org.apache.commons.math3.optimization.GoalType;
import org.apache.commons.math3.optimization.MultivariateDifferentiableOptimizer;
import org.apache.commons.math3.optimization.PointValuePair;
import org.apache.commons.math3.optimization.direct.BaseAbstractMultivariateOptimizer;
/**
* Base class for implementing optimizers for multivariate scalar
* differentiable functions.
* It contains boiler-plate code for dealing with gradient evaluation.
*
* @version $Id$
* @since 3.1
*/
public abstract class AbstractDifferentiableOptimizer
extends BaseAbstractMultivariateOptimizer<MultivariateDifferentiableFunction>
implements MultivariateDifferentiableOptimizer {
/**
* Objective function gradient.
*/
private MultivariateVectorFunction gradient;
/**
* @param checker Convergence checker.
*/
protected AbstractDifferentiableOptimizer(ConvergenceChecker<PointValuePair> checker) {
super(checker);
}
/**
* Compute the gradient vector.
*
* @param evaluationPoint Point at which the gradient must be evaluated.
* @return the gradient at the specified point.
*/
protected double[] computeObjectiveGradient(final double[] evaluationPoint) {
return gradient.value(evaluationPoint);
}
/** {@inheritDoc} */
@Override
public PointValuePair optimize(final int maxEval, final MultivariateDifferentiableFunction f,
final GoalType goalType, final double[] startPoint) {
// store optimization problem characteristics
gradient = new GradientFunction(f);
// perform optimization
return super.optimize(maxEval, f, goalType, startPoint);
}
}

View File

@ -21,6 +21,8 @@ import org.apache.commons.math3.exception.NumberIsTooSmallException;
import org.apache.commons.math3.exception.DimensionMismatchException;
import org.apache.commons.math3.analysis.DifferentiableMultivariateVectorFunction;
import org.apache.commons.math3.analysis.MultivariateMatrixFunction;
import org.apache.commons.math3.analysis.differentiation.JacobianFunction;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableVectorFunction;
import org.apache.commons.math3.exception.util.LocalizedFormats;
import org.apache.commons.math3.linear.QRDecomposition;
import org.apache.commons.math3.linear.DecompositionSolver;
@ -36,11 +38,10 @@ import org.apache.commons.math3.util.FastMath;
* It handles the boilerplate methods associated to thresholds settings,
* jacobian and error estimation.
* <br/>
* This class uses the {@link DifferentiableMultivariateVectorFunction#jacobian()}
* of the function argument in method
* {@link #optimize(int,DifferentiableMultivariateVectorFunction,double[],double[],double[])
* This class uses the {@link JacobianFunction Jacobian} of the function argument in method
* {@link #optimize(int, MultivariateDifferentiableVectorFunction, double[], double[], double[])
* optimize} and assumes that, in the matrix returned by the
* {@link MultivariateMatrixFunction#value(double[]) value} method, the rows
* {@link JacobianFunction#value(double[]) value} method, the rows
* iterate on the model functions while the columns iterate on the parameters; thus,
* the numbers of rows is equal to the length of the {@code target} array while the
* number of columns is equal to the length of the {@code startPoint} array.
@ -292,8 +293,12 @@ public abstract class AbstractLeastSquaresOptimizer
return sig;
}
/** {@inheritDoc} */
/** {@inheritDoc}
* @deprecated as of 3.1 replaced by {@link #optimize(int,
* MultivariateDifferentiableVectorFunction, double[], double[], double[])}
*/
@Override
@Deprecated
public PointVectorValuePair optimize(int maxEval,
final DifferentiableMultivariateVectorFunction f,
final double[] target, final double[] weights,
@ -314,6 +319,51 @@ public abstract class AbstractLeastSquaresOptimizer
cost = Double.POSITIVE_INFINITY;
return super.optimize(maxEval, f, target, weights, startPoint);
return optimizeInternal(maxEval, f, target, weights, startPoint);
}
/**
* Optimize an objective function.
* Optimization is considered to be a weighted least-squares minimization.
* The cost function to be minimized is
* <code>&sum;weight<sub>i</sub>(objective<sub>i</sub> - target<sub>i</sub>)<sup>2</sup></code>
*
* @param f Objective function.
* @param target Target value for the objective functions at optimum.
* @param weight Weights for the least squares cost computation.
* @param startPoint Start point for optimization.
* @return the point/value pair giving the optimal value for objective
* function.
* @param maxEval Maximum number of function evaluations.
* @throws org.apache.commons.math3.exception.DimensionMismatchException
* if the start point dimension is wrong.
* @throws org.apache.commons.math3.exception.TooManyEvaluationsException
* if the maximal number of evaluations is exceeded.
* @throws org.apache.commons.math3.exception.NullArgumentException if
* any argument is {@code null}.
*/
public PointVectorValuePair optimize(final int maxEval,
final MultivariateDifferentiableVectorFunction f,
final double[] target, final double[] weights,
final double[] startPoint) {
// Reset counter.
jacobianEvaluations = 0;
// Store least squares problem characteristics.
jF = new JacobianFunction(f);
// Arrays shared with the other private methods.
point = startPoint.clone();
rows = target.length;
cols = point.length;
weightedResidualJacobian = new double[rows][cols];
this.weightedResiduals = new double[rows];
cost = Double.POSITIVE_INFINITY;
return optimizeInternal(maxEval, f, target, weights, startPoint);
}
}

View File

@ -19,6 +19,8 @@ package org.apache.commons.math3.optimization.general;
import org.apache.commons.math3.analysis.DifferentiableMultivariateFunction;
import org.apache.commons.math3.analysis.MultivariateVectorFunction;
import org.apache.commons.math3.analysis.differentiation.GradientFunction;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableFunction;
import org.apache.commons.math3.optimization.DifferentiableMultivariateOptimizer;
import org.apache.commons.math3.optimization.GoalType;
import org.apache.commons.math3.optimization.ConvergenceChecker;
@ -32,7 +34,9 @@ import org.apache.commons.math3.optimization.direct.BaseAbstractMultivariateOpti
*
* @version $Id$
* @since 2.0
* @deprecated as of 3.1 replaced by {@link AbstractDifferentiableOptimizer}
*/
@Deprecated
public abstract class AbstractScalarDifferentiableOptimizer
extends BaseAbstractMultivariateOptimizer<DifferentiableMultivariateFunction>
implements DifferentiableMultivariateOptimizer {
@ -79,6 +83,33 @@ public abstract class AbstractScalarDifferentiableOptimizer
// Store optimization problem characteristics.
gradient = f.gradient();
return super.optimize(maxEval, f, goalType, startPoint);
return optimizeInternal(maxEval, f, goalType, startPoint);
}
/**
* Optimize an objective function.
*
* @param f Objective function.
* @param goalType Type of optimization goal: either
* {@link GoalType#MAXIMIZE} or {@link GoalType#MINIMIZE}.
* @param startPoint Start point for optimization.
* @param maxEval Maximum number of function evaluations.
* @return the point/value pair giving the optimal value for objective
* function.
* @throws org.apache.commons.math3.exception.DimensionMismatchException
* if the start point dimension is wrong.
* @throws org.apache.commons.math3.exception.TooManyEvaluationsException
* if the maximal number of evaluations is exceeded.
* @throws org.apache.commons.math3.exception.NullArgumentException if
* any argument is {@code null}.
*/
public PointValuePair optimize(final int maxEval,
final MultivariateDifferentiableFunction f,
final GoalType goalType,
final double[] startPoint) {
// Store optimization problem characteristics.
gradient = new GradientFunction(f);
return optimizeInternal(maxEval, f, goalType, startPoint);
}
}

View File

@ -38,13 +38,13 @@
* <li>{@link org.apache.commons.math3.optimization.MultivariateOptimizer
* MultivariateOptimizer} for {@link org.apache.commons.math3.analysis.MultivariateFunction
* multivariate real functions}</li>
* <li>{@link org.apache.commons.math3.optimization.DifferentiableMultivariateOptimizer
* DifferentiableMultivariateOptimizer} for {@link
* org.apache.commons.math3.analysis.DifferentiableMultivariateFunction
* <li>{@link org.apache.commons.math3.optimization.MultivariateDifferentiableOptimizer
* MultivariateDifferentiableOptimizer} for {@link
* org.apache.commons.math3.analysis.MultivariateDifferentiableFunction
* differentiable multivariate real functions}</li>
* <li>{@link org.apache.commons.math3.optimization.DifferentiableMultivariateVectorOptimizer
* DifferentiableMultivariateVectorOptimizer} for {@link
* org.apache.commons.math3.analysis.DifferentiableMultivariateVectorFunction
* <li>{@link org.apache.commons.math3.optimization.MultivariateDifferentiableVectorOptimizer
* MultivariateDifferentiableVectorOptimizer} for {@link
* org.apache.commons.math3.analysis.MultivariateDifferentiableVectorFunction
* differentiable multivariate vectorial functions}</li>
* </ul>
* </p>

View File

@ -1,140 +0,0 @@
/*
* Licensed to the Apache Software Foundation (ASF) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* The ASF licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.apache.commons.math3.optimization;
import java.awt.geom.Point2D;
import java.util.ArrayList;
import org.apache.commons.math3.analysis.DifferentiableMultivariateFunction;
import org.apache.commons.math3.analysis.MultivariateFunction;
import org.apache.commons.math3.analysis.MultivariateVectorFunction;
import org.apache.commons.math3.optimization.general.ConjugateGradientFormula;
import org.apache.commons.math3.optimization.general.NonLinearConjugateGradientOptimizer;
import org.apache.commons.math3.random.GaussianRandomGenerator;
import org.apache.commons.math3.random.JDKRandomGenerator;
import org.apache.commons.math3.random.RandomVectorGenerator;
import org.apache.commons.math3.random.UncorrelatedRandomVectorGenerator;
import org.junit.Assert;
import org.junit.Test;
public class DifferentiableMultivariateMultiStartOptimizerTest {
@Test
public void testCircleFitting() {
Circle circle = new Circle();
circle.addPoint( 30.0, 68.0);
circle.addPoint( 50.0, -6.0);
circle.addPoint(110.0, -20.0);
circle.addPoint( 35.0, 15.0);
circle.addPoint( 45.0, 97.0);
NonLinearConjugateGradientOptimizer underlying =
new NonLinearConjugateGradientOptimizer(ConjugateGradientFormula.POLAK_RIBIERE,
new SimpleValueChecker(1.0e-10, 1.0e-10));
JDKRandomGenerator g = new JDKRandomGenerator();
g.setSeed(753289573253l);
RandomVectorGenerator generator =
new UncorrelatedRandomVectorGenerator(new double[] { 50.0, 50.0 }, new double[] { 10.0, 10.0 },
new GaussianRandomGenerator(g));
DifferentiableMultivariateMultiStartOptimizer optimizer =
new DifferentiableMultivariateMultiStartOptimizer(underlying, 10, generator);
PointValuePair optimum =
optimizer.optimize(200, circle, GoalType.MINIMIZE, new double[] { 98.680, 47.345 });
Assert.assertEquals(200, optimizer.getMaxEvaluations());
PointValuePair[] optima = optimizer.getOptima();
for (PointValuePair o : optima) {
Point2D.Double center = new Point2D.Double(o.getPointRef()[0], o.getPointRef()[1]);
Assert.assertEquals(69.960161753, circle.getRadius(center), 1.0e-8);
Assert.assertEquals(96.075902096, center.x, 1.0e-8);
Assert.assertEquals(48.135167894, center.y, 1.0e-8);
}
Assert.assertTrue(optimizer.getEvaluations() > 70);
Assert.assertTrue(optimizer.getEvaluations() < 90);
Assert.assertEquals(3.1267527, optimum.getValue(), 1.0e-8);
}
private static class Circle implements DifferentiableMultivariateFunction {
private ArrayList<Point2D.Double> points;
public Circle() {
points = new ArrayList<Point2D.Double>();
}
public void addPoint(double px, double py) {
points.add(new Point2D.Double(px, py));
}
public double getRadius(Point2D.Double center) {
double r = 0;
for (Point2D.Double point : points) {
r += point.distance(center);
}
return r / points.size();
}
private double[] gradient(double[] point) {
// optimal radius
Point2D.Double center = new Point2D.Double(point[0], point[1]);
double radius = getRadius(center);
// gradient of the sum of squared residuals
double dJdX = 0;
double dJdY = 0;
for (Point2D.Double pk : points) {
double dk = pk.distance(center);
dJdX += (center.x - pk.x) * (dk - radius) / dk;
dJdY += (center.y - pk.y) * (dk - radius) / dk;
}
dJdX *= 2;
dJdY *= 2;
return new double[] { dJdX, dJdY };
}
public double value(double[] variables) {
Point2D.Double center = new Point2D.Double(variables[0], variables[1]);
double radius = getRadius(center);
double sum = 0;
for (Point2D.Double point : points) {
double di = point.distance(center) - radius;
sum += di * di;
}
return sum;
}
public MultivariateVectorFunction gradient() {
return new MultivariateVectorFunction() {
public double[] value(double[] point) {
return gradient(point);
}
};
}
public MultivariateFunction partialDerivative(final int k) {
return new MultivariateFunction() {
public double value(double[] point) {
return gradient(point)[k];
}
};
}
}
}

View File

@ -0,0 +1,93 @@
/*
* Licensed to the Apache Software Foundation (ASF) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* The ASF licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.apache.commons.math3.optimization;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableFunction;
import org.apache.commons.math3.geometry.euclidean.twod.Vector2D;
import org.apache.commons.math3.optimization.general.CircleScalar;
import org.apache.commons.math3.optimization.general.ConjugateGradientFormula;
import org.apache.commons.math3.optimization.general.NonLinearConjugateGradientOptimizer;
import org.apache.commons.math3.random.GaussianRandomGenerator;
import org.apache.commons.math3.random.JDKRandomGenerator;
import org.apache.commons.math3.random.RandomVectorGenerator;
import org.apache.commons.math3.random.UncorrelatedRandomVectorGenerator;
import org.junit.Assert;
import org.junit.Test;
public class MultivariateDifferentiableMultiStartOptimizerTest {
@Test
public void testCircleFitting() {
CircleScalar circle = new CircleScalar();
circle.addPoint( 30.0, 68.0);
circle.addPoint( 50.0, -6.0);
circle.addPoint(110.0, -20.0);
circle.addPoint( 35.0, 15.0);
circle.addPoint( 45.0, 97.0);
// TODO: the wrapper around NonLinearConjugateGradientOptimizer is a temporary hack for
// version 3.1 of the library. It should be removed when NonLinearConjugateGradientOptimizer
// will officially be declared as implementing MultivariateDifferentiableOptimizer
MultivariateDifferentiableOptimizer underlying =
new MultivariateDifferentiableOptimizer() {
private final NonLinearConjugateGradientOptimizer cg =
new NonLinearConjugateGradientOptimizer(ConjugateGradientFormula.POLAK_RIBIERE,
new SimpleValueChecker(1.0e-10, 1.0e-10));
public PointValuePair optimize(int maxEval,
MultivariateDifferentiableFunction f,
GoalType goalType,
double[] startPoint) {
return cg.optimize(maxEval, f, goalType, startPoint);
}
public int getMaxEvaluations() {
return cg.getMaxEvaluations();
}
public int getEvaluations() {
return cg.getEvaluations();
}
public ConvergenceChecker<PointValuePair> getConvergenceChecker() {
return cg.getConvergenceChecker();
}
};
JDKRandomGenerator g = new JDKRandomGenerator();
g.setSeed(753289573253l);
RandomVectorGenerator generator =
new UncorrelatedRandomVectorGenerator(new double[] { 50.0, 50.0 }, new double[] { 10.0, 10.0 },
new GaussianRandomGenerator(g));
MultivariateDifferentiableMultiStartOptimizer optimizer =
new MultivariateDifferentiableMultiStartOptimizer(underlying, 10, generator);
PointValuePair optimum =
optimizer.optimize(200, circle, GoalType.MINIMIZE, new double[] { 98.680, 47.345 });
Assert.assertEquals(200, optimizer.getMaxEvaluations());
PointValuePair[] optima = optimizer.getOptima();
for (PointValuePair o : optima) {
Vector2D center = new Vector2D(o.getPointRef()[0], o.getPointRef()[1]);
Assert.assertEquals(69.960161753, circle.getRadius(center), 1.0e-8);
Assert.assertEquals(96.075902096, center.getX(), 1.0e-8);
Assert.assertEquals(48.135167894, center.getY(), 1.0e-8);
}
Assert.assertTrue(optimizer.getEvaluations() > 70);
Assert.assertTrue(optimizer.getEvaluations() < 90);
Assert.assertEquals(3.1267527, optimum.getValue(), 1.0e-8);
}
}

View File

@ -18,8 +18,8 @@
package org.apache.commons.math3.optimization;
import org.apache.commons.math3.analysis.DifferentiableMultivariateVectorFunction;
import org.apache.commons.math3.analysis.MultivariateMatrixFunction;
import org.apache.commons.math3.analysis.differentiation.DerivativeStructure;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableVectorFunction;
import org.apache.commons.math3.exception.MathIllegalStateException;
import org.apache.commons.math3.linear.BlockRealMatrix;
import org.apache.commons.math3.linear.RealMatrix;
@ -93,21 +93,47 @@ import org.junit.Test;
* @author Jorge J. More (original fortran minpack tests)
* @author Luc Maisonobe (non-minpack tests and minpack tests Java translation)
*/
public class DifferentiableMultivariateVectorMultiStartOptimizerTest {
public class MultivariateDifferentiableVectorMultiStartOptimizerTest {
@Test
public void testTrivial() {
LinearProblem problem =
new LinearProblem(new double[][] { { 2 } }, new double[] { 3 });
DifferentiableMultivariateVectorOptimizer underlyingOptimizer =
new GaussNewtonOptimizer(true,
new SimpleVectorValueChecker(1.0e-6, 1.0e-6));
// TODO: the wrapper around GaussNewtonOptimizer is a temporary hack for
// version 3.1 of the library. It should be removed when GaussNewtonOptimizer
// will officialy be declared as implementing MultivariateDifferentiableVectorOptimizer
MultivariateDifferentiableVectorOptimizer underlyingOptimizer =
new MultivariateDifferentiableVectorOptimizer() {
private GaussNewtonOptimizer gn =
new GaussNewtonOptimizer(true,
new SimpleVectorValueChecker(1.0e-6, 1.0e-6));
public PointVectorValuePair optimize(int maxEval,
MultivariateDifferentiableVectorFunction f,
double[] target,
double[] weight,
double[] startPoint) {
return gn.optimize(maxEval, f, target, weight, startPoint);
}
public int getMaxEvaluations() {
return gn.getMaxEvaluations();
}
public int getEvaluations() {
return gn.getEvaluations();
}
public ConvergenceChecker<PointVectorValuePair> getConvergenceChecker() {
return gn.getConvergenceChecker();
}
};
JDKRandomGenerator g = new JDKRandomGenerator();
g.setSeed(16069223052l);
RandomVectorGenerator generator =
new UncorrelatedRandomVectorGenerator(1, new GaussianRandomGenerator(g));
DifferentiableMultivariateVectorMultiStartOptimizer optimizer =
new DifferentiableMultivariateVectorMultiStartOptimizer(underlyingOptimizer,
MultivariateDifferentiableVectorMultiStartOptimizer optimizer =
new MultivariateDifferentiableVectorMultiStartOptimizer(underlyingOptimizer,
10, generator);
// no optima before first optimization attempt
@ -134,23 +160,50 @@ public class DifferentiableMultivariateVectorMultiStartOptimizerTest {
@Test(expected=TestException.class)
public void testNoOptimum() {
DifferentiableMultivariateVectorOptimizer underlyingOptimizer =
new GaussNewtonOptimizer(true,
new SimpleVectorValueChecker(1.0e-6, 1.0e-6));
// TODO: the wrapper around GaussNewtonOptimizer is a temporary hack for
// version 3.1 of the library. It should be removed when GaussNewtonOptimizer
// will officialy be declared as implementing MultivariateDifferentiableVectorOptimizer
MultivariateDifferentiableVectorOptimizer underlyingOptimizer =
new MultivariateDifferentiableVectorOptimizer() {
private GaussNewtonOptimizer gn =
new GaussNewtonOptimizer(true,
new SimpleVectorValueChecker(1.0e-6, 1.0e-6));
public PointVectorValuePair optimize(int maxEval,
MultivariateDifferentiableVectorFunction f,
double[] target,
double[] weight,
double[] startPoint) {
return gn.optimize(maxEval, f, target, weight, startPoint);
}
public int getMaxEvaluations() {
return gn.getMaxEvaluations();
}
public int getEvaluations() {
return gn.getEvaluations();
}
public ConvergenceChecker<PointVectorValuePair> getConvergenceChecker() {
return gn.getConvergenceChecker();
}
};
JDKRandomGenerator g = new JDKRandomGenerator();
g.setSeed(12373523445l);
RandomVectorGenerator generator =
new UncorrelatedRandomVectorGenerator(1, new GaussianRandomGenerator(g));
DifferentiableMultivariateVectorMultiStartOptimizer optimizer =
new DifferentiableMultivariateVectorMultiStartOptimizer(underlyingOptimizer,
MultivariateDifferentiableVectorMultiStartOptimizer optimizer =
new MultivariateDifferentiableVectorMultiStartOptimizer(underlyingOptimizer,
10, generator);
optimizer.optimize(100, new DifferentiableMultivariateVectorFunction() {
public MultivariateMatrixFunction jacobian() {
return null;
}
public double[] value(double[] point) {
throw new TestException();
}
optimizer.optimize(100, new MultivariateDifferentiableVectorFunction() {
public double[] value(double[] point) {
throw new TestException();
}
public DerivativeStructure[] value(DerivativeStructure[] point) {
return point;
}
}, new double[] { 2 }, new double[] { 1 }, new double[] { 0 });
}
@ -158,7 +211,7 @@ public class DifferentiableMultivariateVectorMultiStartOptimizerTest {
private static final long serialVersionUID = -7809988995389067683L;
}
private static class LinearProblem implements DifferentiableMultivariateVectorFunction {
private static class LinearProblem implements MultivariateDifferentiableVectorFunction {
final RealMatrix factors;
final double[] target;
@ -171,12 +224,15 @@ public class DifferentiableMultivariateVectorMultiStartOptimizerTest {
return factors.operate(variables);
}
public MultivariateMatrixFunction jacobian() {
return new MultivariateMatrixFunction() {
public double[][] value(double[] point) {
return factors.getData();
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure[] y = new DerivativeStructure[factors.getRowDimension()];
for (int i = 0; i < y.length; ++i) {
y[i] = variables[0].getField().getZero();
for (int j = 0; j < factors.getColumnDimension(); ++j) {
y[i] = y[i].add(variables[j].multiply(factors.getEntry(i, j)));
}
};
}
return y;
}
}

View File

@ -16,16 +16,16 @@
*/
package org.apache.commons.math3.optimization.general;
import java.awt.geom.Point2D;
import java.io.IOException;
import java.io.Serializable;
import java.util.Arrays;
import org.apache.commons.math3.analysis.DifferentiableMultivariateVectorFunction;
import org.apache.commons.math3.analysis.MultivariateMatrixFunction;
import org.apache.commons.math3.analysis.differentiation.DerivativeStructure;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableVectorFunction;
import org.apache.commons.math3.exception.ConvergenceException;
import org.apache.commons.math3.exception.DimensionMismatchException;
import org.apache.commons.math3.exception.NumberIsTooSmallException;
import org.apache.commons.math3.geometry.euclidean.twod.Vector2D;
import org.apache.commons.math3.linear.BlockRealMatrix;
import org.apache.commons.math3.linear.RealMatrix;
import org.apache.commons.math3.optimization.PointVectorValuePair;
@ -360,10 +360,10 @@ public abstract class AbstractLeastSquaresOptimizerAbstractTest {
Assert.assertTrue(optimizer.getJacobianEvaluations() < 10);
double rms = optimizer.getRMS();
Assert.assertEquals(1.768262623567235, FastMath.sqrt(circle.getN()) * rms, 1.0e-10);
Point2D.Double center = new Point2D.Double(optimum.getPointRef()[0], optimum.getPointRef()[1]);
Vector2D center = new Vector2D(optimum.getPointRef()[0], optimum.getPointRef()[1]);
Assert.assertEquals(69.96016176931406, circle.getRadius(center), 1.0e-6);
Assert.assertEquals(96.07590211815305, center.x, 1.0e-6);
Assert.assertEquals(48.13516790438953, center.y, 1.0e-6);
Assert.assertEquals(96.07590211815305, center.getX(), 1.0e-6);
Assert.assertEquals(48.13516790438953, center.getY(), 1.0e-6);
double[][] cov = optimizer.getCovariances();
Assert.assertEquals(1.839, cov[0][0], 0.001);
Assert.assertEquals(0.731, cov[0][1], 0.001);
@ -376,7 +376,7 @@ public abstract class AbstractLeastSquaresOptimizerAbstractTest {
// add perfect measurements and check errors are reduced
double r = circle.getRadius(center);
for (double d= 0; d < 2 * FastMath.PI; d += 0.01) {
circle.addPoint(center.x + r * FastMath.cos(d), center.y + r * FastMath.sin(d));
circle.addPoint(center.getX() + r * FastMath.cos(d), center.getY() + r * FastMath.sin(d));
}
double[] target = new double[circle.getN()];
Arrays.fill(target, 0.0);
@ -407,13 +407,13 @@ public abstract class AbstractLeastSquaresOptimizerAbstractTest {
AbstractLeastSquaresOptimizer optimizer = createOptimizer();
PointVectorValuePair optimum =
optimizer.optimize(100, circle, target, weights, new double[] { -12, -12 });
Point2D.Double center = new Point2D.Double(optimum.getPointRef()[0], optimum.getPointRef()[1]);
Vector2D center = new Vector2D(optimum.getPointRef()[0], optimum.getPointRef()[1]);
Assert.assertTrue(optimizer.getEvaluations() < 25);
Assert.assertTrue(optimizer.getJacobianEvaluations() < 20);
Assert.assertEquals( 0.043, optimizer.getRMS(), 1.0e-3);
Assert.assertEquals( 0.292235, circle.getRadius(center), 1.0e-6);
Assert.assertEquals(-0.151738, center.x, 1.0e-6);
Assert.assertEquals( 0.2075001, center.y, 1.0e-6);
Assert.assertEquals(-0.151738, center.getX(), 1.0e-6);
Assert.assertEquals( 0.2075001, center.getY(), 1.0e-6);
}
@Test
@ -475,7 +475,7 @@ public abstract class AbstractLeastSquaresOptimizerAbstractTest {
final double[][] data = dataset.getData();
final double[] initial = dataset.getStartingPoint(0);
final DifferentiableMultivariateVectorFunction problem;
final MultivariateDifferentiableVectorFunction problem;
problem = dataset.getLeastSquaresProblem();
final PointVectorValuePair optimum;
optimum = optimizer.optimize(100, problem, data[1], w, initial);
@ -504,7 +504,7 @@ public abstract class AbstractLeastSquaresOptimizerAbstractTest {
doTestStRD(StatisticalReferenceDatasetFactory.createHahn1(), 1E-7, 1E-4);
}
static class LinearProblem implements DifferentiableMultivariateVectorFunction, Serializable {
static class LinearProblem implements MultivariateDifferentiableVectorFunction, Serializable {
private static final long serialVersionUID = 703247177355019415L;
final RealMatrix factors;
@ -518,12 +518,17 @@ public abstract class AbstractLeastSquaresOptimizerAbstractTest {
return factors.operate(variables);
}
public MultivariateMatrixFunction jacobian() {
return new MultivariateMatrixFunction() {
public double[][] value(double[] point) {
return factors.getData();
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure[] value = new DerivativeStructure[factors.getRowDimension()];
for (int i = 0; i < value.length; ++i) {
value[i] = variables[0].getField().getZero();
for (int j = 0; j < factors.getColumnDimension(); ++j) {
value[i] = value[i].add(variables[j].multiply(factors.getEntry(i, j)));
}
};
}
return value;
}
}
}

View File

@ -92,7 +92,7 @@ public class AbstractLeastSquaresOptimizerTest {
for (int i = 0; i < sig.length; i++) {
final double actual = FastMath.sqrt(optimizer.getChiSquare()/dof)*sig[i];
Assert.assertEquals(dataset.getName() + ", parameter #" + i,
actual, expected[i], 1E-8 * expected[i]);
expected[i], actual, 1.3e-8 * expected[i]);
}
}
}

View File

@ -18,9 +18,10 @@
package org.apache.commons.math3.optimization.general;
import java.util.ArrayList;
import org.apache.commons.math3.analysis.DifferentiableMultivariateVectorFunction;
import org.apache.commons.math3.analysis.MultivariateMatrixFunction;
import org.apache.commons.math3.util.MathUtils;
import org.apache.commons.math3.analysis.differentiation.DerivativeStructure;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableVectorFunction;
import org.apache.commons.math3.geometry.euclidean.twod.Vector2D;
import org.apache.commons.math3.util.FastMath;
/**
@ -37,31 +38,13 @@ import org.apache.commons.math3.util.FastMath;
* corresponding circle.</li>
* </ul>
*/
class CircleProblem implements DifferentiableMultivariateVectorFunction {
class CircleProblem implements MultivariateDifferentiableVectorFunction {
/** Cloud of points assumed to be fitted by a circle. */
private final ArrayList<double[]> points;
private final ArrayList<Vector2D> points;
/** Error on the x-coordinate of the points. */
private final double xSigma;
/** Error on the y-coordinate of the points. */
private final double ySigma;
/** Number of points on the circumference (when searching which
model point is closest to a given "observation". */
private final int resolution;
/**
* @param xError Assumed error for the x-coordinate of the circle points.
* @param yError Assumed error for the y-coordinate of the circle points.
* @param searchResolution Number of points to try when searching the one
* that is closest to a given "observed" point.
*/
public CircleProblem(double xError,
double yError,
int searchResolution) {
points = new ArrayList<double[]>();
xSigma = xError;
ySigma = yError;
resolution = searchResolution;
}
/**
* @param xError Assumed error for the x-coordinate of the circle points.
@ -69,20 +52,22 @@ class CircleProblem implements DifferentiableMultivariateVectorFunction {
*/
public CircleProblem(double xError,
double yError) {
this(xError, yError, 500);
points = new ArrayList<Vector2D>();
xSigma = xError;
ySigma = yError;
}
public void addPoint(double px, double py) {
points.add(new double[] { px, py });
public void addPoint(Vector2D p) {
points.add(p);
}
public double[] target() {
final double[] t = new double[points.size() * 2];
for (int i = 0; i < points.size(); i++) {
final double[] p = points.get(i);
final Vector2D p = points.get(i);
final int index = i * 2;
t[index] = p[0];
t[index + 1] = p[1];
t[index] = p.getX();
t[index + 1] = p.getY();
}
return t;
@ -108,65 +93,46 @@ class CircleProblem implements DifferentiableMultivariateVectorFunction {
final double[] model = new double[points.size() * 2];
final double deltaTheta = MathUtils.TWO_PI / resolution;
for (int i = 0; i < points.size(); i++) {
final double[] p = points.get(i);
final double px = p[0];
final double py = p[1];
final Vector2D p = points.get(i);
double bestX = 0;
double bestY = 0;
double dMin = Double.POSITIVE_INFINITY;
// Find the circle point closest to the observed point
// (observed points are points add through the addPoint method above)
final double dX = cx - p.getX();
final double dY = cy - p.getY();
final double scaling = r / FastMath.hypot(dX, dY);
final int index = i * 2;
model[index] = cx - scaling * dX;
model[index + 1] = cy - scaling * dY;
// Find the angle for which the circle passes closest to the
// current point (using a resolution of 100 points along the
// circumference).
for (double theta = 0; theta <= MathUtils.TWO_PI; theta += deltaTheta) {
final double currentX = cx + r * FastMath.cos(theta);
final double currentY = cy + r * FastMath.sin(theta);
final double dX = currentX - px;
final double dY = currentY - py;
final double d = dX * dX + dY * dY;
if (d < dMin) {
dMin = d;
bestX = currentX;
bestY = currentY;
}
}
final int index = i * 2;
model[index] = bestX;
model[index + 1] = bestY;
}
return model;
}
public MultivariateMatrixFunction jacobian() {
return new MultivariateMatrixFunction() {
public double[][] value(double[] point) {
return jacobian(point);
}
};
}
public DerivativeStructure[] value(DerivativeStructure[] params) {
final DerivativeStructure cx = params[0];
final DerivativeStructure cy = params[1];
final DerivativeStructure r = params[2];
private double[][] jacobian(double[] params) {
final double[][] jacobian = new double[points.size() * 2][3];
final DerivativeStructure[] model = new DerivativeStructure[points.size() * 2];
for (int i = 0; i < points.size(); i++) {
final int index = i * 2;
// Partial derivative wrt x-coordinate of center.
jacobian[index][0] = 1;
jacobian[index + 1][0] = 0;
// Partial derivative wrt y-coordinate of center.
jacobian[index][1] = 0;
jacobian[index + 1][1] = 1;
// Partial derivative wrt radius.
final double[] p = points.get(i);
jacobian[index][2] = (p[0] - params[0]) / params[2];
jacobian[index + 1][2] = (p[1] - params[1]) / params[2];
final Vector2D p = points.get(i);
// Find the circle point closest to the observed point
// (observed points are points add through the addPoint method above)
final DerivativeStructure dX = cx.subtract(p.getX());
final DerivativeStructure dY = cy.subtract(p.getY());
final DerivativeStructure scaling = r.divide(dX.multiply(dX).add(dY.multiply(dY)).sqrt());
final int index = i * 2;
model[index] = cx.subtract(scaling.multiply(dX));
model[index + 1] = cy.subtract(scaling.multiply(dY));
}
return jacobian;
return model;
}
}

View File

@ -17,59 +17,55 @@
package org.apache.commons.math3.optimization.general;
import java.awt.geom.Point2D;
import java.util.ArrayList;
import org.apache.commons.math3.analysis.DifferentiableMultivariateFunction;
import org.apache.commons.math3.analysis.MultivariateFunction;
import org.apache.commons.math3.analysis.MultivariateVectorFunction;
import org.apache.commons.math3.analysis.differentiation.DerivativeStructure;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableFunction;
import org.apache.commons.math3.geometry.euclidean.twod.Vector2D;
/**
* Class used in the tests.
*/
class CircleScalar implements DifferentiableMultivariateFunction {
private ArrayList<Point2D.Double> points;
public class CircleScalar implements MultivariateDifferentiableFunction {
private ArrayList<Vector2D> points;
public CircleScalar() {
points = new ArrayList<Point2D.Double>();
points = new ArrayList<Vector2D>();
}
public void addPoint(double px, double py) {
points.add(new Point2D.Double(px, py));
points.add(new Vector2D(px, py));
}
public double getRadius(Point2D.Double center) {
public double getRadius(Vector2D center) {
double r = 0;
for (Point2D.Double point : points) {
for (Vector2D point : points) {
r += point.distance(center);
}
return r / points.size();
}
private double[] gradient(double[] point) {
// optimal radius
Point2D.Double center = new Point2D.Double(point[0], point[1]);
double radius = getRadius(center);
private DerivativeStructure distance(Vector2D point,
DerivativeStructure cx, DerivativeStructure cy) {
DerivativeStructure dx = cx.subtract(point.getX());
DerivativeStructure dy = cy.subtract(point.getY());
return dx.multiply(dx).add(dy.multiply(dy)).sqrt();
}
// gradient of the sum of squared residuals
double dJdX = 0;
double dJdY = 0;
for (Point2D.Double pk : points) {
double dk = pk.distance(center);
dJdX += (center.x - pk.x) * (dk - radius) / dk;
dJdY += (center.y - pk.y) * (dk - radius) / dk;
public DerivativeStructure getRadius(DerivativeStructure cx, DerivativeStructure cy) {
DerivativeStructure r = cx.getField().getZero();
for (Vector2D point : points) {
r = r.add(distance(point, cx, cy));
}
dJdX *= 2;
dJdY *= 2;
return new double[] { dJdX, dJdY };
return r.divide(points.size());
}
public double value(double[] variables) {
Point2D.Double center = new Point2D.Double(variables[0], variables[1]);
Vector2D center = new Vector2D(variables[0], variables[1]);
double radius = getRadius(center);
double sum = 0;
for (Point2D.Double point : points) {
for (Vector2D point : points) {
double di = point.distance(center) - radius;
sum += di * di;
}
@ -77,19 +73,16 @@ class CircleScalar implements DifferentiableMultivariateFunction {
return sum;
}
public MultivariateVectorFunction gradient() {
return new MultivariateVectorFunction() {
public double[] value(double[] point) {
return gradient(point);
}
};
public DerivativeStructure value(DerivativeStructure[] variables) {
DerivativeStructure radius = getRadius(variables[0], variables[1]);
DerivativeStructure sum = variables[0].getField().getZero();
for (Vector2D point : points) {
DerivativeStructure di = distance(point, variables[0], variables[1]).subtract(radius);
sum = sum.add(di.multiply(di));
}
return sum;
}
public MultivariateFunction partialDerivative(final int k) {
return new MultivariateFunction() {
public double value(double[] point) {
return gradient(point)[k];
}
};
}
}

View File

@ -17,66 +17,55 @@
package org.apache.commons.math3.optimization.general;
import java.awt.geom.Point2D;
import java.util.ArrayList;
import org.apache.commons.math3.analysis.DifferentiableMultivariateVectorFunction;
import org.apache.commons.math3.analysis.MultivariateMatrixFunction;
import org.apache.commons.math3.analysis.differentiation.DerivativeStructure;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableVectorFunction;
import org.apache.commons.math3.geometry.euclidean.twod.Vector2D;
/**
* Class used in the tests.
*/
class CircleVectorial implements DifferentiableMultivariateVectorFunction {
private ArrayList<Point2D.Double> points;
class CircleVectorial implements MultivariateDifferentiableVectorFunction {
private ArrayList<Vector2D> points;
public CircleVectorial() {
points = new ArrayList<Point2D.Double>();
points = new ArrayList<Vector2D>();
}
public void addPoint(double px, double py) {
points.add(new Point2D.Double(px, py));
points.add(new Vector2D(px, py));
}
public int getN() {
return points.size();
}
public double getRadius(Point2D.Double center) {
public double getRadius(Vector2D center) {
double r = 0;
for (Point2D.Double point : points) {
for (Vector2D point : points) {
r += point.distance(center);
}
return r / points.size();
}
private double[][] jacobian(double[] point) {
int n = points.size();
Point2D.Double center = new Point2D.Double(point[0], point[1]);
private DerivativeStructure distance(Vector2D point,
DerivativeStructure cx, DerivativeStructure cy) {
DerivativeStructure dx = cx.subtract(point.getX());
DerivativeStructure dy = cy.subtract(point.getY());
return dx.multiply(dx).add(dy.multiply(dy)).sqrt();
}
// gradient of the optimal radius
double dRdX = 0;
double dRdY = 0;
for (Point2D.Double pk : points) {
double dk = pk.distance(center);
dRdX += (center.x - pk.x) / dk;
dRdY += (center.y - pk.y) / dk;
public DerivativeStructure getRadius(DerivativeStructure cx, DerivativeStructure cy) {
DerivativeStructure r = cx.getField().getZero();
for (Vector2D point : points) {
r = r.add(distance(point, cx, cy));
}
dRdX /= n;
dRdY /= n;
// jacobian of the radius residuals
double[][] jacobian = new double[n][2];
for (int i = 0; i < n; ++i) {
Point2D.Double pi = points.get(i);
double di = pi.distance(center);
jacobian[i][0] = (center.x - pi.x) / di - dRdX;
jacobian[i][1] = (center.y - pi.y) / di - dRdY;
}
return jacobian;
return r.divide(points.size());
}
public double[] value(double[] variables) {
Point2D.Double center = new Point2D.Double(variables[0], variables[1]);
Vector2D center = new Vector2D(variables[0], variables[1]);
double radius = getRadius(center);
double[] residuals = new double[points.size()];
@ -87,11 +76,15 @@ class CircleVectorial implements DifferentiableMultivariateVectorFunction {
return residuals;
}
public MultivariateMatrixFunction jacobian() {
return new MultivariateMatrixFunction() {
public double[][] value(double[] point) {
return jacobian(point);
}
};
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure radius = getRadius(variables[0], variables[1]);
DerivativeStructure[] residuals = new DerivativeStructure[points.size()];
for (int i = 0; i < residuals.length; ++i) {
residuals[i] = distance(points.get(i), variables[0], variables[1]).subtract(radius);
}
return residuals;
}
}

View File

@ -17,16 +17,16 @@
package org.apache.commons.math3.optimization.general;
import java.awt.geom.Point2D;
import java.io.Serializable;
import java.util.ArrayList;
import java.util.List;
import org.apache.commons.math3.analysis.DifferentiableMultivariateVectorFunction;
import org.apache.commons.math3.analysis.MultivariateMatrixFunction;
import org.apache.commons.math3.analysis.differentiation.DerivativeStructure;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableVectorFunction;
import org.apache.commons.math3.exception.ConvergenceException;
import org.apache.commons.math3.exception.DimensionMismatchException;
import org.apache.commons.math3.exception.TooManyEvaluationsException;
import org.apache.commons.math3.geometry.euclidean.twod.Vector2D;
import org.apache.commons.math3.linear.SingularMatrixException;
import org.apache.commons.math3.optimization.PointVectorValuePair;
import org.apache.commons.math3.util.FastMath;
@ -120,7 +120,7 @@ public class LevenbergMarquardtOptimizerTest extends AbstractLeastSquaresOptimiz
optimizer.optimize(100, problem, problem.target, new double[] { 1, 1, 1 }, new double[] { 0, 0, 0 });
Assert.assertTrue(FastMath.sqrt(problem.target.length) * optimizer.getRMS() > 0.6);
double[][] cov = optimizer.getCovariances(1.5e-14);
optimizer.getCovariances(1.5e-14);
}
@Test
@ -138,7 +138,7 @@ public class LevenbergMarquardtOptimizerTest extends AbstractLeastSquaresOptimiz
checkEstimate(circle, 0.1, 20, 1.0e-18, 1.0e-16, 1.0e-10, true);
}
private void checkEstimate(DifferentiableMultivariateVectorFunction problem,
private void checkEstimate(MultivariateDifferentiableVectorFunction problem,
double initialStepBoundFactor, int maxCostEval,
double costRelativeTolerance, double parRelativeTolerance,
double orthoTolerance, boolean shouldFail) {
@ -225,7 +225,6 @@ public class LevenbergMarquardtOptimizerTest extends AbstractLeastSquaresOptimiz
optimizer.optimize(100, problem, dataPoints[1], weights,
new double[] { 10, 900, 80, 27, 225 });
final double chi2 = optimizer.getChiSquare();
final double[] solution = optimum.getPoint();
final double[] expectedSolution = { 10.4, 958.3, 131.4, 33.9, 205.0 };
@ -274,8 +273,8 @@ public class LevenbergMarquardtOptimizerTest extends AbstractLeastSquaresOptimiz
final CircleProblem circle = new CircleProblem(xSigma, ySigma);
final int numPoints = 10;
for (Point2D.Double p : factory.generate(numPoints)) {
circle.addPoint(p.x, p.y);
for (Vector2D p : factory.generate(numPoints)) {
circle.addPoint(p);
// System.out.println(p.x + " " + p.y);
}
@ -309,7 +308,7 @@ public class LevenbergMarquardtOptimizerTest extends AbstractLeastSquaresOptimiz
Assert.assertEquals(radius, paramFound[2], asymptoticStandardErrorFound[2]);
}
private static class QuadraticProblem implements DifferentiableMultivariateVectorFunction, Serializable {
private static class QuadraticProblem implements MultivariateDifferentiableVectorFunction, Serializable {
private static final long serialVersionUID = 7072187082052755854L;
private List<Double> x;
@ -325,16 +324,6 @@ public class LevenbergMarquardtOptimizerTest extends AbstractLeastSquaresOptimiz
this.y.add(y);
}
private double[][] jacobian(double[] variables) {
double[][] jacobian = new double[x.size()][3];
for (int i = 0; i < jacobian.length; ++i) {
jacobian[i][0] = x.get(i) * x.get(i);
jacobian[i][1] = x.get(i);
jacobian[i][2] = 1.0;
}
return jacobian;
}
public double[] value(double[] variables) {
double[] values = new double[x.size()];
for (int i = 0; i < values.length; ++i) {
@ -343,17 +332,18 @@ public class LevenbergMarquardtOptimizerTest extends AbstractLeastSquaresOptimiz
return values;
}
public MultivariateMatrixFunction jacobian() {
return new MultivariateMatrixFunction() {
public double[][] value(double[] point) {
return jacobian(point);
}
};
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure[] values = new DerivativeStructure[x.size()];
for (int i = 0; i < values.length; ++i) {
values[i] = (variables[0].multiply(x.get(i)).add(variables[1])).multiply(x.get(i)).add(variables[2]);
}
return values;
}
}
private static class BevingtonProblem
implements DifferentiableMultivariateVectorFunction {
implements MultivariateDifferentiableVectorFunction {
private List<Double> time;
private List<Double> count;
@ -367,25 +357,6 @@ public class LevenbergMarquardtOptimizerTest extends AbstractLeastSquaresOptimiz
count.add(c);
}
private double[][] jacobian(double[] params) {
double[][] jacobian = new double[time.size()][5];
for (int i = 0; i < jacobian.length; ++i) {
final double t = time.get(i);
jacobian[i][0] = 1;
final double p3 = params[3];
final double p4 = params[4];
final double tOp3 = t / p3;
final double tOp4 = t / p4;
jacobian[i][1] = Math.exp(-tOp3);
jacobian[i][2] = Math.exp(-tOp4);
jacobian[i][3] = params[1] * Math.exp(-tOp3) * tOp3 / p3;
jacobian[i][4] = params[2] * Math.exp(-tOp4) * tOp4 / p4;
}
return jacobian;
}
public double[] value(double[] params) {
double[] values = new double[time.size()];
for (int i = 0; i < values.length; ++i) {
@ -397,12 +368,16 @@ public class LevenbergMarquardtOptimizerTest extends AbstractLeastSquaresOptimiz
return values;
}
public MultivariateMatrixFunction jacobian() {
return new MultivariateMatrixFunction() {
public double[][] value(double[] point) {
return jacobian(point);
}
};
public DerivativeStructure[] value(DerivativeStructure[] params) {
DerivativeStructure[] values = new DerivativeStructure[time.size()];
for (int i = 0; i < values.length; ++i) {
final double t = time.get(i);
values[i] = params[0].add(
params[1].multiply(params[3].reciprocal().multiply(-t).exp())).add(
params[2].multiply(params[4].reciprocal().multiply(-t).exp()));
}
return values;
}
}
}

View File

@ -22,8 +22,8 @@ import java.util.Arrays;
import org.apache.commons.math3.exception.TooManyEvaluationsException;
import org.apache.commons.math3.analysis.DifferentiableMultivariateVectorFunction;
import org.apache.commons.math3.analysis.MultivariateMatrixFunction;
import org.apache.commons.math3.analysis.differentiation.DerivativeStructure;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableVectorFunction;
import org.apache.commons.math3.optimization.PointVectorValuePair;
import org.apache.commons.math3.util.FastMath;
import org.junit.Assert;
@ -522,7 +522,7 @@ public class MinpackTest {
}
private static abstract class MinpackFunction
implements DifferentiableMultivariateVectorFunction, Serializable {
implements MultivariateDifferentiableVectorFunction, Serializable {
private static final long serialVersionUID = -6209760235478794233L;
protected int n;
@ -590,17 +590,20 @@ public class MinpackTest {
}
}
public MultivariateMatrixFunction jacobian() {
return new MultivariateMatrixFunction() {
public double[][] value(double[] point) {
return jacobian(point);
}
};
public double[] value(double[] variables) {
DerivativeStructure[] dsV = new DerivativeStructure[variables.length];
for (int i = 0; i < variables.length; ++i) {
dsV[i] = new DerivativeStructure(0, 0, variables[i]);
}
DerivativeStructure[] dsY = value(dsV);
double[] y = new double[dsY.length];
for (int i = 0; i < dsY.length; ++i) {
y[i] = dsY[i].getValue();
}
return y;
}
public abstract double[][] jacobian(double[] variables);
public abstract double[] value(double[] variables);
public abstract DerivativeStructure[] value(DerivativeStructure[] variables);
}
@ -616,30 +619,17 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
double t = 2.0 / m;
double[][] jacobian = new double[m][];
for (int i = 0; i < m; ++i) {
jacobian[i] = new double[n];
for (int j = 0; j < n; ++j) {
jacobian[i][j] = (i == j) ? (1 - t) : -t;
}
}
return jacobian;
}
@Override
public double[] value(double[] variables) {
double sum = 0;
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure sum = variables[0].getField().getZero();
for (int i = 0; i < n; ++i) {
sum += variables[i];
sum = sum.add(variables[i]);
}
double t = 1 + 2 * sum / m;
double[] f = new double[m];
DerivativeStructure t = sum.multiply(2.0 / m).add(1);
DerivativeStructure[] f = new DerivativeStructure[m];
for (int i = 0; i < n; ++i) {
f[i] = variables[i] - t;
f[i] = variables[i].subtract(t);
}
Arrays.fill(f, n, m, -t);
Arrays.fill(f, n, m, t.negate());
return f;
}
@ -656,28 +646,16 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
double[][] jacobian = new double[m][];
for (int i = 0; i < m; ++i) {
jacobian[i] = new double[n];
for (int j = 0; j < n; ++j) {
jacobian[i][j] = (i + 1) * (j + 1);
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure[] f = new DerivativeStructure[m];
DerivativeStructure sum = variables[0].getField().getZero();
for (int i = 0; i < n; ++i) {
sum = sum.add(variables[i].multiply(i + 1));
}
}
return jacobian;
}
@Override
public double[] value(double[] variables) {
double[] f = new double[m];
double sum = 0;
for (int i = 0; i < n; ++i) {
sum += (i + 1) * variables[i];
}
for (int i = 0; i < m; ++i) {
f[i] = (i + 1) * sum - 1;
}
return f;
for (int i = 0; i < m; ++i) {
f[i] = sum.multiply(i + 1).subtract(1);
}
return f;
}
}
@ -693,36 +671,16 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
double[][] jacobian = new double[m][];
for (int i = 0; i < m; ++i) {
jacobian[i] = new double[n];
jacobian[i][0] = 0;
for (int j = 1; j < (n - 1); ++j) {
if (i == 0) {
jacobian[i][j] = 0;
} else if (i != (m - 1)) {
jacobian[i][j] = i * (j + 1);
} else {
jacobian[i][j] = 0;
}
}
jacobian[i][n - 1] = 0;
}
return jacobian;
}
@Override
public double[] value(double[] variables) {
double[] f = new double[m];
double sum = 0;
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure[] f = new DerivativeStructure[m];
DerivativeStructure sum = variables[0].getField().getZero();
for (int i = 1; i < (n - 1); ++i) {
sum += (i + 1) * variables[i];
sum = sum.add(variables[i].multiply(i + 1));
}
for (int i = 0; i < (m - 1); ++i) {
f[i] = i * sum - 1;
f[i] = sum.multiply(i).subtract(1);
}
f[m - 1] = -1;
f[m - 1] = variables[0].getField().getOne().negate();
return f;
}
@ -737,16 +695,13 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
double x1 = variables[0];
return new double[][] { { -20 * x1, 10 }, { -1, 0 } };
}
@Override
public double[] value(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
return new double[] { 10 * (x2 - x1 * x1), 1 - x1 };
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure x1 = variables[0];
DerivativeStructure x2 = variables[1];
return new DerivativeStructure[] {
x2.subtract(x1.multiply(x1)).multiply(10),
x1.negate().add(1)
};
}
}
@ -761,39 +716,25 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
double tmpSquare = x1 * x1 + x2 * x2;
double tmp1 = twoPi * tmpSquare;
double tmp2 = FastMath.sqrt(tmpSquare);
return new double[][] {
{ 100 * x2 / tmp1, -100 * x1 / tmp1, 10 },
{ 10 * x1 / tmp2, 10 * x2 / tmp2, 0 },
{ 0, 0, 1 }
};
}
@Override
public double[] value(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
double x3 = variables[2];
double tmp1;
if (x1 == 0) {
tmp1 = (x2 >= 0) ? 0.25 : -0.25;
} else {
tmp1 = FastMath.atan(x2 / x1) / twoPi;
if (x1 < 0) {
tmp1 += 0.5;
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure x1 = variables[0];
DerivativeStructure x2 = variables[1];
DerivativeStructure x3 = variables[2];
DerivativeStructure tmp1 = variables[0].getField().getZero();
if (x1.getValue() == 0) {
tmp1 = tmp1.add((x2.getValue() >= 0) ? 0.25 : -0.25);
} else {
tmp1 = x2.divide(x1).atan().divide(twoPi);
if (x1.getValue() < 0) {
tmp1 = tmp1.add(0.5);
}
}
}
double tmp2 = FastMath.sqrt(x1 * x1 + x2 * x2);
return new double[] {
10.0 * (x3 - 10 * tmp1),
10.0 * (tmp2 - 1),
x3
};
DerivativeStructure tmp2 = x1.multiply(x1).add(x2.multiply(x2)).sqrt();
return new DerivativeStructure[] {
x3.subtract(tmp1.multiply(10)).multiply(10),
tmp2.subtract(1).multiply(10),
x3
};
}
private static final double twoPi = 2.0 * FastMath.PI;
@ -810,30 +751,16 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
double x3 = variables[2];
double x4 = variables[3];
return new double[][] {
{ 1, 10, 0, 0 },
{ 0, 0, sqrt5, -sqrt5 },
{ 0, 2 * (x2 - 2 * x3), -4 * (x2 - 2 * x3), 0 },
{ 2 * sqrt10 * (x1 - x4), 0, 0, -2 * sqrt10 * (x1 - x4) }
};
}
@Override
public double[] value(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
double x3 = variables[2];
double x4 = variables[3];
return new double[] {
x1 + 10 * x2,
sqrt5 * (x3 - x4),
(x2 - 2 * x3) * (x2 - 2 * x3),
sqrt10 * (x1 - x4) * (x1 - x4)
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure x1 = variables[0];
DerivativeStructure x2 = variables[1];
DerivativeStructure x3 = variables[2];
DerivativeStructure x4 = variables[3];
return new DerivativeStructure[] {
x1.add(x2.multiply(10)),
x3.subtract(x4).multiply(sqrt5),
x2.subtract(x3.multiply(2)).multiply(x2.subtract(x3.multiply(2))),
x1.subtract(x4).multiply(x1.subtract(x4)).multiply(sqrt10)
};
}
@ -855,22 +782,13 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
double x2 = variables[1];
return new double[][] {
{ 1, x2 * (10 - 3 * x2) - 2 },
{ 1, x2 * ( 2 + 3 * x2) - 14, }
};
}
@Override
public double[] value(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
return new double[] {
-13.0 + x1 + ((5.0 - x2) * x2 - 2.0) * x2,
-29.0 + x1 + ((1.0 + x2) * x2 - 14.0) * x2
};
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure x1 = variables[0];
DerivativeStructure x2 = variables[1];
return new DerivativeStructure[] {
x1.subtract(13.0).add(x2.negate().add(5.0).multiply(x2).subtract(2).multiply(x2)),
x1.subtract(29.0).add(x2.add(1).multiply(x2).subtract(14).multiply(x2))
};
}
}
@ -888,32 +806,16 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
double x2 = variables[1];
double x3 = variables[2];
double[][] jacobian = new double[m][];
for (int i = 0; i < m; ++i) {
double tmp1 = i + 1;
double tmp2 = 15 - i;
double tmp3 = (i <= 7) ? tmp1 : tmp2;
double tmp4 = x2 * tmp2 + x3 * tmp3;
tmp4 *= tmp4;
jacobian[i] = new double[] { -1, tmp1 * tmp2 / tmp4, tmp1 * tmp3 / tmp4 };
}
return jacobian;
}
@Override
public double[] value(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
double x3 = variables[2];
double[] f = new double[m];
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure x1 = variables[0];
DerivativeStructure x2 = variables[1];
DerivativeStructure x3 = variables[2];
DerivativeStructure[] f = new DerivativeStructure[m];
for (int i = 0; i < m; ++i) {
double tmp1 = i + 1;
double tmp2 = 15 - i;
double tmp3 = (i <= 7) ? tmp1 : tmp2;
f[i] = y[i] - (x1 + tmp1 / (x2 * tmp2 + x3 * tmp3));
f[i] = x1.add(x2.multiply(tmp2).add(x3.multiply(tmp3)).reciprocal().multiply(tmp1)).negate().add(y[i]);
}
return f;
}
@ -943,34 +845,16 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
double x3 = variables[2];
double x4 = variables[3];
double[][] jacobian = new double[m][];
for (int i = 0; i < m; ++i) {
double tmp = v[i] * (v[i] + x3) + x4;
double j1 = -v[i] * (v[i] + x2) / tmp;
double j2 = -v[i] * x1 / tmp;
double j3 = j1 * j2;
double j4 = j3 / v[i];
jacobian[i] = new double[] { j1, j2, j3, j4 };
}
return jacobian;
}
@Override
public double[] value(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
double x3 = variables[2];
double x4 = variables[3];
double[] f = new double[m];
for (int i = 0; i < m; ++i) {
f[i] = y[i] - x1 * (v[i] * (v[i] + x2)) / (v[i] * (v[i] + x3) + x4);
}
return f;
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure x1 = variables[0];
DerivativeStructure x2 = variables[1];
DerivativeStructure x3 = variables[2];
DerivativeStructure x4 = variables[3];
DerivativeStructure[] f = new DerivativeStructure[m];
for (int i = 0; i < m; ++i) {
f[i] = x1.multiply(x2.add(v[i]).multiply(v[i])).divide(x4.add(x3.add(v[i]).multiply(v[i]))).negate().add(y[i]);
}
return f;
}
private static final double[] v = {
@ -1001,29 +885,13 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
double x3 = variables[2];
double[][] jacobian = new double[m][];
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure x1 = variables[0];
DerivativeStructure x2 = variables[1];
DerivativeStructure x3 = variables[2];
DerivativeStructure[] f = new DerivativeStructure[m];
for (int i = 0; i < m; ++i) {
double temp = 5.0 * (i + 1) + 45.0 + x3;
double tmp1 = x2 / temp;
double tmp2 = FastMath.exp(tmp1);
double tmp3 = x1 * tmp2 / temp;
jacobian[i] = new double[] { tmp2, tmp3, -tmp1 * tmp3 };
}
return jacobian;
}
@Override
public double[] value(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
double x3 = variables[2];
double[] f = new double[m];
for (int i = 0; i < m; ++i) {
f[i] = x1 * FastMath.exp(x2 / (5.0 * (i + 1) + 45.0 + x3)) - y[i];
f[i] = x1.multiply(x2.divide(x3.add(5.0 * (i + 1) + 45.0)).exp()).subtract(y[i]);
}
return f;
}
@ -1050,64 +918,31 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
double[][] jacobian = new double[m][];
for (int i = 0; i < (m - 2); ++i) {
double div = (i + 1) / 29.0;
double s2 = 0.0;
double dx = 1.0;
for (int j = 0; j < n; ++j) {
s2 += dx * variables[j];
dx *= div;
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure[] f = new DerivativeStructure[m];
for (int i = 0; i < (m - 2); ++i) {
double div = (i + 1) / 29.0;
DerivativeStructure s1 = variables[0].getField().getZero();
DerivativeStructure dx = variables[0].getField().getOne();
for (int j = 1; j < n; ++j) {
s1 = s1.add(dx.multiply(j).multiply(variables[j]));
dx = dx.multiply(div);
}
DerivativeStructure s2 = variables[0].getField().getZero();
dx = variables[0].getField().getOne();
for (int j = 0; j < n; ++j) {
s2 = s2.add(dx.multiply(variables[j]));
dx = dx.multiply(div);
}
f[i] = s1.subtract(s2.multiply(s2)).subtract(1);
}
double temp= 2 * div * s2;
dx = 1.0 / div;
jacobian[i] = new double[n];
for (int j = 0; j < n; ++j) {
jacobian[i][j] = dx * (j - temp);
dx *= div;
}
}
jacobian[m - 2] = new double[n];
jacobian[m - 2][0] = 1;
DerivativeStructure x1 = variables[0];
DerivativeStructure x2 = variables[1];
f[m - 2] = x1;
f[m - 1] = x2.subtract(x1.multiply(x1)).subtract(1);
jacobian[m - 1] = new double[n];
jacobian[m - 1][0]= -2 * variables[0];
jacobian[m - 1][1]= 1;
return jacobian;
}
@Override
public double[] value(double[] variables) {
double[] f = new double[m];
for (int i = 0; i < (m - 2); ++i) {
double div = (i + 1) / 29.0;
double s1 = 0;
double dx = 1;
for (int j = 1; j < n; ++j) {
s1 += j * dx * variables[j];
dx *= div;
}
double s2 =0;
dx =1;
for (int j = 0; j < n; ++j) {
s2 += dx * variables[j];
dx *= div;
}
f[i] = s1 - s2 * s2 - 1;
}
double x1 = variables[0];
double x2 = variables[1];
f[m - 2] = x1;
f[m - 1] = x2 - x1 * x1 - 1;
return f;
return f;
}
@ -1124,31 +959,15 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
double[][] jacobian = new double[m][];
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure x1 = variables[0];
DerivativeStructure x2 = variables[1];
DerivativeStructure x3 = variables[2];
DerivativeStructure[] f = new DerivativeStructure[m];
for (int i = 0; i < m; ++i) {
double tmp = (i + 1) / 10.0;
jacobian[i] = new double[] {
-tmp * FastMath.exp(-tmp * x1),
tmp * FastMath.exp(-tmp * x2),
FastMath.exp(-i - 1) - FastMath.exp(-tmp)
};
}
return jacobian;
}
@Override
public double[] value(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
double x3 = variables[2];
double[] f = new double[m];
for (int i = 0; i < m; ++i) {
double tmp = (i + 1) / 10.0;
f[i] = FastMath.exp(-tmp * x1) - FastMath.exp(-tmp * x2)
+ (FastMath.exp(-i - 1) - FastMath.exp(-tmp)) * x3;
f[i] = x1.multiply(-tmp).exp().subtract(x2.multiply(-tmp).exp()).add(
x3.multiply(FastMath.exp(-i - 1) - FastMath.exp(-tmp)));
}
return f;
}
@ -1168,27 +987,15 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
double[][] jacobian = new double[m][];
for (int i = 0; i < m; ++i) {
double t = i + 1;
jacobian[i] = new double[] { -t * FastMath.exp(t * x1), -t * FastMath.exp(t * x2) };
}
return jacobian;
}
@Override
public double[] value(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
double[] f = new double[m];
for (int i = 0; i < m; ++i) {
double temp = i + 1;
f[i] = 2 + 2 * temp - FastMath.exp(temp * x1) - FastMath.exp(temp * x2);
}
return f;
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure x1 = variables[0];
DerivativeStructure x2 = variables[1];
DerivativeStructure[] f = new DerivativeStructure[m];
for (int i = 0; i < m; ++i) {
double temp = i + 1;
f[i] = x1.multiply(temp).exp().add(x2.multiply(temp).exp()).subtract(2 + 2 * temp).negate();
}
return f;
}
}
@ -1207,38 +1014,19 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
double x3 = variables[2];
double x4 = variables[3];
double[][] jacobian = new double[m][];
for (int i = 0; i < m; ++i) {
double temp = (i + 1) / 5.0;
double ti = FastMath.sin(temp);
double tmp1 = x1 + temp * x2 - FastMath.exp(temp);
double tmp2 = x3 + ti * x4 - FastMath.cos(temp);
jacobian[i] = new double[] {
2 * tmp1, 2 * temp * tmp1, 2 * tmp2, 2 * ti * tmp2
};
}
return jacobian;
}
@Override
public double[] value(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
double x3 = variables[2];
double x4 = variables[3];
double[] f = new double[m];
for (int i = 0; i < m; ++i) {
double temp = (i + 1) / 5.0;
double tmp1 = x1 + temp * x2 - FastMath.exp(temp);
double tmp2 = x3 + FastMath.sin(temp) * x4 - FastMath.cos(temp);
f[i] = tmp1 * tmp1 + tmp2 * tmp2;
}
return f;
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure x1 = variables[0];
DerivativeStructure x2 = variables[1];
DerivativeStructure x3 = variables[2];
DerivativeStructure x4 = variables[3];
DerivativeStructure[] f = new DerivativeStructure[m];
for (int i = 0; i < m; ++i) {
double temp = (i + 1) / 5.0;
DerivativeStructure tmp1 = x1.add(x2.multiply(temp)).subtract(FastMath.exp(temp));
DerivativeStructure tmp2 = x3.add(x4.multiply(FastMath.sin(temp))).subtract(FastMath.cos(temp));
f[i] = tmp1.multiply(tmp1).add(tmp2.multiply(tmp2));
}
return f;
}
}
@ -1265,63 +1053,34 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
public DerivativeStructure[] value(DerivativeStructure[] variables) {
double[][] jacobian = new double[m][];
for (int i = 0; i < m; ++i) {
jacobian[i] = new double[n];
}
DerivativeStructure[] f = new DerivativeStructure[m];
Arrays.fill(f, variables[0].getField().getZero());
double dx = 1.0 / n;
for (int j = 0; j < n; ++j) {
double tmp1 = 1;
double tmp2 = 2 * variables[j] - 1;
double temp = 2 * tmp2;
double tmp3 = 0;
double tmp4 = 2;
for (int j = 0; j < n; ++j) {
DerivativeStructure tmp1 = variables[0].getField().getOne();
DerivativeStructure tmp2 = variables[j].multiply(2).subtract(1);
DerivativeStructure temp = tmp2.multiply(2);
for (int i = 0; i < m; ++i) {
f[i] = f[i].add(tmp2);
DerivativeStructure ti = temp.multiply(tmp2).subtract(tmp1);
tmp1 = tmp2;
tmp2 = ti;
}
}
double dx = 1.0 / n;
boolean iev = false;
for (int i = 0; i < m; ++i) {
jacobian[i][j] = dx * tmp4;
double ti = 4 * tmp2 + temp * tmp4 - tmp3;
tmp3 = tmp4;
tmp4 = ti;
ti = temp * tmp2 - tmp1;
tmp1 = tmp2;
tmp2 = ti;
f[i] = f[i].multiply(dx);
if (iev) {
f[i] = f[i].add(1.0 / (i * (i + 2)));
}
iev = ! iev;
}
}
return jacobian;
}
@Override
public double[] value(double[] variables) {
double[] f = new double[m];
for (int j = 0; j < n; ++j) {
double tmp1 = 1;
double tmp2 = 2 * variables[j] - 1;
double temp = 2 * tmp2;
for (int i = 0; i < m; ++i) {
f[i] += tmp2;
double ti = temp * tmp2 - tmp1;
tmp1 = tmp2;
tmp2 = ti;
}
}
double dx = 1.0 / n;
boolean iev = false;
for (int i = 0; i < m; ++i) {
f[i] *= dx;
if (iev) {
f[i] += 1.0 / (i * (i + 2));
}
iev = ! iev;
}
return f;
return f;
}
@ -1340,52 +1099,18 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
double[][] jacobian = new double[m][];
for (int i = 0; i < m; ++i) {
jacobian[i] = new double[n];
}
double prod = 1;
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure[] f = new DerivativeStructure[m];
DerivativeStructure sum = variables[0].getField().getZero().subtract(n + 1);
DerivativeStructure prod = variables[0].getField().getOne();
for (int j = 0; j < n; ++j) {
prod *= variables[j];
for (int i = 0; i < n; ++i) {
jacobian[i][j] = 1;
}
jacobian[j][j] = 2;
}
for (int j = 0; j < n; ++j) {
double temp = variables[j];
if (temp == 0) {
temp = 1;
prod = 1;
for (int k = 0; k < n; ++k) {
if (k != j) {
prod *= variables[k];
}
}
}
jacobian[n - 1][j] = prod / temp;
}
return jacobian;
}
@Override
public double[] value(double[] variables) {
double[] f = new double[m];
double sum = -(n + 1);
double prod = 1;
for (int j = 0; j < n; ++j) {
sum += variables[j];
prod *= variables[j];
sum = sum.add(variables[j]);
prod = prod.multiply(variables[j]);
}
for (int i = 0; i < n; ++i) {
f[i] = variables[i] + sum;
f[i] = variables[i].add(sum);
}
f[n - 1] = prod - 1;
f[n - 1] = prod.subtract(1);
return f;
}
@ -1404,36 +1129,18 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
double x2 = variables[1];
double x3 = variables[2];
double x4 = variables[3];
double x5 = variables[4];
double[][] jacobian = new double[m][];
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure x1 = variables[0];
DerivativeStructure x2 = variables[1];
DerivativeStructure x3 = variables[2];
DerivativeStructure x4 = variables[3];
DerivativeStructure x5 = variables[4];
DerivativeStructure[] f = new DerivativeStructure[m];
for (int i = 0; i < m; ++i) {
double temp = 10.0 * i;
double tmp1 = FastMath.exp(-temp * x4);
double tmp2 = FastMath.exp(-temp * x5);
jacobian[i] = new double[] {
-1, -tmp1, -tmp2, temp * x2 * tmp1, temp * x3 * tmp2
};
}
return jacobian;
}
@Override
public double[] value(double[] variables) {
double x1 = variables[0];
double x2 = variables[1];
double x3 = variables[2];
double x4 = variables[3];
double x5 = variables[4];
double[] f = new double[m];
for (int i = 0; i < m; ++i) {
double temp = 10.0 * i;
double tmp1 = FastMath.exp(-temp * x4);
double tmp2 = FastMath.exp(-temp * x5);
f[i] = y[i] - (x1 + x2 * tmp1 + x3 * tmp2);
DerivativeStructure tmp1 = x4.multiply(-temp).exp();
DerivativeStructure tmp2 = x5.multiply(-temp).exp();
f[i] = x1.add(x2.multiply(tmp1)).add(x3.multiply(tmp2)).negate().add(y[i]);
}
return f;
}
@ -1459,65 +1166,28 @@ public class MinpackTest {
}
@Override
public double[][] jacobian(double[] variables) {
double x01 = variables[0];
double x02 = variables[1];
double x03 = variables[2];
double x04 = variables[3];
double x05 = variables[4];
double x06 = variables[5];
double x07 = variables[6];
double x08 = variables[7];
double x09 = variables[8];
double x10 = variables[9];
double x11 = variables[10];
double[][] jacobian = new double[m][];
for (int i = 0; i < m; ++i) {
double temp = i / 10.0;
double tmp1 = FastMath.exp(-x05 * temp);
double tmp2 = FastMath.exp(-x06 * (temp - x09) * (temp - x09));
double tmp3 = FastMath.exp(-x07 * (temp - x10) * (temp - x10));
double tmp4 = FastMath.exp(-x08 * (temp - x11) * (temp - x11));
jacobian[i] = new double[] {
-tmp1,
-tmp2,
-tmp3,
-tmp4,
temp * x01 * tmp1,
x02 * (temp - x09) * (temp - x09) * tmp2,
x03 * (temp - x10) * (temp - x10) * tmp3,
x04 * (temp - x11) * (temp - x11) * tmp4,
-2 * x02 * x06 * (temp - x09) * tmp2,
-2 * x03 * x07 * (temp - x10) * tmp3,
-2 * x04 * x08 * (temp - x11) * tmp4
};
}
return jacobian;
}
@Override
public double[] value(double[] variables) {
double x01 = variables[0];
double x02 = variables[1];
double x03 = variables[2];
double x04 = variables[3];
double x05 = variables[4];
double x06 = variables[5];
double x07 = variables[6];
double x08 = variables[7];
double x09 = variables[8];
double x10 = variables[9];
double x11 = variables[10];
double[] f = new double[m];
for (int i = 0; i < m; ++i) {
double temp = i / 10.0;
double tmp1 = FastMath.exp(-x05 * temp);
double tmp2 = FastMath.exp(-x06 * (temp - x09) * (temp - x09));
double tmp3 = FastMath.exp(-x07 * (temp - x10) * (temp - x10));
double tmp4 = FastMath.exp(-x08 * (temp - x11) * (temp - x11));
f[i] = y[i] - (x01 * tmp1 + x02 * tmp2 + x03 * tmp3 + x04 * tmp4);
}
return f;
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure x01 = variables[0];
DerivativeStructure x02 = variables[1];
DerivativeStructure x03 = variables[2];
DerivativeStructure x04 = variables[3];
DerivativeStructure x05 = variables[4];
DerivativeStructure x06 = variables[5];
DerivativeStructure x07 = variables[6];
DerivativeStructure x08 = variables[7];
DerivativeStructure x09 = variables[8];
DerivativeStructure x10 = variables[9];
DerivativeStructure x11 = variables[10];
DerivativeStructure[] f = new DerivativeStructure[m];
for (int i = 0; i < m; ++i) {
double temp = i / 10.0;
DerivativeStructure tmp1 = x05.multiply(-temp).exp();
DerivativeStructure tmp2 = x06.negate().multiply(x09.subtract(temp).multiply(x09.subtract(temp))).exp();
DerivativeStructure tmp3 = x07.negate().multiply(x10.subtract(temp).multiply(x10.subtract(temp))).exp();
DerivativeStructure tmp4 = x08.negate().multiply(x11.subtract(temp).multiply(x11.subtract(temp))).exp();
f[i] = x01.multiply(tmp1).add(x02.multiply(tmp2)).add(x03.multiply(tmp3)).add(x04.multiply(tmp4)).negate().add(y[i]);
}
return f;
}
private static final double[] y = {

View File

@ -17,12 +17,17 @@
package org.apache.commons.math3.optimization.general;
import java.awt.geom.Point2D;
import java.io.Serializable;
import org.apache.commons.math3.analysis.DifferentiableMultivariateFunction;
import org.apache.commons.math3.analysis.MultivariateFunction;
import org.apache.commons.math3.analysis.MultivariateVectorFunction;
import org.apache.commons.math3.analysis.differentiation.DerivativeStructure;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableFunction;
import org.apache.commons.math3.analysis.solvers.BrentSolver;
import org.apache.commons.math3.exception.DimensionMismatchException;
import org.apache.commons.math3.exception.MathIllegalArgumentException;
import org.apache.commons.math3.geometry.euclidean.twod.Vector2D;
import org.apache.commons.math3.linear.BlockRealMatrix;
import org.apache.commons.math3.linear.RealMatrix;
import org.apache.commons.math3.optimization.GoalType;
@ -337,13 +342,13 @@ public class NonLinearConjugateGradientOptimizerTest {
new BrentSolver(1e-15, 1e-13));
PointValuePair optimum =
optimizer.optimize(100, circle, GoalType.MINIMIZE, new double[] { 98.680, 47.345 });
Point2D.Double center = new Point2D.Double(optimum.getPointRef()[0], optimum.getPointRef()[1]);
Vector2D center = new Vector2D(optimum.getPointRef()[0], optimum.getPointRef()[1]);
Assert.assertEquals(69.960161753, circle.getRadius(center), 1.0e-8);
Assert.assertEquals(96.075902096, center.x, 1.0e-8);
Assert.assertEquals(48.135167894, center.y, 1.0e-8);
Assert.assertEquals(96.075902096, center.getX(), 1.0e-8);
Assert.assertEquals(48.135167894, center.getY(), 1.0e-8);
}
private static class LinearProblem implements DifferentiableMultivariateFunction, Serializable {
private static class LinearProblem implements MultivariateDifferentiableFunction, Serializable {
private static final long serialVersionUID = 703247177355019415L;
final RealMatrix factors;
@ -353,18 +358,6 @@ public class NonLinearConjugateGradientOptimizerTest {
this.target = target;
}
private double[] gradient(double[] point) {
double[] r = factors.operate(point);
for (int i = 0; i < r.length; ++i) {
r[i] -= target[i];
}
double[] p = factors.transpose().operate(r);
for (int i = 0; i < p.length; ++i) {
p[i] *= 2;
}
return p;
}
public double value(double[] variables) {
double[] y = factors.operate(variables);
double sum = 0;
@ -375,20 +368,22 @@ public class NonLinearConjugateGradientOptimizerTest {
return sum;
}
public MultivariateVectorFunction gradient() {
return new MultivariateVectorFunction() {
public double[] value(double[] point) {
return gradient(point);
public DerivativeStructure value(DerivativeStructure[] variables) {
DerivativeStructure[] y = new DerivativeStructure[factors.getRowDimension()];
for (int i = 0; i < y.length; ++i) {
y[i] = variables[0].getField().getZero();
for (int j = 0; j < factors.getColumnDimension(); ++j) {
y[i] = y[i].add(variables[j].multiply(factors.getEntry(i, j)));
}
};
}
DerivativeStructure sum = variables[0].getField().getZero();
for (int i = 0; i < y.length; ++i) {
DerivativeStructure ri = y[i].subtract(target[i]);
sum = sum.add(ri.multiply(ri));
}
return sum;
}
public MultivariateFunction partialDerivative(final int k) {
return new MultivariateFunction() {
public double value(double[] point) {
return gradient(point)[k];
}
};
}
}
}

View File

@ -17,7 +17,6 @@
package org.apache.commons.math3.optimization.general;
import java.awt.geom.Point2D;
import org.apache.commons.math3.random.RandomGenerator;
import org.apache.commons.math3.random.Well44497b;
import org.apache.commons.math3.util.MathUtils;
@ -25,6 +24,7 @@ import org.apache.commons.math3.util.FastMath;
import org.apache.commons.math3.distribution.RealDistribution;
import org.apache.commons.math3.distribution.UniformRealDistribution;
import org.apache.commons.math3.distribution.NormalDistribution;
import org.apache.commons.math3.geometry.euclidean.twod.Vector2D;
/**
* Factory for generating a cloud of points that approximate a circle.
@ -69,8 +69,8 @@ public class RandomCirclePointGenerator {
* @param n Number of points to create.
* @return the cloud of {@code n} points.
*/
public Point2D.Double[] generate(int n) {
final Point2D.Double[] cloud = new Point2D.Double[n];
public Vector2D[] generate(int n) {
final Vector2D[] cloud = new Vector2D[n];
for (int i = 0; i < n; i++) {
cloud[i] = create();
}
@ -82,11 +82,11 @@ public class RandomCirclePointGenerator {
*
* @return a point.
*/
private Point2D.Double create() {
private Vector2D create() {
final double t = tP.sample();
final double pX = cX.sample() + radius * FastMath.cos(t);
final double pY = cY.sample() + radius * FastMath.sin(t);
return new Point2D.Double(pX, pY);
return new Vector2D(pX, pY);
}
}

View File

@ -20,8 +20,8 @@ import java.io.BufferedReader;
import java.io.IOException;
import java.util.ArrayList;
import org.apache.commons.math3.analysis.DifferentiableMultivariateVectorFunction;
import org.apache.commons.math3.analysis.MultivariateMatrixFunction;
import org.apache.commons.math3.analysis.differentiation.DerivativeStructure;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableVectorFunction;
import org.apache.commons.math3.util.MathArrays;
/**
@ -67,7 +67,7 @@ public abstract class StatisticalReferenceDataset {
private double residualSumOfSquares;
/** The least-squares problem. */
private final DifferentiableMultivariateVectorFunction problem;
private final MultivariateDifferentiableVectorFunction problem;
/**
* Creates a new instance of this class from the specified data file. The
@ -150,29 +150,30 @@ public abstract class StatisticalReferenceDataset {
}
this.name = dummyString;
this.problem = new DifferentiableMultivariateVectorFunction() {
this.problem = new MultivariateDifferentiableVectorFunction() {
public double[] value(final double[] a) {
DerivativeStructure[] dsA = new DerivativeStructure[a.length];
for (int i = 0; i < a.length; ++i) {
dsA[i] = new DerivativeStructure(a.length, 0, a[i]);
}
final int n = getNumObservations();
final double[] yhat = new double[n];
for (int i = 0; i < n; i++) {
yhat[i] = getModelValue(getX(i), dsA).getValue();
}
return yhat;
}
public DerivativeStructure[] value(final DerivativeStructure[] a) {
final int n = getNumObservations();
final DerivativeStructure[] yhat = new DerivativeStructure[n];
for (int i = 0; i < n; i++) {
yhat[i] = getModelValue(getX(i), a);
}
return yhat;
}
public MultivariateMatrixFunction jacobian() {
return new MultivariateMatrixFunction() {
public double[][] value(final double[] a)
throws IllegalArgumentException {
final int n = getNumObservations();
final double[][] j = new double[n][];
for (int i = 0; i < n; i++) {
j[i] = getModelDerivatives(getX(i), a);
}
return j;
}
};
}
};
}
@ -309,7 +310,7 @@ public abstract class StatisticalReferenceDataset {
*
* @return the least-squares problem
*/
public DifferentiableMultivariateVectorFunction getLeastSquaresProblem() {
public MultivariateDifferentiableVectorFunction getLeastSquaresProblem() {
return problem;
}
@ -321,18 +322,7 @@ public abstract class StatisticalReferenceDataset {
* @param a the parameters
* @return the value of the model
*/
public abstract double getModelValue(final double x, final double[] a);
/**
* Returns the values of the partial derivatives of the model with respect
* to the parameters.
*
* @param x the predictor variable
* @param a the parameters
* @return the partial derivatives
*/
public abstract double[] getModelDerivatives(final double x,
final double[] a);
public abstract DerivativeStructure getModelValue(final double x, final DerivativeStructure[] a);
/**
* <p>

View File

@ -21,7 +21,7 @@ import java.io.IOException;
import java.io.InputStream;
import java.io.InputStreamReader;
import org.apache.commons.math3.util.FastMath;
import org.apache.commons.math3.analysis.differentiation.DerivativeStructure;
/**
* A factory to create instances of {@link StatisticalReferenceDataset} from
@ -59,25 +59,12 @@ public class StatisticalReferenceDatasetFactory {
dataset = new StatisticalReferenceDataset(in) {
@Override
public double getModelValue(final double x, final double[] a) {
final double p = a[0] + x * (a[1] + x * a[2]);
final double q = 1.0 + x * (a[3] + x * a[4]);
return p / q;
public DerivativeStructure getModelValue(final double x, final DerivativeStructure[] a) {
final DerivativeStructure p = a[0].add(a[1].add(a[2].multiply(x)).multiply(x));
final DerivativeStructure q = a[3].add(a[4].multiply(x)).multiply(x).add(1.0);
return p.divide(q);
}
@Override
public double[] getModelDerivatives(final double x,
final double[] a) {
final double[] dy = new double[5];
final double p = a[0] + x * (a[1] + x * a[2]);
final double q = 1.0 + x * (a[3] + x * a[4]);
dy[0] = 1.0 / q;
dy[1] = x / q;
dy[2] = x * dy[1];
dy[3] = -x * p / (q * q);
dy[4] = x * dy[3];
return dy;
}
};
} finally {
in.close();
@ -93,27 +80,12 @@ public class StatisticalReferenceDatasetFactory {
dataset = new StatisticalReferenceDataset(in) {
@Override
public double getModelValue(final double x, final double[] a) {
final double p = a[0] + x * (a[1] + x * (a[2] + x * a[3]));
final double q = 1.0 + x * (a[4] + x * (a[5] + x * a[6]));
return p / q;
public DerivativeStructure getModelValue(final double x, final DerivativeStructure[] a) {
final DerivativeStructure p = a[0].add(a[1].add(a[2].add(a[3].multiply(x)).multiply(x)).multiply(x));
final DerivativeStructure q = a[4].add(a[5].add(a[6].multiply(x)).multiply(x)).multiply(x).add(1.0);
return p.divide(q);
}
@Override
public double[] getModelDerivatives(final double x,
final double[] a) {
final double[] dy = new double[7];
final double p = a[0] + x * (a[1] + x * (a[2] + x * a[3]));
final double q = 1.0 + x * (a[4] + x * (a[5] + x * a[6]));
dy[0] = 1.0 / q;
dy[1] = x * dy[0];
dy[2] = x * dy[1];
dy[3] = x * dy[2];
dy[4] = -x * p / (q * q);
dy[5] = x * dy[4];
dy[6] = x * dy[5];
return dy;
}
};
} finally {
in.close();
@ -129,22 +101,10 @@ public class StatisticalReferenceDatasetFactory {
dataset = new StatisticalReferenceDataset(in) {
@Override
public double getModelValue(final double x, final double[] a) {
return a[0] + a[1] * FastMath.exp(-a[3] * x) + a[2] *
FastMath.exp(-a[4] * x);
public DerivativeStructure getModelValue(final double x, final DerivativeStructure[] a) {
return a[0].add(a[1].multiply(a[3].multiply(-x).exp())).add(a[2].multiply(a[4].multiply(-x).exp()));
}
@Override
public double[] getModelDerivatives(final double x,
final double[] a) {
final double[] dy = new double[5];
dy[0] = 1.0;
dy[1] = FastMath.exp(-x * a[3]);
dy[2] = FastMath.exp(-x * a[4]);
dy[3] = -x * a[1] * dy[1];
dy[4] = -x * a[2] * dy[2];
return dy;
}
};
} finally {
in.close();
@ -161,25 +121,12 @@ public class StatisticalReferenceDatasetFactory {
dataset = new StatisticalReferenceDataset(in) {
@Override
public double getModelValue(final double x, final double[] a) {
System.out.println(a[0]+", "+a[1]+", "+a[2]+", "+a[3]+", "+a[4]+", "+a[5]);
return a[0] * FastMath.exp(-a[3] * x) +
a[1] * FastMath.exp(-a[4] * x) +
a[2] * FastMath.exp(-a[5] * x);
public DerivativeStructure getModelValue(final double x, final DerivativeStructure[] a) {
return a[0].multiply(a[3].multiply(-x).exp()).add(
a[1].multiply(a[4].multiply(-x).exp())).add(
a[2].multiply(a[5].multiply(-x).exp()));
}
@Override
public double[] getModelDerivatives(final double x,
final double[] a) {
final double[] dy = new double[6];
dy[0] = FastMath.exp(-x * a[3]);
dy[1] = FastMath.exp(-x * a[4]);
dy[2] = FastMath.exp(-x * a[5]);
dy[3] = -x * a[0] * dy[0];
dy[4] = -x * a[1] * dy[1];
dy[5] = -x * a[2] * dy[2];
return dy;
}
};
} finally {
in.close();

View File

@ -18,9 +18,10 @@
package org.apache.commons.math3.optimization.general;
import java.util.ArrayList;
import org.apache.commons.math3.analysis.DifferentiableMultivariateVectorFunction;
import org.apache.commons.math3.analysis.MultivariateMatrixFunction;
import org.apache.commons.math3.analysis.UnivariateFunction;
import org.apache.commons.math3.analysis.differentiation.DerivativeStructure;
import org.apache.commons.math3.analysis.differentiation.MultivariateDifferentiableVectorFunction;
import org.apache.commons.math3.analysis.differentiation.UnivariateDifferentiableFunction;
import org.apache.commons.math3.stat.regression.SimpleRegression;
/**
@ -35,7 +36,7 @@ import org.apache.commons.math3.stat.regression.SimpleRegression;
* <li>for each pair (a, b), the y-coordinate of the line.</li>
* </ul>
*/
class StraightLineProblem implements DifferentiableMultivariateVectorFunction {
class StraightLineProblem implements MultivariateDifferentiableVectorFunction {
/** Cloud of points assumed to be fitted by a straight line. */
private final ArrayList<double[]> points;
/** Error (on the y-coordinate of the points). */
@ -94,7 +95,8 @@ class StraightLineProblem implements DifferentiableMultivariateVectorFunction {
}
public double[] value(double[] params) {
final Model line = new Model(params[0], params[1]);
final Model line = new Model(new DerivativeStructure(0, 0, params[0]),
new DerivativeStructure(0, 0, params[1]));
final double[] model = new double[points.size()];
for (int i = 0; i < points.size(); i++) {
@ -105,12 +107,16 @@ class StraightLineProblem implements DifferentiableMultivariateVectorFunction {
return model;
}
public MultivariateMatrixFunction jacobian() {
return new MultivariateMatrixFunction() {
public double[][] value(double[] point) {
return jacobian(point);
}
};
public DerivativeStructure[] value(DerivativeStructure[] params) {
final Model line = new Model(params[0], params[1]);
final DerivativeStructure[] model = new DerivativeStructure[points.size()];
for (int i = 0; i < points.size(); i++) {
final DerivativeStructure p0 = params[0].getField().getZero().add(points.get(i)[0]);
model[i] = line.value(p0);
}
return model;
}
/**
@ -127,35 +133,26 @@ class StraightLineProblem implements DifferentiableMultivariateVectorFunction {
return result;
}
private double[][] jacobian(double[] params) {
final double[][] jacobian = new double[points.size()][2];
for (int i = 0; i < points.size(); i++) {
final double[] p = points.get(i);
// Partial derivative wrt "a".
jacobian[i][0] = p[0];
// Partial derivative wrt "b".
jacobian[i][1] = 1;
}
return jacobian;
}
/**
* Linear function.
*/
public static class Model implements UnivariateFunction {
final double a;
final double b;
public static class Model implements UnivariateDifferentiableFunction {
final DerivativeStructure a;
final DerivativeStructure b;
public Model(double a,
double b) {
public Model(DerivativeStructure a,
DerivativeStructure b) {
this.a = a;
this.b = b;
}
public double value(double x) {
return a * x + b;
return a.getValue() * x + b.getValue();
}
public DerivativeStructure value(DerivativeStructure x) {
return x.multiply(a).add(b);
}
}
}