Java Code Examples for org.apache.commons.math.util.MathUtils#EPSILON

The following examples show how to use org.apache.commons.math.util.MathUtils#EPSILON . You can vote up the ones you like or vote down the ones you don't like, and go to the original project or source file by following the links above each example. You may check out the related API usage on the sidebar.
Example 1
Source File: 1_EigenDecompositionImpl.java    From SimFix with GNU General Public License v2.0 6 votes vote down vote up
/**
 * Check if a matrix is symmetric.
 * @param matrix matrix to check
 * @return true if matrix is symmetric
 */
private boolean isSymmetric(final RealMatrix matrix) {
    final int rows    = matrix.getRowDimension();
    final int columns = matrix.getColumnDimension();
    final double eps  = 10 * rows * columns * MathUtils.EPSILON;
    for (int i = 0; i < rows; ++i) {
        for (int j = i + 1; j < columns; ++j) {
            final double mij = matrix.getEntry(i, j);
            final double mji = matrix.getEntry(j, i);
            if (Math.abs(mij - mji) > (Math.max(Math.abs(mij), Math.abs(mji)) * eps)) {
                return false;
            }
        }
    }
    return true;
}
 
Example 2
Source File: JGenProg2017_0072_s.java    From coming with MIT License 6 votes vote down vote up
/**
 * Check if a matrix is symmetric.
 * @param matrix matrix to check
 * @return true if matrix is symmetric
 */
private boolean isSymmetric(final RealMatrix matrix) {
    final int rows    = matrix.getRowDimension();
    final int columns = matrix.getColumnDimension();
    final double eps  = 10 * rows * columns * MathUtils.EPSILON;
    for (int i = 0; i < rows; ++i) {
        for (int j = i + 1; j < columns; ++j) {
            final double mij = matrix.getEntry(i, j);
            final double mji = matrix.getEntry(j, i);
            if (Math.abs(mij - mji) > (Math.max(Math.abs(mij), Math.abs(mji)) * eps)) {
                return false;
            }
        }
    }
    return true;
}
 
Example 3
Source File: jMutRepair_0026_t.java    From coming with MIT License 6 votes vote down vote up
/**
 * Check if a matrix is symmetric.
 * @param matrix matrix to check
 * @return true if matrix is symmetric
 */
private boolean isSymmetric(final RealMatrix matrix) {
    final int rows    = matrix.getRowDimension();
    final int columns = matrix.getColumnDimension();
    final double eps  = 10 * rows * columns * MathUtils.EPSILON;
    for (int i = 0; i < rows; ++i) {
        for (int j = i + 1; j < columns; ++j) {
            final double mij = matrix.getEntry(i, j);
            final double mji = matrix.getEntry(j, i);
            if (Math.abs(mij - mji) > (Math.max(Math.abs(mij), Math.abs(mji)) * eps)) {
                return false;
            }
        }
    }
    return true;
}
 
Example 4
Source File: jMutRepair_0050_t.java    From coming with MIT License 6 votes vote down vote up
/**
 * Check if a matrix is symmetric.
 * @param matrix matrix to check
 * @return true if matrix is symmetric
 */
private boolean isSymmetric(final RealMatrix matrix) {
    final int rows    = matrix.getRowDimension();
    final int columns = matrix.getColumnDimension();
    final double eps  = 10 * rows * columns * MathUtils.EPSILON;
    for (int i = 0; i < rows; ++i) {
        for (int j = i + 1; j < columns; ++j) {
            final double mij = matrix.getEntry(i, j);
            final double mji = matrix.getEntry(j, i);
            if (Math.abs(mij - mji) > (Math.max(Math.abs(mij), Math.abs(mji)) * eps)) {
                return false;
            }
        }
    }
    return true;
}
 
Example 5
Source File: Arja_00109_t.java    From coming with MIT License 6 votes vote down vote up
/**
 * Check if a matrix is symmetric.
 * @param matrix matrix to check
 * @return true if matrix is symmetric
 */
