added absolute accuracy handling for integrators

git-svn-id: https://svn.apache.org/repos/asf/commons/proper/math/trunk@735536 13f79535-47bb-0310-9956-ffa450edef68
This commit is contained in:
Luc Maisonobe 2009-01-18 20:45:16 +00:00
parent e2048a1dd3
commit 217fe6de71
7 changed files with 67 additions and 267 deletions

View File

@ -283,6 +283,7 @@ public class MessagesResources_fr
"les valeurs de la fonction aux bornes n''ont pas des signes diff\u00e9rents. Bornes : [{0}, {1}], valeurs : [{2}, {3}]" },
// org.apache.commons.math.analysis.solvers.UnivariateRealSolverImpl
// org.apache.commons.math.analysis.integration.UnivariateRealIntegratorImpl
// org.apache.commons.math.transform.FastFourierTransformer
{ "endpoints do not specify an interval: [{0}, {1}]",
"les extr\u00e9mit\u00e9s ne constituent pas un intervalle : [{0}, {1}]" },
@ -295,6 +296,10 @@ public class MessagesResources_fr
{ "function is not differentiable",
"la fonction n''est pas diff\u00e9rentiable" },
// org.apache.commons.math.analysis.integration.UnivariateRealIntegratorImpl
{ "invalid iteration limits: min={0}, max={1}",
"limites d''it\u00e9rations invalides : min = {0}, max = {1}" },
// org.apache.commons.math.fraction.Fraction
{ "zero denominator in fraction {0}/{1}",
"d\u00e9nominateur null dans le nombre rationnel {0}/{1}" },

View File

