org.ejml.alg.dense.decomposition.TriangularSolver Maven / Gradle / Ivy
Show all versions of ejml Show documentation
/*
* Copyright (c) 2009-2011, Peter Abeles. All Rights Reserved.
*
* This file is part of Efficient Java Matrix Library (EJML).
*
* EJML is free software: you can redistribute it and/or modify
* it under the terms of the GNU Lesser General Public License as
* published by the Free Software Foundation, either version 3
* of the License, or (at your option) any later version.
*
* EJML is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU Lesser General Public License for more details.
*
* You should have received a copy of the GNU Lesser General Public
* License along with EJML. If not, see
* This contains algorithms for solving systems of equations where T is a
* non-singular triangular matrix:
*
* T*x = b
*
* where x and b are vectors, and T is an n by n matrix. T can either be a lower or upper triangular matrix.
*
*
* These functions are designed for use inside of other algorithms. To use them directly
* is dangerous since no sanity checks are performed.
*
*
* @author Peter Abeles
*/
public class TriangularSolver {
/**
*
* Inverts a square lower triangular matrix: L = L-1
*
*
*
* @param L
* @param m
*/
public static void invertLower( double L[] , int m ) {
for( int i = 0; i < m; i++ ) {
double L_ii = L[ i*m + i ];
for( int j = 0; j < i; j++ ) {
double val = 0;
for( int k = j; k < i; k++ ) {
val += L[ i*m + k] * L[ k*m + j ];
}
L[ i*m + j ] = -val / L_ii;
}
L[ i*m + i ] = 1.0 / L_ii;
}
}
public static void invertLower( double L[] , double L_inv[] , int m ) {
for( int i = 0; i < m; i++ ) {
double L_ii = L[ i*m + i ];
for( int j = 0; j < i; j++ ) {
double val = 0;
for( int k = j; k < i; k++ ) {
val -= L[ i*m + k] * L_inv[ k*m + j ];
}
L_inv[ i*m + j ] = val / L_ii;
}
L_inv[ i*m + i ] = 1.0 / L_ii;
}
}
/**
*
* Solves for non-singular lower triangular matrices using forward substitution.
*
* b = L-1b
*
* where b is a vector, L is an n by n matrix.
*
*
* @param L An n by n non-singular lower triangular matrix. Not modified.
* @param b A vector of length n. Modified.
* @param n The size of the matrices.
*/
public static void solveL( double L[] , double []b , int n )
{
// for( int i = 0; i < n; i++ ) {
// double sum = b[i];
// for( int k=0; k
* This is a forward substitution solver for non-singular lower triangular matrices.
*
* b = (LT)-1b
*
* where b is a vector, L is an n by n matrix.
*
*
* L is a lower triangular matrix, but it comes up with a solution as if it was
* an upper triangular matrix that was computed by transposing L.
*
*
* @param L An n by n non-singular lower triangular matrix. Not modified.
* @param b A vector of length n. Modified.
* @param n The size of the matrices.
*/
public static void solveTranL( double L[] , double []b , int n )
{
for( int i =n-1; i>=0; i-- ) {
double sum = b[i];
for( int k = i+1; k
* This is a forward substitution solver for non-singular upper triangular matrices.
*
* b = U-1b
*
* where b is a vector, U is an n by n matrix.
*
*
* @param U An n by n non-singular upper triangular matrix. Not modified.
* @param b A vector of length n. Modified.
* @param n The size of the matrices.
*/
public static void solveU( double U[] , double []b , int n )
{
// for( int i =n-1; i>=0; i-- ) {
// double sum = b[i];
// for( int j = i+1; j =0; i-- ) {
double sum = b[i];
int indexU = i*n+i+1;
for( int j = i+1; j =minRow; i-- ) {
// double sum = b[i];
// for( int j = i+1; j =minRow; i-- ) {
double sum = b[i];
int indexU = i*sideLength+i+1;
for( int j = i+1; j