private boolean isSymmetric(final RealMatrix matrix) {
    final int rows    = matrix.getRowDimension();
    final int columns = matrix.getColumnDimension();
    final double eps  = 10 * rows * columns * MathUtils.EPSILON;
    for (int i = 0; i < rows; ++i) {
        for (int j = i + 1; j < columns; ++j) {
            final double mij = matrix.getEntry(i, j);
            final double mji = matrix.getEntry(j, i);
            if (Math.abs(mij - mji) > (Math.max(Math.abs(mij), Math.abs(mji)) * eps)) {
                return false;
            }
        }
    }
    return true;
}
 
Example 6
Source File: 1_EigenDecompositionImpl.java    From SimFix with GNU General Public License v2.0 6 votes vote down vote up
/**
 * Check if a matrix is symmetric.
 * @param matrix matrix to check
 * @return true if matrix is symmetric
 */
private boolean isSymmetric(final RealMatrix matrix) {
    final int rows    = matrix.getRowDimension();
    final int columns = matrix.getColumnDimension();
    final double eps  = 10 * rows * columns * MathUtils.EPSILON;
    for (int i = 0; i < rows; ++i) {
        for (int j = i + 1; j < columns; ++j) {
            final double mij = matrix.getEntry(i, j);
            final double mji = matrix.getEntry(j, i);
            if (Math.abs(mij - mji) > (Math.max(Math.abs(mij), Math.abs(mji)) * eps)) {
                return false;
            }
        }
    }
    return true;
}
 
Example 7
Source File: EigenDecompositionImpl.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
/**
 * Check if a matrix is symmetric.
 * @param matrix matrix to check
 * @return true if matrix is symmetric
 */
private boolean isSymmetric(final RealMatrix matrix) {
    final int rows    = matrix.getRowDimension();
    final int columns = matrix.getColumnDimension();
    final double eps  = 10 * rows * columns * MathUtils.EPSILON;
    for (int i = 0; i < rows; ++i) {
        for (int j = i + 1; j < columns; ++j) {
            final double mij = matrix.getEntry(i, j);
            final double mji = matrix.getEntry(j, i);
            if (Math.abs(mij - mji) > (Math.max(Math.abs(mij), Math.abs(mji)) * eps)) {
                return false;
            }
        }
    }
    return true;
}
 
Example 8
Source File: Cardumen_00181_t.java    From coming with MIT License 6 votes vote down vote up
/**
 * Check if a matrix is symmetric.
 * @param matrix matrix to check
 * @return true if matrix is symmetric
 */
private boolean isSymmetric(final RealMatrix matrix) {
    final int rows    = matrix.getRowDimension();
    final int columns = matrix.getColumnDimension();
    final double eps  = 10 * rows * columns * MathUtils.EPSILON;
    for (int i = 0; i < rows; ++i) {
        for (int j = i + 1; j < columns; ++j) {
            final double mij = matrix.getEntry(i, j);
            final double mji = matrix.getEntry(j, i);
            if (Math.abs(mij - mji) > (Math.max(Math.abs(mij), Math.abs(mji)) * eps)) {
                return false;
            }
        }
    }
    return true;
}
 
Example 9
Source File: jKali_0027_t.java    From coming with MIT License 6 votes vote down vote up
/**
 * Check if a matrix is symmetric.
 * @param matrix matrix to check
 * @return true if matrix is symmetric
 */
private boolean isSymmetric(final RealMatrix matrix) {
    final int rows    = matrix.getRowDimension();
    final int columns = matrix.getColumnDimension();
    final double eps  = 10 * rows * columns * MathUtils.EPSILON;
    for (int i = 0; i < rows; ++i) {
        for (int j = i + 1; j < columns; ++j) {
            final double mij = matrix.getEntry(i, j);
            final double mji = matrix.getEntry(j, i);
            if (Math.abs(mij - mji) > (Math.max(Math.abs(mij), Math.abs(mji)) * eps)) {
                return false;
            }
        }
    }
    return true;
}
 
Example 10
Source File: EigenDecompositionImpl.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
/**
 * Check if a matrix is symmetric.
 * @param matrix
 *            matrix to check
 * @return true if matrix is symmetric
 */
