Package mikera.matrixx.decompose.impl.hessenberg

Source Code of mikera.matrixx.decompose.impl.hessenberg.TestHessenbergSimilarDecomposition

/*
* Copyright (c) 2009-2013, 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 mikera.matrixx.decompose.impl.hessenberg;

import mikera.matrixx.Matrix;
import mikera.matrixx.algo.Multiplications;
import mikera.vectorz.Vector;

import org.junit.Test;

import static org.junit.Assert.assertTrue;

public class TestHessenbergSimilarDecomposition {

  /**
     * Default tolerance.
     */
    private static double TOLERANCE = 1e-8;


    /**
     * Decomposes the matrix, extracts H and Q, then sees if it can recompute A using similar matrix stuff.
     */
    @Test
    public void testItAllTogether() {
        Matrix A = Matrix.createRandom(5,5);

        checkItAll(A);
    }

    private void checkItAll(Matrix A) {
        HessenbergResult result = HessenbergSimilarDecomposition.decompose(A);
        Matrix Q = result.getQ().toMatrix();
        Matrix H = result.getH().toMatrix();
//        System.out.println("-------- H ---------");
//        UtilEjml.print(H,"%8.2e");
//        System.out.println("-------- Q ---------");
//        UtilEjml.print(Q,"%8.2e");
        assertTrue(isOrthogonal(Q, TOLERANCE));

        H = Multiplications.multiply(Q, Multiplications.multiply(H, Q.getTranspose()));

//        System.out.println("------- A ----------");
//        UtilEjml.print(A,"%8.2e");
//        System.out.println("----- Found A ------");
//        UtilEjml.print(H,"%8.2e");

        assertTrue(!H.hasUncountable());

        assertTrue(A.epsilonEquals(H,TOLERANCE));
    }

    /**
     * Make sure it doesn't change the input
     */
    @Test
    public void testInputUnmodified() {
        Matrix A = Matrix.createRandom(4,4);
        Matrix B = A.copy();
        HessenbergSimilarDecomposition.decompose(A);
        assertTrue(A.equals(B));
    }

    /**
     * Give it a matrix that is already a Hessenberg matrix and see if its comes out the same.
     */
//    @Test
//    public void testNoChange() {
//        DenseMatrix64F A = RandomMatrices.createUpperTriangle(4,1,-1,1,rand);
//
//        HessenbergSimilarDecomposition decomp = new HessenbergSimilarDecomposition(A.numRows);
//
//        assertTrue(decomp.decompose(A));
//
//        DenseMatrix64F H = decomp.getH(null);
//
//        assertTrue(MatrixFeatures.isIdentical(A,H,0));
//    }

//    /**
//     * This checks to see the gammas and if the householder vectors stored in QH are correct. This
//     * is done by extracting the vectors, computing reflectors, and multipling them by A and seeing
//     * if it has the expected response.
//     */
//    @Test
//    public void testHouseholderVectors()
//    {
//        int N = 5;
//        Matrix A = Matrix.createRandom(N,N);
//        Matrix B = Matrix.create(N,N);
//
//        HessenbergSimilarDecomposition decomp = new HessenbergSimilarDecomposition();
//        HessenbergResult result = decomp.decompose(A);
//       
//       
//        Matrix QH = decomp.getQH().toMatrix();
////        System.out.println("------------ QH -----------");
////        UtilEjml.print(QH);
//
//        double gammas[] = decomp.getGammas();
//
//        DenseMatrix64F u = new DenseMatrix64F(N,1);
//
////        UtilEjml.print(A);
////        System.out.println("-------------------");
//
//        for( int i = 0; i < N-1; i++ ) {
//            u.zero();
//            u.data[i+1] = 1;
//            for( int j = i+2; j < N; j++ ) {
//                u.data[j] = QH.get(j,i);
//            }
//
//            DenseMatrix64F Q = SpecializedOps.createReflector(u,gammas[i]);
//            CommonOps.mult(Q,A,B);
////            System.out.println("----- u ------");
////            UtilEjml.print(u);
////            System.out.println("----- Q ------");
////            UtilEjml.print(Q);
////            System.out.println("----- B ------");
////            UtilEjml.print(B);
//
//            for( int j = 0; j < i+2; j++ ) {
//                assertTrue(Math.abs(B.get(j,i))>UtilEjml.TOLERANCE);
//            }
//            for( int j = i+2; j < N; j++ ) {
//                assertEquals(0,B.get(j,i),UtilEjml.TOLERANCE);
//            }
//            CommonOps.mult(B,Q,A);
//
////            System.out.println("-------------------");
////            UtilEjml.print(A);
////            System.out.println("-------------------");
//        }
//    }

