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 2008/03/23 14:36:06 UTC
svn commit: r640204 - in /commons/proper/math/trunk/src:
java/org/apache/commons/math/estimation/LevenbergMarquardtEstimator.java
site/xdoc/changes.xml
test/org/apache/commons/math/estimation/LevenbergMarquardtEstimatorTest.java
Author: luc
Date: Sun Mar 23 06:36:03 2008
New Revision: 640204
URL: http://svn.apache.org/viewvc?rev=640204&view=rev
Log:
detect numerical problems in Q.R decomposition for Levenberg-Marquardt estimator
and report them appropriately
JIRA: MATH-199
Modified:
commons/proper/math/trunk/src/java/org/apache/commons/math/estimation/LevenbergMarquardtEstimator.java
commons/proper/math/trunk/src/site/xdoc/changes.xml
commons/proper/math/trunk/src/test/org/apache/commons/math/estimation/LevenbergMarquardtEstimatorTest.java
Modified: commons/proper/math/trunk/src/java/org/apache/commons/math/estimation/LevenbergMarquardtEstimator.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/java/org/apache/commons/math/estimation/LevenbergMarquardtEstimator.java?rev=640204&r1=640203&r2=640204&view=diff
==============================================================================
--- commons/proper/math/trunk/src/java/org/apache/commons/math/estimation/LevenbergMarquardtEstimator.java (original)
+++ commons/proper/math/trunk/src/java/org/apache/commons/math/estimation/LevenbergMarquardtEstimator.java Sun Mar 23 06:36:03 2008
@@ -733,8 +733,9 @@
* are performed in non-increasing columns norms order thanks to columns
* pivoting. The diagonal elements of the R matrix are therefore also in
* non-increasing absolute values order.</p>
+ * @exception EstimationException if the decomposition cannot be performed
*/
- private void qrDecomposition() {
+ private void qrDecomposition() throws EstimationException {
// initializations
for (int k = 0; k < cols; ++k) {
@@ -759,6 +760,10 @@
for (int index = iDiag; index < jacobian.length; index += cols) {
double aki = jacobian[index];
norm2 += aki * aki;
+ }
+ if (Double.isInfinite(norm2) || Double.isNaN(norm2)) {
+ throw new EstimationException("unable to perform Q.R decomposition on the {0}x{1} jacobian matrix",
+ new Object[] { new Integer(rows), new Integer(cols) });
}
if (norm2 > ak2) {
nextColumn = i;
Modified: commons/proper/math/trunk/src/site/xdoc/changes.xml
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/site/xdoc/changes.xml?rev=640204&r1=640203&r2=640204&view=diff
==============================================================================
--- commons/proper/math/trunk/src/site/xdoc/changes.xml (original)
+++ commons/proper/math/trunk/src/site/xdoc/changes.xml Sun Mar 23 06:36:03 2008
@@ -47,6 +47,10 @@
<action dev="luc" type="fix" issue="MATH-198" due-to="Frederick Salardi">
added an error detection for missing imaginary character while parsing complex string
</action>
+ <action dev="luc" type="fix" issue="MATH-199" due-to="Mick">
+ detect numerical problems in Q.R decomposition for Levenberg-Marquardt estimator
+ and report them appropriately
+ </action>
</release>
<release version="1.2" date="2008-02-24"
description="This release combines bug fixes and new features. Most notable
Modified: commons/proper/math/trunk/src/test/org/apache/commons/math/estimation/LevenbergMarquardtEstimatorTest.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/test/org/apache/commons/math/estimation/LevenbergMarquardtEstimatorTest.java?rev=640204&r1=640203&r2=640204&view=diff
==============================================================================
--- commons/proper/math/trunk/src/test/org/apache/commons/math/estimation/LevenbergMarquardtEstimatorTest.java (original)
+++ commons/proper/math/trunk/src/test/org/apache/commons/math/estimation/LevenbergMarquardtEstimatorTest.java Sun Mar 23 06:36:03 2008
@@ -20,6 +20,7 @@
import java.util.ArrayList;
import java.util.HashMap;
import java.util.Iterator;
+import java.util.Locale;
import java.util.Set;
import org.apache.commons.math.estimation.EstimatedParameter;
@@ -587,6 +588,22 @@
assertEquals( 0.20750021499570379, circle.getY(), 1.0e-8);
}
+ public void testMath199() {
+ try {
+ QuadraticProblem problem = new QuadraticProblem();
+ problem.addPoint (0, -3.182591015485607, 0.0);
+ problem.addPoint (1, -2.5581184967730577, 4.4E-323);
+ problem.addPoint (2, -2.1488478161387325, 1.0);
+ problem.addPoint (3, -1.9122489313410047, 4.4E-323);
+ problem.addPoint (4, 1.7785661310051026, 0.0);
+ new LevenbergMarquardtEstimator().estimate(problem);
+ fail("an exception should have been thrown");
+ } catch (EstimationException ee) {
+ // expected behavior
+ }
+
+ }
+
private static class LinearProblem implements EstimationProblem {
public LinearProblem(LinearMeasurement[] measurements) {
@@ -758,6 +775,71 @@
private EstimatedParameter cy;
private ArrayList points;
+ }
+ public class QuadraticProblem extends SimpleEstimationProblem {
+
+ private EstimatedParameter a;
+ private EstimatedParameter b;
+ private EstimatedParameter c;
+
+ public QuadraticProblem() {
+ a = new EstimatedParameter("a", 0.0);
+ b = new EstimatedParameter("b", 0.0);
+ c = new EstimatedParameter("c", 0.0);
+ addParameter(a);
+ addParameter(b);
+ addParameter(c);
+ }
+
+ public void addPoint(double x, double y, double w) {
+ addMeasurement(new LocalMeasurement(x, y, w));
+ }
+
+ public double getA() {
+ return a.getEstimate();
+ }
+
+ public double getB() {
+ return b.getEstimate();
+ }
+
+ public double getC() {
+ return c.getEstimate();
+ }
+
+ public double theoreticalValue(double x) {
+ return ( (a.getEstimate() * x + b.getEstimate() ) * x + c.getEstimate());
+ }
+
+ private double partial(double x, EstimatedParameter parameter) {
+ if (parameter == a) {
+ return x * x;
+ } else if (parameter == b) {
+ return x;
+ } else {
+ return 1.0;
+ }
+ }
+
+ private class LocalMeasurement extends WeightedMeasurement {
+
+ private final double x;
+
+ // constructor
+ public LocalMeasurement(double x, double y, double w) {
+ super(w, y);
+ this.x = x;
+ }
+
+ public double getTheoreticalValue() {
+ return theoreticalValue(x);
+ }
+
+ public double getPartial(EstimatedParameter parameter) {
+ return partial(x, parameter);
+ }
+
+ }
}
public static Test suite() {