private boolean isSymmetric(final RealMatrix matrix) {
    final int rows = matrix.getRowDimension();
    final int columns = matrix.getColumnDimension();
    final double eps = 10 * rows * columns * MathUtils.EPSILON;
    for (int i = 0; i < rows; ++i) {
        for (int j = i + 1; j < columns; ++j) {
            final double mij = matrix.getEntry(i, j);
            final double mji = matrix.getEntry(j, i);
            if (FastMath.abs(mij - mji) > (FastMath.max(FastMath.abs(mij), Math
                    .abs(mji)) * eps)) {
                return false;
            }
        }
    }
    return true;
}
 
Example 11
Source File: patch3-Math-81-jMutRepair_patch3-Math-81-jMutRepair_t.java    From coming with MIT License 4 votes vote down vote up
/**
 * Find the realEigenvalues.
 * @exception InvalidMatrixException if a block cannot be diagonalized
 */
private void findEigenvalues()
    throws InvalidMatrixException {

    // compute splitting points
    List<Integer> splitIndices = computeSplits();

    // find realEigenvalues in each block
    realEigenvalues = new double[main.length];
    imagEigenvalues = new double[main.length];
    int begin = 0;
    for (final int end : splitIndices) {
        final int n = end - begin;
        switch (n) {

        case 1:
            // apply dedicated method for dimension 1
            process1RowBlock(begin);
            break;

        case 2:
            // apply dedicated method for dimension 2
            process2RowsBlock(begin);
            break;

        case 3:
            // apply dedicated method for dimension 3
            process3RowsBlock(begin);
            break;

        default:

            // choose an initial shift for LDL<sup>T</sup> decomposition
            final double[] range       = eigenvaluesRange(begin, n);
            final double oneFourth     = 0.25 * (3 * range[0] + range[1]);
            final int oneFourthCount   = countEigenValues(oneFourth, begin, n);
            final double threeFourth   = 0.25 * (range[0] + 3 * range[1]);
            final int threeFourthCount = countEigenValues(threeFourth, begin, n);
            final boolean chooseLeft   = (oneFourthCount - 1) >= (n - threeFourthCount);
            final double lambda        = chooseLeft ? range[0] : range[1];

            tau = (range[1] - range[0]) * MathUtils.EPSILON * n + 2 * minPivot;

            // decompose T-&lambda;I as LDL<sup>T</sup>
            ldlTDecomposition(lambda, begin, n);

            // apply general dqd/dqds method
            processGeneralBlock(n);

            // extract realEigenvalues
            if (chooseLeft) {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda + work[4 * i];
                }
            } else {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda - work[4 * i];
                }
            }

        }
        begin = end;
    }

    // sort the realEigenvalues in decreasing order
    Arrays.sort(realEigenvalues);
    int j = realEigenvalues.length - 1;
    for (int i = 0; i < j; ++i) {
        final double tmp = realEigenvalues[i];
        realEigenvalues[i] = realEigenvalues[j];
        realEigenvalues[j] = tmp;
        --j;
    }

}
 
Example 12
Source File: jKali_0013_s.java    From coming with MIT License 4 votes vote down vote up
/**
 * Find the realEigenvalues.
 * @exception InvalidMatrixException if a block cannot be diagonalized
 */
