Example usage for org.apache.commons.math.linear LUDecompositionImpl LUDecompositionImpl

List of usage examples for org.apache.commons.math.linear LUDecompositionImpl LUDecompositionImpl

Introduction

In this page you can find the example usage for org.apache.commons.math.linear LUDecompositionImpl LUDecompositionImpl.

Prototype

public LUDecompositionImpl(RealMatrix matrix) throws InvalidMatrixException 

Source Link

Document

Calculates the LU-decomposition of the given matrix.

Usage

From source file:com.zinnia.nectar.util.math.MatrixSolver.java

public static double[] solveMatrix(double[][] coefficientMatrix, double[] rhsVector) {
    RealVector x = null;//from   w w  w .  j  a  va2 s. c o m
    try {
        RealMatrix a = new Array2DRowRealMatrix(coefficientMatrix);
        RealVector b = new ArrayRealVector(rhsVector);
        DecompositionSolver solver = new LUDecompositionImpl(a).getSolver();
        x = solver.solve(b);
    } catch (Exception e) {
        e.printStackTrace();
    }

    double[] result;
    result = x.toArray();

    return result;
}

From source file:com.opengamma.analytics.math.linearalgebra.LUDecompositionCommons.java

/**
 * {@inheritDoc}//ww  w.  j ava2 s.com
 */
@Override
public LUDecompositionResult evaluate(final DoubleMatrix2D x) {
    Validate.notNull(x);
    final RealMatrix temp = CommonsMathWrapper.wrap(x);
    final LUDecomposition lu = new LUDecompositionImpl(temp);
    return new LUDecompositionCommonsResult(lu);
}

From source file:com.opengamma.analytics.math.matrix.CommonsMatrixAlgebra.java

/**
 * {@inheritDoc}/*from  ww  w  .j  ava2 s  .  c om*/
 */
@Override
public double getDeterminant(final Matrix<?> m) {
    Validate.notNull(m, "m");
    if (m instanceof DoubleMatrix2D) {
        final RealMatrix temp = CommonsMathWrapper.wrap((DoubleMatrix2D) m);
        final LUDecomposition lud = new LUDecompositionImpl(temp);
        return lud.getDeterminant();
    }
    throw new IllegalArgumentException("Can only find determinant of DoubleMatrix2D; have " + m.getClass());
}

From source file:adams.data.utils.SavitzkyGolay.java

/**
 * Determines the coefficients for the smoothing, with optional debugging
 * output./*from  www.j  a v a 2 s .  c o  m*/
 *
 * @param numLeft   the number of points to the left
 * @param numRight   the number of points to the right
 * @param polyOrder   the polynomial order
 * @param derOrder   the derivative order
 * @param debug   whether to output debugging information
 * @return      the coefficients
 */
public static double[] determineCoefficients(int numLeft, int numRight, int polyOrder, int derOrder,
        boolean debug) {
    double[] result;
    RealMatrix A;
    int i;
    int j;
    int k;
    float sum;
    RealMatrix b;
    LUDecomposition lu;
    RealMatrix solution;

    result = new double[numLeft + numRight + 1];

    // no window?
    if (result.length == 1) {
        result[0] = 1.0;
        return result;
    }

    // Note: "^" = superscript, "." = subscript

    // {A^T*A}.ij = Sum[k:-nl..nr](k^(i+j))
    A = new Array2DRowRealMatrix(polyOrder + 1, polyOrder + 1);
    for (i = 0; i < A.getRowDimension(); i++) {
        for (j = 0; j < A.getColumnDimension(); j++) {
            sum = 0;
            for (k = -numLeft; k <= numRight; k++)
                sum += Math.pow(k, i + j);
            A.setEntry(i, j, sum);
        }
    }
    if (debug)
        System.out.println("A:\n" + A);

    // LU decomp for inverse matrix
    b = new Array2DRowRealMatrix(polyOrder + 1, 1);
    b.setEntry(derOrder, 0, 1.0);
    if (debug)
        System.out.println("b:\n" + b);

    lu = new LUDecompositionImpl(A);
    solution = lu.getSolver().solve(b);
    if (debug)
        System.out.println("LU decomp. - solution:\n" + solution);

    // coefficients: c.n = Sum[m:0..M]((A^T*A)^-1).0m * n^m with n=-nl..nr
    for (i = -numLeft; i <= numRight; i++) {
        sum = 0;
        for (j = 0; j <= polyOrder; j++)
            sum += solution.getEntry(j, 0) * Math.pow(i, j);
        result[i + numLeft] = sum;
    }
    if (debug)
        System.out.println("Coefficients:\n" + Utils.arrayToString(result));

    return result;
}

From source file:edu.cudenver.bios.matrixsvc.resource.MatrixInversionResource.java

/**
 * Handle POST request to invert a matrix
 * @param entity XMl formatted matrix //from   www  .  j  a  v a  2  s  .  c  o m
 * @return XML representation of the matrix inverse
 * @throws ResourceException
 */
