/*
 * Cholesky.java by Richard J. Davies
 * from `Introductory Java for Scientists and Engineers'
 * chapter: `JNL, a Numerical Library for Java'
 * section: `Example - Using DoubleCholesky'
 *
 * We use the class DoubleCholesky to perform a Cholesky
 * decomposition of a symmetric positive definite matrix.
 */
import VisualNumerics.math.*;


public class Cholesky
{
  // A utility method to print a matrix
  public static void printMatrix(double[][] m)
  {
    for (int i=0; i<3; i++)
    {
      // Print each row, elements separated by tabs
      for (int j=0; j<3; j++)
        System.out.print(m[i][j] + "\t");

      // Start a new line at the row's end
      System.out.println();
    }

    // Leave a gap after the entire matrix
    System.out.println();
  }


  // The `main' method contains the body of the program
  public static void main(String[] argv)
  {
    // Define the matrix that we are going to operate on.
    double[][] A = { { 1.3, 1.2, 0.1 },
                     { 1.2, 3.4, 1.1 },
                     { 0.1, 1.1, 3.1} };

    // Print it out
    System.out.println("A is the matrix");
    printMatrix(A);

    try
    {
      // Invert the matrix and print the inverse.
      System.out.println("The inverse of A is");
      printMatrix(DoubleMatrix.inverse(A));;
    }
    catch(MathException e)
    {
      // Error message if inverse fails.
      System.out.println("A is singular");
    }

    try
    {
      // Perform the Cholesky decomposition
      DoubleCholesky decomp = new DoubleCholesky(A);

      // Print the Cholesky factor
      System.out.println("The Cholesky factor of A is");
      printMatrix(decomp.R());

      // Invert the matrix using the decomposition
      System.out.println("The inverse of A via Cholesky is");
      printMatrix(decomp.inverse());
    }
    catch(MathException e)
    {
      // Error message if decomposition fails
      System.out.println("That matrix is not symmetric "
                         + "positive definite");
    }
  }
}