private void findEigenvalues()
    throws InvalidMatrixException {

    // compute splitting points
    List<Integer> splitIndices = computeSplits();

    // find realEigenvalues in each block
    realEigenvalues = new double[main.length];
    imagEigenvalues = new double[main.length];
    int begin = 0;
    for (final int end : splitIndices) {
        final int n = end - begin;
        switch (n) {

        case 1:
            // apply dedicated method for dimension 1
            process1RowBlock(begin);
            break;

        case 2:
            // apply dedicated method for dimension 2
            process2RowsBlock(begin);
            break;

        case 3:
            // apply dedicated method for dimension 3
            process3RowsBlock(begin);
            break;

        default:

            // choose an initial shift for LDL<sup>T</sup> decomposition
            final double[] range       = eigenvaluesRange(begin, n);
            final double oneFourth     = 0.25 * (3 * range[0] + range[1]);
            final int oneFourthCount   = countEigenValues(oneFourth, begin, n);
            final double threeFourth   = 0.25 * (range[0] + 3 * range[1]);
            final int threeFourthCount = countEigenValues(threeFourth, begin, n);
            final boolean chooseLeft   = (oneFourthCount - 1) >= (n - threeFourthCount);
            final double lambda        = chooseLeft ? range[0] : range[1];

            tau = (range[1] - range[0]) * MathUtils.EPSILON * n + 2 * minPivot;

            // decompose T-&lambda;I as LDL<sup>T</sup>
            ldlTDecomposition(lambda, begin, n);

            // apply general dqd/dqds method
            processGeneralBlock(n);

            // extract realEigenvalues
            if (chooseLeft) {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda + work[4 * i];
                }
            } else {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda - work[4 * i];
                }
            }

        }
        begin = end;
    }

    // sort the realEigenvalues in decreasing order
    Arrays.sort(realEigenvalues);
    int j = realEigenvalues.length - 1;
    for (int i = 0; i < j; ++i) {
        final double tmp = realEigenvalues[i];
        realEigenvalues[i] = realEigenvalues[j];
        realEigenvalues[j] = tmp;
        --j;
    }

}
 
Example 13
Source File: Math_80_EigenDecompositionImpl_t.java    From coming with MIT License 4 votes vote down vote up
/**
 * Find the realEigenvalues.
 * @exception InvalidMatrixException if a block cannot be diagonalized
 */
private void findEigenvalues()
    throws InvalidMatrixException {

    // compute splitting points
    List<Integer> splitIndices = computeSplits();

    // find realEigenvalues in each block
    realEigenvalues = new double[main.length];
    imagEigenvalues = new double[main.length];
    int begin = 0;
    for (final int end : splitIndices) {
        final int n = end - begin;
        switch (n) {

        case 1:
            // apply dedicated method for dimension 1
            process1RowBlock(begin);
            break;

        case 2:
            // apply dedicated method for dimension 2
            process2RowsBlock(begin);
            break;

        case 3:
            // apply dedicated method for dimension 3
            process3RowsBlock(begin);
            break;

        default:

            // choose an initial shift for LDL<sup>T</sup> decomposition
            final double[] range       = eigenvaluesRange(begin, n);
            final double oneFourth     = 0.25 * (3 * range[0] + range[1]);
            final int oneFourthCount   = countEigenValues(oneFourth, begin, n);
            final double threeFourth   = 0.25 * (range[0] + 3 * range[1]);
            final int threeFourthCount = countEigenValues(threeFourth, begin, n);
            final boolean chooseLeft   = (oneFourthCount - 1) >= (n - threeFourthCount);
            final double lambda        = chooseLeft ? range[0] : range[1];

            tau = (range[1] - range[0]) * MathUtils.EPSILON * n + 2 * minPivot;

            // decompose T-&lambda;I as LDL<sup>T</sup>
            ldlTDecomposition(lambda, begin, n);

            // apply general dqd/dqds method
            processGeneralBlock(n);

            // extract realEigenvalues
            if (chooseLeft) {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda + work[4 * i];
                }
            } else {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda - work[4 * i];
                }
            }

        }
        begin = end;
    }

    // sort the realEigenvalues in decreasing order
    Arrays.sort(realEigenvalues);
    int j = realEigenvalues.length - 1;
    for (int i = 0; i < j; ++i) {
        final double tmp = realEigenvalues[i];
        realEigenvalues[i] = realEigenvalues[j];
        realEigenvalues[j] = tmp;
        --j;
    }

}
 
Example 14
Source File: jKali_0014_s.java    From coming with MIT License 4 votes vote down vote up
/**
 * Find the realEigenvalues.
 * @exception InvalidMatrixException if a block cannot be diagonalized
 */
