Author: psteitz
Date: Fri Aug 12 05:39:55 2011
New Revision: 1156968

URL: http://svn.apache.org/viewvc?rev=1156968&view=rev
Log:
Added methods to solve upper and lower triangular systems.  JIRA: MATH-624.  
Contributed by Greg Sterijevski.

Modified:
    
commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/MatrixUtils.java
    commons/proper/math/trunk/src/site/xdoc/changes.xml
    
commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/MatrixUtilsTest.java

Modified: 
commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/MatrixUtils.java
URL: 
http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/MatrixUtils.java?rev=1156968&r1=1156967&r2=1156968&view=diff
==============================================================================
--- 
commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/MatrixUtils.java
 (original)
+++ 
commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/MatrixUtils.java
 Fri Aug 12 05:39:55 2011
@@ -25,6 +25,8 @@ import java.util.Arrays;
 
 import org.apache.commons.math.Field;
 import org.apache.commons.math.FieldElement;
+import org.apache.commons.math.exception.MathArithmeticException;
+import org.apache.commons.math.exception.MathIllegalArgumentException;
 import org.apache.commons.math.exception.OutOfRangeException;
 import org.apache.commons.math.exception.NoDataException;
 import org.apache.commons.math.exception.NumberIsTooSmallException;
@@ -34,6 +36,8 @@ import org.apache.commons.math.exception
 import org.apache.commons.math.exception.util.LocalizedFormats;
 import org.apache.commons.math.fraction.BigFraction;
 import org.apache.commons.math.fraction.Fraction;