@ -17,6 +17,7 @@
package org.apache.commons.math.analysis.integration;
import org.apache.commons.math.FunctionEvaluationException;
import org.apache.commons.math.MathRuntimeException;
import org.apache.commons.math.MaxIterationsExceededException;
import org.apache.commons.math.analysis.UnivariateRealFunction;
@ -26,7 +27,7 @@ import org.apache.commons.math.analysis.UnivariateRealFunction;
* reference, see <b>Introduction to Numerical Analysis</b>, ISBN 038795452X,
* chapter 3.
* <p>
* Romberg integration employs k successvie refinements of the trapezoid
* Romberg integration employs k successive refinements of the trapezoid
* rule to remove error terms less than order O(N^(-2k)). Simpson's rule
* is a special case of k = 2.</p>
*
@ -47,23 +48,12 @@ public class RombergIntegrator extends UnivariateRealIntegratorImpl {
super(f, 32);
}
/**
* Integrate the function in the given interval.
*
* @param min the lower bound for the interval
* @param max the upper bound for the interval
* @return the value of integral
* @throws MaxIterationsExceededException if the maximum iteration count is exceeded
* or the integrator detects convergence problems otherwise
* @throws FunctionEvaluationException if an error occurs evaluating the
* function
* @throws IllegalArgumentException if any parameters are invalid
*/
public double integrate(double min, double max) throws MaxIterationsExceededException,
/** {@inheritDoc} */
public double integrate(double min, double max) throws MaxIterationsExceededException,
FunctionEvaluationException, IllegalArgumentException {
int i = 1, j, m = maximalIterationCount + 1;
// Array strcture here can be improved for better space
// Array structure here can be improved for better space
// efficiency because only the lower triangle is used.
double r, t[][] = new double[m][m], s, olds;
@ -83,7 +73,10 @@ public class RombergIntegrator extends UnivariateRealIntegratorImpl {
}
s = t[i][i];
if (i >= minimalIterationCount) {
if (Math.abs(s - olds) <= Math.abs(relativeAccuracy * olds)) {
final double delta = Math.abs(s - olds);
final double rLimit =
relativeAccuracy * (Math.abs(olds) + Math.abs(s)) * 0.5;
if ((delta <= rLimit) || (delta <= absoluteAccuracy)) {
setResult(s, i);
return result;
}
@ -94,18 +87,14 @@ public class RombergIntegrator extends UnivariateRealIntegratorImpl {
throw new MaxIterationsExceededException(maximalIterationCount);
}
/**
* Verifies that the iteration limits are valid and within the range.
*
* @throws IllegalArgumentException if not
*/
/** {@inheritDoc} */
protected void verifyIterationCount() throws IllegalArgumentException {
super.verifyIterationCount();
// at most 32 bisection refinements due to higher order divider
if (maximalIterationCount > 32) {
throw new IllegalArgumentException
("Iteration upper limit out of [0, 32] range: " +
maximalIterationCount);
throw MathRuntimeException.createIllegalArgumentException(
"invalid iteration limits: min={0}, max={1}",
new Object[] { 0, 32 });
}
}
}

View File

@ -17,6 +17,7 @@
package org.apache.commons.math.analysis.integration;
import org.apache.commons.math.FunctionEvaluationException;
import org.apache.commons.math.MathRuntimeException;
import org.apache.commons.math.MaxIterationsExceededException;
import org.apache.commons.math.analysis.UnivariateRealFunction;
@ -46,18 +47,7 @@ public class SimpsonIntegrator extends UnivariateRealIntegratorImpl {
super(f, 64);
}
/**
* Integrate the function in the given interval.
*
* @param min the lower bound for the interval
* @param max the upper bound for the interval
* @return the value of integral
* @throws MaxIterationsExceededException if the maximum iteration count is exceeded
* or the integrator detects convergence problems otherwise
* @throws FunctionEvaluationException if an error occurs evaluating the
* function
* @throws IllegalArgumentException if any parameters are invalid
*/
/** {@inheritDoc} */
public double integrate(double min, double max) throws MaxIterationsExceededException,
FunctionEvaluationException, IllegalArgumentException {
@ -81,7 +71,10 @@ public class SimpsonIntegrator extends UnivariateRealIntegratorImpl {
t = qtrap.stage(min, max, i);
s = (4 * t - oldt) / 3.0;
if (i >= minimalIterationCount) {
if (Math.abs(s - olds) <= Math.abs(relativeAccuracy * olds)) {
final double delta = Math.abs(s - olds);
final double rLimit =
relativeAccuracy * (Math.abs(olds) + Math.abs(s)) * 0.5;
if ((delta <= rLimit) || (delta <= absoluteAccuracy)) {
setResult(s, i);
return result;
}
@ -93,18 +86,14 @@ public class SimpsonIntegrator extends UnivariateRealIntegratorImpl {
throw new MaxIterationsExceededException(maximalIterationCount);
}
/**
* Verifies that the iteration limits are valid and within the range.
*
* @throws IllegalArgumentException if not
*/
/** {@inheritDoc} */
protected void verifyIterationCount() throws IllegalArgumentException {
super.verifyIterationCount();
// at most 64 bisection refinements
if (maximalIterationCount > 64) {
throw new IllegalArgumentException
("Iteration upper limit out of [0, 64] range: " +
maximalIterationCount);
throw MathRuntimeException.createIllegalArgumentException(
"invalid iteration limits: min={0}, max={1}",
new Object[] { 0, 64 });
}
}
}

View File

@ -17,6 +17,7 @@
package org.apache.commons.math.analysis.integration;
import org.apache.commons.math.FunctionEvaluationException;
import org.apache.commons.math.MathRuntimeException;
import org.apache.commons.math.MaxIterationsExceededException;
import org.apache.commons.math.analysis.UnivariateRealFunction;
@ -87,18 +88,7 @@ public class TrapezoidIntegrator extends UnivariateRealIntegratorImpl {
}
}
/**
* Integrate the function in the given interval.
*
* @param min the lower bound for the interval
* @param max the upper bound for the interval
* @return the value of integral
* @throws MaxIterationsExceededException if the maximum iteration count is exceeded
* or the integrator detects convergence problems otherwise
* @throws FunctionEvaluationException if an error occurs evaluating the
* function
* @throws IllegalArgumentException if any parameters are invalid
*/
/** {@inheritDoc} */
public double integrate(double min, double max) throws MaxIterationsExceededException,
FunctionEvaluationException, IllegalArgumentException {
@ -113,7 +103,10 @@ public class TrapezoidIntegrator extends UnivariateRealIntegratorImpl {
while (i <= maximalIterationCount) {
t = stage(min, max, i);
if (i >= minimalIterationCount) {
if (Math.abs(t - oldt) <= Math.abs(relativeAccuracy * oldt)) {
final double delta = Math.abs(t - oldt);
final double rLimit =
relativeAccuracy * (Math.abs(oldt) + Math.abs(t)) * 0.5;
if ((delta <= rLimit) || (delta <= absoluteAccuracy)) {
setResult(t, i);
return result;
}
@ -124,18 +117,14 @@ public class TrapezoidIntegrator extends UnivariateRealIntegratorImpl {
throw new MaxIterationsExceededException(maximalIterationCount);
}
/**
* Verifies that the iteration limits are valid and within the range.
*
* @throws IllegalArgumentException if not
*/
/** {@inheritDoc} */
protected void verifyIterationCount() throws IllegalArgumentException {
super.verifyIterationCount();
// at most 64 bisection refinements
if (maximalIterationCount > 64) {
throw new IllegalArgumentException
("Iteration upper limit out of [0, 64] range: " +
maximalIterationCount);
throw MathRuntimeException.createIllegalArgumentException(
"invalid iteration limits: min={0}, max={1}",
new Object[] { 0, 64 });
}
}
}

View File

@ -17,6 +17,7 @@
package org.apache.commons.math.analysis.integration;
import org.apache.commons.math.ConvergenceException;
import org.apache.commons.math.ConvergingAlgorithm;
import org.apache.commons.math.FunctionEvaluationException;
/**
@ -25,39 +26,9 @@ import org.apache.commons.math.FunctionEvaluationException;
* @version $Revision$ $Date$
* @since 1.2
*/
public interface UnivariateRealIntegrator {
public interface UnivariateRealIntegrator extends ConvergingAlgorithm {
/**
* Set the upper limit for the number of iterations.
* <p>
* Usually a high iteration count indicates convergence problem. However,
* the "reasonable value" varies widely for different cases. Users are
* advised to use the default value.</p>
* <p>
* A <code>ConvergenceException</code> will be thrown if this number
* is exceeded.</p>
*
* @param count maximum number of iterations
*/
void setMaximalIterationCount(int count);
/**
* Get the upper limit for the number of iterations.
*
* @return the actual upper limit
*/
int getMaximalIterationCount();
/**
* Reset the upper limit for the number of iterations to the default.
* <p>
* The default value is supplied by the implementation.</p>
*
* @see #setMaximalIterationCount(int)
*/
void resetMaximalIterationCount();
/**
/**
* Set the lower limit for the number of iterations.
* <p>
* Minimal iteration is needed to avoid false early convergence, e.g.
@ -87,33 +58,6 @@ public interface UnivariateRealIntegrator {
*/
void resetMinimalIterationCount();
/**
* Set the relative accuracy.
* <p>
* This is used to stop iterations.</p>
*
* @param accuracy the relative accuracy
* @throws IllegalArgumentException if the accuracy can't be achieved
* or is otherwise deemed unreasonable
*/
void setRelativeAccuracy(double accuracy);
/**
* Get the actual relative accuracy.
*
* @return the accuracy
*/
double getRelativeAccuracy();
/**
* Reset the relative accuracy to the default.
* <p>
* The default value is provided by the implementation.</p>
*
* @see #setRelativeAccuracy(double)
*/
void resetRelativeAccuracy();
/**
* Integrate the function in the given interval.
*
@ -139,18 +83,4 @@ public interface UnivariateRealIntegrator {
*/
double getResult() throws IllegalStateException;
/**
* Get the number of iterations in the last run of the integrator.
* <p>
* This is mainly meant for testing purposes. It may occasionally
* help track down performance problems: if the iteration count
* is notoriously high, check whether the function is evaluated
* properly, and whether another integrator is more amenable to the
* problem.</p>
*
* @return the last iteration count
* @throws IllegalStateException if there is no result available, either
* because no result was yet computed or the last attempt failed
*/
int getIterationCount() throws IllegalStateException;
}

View File

@ -16,8 +16,7 @@
*/
package org.apache.commons.math.analysis.integration;
import java.io.Serializable;
import org.apache.commons.math.ConvergingAlgorithmImpl;
import org.apache.commons.math.MathRuntimeException;
import org.apache.commons.math.analysis.UnivariateRealFunction;
@ -27,27 +26,15 @@ import org.apache.commons.math.analysis.UnivariateRealFunction;
* @version $Revision$ $Date$
* @since 1.2
*/
public abstract class UnivariateRealIntegratorImpl implements
UnivariateRealIntegrator, Serializable {
public abstract class UnivariateRealIntegratorImpl
extends ConvergingAlgorithmImpl implements UnivariateRealIntegrator {
/** serializable version identifier */
static final long serialVersionUID = -3365294665201465048L;
/** maximum relative error */
protected double relativeAccuracy;
/** maximum number of iterations */
protected int maximalIterationCount;
/** Serializable version identifier. */
private static final long serialVersionUID = 6248808456637441533L;
/** minimum number of iterations */
protected int minimalIterationCount;
/** default maximum relative error */
protected double defaultRelativeAccuracy;
/** default maximum number of iterations */
protected int defaultMaximalIterationCount;
/** default minimum number of iterations */
protected int defaultMinimalIterationCount;
@ -57,9 +44,6 @@ public abstract class UnivariateRealIntegratorImpl implements
/** the last computed integral */
protected double result;
/** the last iteration count */
protected int iterationCount;
/** the integrand function */
protected UnivariateRealFunction f;
@ -71,21 +55,18 @@ public abstract class UnivariateRealIntegratorImpl implements
* @throws IllegalArgumentException if f is null or the iteration
* limits are not valid
*/
protected UnivariateRealIntegratorImpl(
UnivariateRealFunction f,
int defaultMaximalIterationCount) throws IllegalArgumentException {
protected UnivariateRealIntegratorImpl(final UnivariateRealFunction f,
final int defaultMaximalIterationCount)
throws IllegalArgumentException {
super(defaultMaximalIterationCount, 1.0e-15);
if (f == null) {
throw new IllegalArgumentException("Function can not be null.");
}
this.f = f;
// parameters that may depend on algorithm
this.defaultMaximalIterationCount = defaultMaximalIterationCount;
this.maximalIterationCount = defaultMaximalIterationCount;
// parameters that are problem specific
this.defaultRelativeAccuracy = 1E-6;
this.relativeAccuracy = defaultRelativeAccuracy;
setRelativeAccuracy(1.0e-6);
this.defaultMinimalIterationCount = 3;
this.minimalIterationCount = defaultMinimalIterationCount;
@ -106,20 +87,6 @@ public abstract class UnivariateRealIntegratorImpl implements
}
}
/**
* Access the last iteration count.
*
* @return the last iteration count
* @throws IllegalStateException if no integral has been computed
*/
public int getIterationCount() throws IllegalStateException {
if (resultComputed) {
return iterationCount;
} else {
throw MathRuntimeException.createIllegalStateException("no result available", null);
}
}
/**
* Convenience function for implementations.
*
@ -136,98 +103,25 @@ public abstract class UnivariateRealIntegratorImpl implements
* Convenience function for implementations.
*/
protected final void clearResult() {
this.iterationCount = 0;
this.resultComputed = false;
}
/**
* Set the upper limit for the number of iterations.
*
* @param count maximum number of iterations
*/
public void setMaximalIterationCount(int count) {
maximalIterationCount = count;
}
/**
* Get the upper limit for the number of iterations.
*
* @return the actual upper limit
*/
public int getMaximalIterationCount() {
return maximalIterationCount;
}
/**
* Reset the upper limit for the number of iterations to the default.
*/
public void resetMaximalIterationCount() {
maximalIterationCount = defaultMaximalIterationCount;
}
/**
* Set the lower limit for the number of iterations.
*
* @param count minimum number of iterations
*/
/** {@inheritDoc} */
public void setMinimalIterationCount(int count) {
minimalIterationCount = count;
}
/**
* Get the lower limit for the number of iterations.
*
* @return the actual lower limit
*/
/** {@inheritDoc} */
public int getMinimalIterationCount() {
return minimalIterationCount;
}
/**
* Reset the lower limit for the number of iterations to the default.
*/
/** {@inheritDoc} */
public void resetMinimalIterationCount() {
minimalIterationCount = defaultMinimalIterationCount;
}
/**
* Set the relative accuracy.
*
* @param accuracy the relative accuracy
* @throws IllegalArgumentException if the accuracy can't be achieved by
* the integrator or is otherwise deemed unreasonable
*/
public void setRelativeAccuracy(double accuracy) {
relativeAccuracy = accuracy;
}
/**
* Get the actual relative accuracy.
*
* @return the accuracy
*/
public double getRelativeAccuracy() {
return relativeAccuracy;
}
/**
* Reset the relative accuracy to the default.
*/
public void resetRelativeAccuracy() {
relativeAccuracy = defaultRelativeAccuracy;
}
/**
* Returns true if the arguments form a (strictly) increasing sequence
*
* @param start first number
* @param mid second number
* @param end third number
* @return true if the arguments form an increasing sequence
*/
protected boolean isSequence(double start, double mid, double end) {
return (start < mid) && (mid < end);
}
/**
* Verifies that the endpoints specify an interval.
*
@ -238,9 +132,9 @@ public abstract class UnivariateRealIntegratorImpl implements
protected void verifyInterval(double lower, double upper) throws
IllegalArgumentException {
if (lower >= upper) {
throw new IllegalArgumentException
("Endpoints do not specify an interval: [" + lower +
", " + upper + "]");
throw MathRuntimeException.createIllegalArgumentException(
"endpoints do not specify an interval: [{0}, {1}]",
new Object[] { lower, upper });
}
}
@ -250,10 +144,10 @@ public abstract class UnivariateRealIntegratorImpl implements
* @throws IllegalArgumentException if not valid
*/
protected void verifyIterationCount() throws IllegalArgumentException {
if (!isSequence(0, minimalIterationCount, maximalIterationCount+1)) {
throw new IllegalArgumentException
("Invalid iteration limits: min=" + minimalIterationCount +
" max=" + maximalIterationCount);
if ((minimalIterationCount <= 0) || (maximalIterationCount <= minimalIterationCount)) {
throw MathRuntimeException.createIllegalArgumentException(
"invalid iteration limits: min={0}, max={1}",
new Object[] { minimalIterationCount, maximalIterationCount });
}
}
}

View File

@ -39,6 +39,10 @@ The <action> type attribute can be add,update,fix,remove.
</properties>
<body>
<release version="2.0" date="TBD" description="TBD">
<action dev="luc" type="update" >
Integration algorithms now can have both relative and absolute
accuracy settings.
</action>
<action dev="luc" type="add" issue="MATH-177" due-to="Gilles Sadowski">
Added a new minimization sub-package below the analysis package.
</action>