private void findEigenvalues()
    throws InvalidMatrixException {

    // compute splitting points
    List<Integer> splitIndices = computeSplits();

    // find realEigenvalues in each block
    realEigenvalues = new double[main.length];
    imagEigenvalues = new double[main.length];
    int begin = 0;
    for (final int end : splitIndices) {
        final int n = end - begin;
        switch (n) {

        case 1:
            // apply dedicated method for dimension 1
            process1RowBlock(begin);
            break;

        case 2:
            // apply dedicated method for dimension 2
            process2RowsBlock(begin);
            break;

        case 3:
            // apply dedicated method for dimension 3
            process3RowsBlock(begin);
            break;

        default:

            // choose an initial shift for LDL<sup>T</sup> decomposition
            final double[] range       = eigenvaluesRange(begin, n);
            final double oneFourth     = 0.25 * (3 * range[0] + range[1]);
            final int oneFourthCount   = countEigenValues(oneFourth, begin, n);
            final double threeFourth   = 0.25 * (range[0] + 3 * range[1]);
            final int threeFourthCount = countEigenValues(threeFourth, begin, n);
            final boolean chooseLeft   = (oneFourthCount - 1) >= (n - threeFourthCount);
            final double lambda        = chooseLeft ? range[0] : range[1];

            tau = (range[1] - range[0]) * MathUtils.EPSILON * n + 2 * minPivot;

            // decompose T-&lambda;I as LDL<sup>T</sup>
            ldlTDecomposition(lambda, begin, n);

            // apply general dqd/dqds method
            processGeneralBlock(n);

            // extract realEigenvalues
            if (chooseLeft) {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda + work[4 * i];
                }
            } else {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda - work[4 * i];
                }
            }

        }
        begin = end;
    }

    // sort the realEigenvalues in decreasing order
    Arrays.sort(realEigenvalues);
    int j = realEigenvalues.length - 1;
    for (int i = 0; i < j; ++i) {
        final double tmp = realEigenvalues[i];
        realEigenvalues[i] = realEigenvalues[j];
        realEigenvalues[j] = tmp;
        --j;
    }

}
 
Example 15
Source File: Arja_00145_t.java    From coming with MIT License 4 votes vote down vote up
/**
 * Find the realEigenvalues.
 * @exception InvalidMatrixException if a block cannot be diagonalized
 */
private void findEigenvalues()
    throws InvalidMatrixException {

    // compute splitting points
    List<Integer> splitIndices = computeSplits();

    // find realEigenvalues in each block
    realEigenvalues = new double[main.length];
    imagEigenvalues = new double[main.length];
    int begin = 0;
    for (final int end : splitIndices) {
        final int n = end - begin;
        switch (n) {

        case 1:
            // apply dedicated method for dimension 1
            process1RowBlock(begin);
            break;

        case 2:
            // apply dedicated method for dimension 2
            process2RowsBlock(begin);
            break;

        case 3:
            // apply dedicated method for dimension 3
            process3RowsBlock(begin);
            break;

        default:

            // choose an initial shift for LDL<sup>T</sup> decomposition
            final double[] range       = eigenvaluesRange(begin, n);
            final double oneFourth     = 0.25 * (3 * range[0] + range[1]);
            final int oneFourthCount   = countEigenValues(oneFourth, begin, n);
            final double threeFourth   = 0.25 * (range[0] + 3 * range[1]);
            final int threeFourthCount = countEigenValues(threeFourth, begin, n);
            final boolean chooseLeft   = (oneFourthCount - 1) >= (n - threeFourthCount);
            final double lambda        = chooseLeft ? range[0] : range[1];

            tau = (range[1] - range[0]) * MathUtils.EPSILON * n + 2 * minPivot;

            // decompose T-&lambda;I as LDL<sup>T</sup>
            ldlTDecomposition(lambda, begin, n);

            // apply general dqd/dqds method
            processGeneralBlock(n);

            // extract realEigenvalues
            if (chooseLeft) {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda + work[4 * i];
                }
            } else {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda - work[4 * i];
                }
            }

        }
        begin = end;
    }

    // sort the realEigenvalues in decreasing order
    Arrays.sort(realEigenvalues);
    int j = realEigenvalues.length - 1;
    for (int i = 0; i < j; ++i) {
        final double tmp = realEigenvalues[i];
        realEigenvalues[i] = realEigenvalues[j];
        realEigenvalues[j] = tmp;
        --j;
    }

}
 