+import org.apache.commons.math.util.FastMath;
+import org.apache.commons.math.util.MathUtils;
 
 /**
  * A collection of static methods that operate on or return matrices.
@@ -803,4 +807,84 @@ public class MatrixUtils {
             throw ioe;
         }
     }
+
+    /**Solve  a  system of composed of a Lower Triangular Matrix
+     * {@link RealMatrix}.
+     * <p>
+     * This method is called to solve systems of equations which are
+     * of the lower triangular form. The matrix {@link RealMatrix}
+     * is assumed, though not checked, to be in lower triangular form.
+     * The vector {@link RealVector} is overwritten with the solution.
+     * The matrix is checked that it is square and its dimensions match
+     * the length of the vector.
+     * </p>
+     * @param rm RealMatrix which is lower triangular
+     * @param b  RealVector this is overwritten
+     * @exception IllegalArgumentException if the matrix and vector are not 
conformable
+     * @exception ArithmeticException there is a zero or near zero on the 
diagonal of rm
+     */
+    public static void solveLowerTriangularSystem( RealMatrix rm, RealVector 
b){
+        if ((rm == null) || (b == null) || ( rm.getRowDimension() != 
b.getDimension())) {
+            throw new 
MathIllegalArgumentException(LocalizedFormats.DIMENSIONS_MISMATCH_SIMPLE,
+                    (rm == null) ? 0 : rm.getRowDimension(),
+                    (b == null) ? 0 : b.getDimension());
+        }
+        if( rm.getColumnDimension() != rm.getRowDimension() ){
+            throw new 
MathIllegalArgumentException(LocalizedFormats.DIMENSIONS_MISMATCH_2x2,
+                    rm.getRowDimension(),rm.getRowDimension(),
+                    rm.getRowDimension(),rm.getColumnDimension());
+        }
+        int rows = rm.getRowDimension();
+        for( int i = 0 ; i < rows ; i++ ){
+            double diag = rm.getEntry(i, i);
+            if( FastMath.abs(diag) < MathUtils.SAFE_MIN ){
+                throw new 
MathArithmeticException(LocalizedFormats.ZERO_DENOMINATOR);
+            }
+            double bi = b.getEntry(i)/diag;
+            b.setEntry(i,  bi );
+            for( int j = i+1; j< rows; j++ ){
+                b.setEntry(j, b.getEntry(j)-bi*rm.getEntry(j,i)  );
+            }
+        }
+    }
+
+    /** Solver a  system composed  of an Upper Triangular Matrix
+     * {@link RealMatrix}.
+     * <p>
+     * This method is called to solve systems of equations which are
+     * of the lower triangular form. The matrix {@link RealMatrix}
+     * is assumed, though not checked, to be in upper triangular form.
+     * The vector {@link RealVector} is overwritten with the solution.
+     * The matrix is checked that it is square and its dimensions match
+     * the length of the vector.
+     * </p>
+     * @param rm RealMatrix which is upper triangular
+     * @param b  RealVector this is overwritten
+     * @exception IllegalArgumentException if the matrix and vector are not 
conformable
+     * @exception ArithmeticException there is a zero or near zero on the 
diagonal of rm
+     */
+    public static void solveUpperTriangularSystem( RealMatrix rm, RealVector 
b){
+        if ((rm == null) || (b == null) || ( rm.getRowDimension() != 
b.getDimension())) {
+            throw new 
MathIllegalArgumentException(LocalizedFormats.DIMENSIONS_MISMATCH_SIMPLE,
+                    (rm == null) ? 0 : rm.getRowDimension(),
+                    (b == null) ? 0 : b.getDimension());
+        }
+        if( rm.getColumnDimension() != rm.getRowDimension() ){
+            throw new 
MathIllegalArgumentException(LocalizedFormats.DIMENSIONS_MISMATCH_2x2,
+                    rm.getRowDimension(),rm.getRowDimension(),
+                    rm.getRowDimension(),rm.getColumnDimension());
+        }
+        int rows = rm.getRowDimension();
+        for( int i = rows-1 ; i >-1 ; i-- ){
+            double diag = rm.getEntry(i, i);
+            if( FastMath.abs(diag) < MathUtils.SAFE_MIN ){
+                throw new 
MathArithmeticException(LocalizedFormats.ZERO_DENOMINATOR);
+            }
+            double bi = b.getEntry(i)/diag;
+            b.setEntry(i,  bi );
+            for( int j = i-1; j>-1; j-- ){
+                b.setEntry(j, b.getEntry(j)-bi*rm.getEntry(j,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=1156968&r1=1156967&r2=1156968&view=diff
==============================================================================
--- commons/proper/math/trunk/src/site/xdoc/changes.xml (original)
+++ commons/proper/math/trunk/src/site/xdoc/changes.xml Fri Aug 12 05:39:55 2011
@@ -52,6 +52,9 @@ The <action> type attribute can be add,u
     If the output is not quite correct, check for invisible trailing spaces!
      -->
     <release version="3.0" date="TBD" description="TBD">
+      <action dev="psteitz" type="update" issue="MATH-624" due-to="Greg 
Sterijevski">
+        Added methods to solve upper and lower triangular systems to 
MatrixUtils.
+      </action>
       <action dev="psteitz" type="update" issue="MATH-642">
         Improved performance of nextInt(int) in BitsStreamGenerator.
       </action>

Modified: 
commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/MatrixUtilsTest.java
URL: 
http://svn.apache.org/viewvc/commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/MatrixUtilsTest.java?rev=1156968&r1=1156967&r2=1156968&view=diff
==============================================================================
--- 
commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/MatrixUtilsTest.java
 (original)
+++ 
commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/MatrixUtilsTest.java
 Fri Aug 12 05:39:55 2011
@@ -17,8 +17,7 @@
 package org.apache.commons.math.linear;
 
 import java.math.BigDecimal;
-
-
+import org.apache.commons.math.TestUtils;
 import org.apache.commons.math.fraction.BigFraction;
 import org.apache.commons.math.fraction.Fraction;
 import org.apache.commons.math.fraction.FractionConversionException;
@@ -302,5 +301,29 @@ public final class MatrixUtilsTest {
         }
         return d;
     }
+  
+    @Test 
+    public void testSolveLowerTriangularSystem(){
+        RealMatrix rm = new Array2DRowRealMatrix(
+                new double[][] { {2,0,0,0 }, { 1,1,0,0 }, { 3,3,3,0 }, { 
3,3,3,4 } },
+                       false);
+        RealVector b = new ArrayRealVector(new double[] { 2,3,4,8 }, false);
+        MatrixUtils.solveLowerTriangularSystem(rm, b);
+        TestUtils.assertEquals( new double[]{1,2,-1.66666666666667, 1.0}  , 
b.getData() , 1.0e-12);
+    }
+    
+     
+    /*
+     * Taken from R manual 
http://stat.ethz.ch/R-manual/R-patched/library/base/html/backsolve.html
+     */
+    @Test
+    public void testSolveUpperTriangularSystem(){
+        RealMatrix rm = new Array2DRowRealMatrix(
+                new double[][] { {1,2,3 }, { 0,1,1 }, { 0,0,2 } },
+                       false);
+        RealVector b = new ArrayRealVector(new double[] { 8,4,2 }, false);
+        MatrixUtils.solveUpperTriangularSystem(rm, b);
+        TestUtils.assertEquals( new double[]{-1,3,1}  , b.getData() , 1.0e-12);
+    }
 }
 


Reply via email to