@Post
public Representation matrixInversion(Representation entity) throws ResourceException {
    if (entity == null) {
        throw new ResourceException(Status.CLIENT_ERROR_BAD_REQUEST, "No Input- Matrix Inversion not possible");
    }
    DomRepresentation rep = new DomRepresentation(entity);
    NamedRealMatrix matrixInput = null;
    NamedRealMatrix inverseMatrix = null;
    try {
        // parse the parameters from the entity body
        matrixInput = MatrixParamParser
                .getMatrixInversionParamsFromDomNode(rep.getDocument().getDocumentElement());

        if (matrixInput == null || !matrixInput.getName().equalsIgnoreCase("A")) {
            throw new ResourceException(Status.CLIENT_ERROR_BAD_REQUEST,
                    "Couldn't retrieve the matrix for inversion.");
        }

        // perform inversion
        inverseMatrix = new NamedRealMatrix(new LUDecompositionImpl(matrixInput).getSolver().getInverse());

        // set name
        inverseMatrix.setName(MatrixConstants.INVERSION_MATRIX_RETURN_NAME);

        //create our response representation
        MatrixXmlRepresentation response = new MatrixXmlRepresentation(inverseMatrix);
        return response;
    } catch (IllegalArgumentException iae) {
        MatrixLogger.getInstance().error(iae.getMessage());
        throw new ResourceException(Status.CLIENT_ERROR_BAD_REQUEST, iae.getMessage());

    } catch (IOException ioe) {
        MatrixLogger.getInstance().error(ioe.getMessage());
        throw new ResourceException(Status.CLIENT_ERROR_BAD_REQUEST, ioe.getMessage());

    }
}

From source file:eu.amidst.core.exponentialfamily.EF_Normal_NormalParents2.java

/**
 * {@inheritDoc}/*from  w  w  w  .  j  av a2 s .c o m*/
 */
@Override
public void updateNaturalFromMomentParameters() {

    /*
     * First step: means and convariances
     */
    CompoundVector globalMomentParam = (CompoundVector) this.momentParameters;
    double mean_X = globalMomentParam.getXYbaseMatrix().getEntry(0);
    RealVector mean_Y = globalMomentParam.getTheta_beta0BetaRV();

    double cov_XX = globalMomentParam.getcovbaseMatrix().getEntry(0, 0) - mean_X * mean_X;
    RealMatrix cov_YY = globalMomentParam.getcovbaseMatrix().getSubMatrix(1, nOfParents, 1, nOfParents)
            .subtract(mean_Y.outerProduct(mean_Y));
    RealVector cov_XY = globalMomentParam.getcovbaseMatrix().getSubMatrix(0, 0, 1, nOfParents).getRowVector(0)
            .subtract(mean_Y.mapMultiply(mean_X));
    //RealVector cov_YX = cov_XY; //outerProduct transposes the vector automatically

    /*
     * Second step: betas and variance
     */
    RealMatrix cov_YYInverse = new LUDecompositionImpl(cov_YY).getSolver().getInverse();
    RealVector beta = cov_YYInverse.preMultiply(cov_XY);

    double beta_0 = mean_X - beta.dotProduct(mean_Y);
    double variance = cov_XX - beta.dotProduct(cov_XY);

    /*
     * Third step: natural parameters (5 in total)
     */

    /*
     * 1) theta_0
     */
    double theta_0 = beta_0 / variance;
    double[] theta_0array = { theta_0 };

    /*
     * 2) theta_0Theta
     */
    double variance2Inv = 1.0 / (2 * variance);
    RealVector theta_0Theta = beta.mapMultiply(-beta_0 / variance);
    ((CompoundVector) this.naturalParameters)
            .setXYbaseVector(new ArrayRealVector(theta_0array, theta_0Theta.getData()));

    /*
     * 3) theta_Minus1
     */
    double theta_Minus1 = -variance2Inv;

    /*
     * 4) theta_beta
     */
    RealVector theta_beta = beta.mapMultiply(variance2Inv);

    /*
     * 5) theta_betaBeta
     */
    RealMatrix theta_betaBeta = beta.outerProduct(beta).scalarMultiply(-variance2Inv * 2);

    /*
     * Store natural parameters
     */
    RealMatrix natural_XY = new Array2DRowRealMatrix(nOfParents + 1, nOfParents + 1);
    double[] theta_Minus1array = { theta_Minus1 };
    RealVector covXY = new ArrayRealVector(theta_Minus1array, theta_beta.getData());
    natural_XY.setColumnVector(0, covXY);
    natural_XY.setRowVector(0, covXY);
    natural_XY.setSubMatrix(theta_betaBeta.getData(), 1, 1);
    ((CompoundVector) this.naturalParameters).setcovbaseVector(natural_XY);

}

From source file:eu.amidst.core.exponentialfamily.EF_Normal_NormalParents.java

/**
 * {@inheritDoc}//from w w  w. j  av a 2  s  . co m
 */