Example 16
Source File: jKali_0027_s.java    From coming with MIT License 4 votes vote down vote up
/**
 * Find the realEigenvalues.
 * @exception InvalidMatrixException if a block cannot be diagonalized
 */
private void findEigenvalues()
    throws InvalidMatrixException {

    // compute splitting points
    List<Integer> splitIndices = computeSplits();

    // find realEigenvalues in each block
    realEigenvalues = new double[main.length];
    imagEigenvalues = new double[main.length];
    int begin = 0;
    for (final int end : splitIndices) {
        final int n = end - begin;
        switch (n) {

        case 1:
            // apply dedicated method for dimension 1
            process1RowBlock(begin);
            break;

        case 2:
            // apply dedicated method for dimension 2
            process2RowsBlock(begin);
            break;

        case 3:
            // apply dedicated method for dimension 3
            process3RowsBlock(begin);
            break;

        default:

            // choose an initial shift for LDL<sup>T</sup> decomposition
            final double[] range       = eigenvaluesRange(begin, n);
            final double oneFourth     = 0.25 * (3 * range[0] + range[1]);
            final int oneFourthCount   = countEigenValues(oneFourth, begin, n);
            final double threeFourth   = 0.25 * (range[0] + 3 * range[1]);
            final int threeFourthCount = countEigenValues(threeFourth, begin, n);
            final boolean chooseLeft   = (oneFourthCount - 1) >= (n - threeFourthCount);
            final double lambda        = chooseLeft ? range[0] : range[1];

            tau = (range[1] - range[0]) * MathUtils.EPSILON * n + 2 * minPivot;

            // decompose T-&lambda;I as LDL<sup>T</sup>
            ldlTDecomposition(lambda, begin, n);

            // apply general dqd/dqds method
            processGeneralBlock(n);

            // extract realEigenvalues
            if (chooseLeft) {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda + work[4 * i];
                }
            } else {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda - work[4 * i];
                }
            }

        }
        begin = end;
    }

    // sort the realEigenvalues in decreasing order
    Arrays.sort(realEigenvalues);
    int j = realEigenvalues.length - 1;
    for (int i = 0; i < j; ++i) {
        final double tmp = realEigenvalues[i];
        realEigenvalues[i] = realEigenvalues[j];
        realEigenvalues[j] = tmp;
        --j;
    }

}
 
Example 17
Source File: Arja_00159_t.java    From coming with MIT License 4 votes vote down vote up
/**
 * Find the realEigenvalues.
 * @exception InvalidMatrixException if a block cannot be diagonalized
 */
private void findEigenvalues()
    throws InvalidMatrixException {

    // compute splitting points
    List<Integer> splitIndices = computeSplits();

    // find realEigenvalues in each block
    realEigenvalues = new double[main.length];
    imagEigenvalues = new double[main.length];
    int begin = 0;
    for (final int end : splitIndices) {
        final int n = end - begin;
        switch (n) {

        case 1:
            // apply dedicated method for dimension 1
            process1RowBlock(begin);
            break;

        case 2:
            // apply dedicated method for dimension 2
            process2RowsBlock(begin);
            break;

        case 3:
            // apply dedicated method for dimension 3
            process3RowsBlock(begin);
            break;

        default:

            // choose an initial shift for LDL<sup>T</sup> decomposition
            final double[] range       = eigenvaluesRange(begin, n);
            final double oneFourth     = 0.25 * (3 * range[0] + range[1]);
            final int oneFourthCount   = countEigenValues(oneFourth, begin, n);
            final double threeFourth   = 0.25 * (range[0] + 3 * range[1]);
            final int threeFourthCount = countEigenValues(threeFourth, begin, n);
            final boolean chooseLeft   = (oneFourthCount - 1) >= (n - threeFourthCount);
            final double lambda        = chooseLeft ? range[0] : range[1];

            tau = (range[1] - range[0]) * MathUtils.EPSILON * n + 2 * minPivot;

            // decompose T-&lambda;I as LDL<sup>T</sup>
            ldlTDecomposition(lambda, begin, n);

            // apply general dqd/dqds method
            processGeneralBlock(n);

            // extract realEigenvalues
            if (chooseLeft) {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda + work[4 * i];
                }
            } else {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda - work[4 * i];
                }
            }

        }
        begin = end;
    }

    // sort the realEigenvalues in decreasing order
    Arrays.sort(realEigenvalues);
    int j = realEigenvalues.length - 1;
    for (int i = 0; i < j; ++i) {
        final double tmp = realEigenvalues[i];
        realEigenvalues[i] = realEigenvalues[j];
        realEigenvalues[j] = tmp;
        --j;
    }

}
 