//    /**
//     * Compute the overall Q matrix from the stored u vectors.  See if the extract H is the same as the expected H.
//     */
//    @Test
//    public void testH() {
//        int N = 5;
//        DenseMatrix64F A = RandomMatrices.createRandom(N,N,rand);
//
//        HessenbergSimilarDecomposition decomp = new HessenbergSimilarDecomposition(A.numRows);
//
//        assertTrue(safeDecomposition(decomp,A));
//
//        DenseMatrix64F QH = decomp.getQH();
//
//        double gammas[] = decomp.getGammas();
//
//        DenseMatrix64F u = new DenseMatrix64F(N,1);
//
//
//        DenseMatrix64F Q = CommonOps.identity(N);
//        DenseMatrix64F temp = new DenseMatrix64F(N,N);
//
//        for( int i = N-2; i >= 0; i-- ) {
//            u.zero();
//            u.data[i+1] = 1;
//            for( int j = i+2; j < N; j++ ) {
//                u.data[j] = QH.get(j,i);
//            }
//
//            DenseMatrix64F Qi = SpecializedOps.createReflector(u,gammas[i]);
//
//            CommonOps.mult(Qi,Q,temp);
//            Q.set(temp);
//        }
//        DenseMatrix64F expectedH = new DenseMatrix64F(N,N);
//
//        CommonOps.multTransA(Q,A,temp);
//        CommonOps.mult(temp,Q,expectedH);
//
////        UtilEjml.print(expectedH);
//
//        DenseMatrix64F foundH = decomp.getH(null);
//
////        UtilEjml.print(foundH);
//
//        assertTrue(MatrixFeatures.isIdentical(expectedH,foundH,UtilEjml.TOLERANCE));
//
//        System.out.println();
//    }
   
    /**
     * <p>
     * Checks to see if a matrix is orthogonal or isometric.
     * </p>
     *
     * @param Q The matrix being tested. Not modified.
     * @param tol Tolerance.
     * @return True if it passes the test.
     */
    private boolean isOrthogonal( Matrix Q , double tol ) {
       if( Q.rowCount() < Q.columnCount() ) {
            throw new IllegalArgumentException("The number of rows must be more than or equal to the number of columns");
        }

        Vector u[] = columnsToVector(Q);

        for( int i = 0; i < u.length; i++ ) {
            Vector a = u[i];

            for( int j = i+1; j < u.length; j++ ) {
//                double val = VectorVectorMult.innerProd(a,u[j]);
                double val = a.innerProduct(u[j]).get();
                if( !(Math.abs(val) <= tol))
                    return false;
            }
        }

        return true;
    }
   
    /**
     * Converts the columns in a matrix into a set of vectors.
     *
     * @param A Matrix.  Not modified.
     * @param v
     * @return An array of vectors.
     */
    public Vector[] columnsToVector(Matrix A) {
      Vector[] ret = new Vector[A.columnCount()];
        for( int i = 0; i < ret.length; i++ ) {
          ret[i] = Vector.createLength(A.rowCount());

            Vector u = ret[i];

            for( int j = 0; j < A.rowCount(); j++ ) {
                u.set(j, A.get(j,i));
            }
        }

        return ret;
    }
}
TOP

Related Classes of mikera.matrixx.decompose.impl.hessenberg.TestHessenbergSimilarDecomposition

TOP
Copyright © 2018 www.massapi.com. All rights reserved.
All source code are property of their respective owners. Java is a trademark of Sun Microsystems, Inc and owned by ORACLE Inc. Contact coftware#gmail.com.