@Override
public void updateNaturalFromMomentParameters() {

    /*
     * First step: means and convariances
     */
    CompoundVector globalMomentParam = (CompoundVector) this.momentParameters;
    double mean_X = globalMomentParam.getXYbaseMatrix().getEntry(0);
    RealVector mean_Y = globalMomentParam.getTheta_beta0BetaRV();

    double cov_XX = globalMomentParam.getcovbaseMatrix().getEntry(0, 0) - mean_X * mean_X;
    RealMatrix cov_YY = globalMomentParam.getcovbaseMatrix().getSubMatrix(1, nOfParents, 1, nOfParents)
            .subtract(mean_Y.outerProduct(mean_Y));
    RealVector cov_XY = globalMomentParam.getcovbaseMatrix().getSubMatrix(0, 0, 1, nOfParents).getRowVector(0)
            .subtract(mean_Y.mapMultiply(mean_X));
    //RealVector cov_YX = cov_XY; //outerProduct transposes the vector automatically

    /*
     * Second step: betas and variance
     */
    RealMatrix cov_YYInverse = new LUDecompositionImpl(cov_YY).getSolver().getInverse();
    RealVector beta = cov_YYInverse.preMultiply(cov_XY);

    this.betas = new double[beta.getDimension()];
    for (int i = 0; i < beta.getDimension(); i++) {
        this.betas[i] = beta.getEntry(i);
    }
    this.beta0 = mean_X - beta.dotProduct(mean_Y);
    this.variance = cov_XX - beta.dotProduct(cov_XY);

}

From source file:edu.valelab.GaussianFit.CoordinateMapper.java

public static PolynomialCoefficients fitPolynomial(ExponentPairs exponentPairs,
        Map<Point2D.Double, Point2D.Double> pointPairs) {
    final List<Point2D.Double> srcPoints = new ArrayList(pointPairs.keySet());
    final RealMatrix matrix = new Array2DRowRealMatrix(srcPoints.size(), exponentPairs.size());
    for (int i = 0; i < srcPoints.size(); ++i) {
        matrix.setRow(i, powerTerms(srcPoints.get(i).x, srcPoints.get(i).y, exponentPairs));
    }/*w w  w .j  a va2  s  . c om*/
    final DecompositionSolver solver = new LUDecompositionImpl(matrix).getSolver();
    final double[] destX = new double[srcPoints.size()];
    final double[] destY = new double[srcPoints.size()];
    for (int i = 0; i < srcPoints.size(); ++i) {
        final Point2D.Double destPoint = pointPairs.get(srcPoints.get(i));
        destX[i] = destPoint.x;
        destY[i] = destPoint.y;
    }
    final PolynomialCoefficients polys = new PolynomialCoefficients();
    polys.polyX = solver.solve(destX);
    polys.polyY = solver.solve(destY);
    return polys;
}

From source file:com.stromberglabs.jopensurf.FastHessian.java

private double[] interpolateStep(int r, int c, ResponseLayer t, ResponseLayer m, ResponseLayer b) {
    double[] values = new double[3];

    RealMatrix partialDerivs = getPartialDerivativeMatrix(r, c, t, m, b);
    RealMatrix hessian3D = getHessian3DMatrix(r, c, t, m, b);

    DecompositionSolver solver = new LUDecompositionImpl(hessian3D).getSolver();
    RealMatrix X = solver.getInverse().multiply(partialDerivs);

    //      System.out.println("X = " + X.getColumnDimension() + " col x " + X.getRowDimension() + " rows");
    //      for ( int i = 0; i < X.getRowDimension(); i++ ){
    //         for ( int j = 0; j < X.getColumnDimension(); j++ ){
    //            System.out.print(X.getEntry(i,j) + (j != X.getColumnDimension()-1 ? " - " : ""));
    //         }/* w  w w .j  a v a 2s .  co  m*/
    //         System.out.println();
    //      }
    //      System.out.println();
    //      
    //values of them are used
    //xi
    values[0] = -X.getEntry(2, 0);
    //xr
    values[1] = -X.getEntry(1, 0);
    //xc
    values[2] = -X.getEntry(0, 0);

    return values;
}

From source file:edu.valelab.gaussianfit.datasettransformations.CoordinateMapper.java

public static PolynomialCoefficients fitPolynomial(ExponentPairs exponentPairs,
        Map<Point2D.Double, Point2D.Double> pointPairs) {
    final List<Point2D.Double> srcPoints = new ArrayList<Point2D.Double>(pointPairs.keySet());
    final RealMatrix matrix = new Array2DRowRealMatrix(srcPoints.size(), exponentPairs.size());
    for (int i = 0; i < srcPoints.size(); ++i) {
        matrix.setRow(i, powerTerms(srcPoints.get(i).x, srcPoints.get(i).y, exponentPairs));
    }//from w ww  . jav a 2 s. c om
    final DecompositionSolver solver = new LUDecompositionImpl(matrix).getSolver();
    final double[] destX = new double[srcPoints.size()];
    final double[] destY = new double[srcPoints.size()];
    for (int i = 0; i < srcPoints.size(); ++i) {
        final Point2D.Double destPoint = pointPairs.get(srcPoints.get(i));
        destX[i] = destPoint.x;
        destY[i] = destPoint.y;
    }
    final PolynomialCoefficients polys = new PolynomialCoefficients();
    polys.polyX = solver.solve(destX);
    polys.polyY = solver.solve(destY);
    return polys;
}