Example 18
Source File: 1_EigenDecompositionImpl.java    From SimFix with GNU General Public License v2.0 4 votes vote down vote up
/**
 * Find the realEigenvalues.
 * @exception InvalidMatrixException if a block cannot be diagonalized
 */
private void findEigenvalues()
    throws InvalidMatrixException {

    // compute splitting points
    List<Integer> splitIndices = computeSplits();

    // find realEigenvalues in each block
    realEigenvalues = new double[main.length];
    imagEigenvalues = new double[main.length];
    int begin = 0;
    for (final int end : splitIndices) {
        final int n = end - begin;
        switch (n) {

        case 1:
            // apply dedicated method for dimension 1
            process1RowBlock(begin);
            break;

        case 2:
            // apply dedicated method for dimension 2
            process2RowsBlock(begin);
            break;

        case 3:
            // apply dedicated method for dimension 3
            process3RowsBlock(begin);
            break;

        default:

            // choose an initial shift for LDL<sup>T</sup> decomposition
            final double[] range       = eigenvaluesRange(begin, n);
            final double oneFourth     = 0.25 * (3 * range[0] + range[1]);
            final int oneFourthCount   = countEigenValues(oneFourth, begin, n);
            final double threeFourth   = 0.25 * (range[0] + 3 * range[1]);
            final int threeFourthCount = countEigenValues(threeFourth, begin, n);
            final boolean chooseLeft   = (oneFourthCount - 1) >= (n - threeFourthCount);
            final double lambda        = chooseLeft ? range[0] : range[1];

            tau = (range[1] - range[0]) * MathUtils.EPSILON * n + 2 * minPivot;

            // decompose T-&lambda;I as LDL<sup>T</sup>
            ldlTDecomposition(lambda, begin, n);

            // apply general dqd/dqds method
            processGeneralBlock(n);

            // extract realEigenvalues
            if (chooseLeft) {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda + work[4 * i];
                }
            } else {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda - work[4 * i];
                }
            }

        }
        begin = end;
    }

    // sort the realEigenvalues in decreasing order
    Arrays.sort(realEigenvalues);
    int j = realEigenvalues.length - 1;
    for (int i = 0; i < j; ++i) {
        final double tmp = realEigenvalues[i];
        realEigenvalues[i] = realEigenvalues[j];
        realEigenvalues[j] = tmp;
        --j;
    }

}
 
Example 19
Source File: MillerUpdatingRegression.java    From astor with GNU General Public License v2.0 4 votes vote down vote up
/**
 * The include method is where the QR decomposition occurs. This statement forms all
 * intermediate data which will be used for all derivative measures.
 * According to the miller paper, note that in the original implementation the x vector
 * is overwritten. In this implementation, the include method is passed a copy of the
 * original data vector so that there is no contamination of the data. Additionally,
 * this method differs slightly from Gentleman's method, in that the assumption is
 * of dense design matrices, there is some advantage in using the original gentleman algorithm
 * on sparse matrices.
 *
 * @param x observations on the regressors
 * @param wi weight of the this observation (-1,1)
 * @param yi observation on the regressand
 */
