Author: celestin
Date: Fri Sep 9 02:13:03 2011
New Revision: 1166963
URL: http://svn.apache.org/viewvc?rev=1166963&view=rev
Log:
Removed double[][] solve(double[][]) from LUDecompositionImpl.Solver
Modified:
commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/LUDecompositionImpl.java
Modified: commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/LUDecompositionImpl.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/LUDecompositionImpl.java?rev=1166963&r1=1166962&r2=1166963&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/LUDecompositionImpl.java
(original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/LUDecompositionImpl.java
Fri Sep 9 02:13:03 2011
@@ -340,7 +340,56 @@ public class LUDecompositionImpl impleme
/** {@inheritDoc} */
public RealMatrix solve(RealMatrix b) {
- return new Array2DRowRealMatrix(solve(b.getData()), false);
+
+ final int m = pivot.length;
+ if (b.getRowDimension() != m) {
+ throw new DimensionMismatchException(b.getRowDimension(), m);
+ }
+ if (singular) {
+ throw new SingularMatrixException();
+ }
+
+ final int nColB = b.getColumnDimension();
+
+ // Apply permutations to b
+ final double[][] bp = new double[m][nColB];
+ for (int row = 0; row < m; row++) {
+ final double[] bpRow = bp[row];
+ final int pRow = pivot[row];
+ for (int col = 0; col < nColB; col++) {
+ bpRow[col] = b.getEntry(pRow, col);
+ }
+ }
+
+ // Solve LY = b
+ for (int col = 0; col < m; col++) {
+ final double[] bpCol = bp[col];
+ for (int i = col + 1; i < m; i++) {
+ final double[] bpI = bp[i];
+ final double luICol = lu[i][col];
+ for (int j = 0; j < nColB; j++) {
+ bpI[j] -= bpCol[j] * luICol;
+ }
+ }
+ }
+
+ // Solve UX = Y
+ for (int col = m - 1; col >= 0; col--) {
+ final double[] bpCol = bp[col];
+ final double luDiag = lu[col][col];
+ for (int j = 0; j < nColB; j++) {
+ bpCol[j] /= luDiag;
+ }
+ for (int i = 0; i < col; i++) {
+ final double[] bpI = bp[i];
+ final double luICol = lu[i][col];
+ for (int j = 0; j < nColB; j++) {
+ bpI[j] -= bpCol[j] * luICol;
+ }
+ }
+ }
+
+ return new Array2DRowRealMatrix(bp, false);
}
/** {@inheritDoc} */
|