org.ejml.dense.row.decomposition.TriangularSolver_FDRM Maven / Gradle / Ivy
Show all versions of ejml-fdense Show documentation
/*
* Copyright (c) 2022, Peter Abeles. All Rights Reserved.
*
* This file is part of Efficient Java Matrix Library (EJML).
*
* Licensed under the Apache License, Version 2.0 (the "License");
* you may not use this file except in compliance with the License.
* You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.ejml.dense.row.decomposition;
import javax.annotation.Generated;
/**
*
* 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
*/
@Generated("org.ejml.dense.row.decomposition.TriangularSolver_DDRM")
public class TriangularSolver_FDRM {
/**
*
* Inverts a square lower triangular matrix: L = L-1
*
*/
public static void invertLower( float[] L, int m ) {
for (int i = 0; i < m; i++) {
float L_ii = L[i*m + i];
for (int j = 0; j < i; j++) {
float 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.0f/L_ii;
}
}
public static void invertLower( float[] L, float[] L_inv, int m ) {
for (int i = 0; i < m; i++) {
float L_ii = L[i*m + i];
for (int j = 0; j < i; j++) {
float 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.0f/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( float[] L, float[] b, int n ) {
// for( int i = 0; i < n; i++ ) {
// float 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( float[] L, float[] b, int n ) {
for (int i = n - 1; i >= 0; i--) {
float sum = b[i];
for (int k = i + 1; k < n; k++) {
sum -= L[k*n + i]*b[k];
}
b[i] = sum/L[i*n + i];
}
}
/**
*
* 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( float[] U, float[] b, int n ) {
// for( int i =n-1; i>=0; i-- ) {
// float sum = b[i];
// for( int j = i+1; j = 0; i--) {
float sum = b[i];
int indexU = i*n + i + 1;
for (int j = i + 1; j < n; j++) {
sum -= U[indexU++]*b[j];
}
b[i] = sum/U[i*n + i];
}
}
public static void solveU( float[] U, float[] b, int sideLength, int minRow, int maxRow ) {
// for( int i =maxRow-1; i>=minRow; i-- ) {
// float sum = b[i];
// for( int j = i+1; j = minRow; i--) {
float sum = b[i];
int indexU = i*sideLength + i + 1;
for (int j = i + 1; j < maxRow; j++) {
sum -= U[indexU++]*b[j];
}
b[i] = sum/U[i*sideLength + i];
}
}
/**
*
* This is a forward substitution solver for non-singular upper triangular matrices which are
* a sub-matrix inside a larger. The columns of 'b' are solved for individually
*
* b = U-1b
*
* where b is a matrix, U is an n by n matrix.
*
*
* @param U Matrix containing the upper triangle system
* @param startU Index of the first element in U
* @param strideU stride between rows
* @param widthU How wide the square matrix is
* @param b Matrix containing the solution to the system. Overwritten with the solution.
* @param startB Index of the first element in B
* @param strideB stride between rows
* @param widthB How wide the matrix is. Length is the same as U's width
*/
public static void solveU( float[] U, int startU, int strideU, int widthU,
float[] b, int startB, int strideB, int widthB ) {
for (int colB = 0; colB < widthB; colB++) {
for (int i = widthU - 1; i >= 0; i--) {
float sum = b[startB + i*strideB + colB];
for (int j = i + 1; j < widthU; j++) {
sum -= U[startU + i*strideU + j]*b[startB + j*strideB + colB];
}
b[startB + i*strideB + colB] = sum/U[startU + i*strideU + i];
}
}
// todo comment out the above and optimize it
}
}