Apache
/
commons-math
mirror of https://github.com/apache/commons-math.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Added methods to solve upper and lower triangular systems. JIRA: MATH-624. Contributed by Greg Sterijevski. 

git-svn-id: https://svn.apache.org/repos/asf/commons/proper/math/trunk@1156968 13f79535-47bb-0310-9956-ffa450edef68
This commit is contained in:
Phil Steitz 2011-08-12 05:39:55 +00:00
parent 2f2bb7783c
commit 4a65b71c2a
3 changed files with 112 additions and 2 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

84
src/main/java/org/apache/commons/math/linear/MatrixUtils.java
View File

@ -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.ZeroException;
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) );
}
}
}
}

3
src/site/xdoc/changes.xml
View File

@ -52,6 +52,9 @@ The <action> type attribute can be add,update,fix,remove.
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>

27
src/test/java/org/apache/commons/math/linear/MatrixUtilsTest.java
View File

@ -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);
}
}