NUMBERS-30: "safe norm" ported from "Commons Math".
Project: http://git-wip-us.apache.org/repos/asf/commons-numbers/repo Commit: http://git-wip-us.apache.org/repos/asf/commons-numbers/commit/9f1d9ef6 Tree: http://git-wip-us.apache.org/repos/asf/commons-numbers/tree/9f1d9ef6 Diff: http://git-wip-us.apache.org/repos/asf/commons-numbers/diff/9f1d9ef6 Branch: refs/heads/master Commit: 9f1d9ef6f79adcbf56d8b71663edbf525a4d5be9 Parents: 7a89652 Author: Gilles Sadowski <[email protected]> Authored: Wed May 24 17:41:06 2017 +0200 Committer: Gilles Sadowski <[email protected]> Committed: Wed May 24 17:41:06 2017 +0200 ---------------------------------------------------------------------- commons-numbers-core/LICENSE.txt | 65 +++++++++++++++ .../apache/commons/numbers/core/SafeNorm.java | 86 ++++++++++++++++++++ .../commons/numbers/core/SafeNormTest.java | 62 ++++++++++++++ 3 files changed, 213 insertions(+) ---------------------------------------------------------------------- http://git-wip-us.apache.org/repos/asf/commons-numbers/blob/9f1d9ef6/commons-numbers-core/LICENSE.txt ---------------------------------------------------------------------- diff --git a/commons-numbers-core/LICENSE.txt b/commons-numbers-core/LICENSE.txt index 261eeb9..67c36eb 100644 --- a/commons-numbers-core/LICENSE.txt +++ b/commons-numbers-core/LICENSE.txt @@ -199,3 +199,68 @@ 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. + +=============================================================================== + +Apache "Commons Numbers" derivative works: + +The "commons-numbers-core" library includes a number of subcomponents +whose implementation is derived from original sources written in C or +Fortran. License terms of the original sources are reproduced below. + +=============================================================================== +Code for "SafeNorm" comes from the MINPACK library. + +Original source copyright and license statement: + +Minpack Copyright Notice (1999) University of Chicago. All rights reserved + +Redistribution and use in source and binary forms, with or +without modification, are permitted provided that the +following conditions are met: + +1. Redistributions of source code must retain the above +copyright notice, this list of conditions and the following +disclaimer. + +2. Redistributions in binary form must reproduce the above +copyright notice, this list of conditions and the following +disclaimer in the documentation and/or other materials +provided with the distribution. + +3. The end-user documentation included with the +redistribution, if any, must include the following +acknowledgment: + + "This product includes software developed by the + University of Chicago, as Operator of Argonne National + Laboratory. + +Alternately, this acknowledgment may appear in the software +itself, if and wherever such third-party acknowledgments +normally appear. + +4. WARRANTY DISCLAIMER. THE SOFTWARE IS SUPPLIED "AS IS" +WITHOUT WARRANTY OF ANY KIND. THE COPYRIGHT HOLDER, THE +UNITED STATES, THE UNITED STATES DEPARTMENT OF ENERGY, AND +THEIR EMPLOYEES: (1) DISCLAIM ANY WARRANTIES, EXPRESS OR +IMPLIED, INCLUDING BUT NOT LIMITED TO ANY IMPLIED WARRANTIES +OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE, TITLE +OR NON-INFRINGEMENT, (2) DO NOT ASSUME ANY LEGAL LIABILITY +OR RESPONSIBILITY FOR THE ACCURACY, COMPLETENESS, OR +USEFULNESS OF THE SOFTWARE, (3) DO NOT REPRESENT THAT USE OF +THE SOFTWARE WOULD NOT INFRINGE PRIVATELY OWNED RIGHTS, (4) +DO NOT WARRANT THAT THE SOFTWARE WILL FUNCTION +UNINTERRUPTED, THAT IT IS ERROR-FREE OR THAT ANY ERRORS WILL +BE CORRECTED. + +5. LIMITATION OF LIABILITY. IN NO EVENT WILL THE COPYRIGHT +HOLDER, THE UNITED STATES, THE UNITED STATES DEPARTMENT OF +ENERGY, OR THEIR EMPLOYEES: BE LIABLE FOR ANY INDIRECT, +INCIDENTAL, CONSEQUENTIAL, SPECIAL OR PUNITIVE DAMAGES OF +ANY KIND OR NATURE, INCLUDING BUT NOT LIMITED TO LOSS OF +PROFITS OR LOSS OF DATA, FOR ANY REASON WHATSOEVER, WHETHER +SUCH LIABILITY IS ASSERTED ON THE BASIS OF CONTRACT, TORT +(INCLUDING NEGLIGENCE OR STRICT LIABILITY), OR OTHERWISE, +EVEN IF ANY OF SAID PARTIES HAS BEEN WARNED OF THE +POSSIBILITY OF SUCH LOSS OR DAMAGES. http://git-wip-us.apache.org/repos/asf/commons-numbers/blob/9f1d9ef6/commons-numbers-core/src/main/java/org/apache/commons/numbers/core/SafeNorm.java ---------------------------------------------------------------------- diff --git a/commons-numbers-core/src/main/java/org/apache/commons/numbers/core/SafeNorm.java b/commons-numbers-core/src/main/java/org/apache/commons/numbers/core/SafeNorm.java new file mode 100644 index 0000000..5d077bc --- /dev/null +++ b/commons-numbers-core/src/main/java/org/apache/commons/numbers/core/SafeNorm.java @@ -0,0 +1,86 @@ +/* + * 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.numbers.core; + +/** + * Computes the Cartesian norm (2-norm), handling both overflow and underflow. + * Translation of the <a href="http://www.netlib.org/minpack">minpack</a> + * "enorm" subroutine. + */ +public class SafeNorm { + /** Constant. */ + private static final double R_DWARF = 3.834e-20; + /** Constant. */ + private static final double R_GIANT = 1.304e+19; + + /** + * @param v Cartesian coordinates. + * @return the 2-norm of the vector. + */ + public static double value(double[] v) { + double s1 = 0; + double s2 = 0; + double s3 = 0; + double x1max = 0; + double x3max = 0; + double floatn = v.length; + double agiant = R_GIANT / floatn; + for (int i = 0; i < v.length; i++) { + double xabs = Math.abs(v[i]); + if (xabs < R_DWARF || xabs > agiant) { + if (xabs > R_DWARF) { + if (xabs > x1max) { + double r = x1max / xabs; + s1 = 1 + s1 * r * r; + x1max = xabs; + } else { + double r = xabs / x1max; + s1 += r * r; + } + } else { + if (xabs > x3max) { + double r = x3max / xabs; + s3 = 1 + s3 * r * r; + x3max = xabs; + } else { + if (xabs != 0) { + double r = xabs / x3max; + s3 += r * r; + } + } + } + } else { + s2 += xabs * xabs; + } + } + double norm; + if (s1 != 0) { + norm = x1max * Math.sqrt(s1 + (s2 / x1max) / x1max); + } else { + if (s2 == 0) { + norm = x3max * Math.sqrt(s3); + } else { + if (s2 >= x3max) { + norm = Math.sqrt(s2 * (1 + (x3max / s2) * (x3max * s3))); + } else { + norm = Math.sqrt(x3max * ((s2 / x3max) + (x3max * s3))); + } + } + } + return norm; + } +} http://git-wip-us.apache.org/repos/asf/commons-numbers/blob/9f1d9ef6/commons-numbers-core/src/test/java/org/apache/commons/numbers/core/SafeNormTest.java ---------------------------------------------------------------------- diff --git a/commons-numbers-core/src/test/java/org/apache/commons/numbers/core/SafeNormTest.java b/commons-numbers-core/src/test/java/org/apache/commons/numbers/core/SafeNormTest.java new file mode 100644 index 0000000..0c4e9bc --- /dev/null +++ b/commons-numbers-core/src/test/java/org/apache/commons/numbers/core/SafeNormTest.java @@ -0,0 +1,62 @@ +/* + * Licensed to the Apache Software Foundation (ASF) under one or more + * contributor license agreements. See the NOTICE file distributed with this + * work for additional information regarding copyright ownership. The ASF + * licenses this file to You under the Apache License, Version 2.0 (the + * "License"); you may not use this file except in compliance with the License. + * You may obtain a copy of the License at + * http://www.apache.org/licenses/LICENSE-2.0 Unless required by applicable law + * or agreed to in writing, software distributed under the License is + * distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY + * KIND, either express or implied. See the License for the specific language + * governing permissions and limitations under the License. + */ +package org.apache.commons.numbers.core; + +import org.junit.Assert; +import org.junit.Test; + +/** + * Test cases for the {@link SafeNorm} class. + */ +public class SafeNormTest { + + @Test + public void testTiny() { + final double s = 1e-320; + final double[] v = new double[] { s, s }; + Assert.assertEquals(Math.sqrt(2) * s, SafeNorm.value(v), 0d); + } + + @Test + public void testBig() { + final double s = 1e300; + final double[] v = new double[] { s, s }; + Assert.assertEquals(Math.sqrt(2) * s, SafeNorm.value(v), 0d); + } + + @Test + public void testOne3D() { + final double s = 1; + final double[] v = new double[] { s, s, s }; + Assert.assertEquals(Math.sqrt(3), SafeNorm.value(v), 0d); + } + + @Test + public void testUnit3D() { + Assert.assertEquals(1, SafeNorm.value(new double[] { 1, 0, 0 }), 0d); + Assert.assertEquals(1, SafeNorm.value(new double[] { 0, 1, 0 }), 0d); + Assert.assertEquals(1, SafeNorm.value(new double[] { 0, 0, 1 }), 0d); + } + + @Test + public void testSimple() { + final double[] v = new double[] { -0.9, 8.7, -6.5, -4.3, -2.1, 0, 1.2, 3.4, -5.6, 7.8, 9.0 }; + double n = 0; + for (int i = 0; i < v.length; i++) { + n += v[i] * v[i]; + } + final double expected = Math.sqrt(n); + Assert.assertEquals(expected, SafeNorm.value(v), 0d); + } +}
