You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@commons.apache.org by lu...@apache.org on 2013/06/03 11:04:40 UTC
svn commit: r1488914 - in /commons/proper/math/trunk: ./ src/changes/
src/main/java/org/apache/commons/math3/analysis/integration/
src/test/java/org/apache/commons/math3/analysis/integration/
Author: luc
Date: Mon Jun 3 09:04:40 2013
New Revision: 1488914
URL: http://svn.apache.org/r1488914
Log:
Added midpoint integration method.
Patch contributed by Oleksandr Kornieiev.
JIRA: MATH-967
Added:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/integration/MidPointIntegrator.java (with props)
commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/integration/MidPointIntegratorTest.java (with props)
Modified:
commons/proper/math/trunk/pom.xml
commons/proper/math/trunk/src/changes/changes.xml
Modified: commons/proper/math/trunk/pom.xml
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/pom.xml?rev=1488914&r1=1488913&r2=1488914&view=diff
==============================================================================
--- commons/proper/math/trunk/pom.xml (original)
+++ commons/proper/math/trunk/pom.xml Mon Jun 3 09:04:40 2013
@@ -217,6 +217,9 @@
<name>Eugene Kirpichov</name>
</contributor>
<contributor>
+ <name>Oleksandr Kornieiev</name>
+ </contributor>
+ <contributor>
<name>Piotr Kochanski</name>
</contributor>
<contributor>
Modified: commons/proper/math/trunk/src/changes/changes.xml
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/changes/changes.xml?rev=1488914&r1=1488913&r2=1488914&view=diff
==============================================================================
--- commons/proper/math/trunk/src/changes/changes.xml (original)
+++ commons/proper/math/trunk/src/changes/changes.xml Mon Jun 3 09:04:40 2013
@@ -51,6 +51,9 @@ If the output is not quite correct, chec
</properties>
<body>
<release version="x.y" date="TBD" description="TBD">
+ <action dev="luc" type="add" issue="MATH-967" due-to="Oleksandr Kornieiev">
+ Added midpoint integration method.
+ </action>
<action dev="luc" type="fix" issue="MATH-988" due-to="Andreas Huber">
Fixed NullPointerException in 2D and 3D sub-line intersections.
</action>
Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/integration/MidPointIntegrator.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/integration/MidPointIntegrator.java?rev=1488914&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/integration/MidPointIntegrator.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/integration/MidPointIntegrator.java Mon Jun 3 09:04:40 2013
@@ -0,0 +1,165 @@
+/*
+ * 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.analysis.integration;
+
+import org.apache.commons.math3.exception.MathIllegalArgumentException;
+import org.apache.commons.math3.exception.MaxCountExceededException;
+import org.apache.commons.math3.exception.NotStrictlyPositiveException;
+import org.apache.commons.math3.exception.NumberIsTooLargeException;
+import org.apache.commons.math3.exception.NumberIsTooSmallException;
+import org.apache.commons.math3.exception.TooManyEvaluationsException;
+import org.apache.commons.math3.util.FastMath;
+
+/**
+ * Implements the <a href="http://en.wikipedia.org/wiki/Midpoint_method">
+ * Midpoint Rule</a> for integration of real univariate functions. For
+ * reference, see <b>Numerical Mathematics</b>, ISBN 0387989595,
+ * chapter 9.2.
+ * <p>
+ * The function should be integrable.</p>
+ *
+ * @version $Id$
+ * @since 3.3
+ */
+public class MidPointIntegrator extends BaseAbstractUnivariateIntegrator {
+
+ /** Maximum number of iterations for midpoint. */
+ public static final int MIDPOINT_MAX_ITERATIONS_COUNT = 64;
+
+ /** Intermediate result. */
+ private double s;
+
+ /**
+ * Build a midpoint integrator with given accuracies and iterations counts.
+ * @param relativeAccuracy relative accuracy of the result
+ * @param absoluteAccuracy absolute accuracy of the result
+ * @param minimalIterationCount minimum number of iterations
+ * @param maximalIterationCount maximum number of iterations
+ * (must be less than or equal to {@link #MIDPOINT_MAX_ITERATIONS_COUNT}
+ * @exception NotStrictlyPositiveException if minimal number of iterations
+ * is not strictly positive
+ * @exception NumberIsTooSmallException if maximal number of iterations
+ * is lesser than or equal to the minimal number of iterations
+ * @exception NumberIsTooLargeException if maximal number of iterations
+ * is greater than {@link #MIDPOINT_MAX_ITERATIONS_COUNT}
+ */
+ public MidPointIntegrator(final double relativeAccuracy,
+ final double absoluteAccuracy,
+ final int minimalIterationCount,
+ final int maximalIterationCount)
+ throws NotStrictlyPositiveException, NumberIsTooSmallException, NumberIsTooLargeException {
+ super(relativeAccuracy, absoluteAccuracy, minimalIterationCount, maximalIterationCount);
+ if (maximalIterationCount > MIDPOINT_MAX_ITERATIONS_COUNT) {
+ throw new NumberIsTooLargeException(maximalIterationCount,
+ MIDPOINT_MAX_ITERATIONS_COUNT, false);
+ }
+ }
+
+ /**
+ * Build a midpoint integrator with given iteration counts.
+ * @param minimalIterationCount minimum number of iterations
+ * @param maximalIterationCount maximum number of iterations
+ * (must be less than or equal to {@link #MIDPOINT_MAX_ITERATIONS_COUNT}
+ * @exception NotStrictlyPositiveException if minimal number of iterations
+ * is not strictly positive
+ * @exception NumberIsTooSmallException if maximal number of iterations
+ * is lesser than or equal to the minimal number of iterations
+ * @exception NumberIsTooLargeException if maximal number of iterations
+ * is greater than {@link #MIDPOINT_MAX_ITERATIONS_COUNT}
+ */
+ public MidPointIntegrator(final int minimalIterationCount,
+ final int maximalIterationCount)
+ throws NotStrictlyPositiveException, NumberIsTooSmallException, NumberIsTooLargeException {
+ super(minimalIterationCount, maximalIterationCount);
+ if (maximalIterationCount > MIDPOINT_MAX_ITERATIONS_COUNT) {
+ throw new NumberIsTooLargeException(maximalIterationCount,
+ MIDPOINT_MAX_ITERATIONS_COUNT, false);
+ }
+ }
+
+ /**
+ * Construct a midpoint integrator with default settings.
+ * (max iteration count set to {@link #MIDPOINT_MAX_ITERATIONS_COUNT})
+ */
+ public MidPointIntegrator() {
+ super(DEFAULT_MIN_ITERATIONS_COUNT, MIDPOINT_MAX_ITERATIONS_COUNT);
+ }
+
+ /**
+ * Compute the n-th stage integral of midpoint rule.
+ * This function should only be called by API <code>integrate()</code> in the package.
+ * To save time it does not verify arguments - caller does.
+ * <p>
+ * The interval is divided equally into 2^n sections rather than an
+ * arbitrary m sections because this configuration can best utilize the
+ * already computed values.</p>
+ *
+ * @param n the stage of 1/2 refinement, n = 0 is no refinement
+ * @return the value of n-th stage integral
+ * @throws TooManyEvaluationsException if the maximal number of evaluations
+ * is exceeded.
+ */
+ private double stage(final int n)
+ throws TooManyEvaluationsException {
+
+ final double max = getMax();
+ final double min = getMin();
+
+ if (n == 0) {
+ final double midPoint = 0.5 * (max - min);
+ s = (max - min) * computeObjectiveValue(midPoint);
+ return s;
+ } else {
+ final long np = 1L << (n - 1); // number of new points in this stage
+ double sum = 0;
+ // spacing between adjacent new points
+ final double spacing = (max - min) / np;
+ double x = min + 0.5 * spacing; // the first new point
+ for (long i = 0; i < np; i++) {
+ sum += computeObjectiveValue(x);
+ x += spacing;
+ }
+ // add the new sum to previously calculated result
+ s = 0.5 * (s + sum * spacing);
+ return s;
+ }
+ }
+
+ /** {@inheritDoc} */
+ protected double doIntegrate()
+ throws MathIllegalArgumentException, TooManyEvaluationsException, MaxCountExceededException {
+
+ double oldt = stage(0);
+ iterations.incrementCount();
+ while (true) {
+ final int i = iterations.getCount();
+ final double t = stage(i);
+ if (i >= getMinimalIterationCount()) {
+ final double delta = FastMath.abs(t - oldt);
+ final double rLimit =
+ getRelativeAccuracy() * (FastMath.abs(oldt) + FastMath.abs(t)) * 0.5;
+ if ((delta <= rLimit) || (delta <= getAbsoluteAccuracy())) {
+ return t;
+ }
+ }
+ oldt = t;
+ iterations.incrementCount();
+ }
+
+ }
+
+}
Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/integration/MidPointIntegrator.java
------------------------------------------------------------------------------
svn:eol-style = native
Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/integration/MidPointIntegrator.java
------------------------------------------------------------------------------
svn:keywords = "Author Date Id Revision"
Added: commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/integration/MidPointIntegratorTest.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/integration/MidPointIntegratorTest.java?rev=1488914&view=auto
==============================================================================
--- commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/integration/MidPointIntegratorTest.java (added)
+++ commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/integration/MidPointIntegratorTest.java Mon Jun 3 09:04:40 2013
@@ -0,0 +1,132 @@
+/*
+ * 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.analysis.integration;
+
+import org.apache.commons.math3.analysis.QuinticFunction;
+import org.apache.commons.math3.analysis.UnivariateFunction;
+import org.apache.commons.math3.analysis.function.Sin;
+import org.apache.commons.math3.exception.NumberIsTooLargeException;
+import org.apache.commons.math3.exception.NumberIsTooSmallException;
+import org.apache.commons.math3.util.FastMath;
+import org.junit.Assert;
+import org.junit.Test;
+
+/**
+ * Test case for midpoint integrator.
+ * <p>
+ * Test runs show that for a default relative accuracy of 1E-6, it generally
+ * takes 10 to 15 iterations for the integral to converge.
+ *
+ * @version $Id: MidPointIntegratorTest.java 1374632 2012-08-18 18:11:11Z luc $
+ */
+public final class MidPointIntegratorTest {
+
+ /**
+ * Test of integrator for the sine function.
+ */
+ @Test
+ public void testSinFunction() {
+ UnivariateFunction f = new Sin();
+ UnivariateIntegrator integrator = new MidPointIntegrator();
+
+ double min = 0;
+ double max = FastMath.PI;
+ double expected = 2;
+ double tolerance = FastMath.abs(expected * integrator.getRelativeAccuracy());
+ double result = integrator.integrate(Integer.MAX_VALUE, f, min, max);
+ Assert.assertTrue(integrator.getEvaluations() < Integer.MAX_VALUE / 2);
+ Assert.assertTrue(integrator.getIterations() < MidPointIntegrator.MIDPOINT_MAX_ITERATIONS_COUNT / 2);
+ Assert.assertEquals(expected, result, tolerance);
+
+ min = -FastMath.PI/3;
+ max = 0;
+ expected = -0.5;
+ tolerance = FastMath.abs(expected * integrator.getRelativeAccuracy());
+ result = integrator.integrate(Integer.MAX_VALUE, f, min, max);
+ Assert.assertTrue(integrator.getEvaluations() < Integer.MAX_VALUE / 2);
+ Assert.assertTrue(integrator.getIterations() < MidPointIntegrator.MIDPOINT_MAX_ITERATIONS_COUNT / 2);
+ Assert.assertEquals(expected, result, tolerance);
+
+ }
+
+ /**
+ * Test of integrator for the quintic function.
+ */
+ @Test
+ public void testQuinticFunction() {
+ UnivariateFunction f = new QuinticFunction();
+ UnivariateIntegrator integrator = new MidPointIntegrator();
+
+ double min = 0;
+ double max = 1;
+ double expected = -1.0 / 48;
+ double tolerance = FastMath.abs(expected * integrator.getRelativeAccuracy());
+ double result = integrator.integrate(Integer.MAX_VALUE, f, min, max);
+ Assert.assertTrue(integrator.getEvaluations() < Integer.MAX_VALUE / 2);
+ Assert.assertTrue(integrator.getIterations() < MidPointIntegrator.MIDPOINT_MAX_ITERATIONS_COUNT / 2);
+ Assert.assertEquals(expected, result, tolerance);
+
+ min = 0;
+ max = 0.5;
+ expected = 11.0 / 768;
+ tolerance = FastMath.abs(expected * integrator.getRelativeAccuracy());
+ result = integrator.integrate(Integer.MAX_VALUE, f, min, max);
+ Assert.assertTrue(integrator.getEvaluations() < Integer.MAX_VALUE / 2);
+ Assert.assertTrue(integrator.getIterations() < MidPointIntegrator.MIDPOINT_MAX_ITERATIONS_COUNT / 2);
+ Assert.assertEquals(expected, result, tolerance);
+
+ min = -1;
+ max = 4;
+ expected = 2048 / 3.0 - 78 + 1.0 / 48;
+ tolerance = FastMath.abs(expected * integrator.getRelativeAccuracy());
+ result = integrator.integrate(Integer.MAX_VALUE, f, min, max);
+ Assert.assertTrue(integrator.getEvaluations() < Integer.MAX_VALUE / 2);
+ Assert.assertTrue(integrator.getIterations() < MidPointIntegrator.MIDPOINT_MAX_ITERATIONS_COUNT / 2);
+ Assert.assertEquals(expected, result, tolerance);
+
+ }
+
+ /**
+ * Test of parameters for the integrator.
+ */
+ @Test
+ public void testParameters() {
+ UnivariateFunction f = new Sin();
+
+ try {
+ // bad interval
+ new MidPointIntegrator().integrate(1000, f, 1, -1);
+ Assert.fail("Expecting NumberIsTooLargeException - bad interval");
+ } catch (NumberIsTooLargeException ex) {
+ // expected
+ }
+ try {
+ // bad iteration limits
+ new MidPointIntegrator(5, 4);
+ Assert.fail("Expecting NumberIsTooSmallException - bad iteration limits");
+ } catch (NumberIsTooSmallException ex) {
+ // expected
+ }
+ try {
+ // bad iteration limits
+ new MidPointIntegrator(10, 99);
+ Assert.fail("Expecting NumberIsTooLargeException - bad iteration limits");
+ } catch (NumberIsTooLargeException ex) {
+ // expected
+ }
+ }
+}
Propchange: commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/integration/MidPointIntegratorTest.java
------------------------------------------------------------------------------
svn:eol-style = native
Propchange: commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/integration/MidPointIntegratorTest.java
------------------------------------------------------------------------------
svn:keywords = "Author Date Id Revision"
Re: svn commit: r1488914 - in /commons/proper/math/trunk: ./ src/changes/
src/main/java/org/apache/commons/math3/analysis/integration/ src/test/java/org/apache/commons/math3/analysis/integration/
Posted by Luc Maisonobe <Lu...@free.fr>.
Le 03/06/2013 12:19, Gilles a écrit :
> Hello.
Hi Gilles,
>
> On Mon, 03 Jun 2013 09:04:40 -0000, luc@apache.org wrote:
>> Author: luc
>> Date: Mon Jun 3 09:04:40 2013
>> New Revision: 1488914
>>
>> URL: http://svn.apache.org/r1488914
>> Log:
>> Added midpoint integration method.
>>
>> Patch contributed by Oleksandr Kornieiev.
>>
>> JIRA: MATH-967
>
>> [...]
>
>> + private double stage(final int n)
>> + throws TooManyEvaluationsException {
>> +
>> + final double max = getMax();
>> + final double min = getMin();
>> +
>> + if (n == 0) {
>> + final double midPoint = 0.5 * (max - min);
> ^^^^^^^^^^^^^^^^^^^^^^^^^^^^
> Is this correct?
No, it was wrong!
The fact it the problem was not detected by the tests cases is because
the contribution of this single wrong initial point is completely
negectible after just a few iteration (each iteration adds 2^n new points).
Solving the error was trivial, but finding a test case that triggers it
was difficult! It involved setting the integration parameters to very
low accuracy, so very few iterations are performed.
Thanks for spotting this.
Luc
>
>> + s = (max - min) * computeObjectiveValue(midPoint);
>> + return s;
>> + } else {
>> + final long np = 1L << (n - 1); // number of
>> new points in this stage
>> + double sum = 0;
>> + // spacing between adjacent new points
>> + final double spacing = (max - min) / np;
>> + double x = min + 0.5 * spacing; // the first new point
>> + for (long i = 0; i < np; i++) {
>> + sum += computeObjectiveValue(x);
>> + x += spacing;
>> + }
>> + // add the new sum to previously calculated result
>> + s = 0.5 * (s + sum * spacing);
>> + return s;
>> + }
>> + }
>> +
>> + /** {@inheritDoc} */
>> + protected double doIntegrate()
>> + throws MathIllegalArgumentException,
>> TooManyEvaluationsException, MaxCountExceededException {
>> +
>> + double oldt = stage(0);
>> + iterations.incrementCount();
>> + while (true) {
>> + final int i = iterations.getCount();
>> + final double t = stage(i);
>> + if (i >= getMinimalIterationCount()) {
>> + final double delta = FastMath.abs(t - oldt);
>> + final double rLimit =
>> + getRelativeAccuracy() * (FastMath.abs(oldt)
>> + FastMath.abs(t)) * 0.5;
>> + if ((delta <= rLimit) || (delta <=
>> getAbsoluteAccuracy())) {
>> + return t;
>> + }
>> + }
>> + oldt = t;
>> + iterations.incrementCount();
>> + }
>> +
>> + }
>> +
>> +}
>>
>> [...]
>
>
> ---------------------------------------------------------------------
> To unsubscribe, e-mail: dev-unsubscribe@commons.apache.org
> For additional commands, e-mail: dev-help@commons.apache.org
>
>
---------------------------------------------------------------------
To unsubscribe, e-mail: dev-unsubscribe@commons.apache.org
For additional commands, e-mail: dev-help@commons.apache.org
Re: svn commit: r1488914 - in /commons/proper/math/trunk: ./
src/changes/ src/main/java/org/apache/commons/math3/analysis/integration/
src/test/java/org/apache/commons/math3/analysis/integration/
Posted by Gilles <gi...@harfang.homelinux.org>.
Hello.
On Mon, 03 Jun 2013 09:04:40 -0000, luc@apache.org wrote:
> Author: luc
> Date: Mon Jun 3 09:04:40 2013
> New Revision: 1488914
>
> URL: http://svn.apache.org/r1488914
> Log:
> Added midpoint integration method.
>
> Patch contributed by Oleksandr Kornieiev.
>
> JIRA: MATH-967
> [...]
> + private double stage(final int n)
> + throws TooManyEvaluationsException {
> +
> + final double max = getMax();
> + final double min = getMin();
> +
> + if (n == 0) {
> + final double midPoint = 0.5 * (max - min);
^^^^^^^^^^^^^^^^^^^^^^^^^^^^
Is this correct?
> + s = (max - min) * computeObjectiveValue(midPoint);
> + return s;
> + } else {
> + final long np = 1L << (n - 1); // number of
> new points in this stage
> + double sum = 0;
> + // spacing between adjacent new points
> + final double spacing = (max - min) / np;
> + double x = min + 0.5 * spacing; // the first new
> point
> + for (long i = 0; i < np; i++) {
> + sum += computeObjectiveValue(x);
> + x += spacing;
> + }
> + // add the new sum to previously calculated result
> + s = 0.5 * (s + sum * spacing);
> + return s;
> + }
> + }
> +
> + /** {@inheritDoc} */
> + protected double doIntegrate()
> + throws MathIllegalArgumentException,
> TooManyEvaluationsException, MaxCountExceededException {
> +
> + double oldt = stage(0);
> + iterations.incrementCount();
> + while (true) {
> + final int i = iterations.getCount();
> + final double t = stage(i);
> + if (i >= getMinimalIterationCount()) {
> + final double delta = FastMath.abs(t - oldt);
> + final double rLimit =
> + getRelativeAccuracy() * (FastMath.abs(oldt)
> + FastMath.abs(t)) * 0.5;
> + if ((delta <= rLimit) || (delta <=
> getAbsoluteAccuracy())) {
> + return t;
> + }
> + }
> + oldt = t;
> + iterations.incrementCount();
> + }
> +
> + }
> +
> +}
>
> [...]
---------------------------------------------------------------------
To unsubscribe, e-mail: dev-unsubscribe@commons.apache.org
For additional commands, e-mail: dev-help@commons.apache.org