diff --git a/src/main/java/com/thealgorithms/matrix/LUDecomposition.java b/src/main/java/com/thealgorithms/matrix/LUDecomposition.java index e41aaa201338..9e7838648422 100644 --- a/src/main/java/com/thealgorithms/matrix/LUDecomposition.java +++ b/src/main/java/com/thealgorithms/matrix/LUDecomposition.java @@ -1,88 +1,186 @@ package com.thealgorithms.matrix; /** - * LU Decomposition algorithm - * -------------------------- - * Decomposes a square matrix a into a product of two matrices: - * a = l * u - * where: - * - l is a lower triangular matrix with 1s on its diagonal - * - u is an upper triangular matrix + * LU Decomposition algorithm for square matrices. + * Decomposes a matrix A into L (lower triangular) and U (upper triangular) + * such that A = L * U * - * Reference: - * https://en.wikipedia.org/wiki/lu_decomposition + *
Time Complexity: O(n^3) + *
Space Complexity: O(n^2) + * + * @author Raghu0703 + * @see LU Decomposition */ public final class LUDecomposition { - private LUDecomposition() { } /** - * A helper class to store both l and u matrices + * Performs LU decomposition on a square matrix using Doolittle's method. + * + * @param matrix The input square matrix + * @return A Result object containing L and U matrices + * @throws IllegalArgumentException if matrix is not square or singular */ - public static class LU { - double[][] l; - double[][] u; + public static Result decompose(double[][] matrix) { + int n = matrix.length; - LU(double[][] l, double[][] u) { - this.l = l; - this.u = u; + // Validate input + if (n == 0) { + throw new IllegalArgumentException("Matrix cannot be empty"); + } + for (double[] row : matrix) { + if (row.length != n) { + throw new IllegalArgumentException("Matrix must be square"); + } } - } - /** - * Performs LU Decomposition on a square matrix a - * - * @param a input square matrix - * @return LU object containing l and u matrices - */ - public static LU decompose(double[][] a) { - int n = a.length; double[][] l = new double[n][n]; double[][] u = new double[n][n]; + // Initialize L with identity matrix for (int i = 0; i < n; i++) { - // upper triangular matrix - for (int k = i; k < n; k++) { - double sum = 0; - for (int j = 0; j < i; j++) { - sum += l[i][j] * u[j][k]; + l[i][i] = 1.0; + } + + // Perform LU decomposition using Doolittle's method + for (int j = 0; j < n; j++) { + // Calculate U matrix elements + for (int i = 0; i <= j; i++) { + double sum = 0.0; + for (int k = 0; k < i; k++) { + sum += l[i][k] * u[k][j]; } - u[i][k] = a[i][k] - sum; + u[i][j] = matrix[i][j] - sum; } - // lower triangular matrix - for (int k = i; k < n; k++) { - if (i == k) { - l[i][i] = 1; // diagonal as 1 - } else { - double sum = 0; - for (int j = 0; j < i; j++) { - sum += l[k][j] * u[j][i]; - } - l[k][i] = (a[k][i] - sum) / u[i][i]; + // Calculate L matrix elements + for (int i = j + 1; i < n; i++) { + double sum = 0.0; + for (int k = 0; k < j; k++) { + sum += l[i][k] * u[k][j]; } + + if (Math.abs(u[j][j]) < 1e-10) { + throw new IllegalArgumentException( + "Matrix is singular or nearly singular" + ); + } + + l[i][j] = (matrix[i][j] - sum) / u[j][j]; } } - return new LU(l, u); + return new Result(l, u); + } + + /** + * Main method for demonstration. + * + * @param args command line arguments (not used) + */ + public static void main(String[] args) { + // Example from the issue + double[][] matrix = { + {2, -1, -2}, + {-4, 6, 3}, + {-4, -2, 8} + }; + + System.out.println("LU Decomposition Example"); + System.out.println("========================\n"); + + Result result = decompose(matrix); + + System.out.println("Original Matrix A:"); + printMatrix(matrix); + + System.out.println("\nLower Triangular Matrix L:"); + printMatrix(result.getL()); + + System.out.println("\nUpper Triangular Matrix U:"); + printMatrix(result.getU()); + + // Verify that L * U = A + System.out.println("\nVerification: L * U = A"); + double[][] product = multiplyMatrices(result.getL(), result.getU()); + printMatrix(product); } /** - * Utility function to print a matrix + * Helper method to print a matrix. * - * @param m matrix to print + * @param matrix the matrix to print */ - public static void printMatrix(double[][] m) { - for (double[] row : m) { - System.out.print("["); - for (int j = 0; j < row.length; j++) { - System.out.printf("%7.3f", row[j]); - if (j < row.length - 1) { + private static void printMatrix(double[][] matrix) { + for (double[] row : matrix) { + System.out.print("{ "); + for (int i = 0; i < row.length; i++) { + System.out.printf("%.3f", row[i]); + if (i < row.length - 1) { System.out.print(", "); } } - System.out.println("]"); + System.out.println(" }"); + } + } + + /** + * Helper method to multiply two matrices. + * + * @param a first matrix + * @param b second matrix + * @return the product matrix + */ + private static double[][] multiplyMatrices(double[][] a, double[][] b) { + int n = a.length; + double[][] result = new double[n][n]; + + for (int i = 0; i < n; i++) { + for (int j = 0; j < n; j++) { + for (int k = 0; k < n; k++) { + result[i][j] += a[i][k] * b[k][j]; + } + } + } + + return result; + } + + /** + * Result class to hold L and U matrices from decomposition. + */ + public static class Result { + private final double[][] lMatrix; + private final double[][] uMatrix; + + /** + * Constructor for Result. + * + * @param l Lower triangular matrix + * @param u Upper triangular matrix + */ + public Result(double[][] l, double[][] u) { + this.lMatrix = l; + this.uMatrix = u; + } + + /** + * Gets the lower triangular matrix. + * + * @return L matrix + */ + public double[][] getL() { + return lMatrix; + } + + /** + * Gets the upper triangular matrix. + * + * @return U matrix + */ + public double[][] getU() { + return uMatrix; } } } diff --git a/src/test/java/com/thealgorithms/matrix/LUDecompositionTest.java b/src/test/java/com/thealgorithms/matrix/LUDecompositionTest.java index d3cc6d64bf42..4f9f4dd31622 100644 --- a/src/test/java/com/thealgorithms/matrix/LUDecompositionTest.java +++ b/src/test/java/com/thealgorithms/matrix/LUDecompositionTest.java @@ -1,40 +1,107 @@ package com.thealgorithms.matrix; import static org.junit.jupiter.api.Assertions.assertArrayEquals; - +import static org.junit.jupiter.api.Assertions.assertThrows; import org.junit.jupiter.api.Test; -public class LUDecompositionTest { +class LUDecompositionTest { + + private static final double EPSILON = 1e-6; @Test - public void testLUDecomposition() { - double[][] a = {{4, 3}, {6, 3}}; + void testBasicLUDecomposition() { + double[][] matrix = { + {2, -1, -2}, + {-4, 6, 3}, + {-4, -2, 8} + }; - // Perform LU decomposition - LUDecomposition.LU lu = LUDecomposition.decompose(a); - double[][] l = lu.l; - double[][] u = lu.u; + LUDecomposition.Result result = LUDecomposition.decompose(matrix); - // Reconstruct a from l and u - double[][] reconstructed = multiplyMatrices(l, u); + double[][] expectedL = { + {1.0, 0.0, 0.0}, + {-2.0, 1.0, 0.0}, + {-2.0, -1.0, 1.0} + }; - // Assert that reconstructed matrix matches original a - for (int i = 0; i < a.length; i++) { - assertArrayEquals(a[i], reconstructed[i], 1e-9); - } + double[][] expectedU = { + {2.0, -1.0, -2.0}, + {0.0, 4.0, -1.0}, + {0.0, 0.0, 3.0} + }; + + assertMatrixEquals(expectedL, result.getL()); + assertMatrixEquals(expectedU, result.getU()); + } + + @Test + void testIdentityMatrix() { + double[][] identity = { + {1, 0, 0}, + {0, 1, 0}, + {0, 0, 1} + }; + + LUDecomposition.Result result = LUDecomposition.decompose(identity); + + assertMatrixEquals(identity, result.getL()); + assertMatrixEquals(identity, result.getU()); + } + + @Test + void testTwoByTwoMatrix() { + double[][] matrix = { + {4, 3}, + {6, 3} + }; + + LUDecomposition.Result result = LUDecomposition.decompose(matrix); + + double[][] expectedL = { + {1.0, 0.0}, + {1.5, 1.0} + }; + + double[][] expectedU = { + {4.0, 3.0}, + {0.0, -1.5} + }; + + assertMatrixEquals(expectedL, result.getL()); + assertMatrixEquals(expectedU, result.getU()); + } + + @Test + void testNonSquareMatrix() { + double[][] nonSquare = { + {1, 2, 3}, + {4, 5, 6} + }; + + assertThrows(IllegalArgumentException.class, () -> LUDecomposition.decompose(nonSquare)); + } + + @Test + void testEmptyMatrix() { + double[][] empty = {}; + + assertThrows(IllegalArgumentException.class, () -> LUDecomposition.decompose(empty)); + } + + @Test + void testSingularMatrix() { + double[][] singular = { + {1, 2, 3}, + {2, 4, 6}, + {3, 6, 9} + }; + + assertThrows(IllegalArgumentException.class, () -> LUDecomposition.decompose(singular)); } - // Helper method to multiply two matrices - private double[][] multiplyMatrices(double[][] a, double[][] b) { - int n = a.length; - double[][] c = new double[n][n]; - for (int i = 0; i < n; i++) { - for (int j = 0; j < n; j++) { - for (int k = 0; k < n; k++) { - c[i][j] += a[i][k] * b[k][j]; - } - } + private void assertMatrixEquals(double[][] expected, double[][] actual) { + for (int i = 0; i < expected.length; i++) { + assertArrayEquals(expected[i], actual[i], EPSILON); } - return c; } }