private void include(final double[] x, final double wi, final double yi) {
    int nextr = 0;
    double w = wi;
    double y = yi;
    double xi;
    double di;
    double wxi;
    double dpi;
    double xk;
    double _w;
    this.rss_set = false;
    sumy = smartAdd(yi, sumy);
    sumsqy = smartAdd(sumsqy, yi * yi);
    for (int i = 0; i < x.length; i++) {
        if (w == 0.0) {
            return;
        }
        xi = x[i];

        if (xi == 0.0) {
            nextr += nvars - i - 1;
            continue;
        }
        di = d[i];
        wxi = w * xi;
        _w = w;
        if (di != 0.0) {
            dpi = smartAdd(di, wxi * xi);
            double tmp = wxi * xi / di;
            if (FastMath.abs(tmp) > MathUtils.EPSILON) {
                w = (di * w) / dpi;
            }
        } else {
            dpi = wxi * xi;
            w = 0.0;
        }
        d[i] = dpi;
        for (int k = i + 1; k < nvars; k++) {
            xk = x[k];
            x[k] = smartAdd(xk, -xi * r[nextr]);
            if (di != 0.0) {
                r[nextr] = smartAdd(di * r[nextr], (_w * xi) * xk) / dpi;
            } else {
                r[nextr] = xk / xi;
            }
            ++nextr;
        }
        xk = y;
        y = smartAdd(xk, -xi * rhs[i]);
        if (di != 0.0) {
            rhs[i] = smartAdd(di * rhs[i], wxi * xk) / dpi;
        } else {
            rhs[i] = xk / xi;
        }
    }
    sserr = smartAdd(sserr, w * y * y);
    return;
}
 
Example 20
Source File: JGenProg2017_0097_t.java    From coming with MIT License 4 votes vote down vote up
/**
 * Find the realEigenvalues.
 * @exception InvalidMatrixException if a block cannot be diagonalized
 */
private void findEigenvalues()
    throws InvalidMatrixException {

    // compute splitting points
    List<Integer> splitIndices = computeSplits();

    // find realEigenvalues in each block
    realEigenvalues = new double[main.length];
    imagEigenvalues = new double[main.length];
    int begin = 0;
    for (final int end : splitIndices) {
        final int n = end - begin;
        switch (n) {

        case 1:
            // apply dedicated method for dimension 1
            process1RowBlock(begin);
            break;

        case 2:
            // apply dedicated method for dimension 2
            process2RowsBlock(begin);
            break;

        case 3:
            // apply dedicated method for dimension 3
            process3RowsBlock(begin);
            break;

        default:

            // choose an initial shift for LDL<sup>T</sup> decomposition
            final double[] range       = eigenvaluesRange(begin, n);
            final double oneFourth     = 0.25 * (3 * range[0] + range[1]);
            final int oneFourthCount   = countEigenValues(oneFourth, begin, n);
            final double threeFourth   = 0.25 * (range[0] + 3 * range[1]);
            final int threeFourthCount = countEigenValues(threeFourth, begin, n);
            final boolean chooseLeft   = (oneFourthCount - 1) >= (n - threeFourthCount);
            final double lambda        = chooseLeft ? range[0] : range[1];

            tau = (range[1] - range[0]) * MathUtils.EPSILON * n + 2 * minPivot;

            // decompose T-&lambda;I as LDL<sup>T</sup>
            ldlTDecomposition(lambda, begin, n);

            // apply general dqd/dqds method
            processGeneralBlock(n);

            // extract realEigenvalues
            if (chooseLeft) {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda + work[4 * i];
                }
            } else {
                for (int i = 0; i < n; ++i) {
                    realEigenvalues[begin + i] = lambda - work[4 * i];
                }
            }

        }
        begin = end;
    }

    // sort the realEigenvalues in decreasing order
    Arrays.sort(realEigenvalues);
    int j = realEigenvalues.length - 1;
    for (int i = 0; i < j; ++i) {
        final double tmp = realEigenvalues[i];
        realEigenvalues[i] = realEigenvalues[j];
        realEigenvalues[j] = tmp;
        --j;
    }

}