Phil Steitz a écrit : > Any ideas on where this stuff can go??? I am beginning to think we should have a few different matrices shapes: upper triangular, lower triangular, diagonal, tri-diagonal and perhaps sparse.
In this case, this could be put as an specialized implementation of solve in a triangular matrix class. Luc >> - >> + + /** TODO: Find a home for the following methods in the >> linear package */ + + /** >> + * <p>Uses back substitution to solve the system</p> >> + * + * <p>coefficients X = constants</p> >> + * + * <p>coefficients must upper-triangular and constants >> must be a column + * matrix. The solution is returned as a column >> matrix.</p> >> + * + * <p>The number of columns in coefficients determines >> the length >> + * of the returned solution vector (column matrix). If constants >> + * has more rows than coefficients has columns, excess rows are >> ignored. >> + * Similarly, extra (zero) rows in coefficients are ignored</p> >> + * + * @param coefficients upper-triangular coefficients matrix >> + * @param constants column RHS constants matrix >> + * @return solution matrix as a column matrix >> + * + */ >> + private static RealMatrix solveUpperTriangular(RealMatrixImpl >> coefficients, >> + RealMatrixImpl constants) { >> + if (!isUpperTriangular(coefficients, 1E-12)) { >> + throw new IllegalArgumentException( >> + "Coefficients is not upper-triangular"); >> + } >> + if (constants.getColumnDimension() != 1) { >> + throw new IllegalArgumentException( >> + "Constants not a column matrix."); >> + } >> + int length = coefficients.getColumnDimension(); >> + double[][] cons = constants.getDataRef(); >> + double[][] coef = coefficients.getDataRef(); >> + double x[] = new double[length]; >> + for (int i = 0; i < length; i++) { >> + int index = length - 1 - i; >> + double sum = 0; >> + for (int j = index + 1; j < length; j++) { >> + sum += coef[index][j] * x[j]; >> + } >> + x[index] = (cons[index][0] - sum) / coef[index][index]; >> + } + return new RealMatrixImpl(x); >> + } >> + + /** >> + * <p>Returns true iff m is an upper-triangular matrix.</p> >> + * + * <p>Makes sure all below-diagonal elements are within >> epsilon of 0.</p> >> + * + * @param m matrix to check >> + * @param epsilon maximum allowable absolute value for elements >> below >> + * the main diagonal >> + * + * @return true if m is upper-triangular; false otherwise >> + * @throws NullPointerException if m is null >> + */ >> + private static boolean isUpperTriangular(RealMatrixImpl m, double >> epsilon) { >> + double[][] data = m.getDataRef(); >> + int nCols = m.getColumnDimension(); >> + int nRows = m.getRowDimension(); >> + for (int r = 0; r < nRows; r++) { >> + int bound = Math.min(r, nCols); >> + for (int c = 0; c < bound; c++) { >> + if (Math.abs(data[r][c]) > epsilon) { >> + return false; >> + } >> + } >> + } >> + return true; >> + } >> } >> >> > > > --------------------------------------------------------------------- > To unsubscribe, e-mail: [EMAIL PROTECTED] > For additional commands, e-mail: [EMAIL PROTECTED] > > --------------------------------------------------------------------- To unsubscribe, e-mail: [EMAIL PROTECTED] For additional commands, e-mail: [EMAIL PROTECTED]