Example usage for org.apache.commons.math3.linear LUDecomposition getDeterminant

List of usage examples for org.apache.commons.math3.linear LUDecomposition getDeterminant

Introduction

In this page you can find the example usage for org.apache.commons.math3.linear LUDecomposition getDeterminant.

Prototype

public double getDeterminant() 

Source Link

Document

Return the determinant of the matrix

Usage

From source file:Matrix_Operations.java

public static void main(String[] args) {

    Scanner keyboard = new Scanner(System.in);

    //Allow user to enter number of columns
    int rows = getRowsNumberFromUser();
    int columns = getColumnNumberFromUser();
    double matrixA[][] = new double[rows][columns];

    //Enter values for matrix
    System.out.println("Enter values for each position in the matrix below. ");
    for (int i = 0; i < rows; i++) {
        for (int j = 0; j < columns; j++) {
            matrixA[i][j] = keyboard.nextDouble();

        }//  www .  j a  v a  2 s. co  m
        System.out.println("");
    }

    showMatrix(matrixA);
    System.out.println("");
    RealMatrix A = MatrixUtils.createRealMatrix(matrixA);

    LUDecomposition lu = new LUDecomposition(A);

    showMatrix(matrixA);

    System.out.println(lu.getDeterminant());

}

From source file:hivemall.utils.math.StatsUtils.java

/**
 * pdf(x, x_hat) = exp(-0.5 * (x-x_hat) * inv() * (x-x_hat)T) / ( 2^0.5d * det()^0.5)
 * /* w  w  w .  j  av a  2s .  c  o m*/
 * @return value of probabilistic density function
 * @link https://en.wikipedia.org/wiki/Multivariate_normal_distribution#Density_function
 */
public static double pdf(@Nonnull final RealVector x, @Nonnull final RealVector x_hat,
        @Nonnull final RealMatrix sigma) {
    final int dim = x.getDimension();
    Preconditions.checkArgument(x_hat.getDimension() == dim,
            "|x| != |x_hat|, |x|=" + dim + ", |x_hat|=" + x_hat.getDimension());
    Preconditions.checkArgument(sigma.getRowDimension() == dim,
            "|x| != |sigma|, |x|=" + dim + ", |sigma|=" + sigma.getRowDimension());
    Preconditions.checkArgument(sigma.isSquare(), "Sigma is not square matrix");

    LUDecomposition LU = new LUDecomposition(sigma);
    final double detSigma = LU.getDeterminant();
    double denominator = Math.pow(2.d * Math.PI, 0.5d * dim) * Math.pow(detSigma, 0.5d);
    if (denominator == 0.d) { // avoid divide by zero
        return 0.d;
    }

    final RealMatrix invSigma;
    DecompositionSolver solver = LU.getSolver();
    if (solver.isNonSingular() == false) {
        SingularValueDecomposition svd = new SingularValueDecomposition(sigma);
        invSigma = svd.getSolver().getInverse(); // least square solution
    } else {
        invSigma = solver.getInverse();
    }
    //EigenDecomposition eigen = new EigenDecomposition(sigma);
    //double detSigma = eigen.getDeterminant();
    //RealMatrix invSigma = eigen.getSolver().getInverse();

    RealVector diff = x.subtract(x_hat);
    RealVector premultiplied = invSigma.preMultiply(diff);
    double sum = premultiplied.dotProduct(diff);
    double numerator = Math.exp(-0.5d * sum);

    return numerator / denominator;
}

From source file:hivemall.utils.math.StatsUtils.java

/**
 * @param mu1 mean vector of the first normal distribution
 * @param sigma1 covariance matrix of the first normal distribution
 * @param mu2 mean vector of the second normal distribution
 * @param sigma2 covariance matrix of the second normal distribution
 * @return the Hellinger distance between two multivariate normal distributions
 * @link https://en.wikipedia.org/wiki/Hellinger_distance#Examples
 *///from  ww  w  .  j  a v a2  s . c  om
public static double hellingerDistance(@Nonnull final RealVector mu1, @Nonnull final RealMatrix sigma1,
        @Nonnull final RealVector mu2, @Nonnull final RealMatrix sigma2) {
    RealVector muSub = mu1.subtract(mu2);
    RealMatrix sigmaMean = sigma1.add(sigma2).scalarMultiply(0.5d);
    LUDecomposition LUsigmaMean = new LUDecomposition(sigmaMean);
    double denominator = Math.sqrt(LUsigmaMean.getDeterminant());
    if (denominator == 0.d) {
        return 1.d; // avoid divide by zero
    }
    RealMatrix sigmaMeanInv = LUsigmaMean.getSolver().getInverse(); // has inverse iff det != 0
    double sigma1Det = MatrixUtils.det(sigma1);
    double sigma2Det = MatrixUtils.det(sigma2);

    double numerator = Math.pow(sigma1Det, 0.25d) * Math.pow(sigma2Det, 0.25d)
            * Math.exp(-0.125d * sigmaMeanInv.preMultiply(muSub).dotProduct(muSub));
    return 1.d - numerator / denominator;
}

From source file:hivemall.utils.math.MatrixUtils.java

public static double det(@Nonnull final RealMatrix m) {
    LUDecomposition LU = new LUDecomposition(m);
    return LU.getDeterminant();
}

From source file:com.itemanalysis.psychometrics.cfa.MaximumLikelihoodEstimation.java

public double value(double[] argument) {
    model.setParameters(argument);//ww w. j  a v  a2s.c o m

    //            Linesearch method in QNMinimizer is causing NaN values after repeated calls to here
    //            Next libe is for monitoring values when called from line search
    //            No problem occurs with CGMinimizer
    //            System.out.println("valueAt: " + argument[0] + " " + argument[1]);
    SIGMA = model.getImpliedCovariance(argument);

    //compute determinant of SIGMA
    LUDecomposition SLUD = new LUDecomposition(SIGMA);
    double detSig = SLUD.getDeterminant();

    //compute inverse of SIGMA
    RealMatrix SIGMAinv = SLUD.getSolver().getInverse();
    RealMatrix VC_SIGMA_INV = varcov.multiply(SIGMAinv);
    double trace = VC_SIGMA_INV.getTrace();

    //convert number of items to double
    double p = Double.valueOf(model.getNumberOfItems()).doubleValue();

    //compute objective function
    F = Math.log(detSig) + trace - Math.log(detVc) - p;
    return F;
}

From source file:com.opengamma.strata.math.impl.matrix.CommonsMatrixAlgebra.java

@Override
public double getDeterminant(Matrix m) {
    ArgChecker.notNull(m, "m");
    if (m instanceof DoubleMatrix) {
        RealMatrix temp = CommonsMathWrapper.wrap((DoubleMatrix) m);
        LUDecomposition lud = new LUDecomposition(temp);
        return lud.getDeterminant();
    }//from   ww w  . j  a v  a2 s. c o m
    throw new IllegalArgumentException("Can only find determinant of DoubleMatrix; have " + m.getClass());
}

From source file:com.opengamma.strata.math.impl.linearalgebra.LUDecompositionCommonsResult.java

/**
 * @param lu The result of the LU decomposition, not null. $\mathbf{L}$ cannot be singular.
 *///from  w  w  w  .  j  a v a  2 s . c  o  m
public LUDecompositionCommonsResult(LUDecomposition lu) {
    ArgChecker.notNull(lu, "LU decomposition");
    ArgChecker.notNull(lu.getL(), "Matrix is singular; could not perform LU decomposition");
    _determinant = lu.getDeterminant();
    _l = CommonsMathWrapper.unwrap(lu.getL());
    _p = CommonsMathWrapper.unwrap(lu.getP());
    _pivot = lu.getPivot();
    _solver = lu.getSolver();
    _u = CommonsMathWrapper.unwrap(lu.getU());
}

From source file:lanchester.MultiArena.java

public void step() {
    int numFoes = forces.size();
    if (isMyTurn == null) {
        isMyTurn = new boolean[numFoes][numFoes];
        stanceArray = new int[numFoes][numFoes];
        currentFloor = new double[numFoes][numFoes];
        for (int i1 = 0; i1 < numFoes; i1++) {
            int ind1 = forceMap.get(forces.get(i1));
            for (int i2 = 0; i2 < numFoes; i2++) {
                int ind2 = forceMap.get(forces.get(i2));
                isMyTurn[i1][i2] = true;
                if (i1 == i2) {
                    stanceArray[i1][i2] = AthenaConstants.ALLIED_POSTURE;
                    currentFloor[i1][i2] = 0.;
                } else {
                    stanceArray[i1][i2] = initializeStance(forces.get(i1), forces.get(i2));
                    setFloor(i1, i2);//from w  w  w  .  ja va  2 s . c  o  m
                }
            }
        }
    }
    Array2DRowRealMatrix mat = getMat();
    EigenDecomposition eigen = new EigenDecomposition(mat);
    double det = eigen.getDeterminant();
    double[] eVals = eigen.getRealEigenvalues();
    //        for(int i1=0;i1<eVals.length;i1++){
    //            System.out.println("eVals["+i1+"] = "+eVals[i1]);
    //        }
    if (eigen.hasComplexEigenvalues()) {
        System.out.println("Complex eigenvalues");
        for (int i1 = 0; i1 < forces.size(); i1++) {
            AthenaForce f = forces.get(i1);
            System.out.println(f.getName() + " has " + f.getForceSize() + " forces remaining");
        }
    }
    double[] initialNums = getInitialNumbers(forces);
    Array2DRowRealMatrix eVectors = (Array2DRowRealMatrix) eigen.getV();
    LUDecomposition lu = new LUDecomposition(eVectors);
    double det2 = lu.getDeterminant();
    double[] coeffs = new double[numFoes];
    for (int i1 = 0; i1 < numFoes; i1++) {
        Array2DRowRealMatrix tmpMat = (Array2DRowRealMatrix) eVectors.copy();
        tmpMat.setColumn(i1, initialNums);
        LUDecomposition tmpLU = new LUDecomposition(tmpMat);
        double tmpDet = tmpLU.getDeterminant();
        coeffs[i1] = tmpDet / det2;
    }
    MultiTimeStep currentStep = new MultiTimeStep(numFoes);
    currentTime += timeStep;
    currentStep.setTime(currentTime);
    for (int i1 = 0; i1 < numFoes; i1++) {
        double updatedForce = 0.;
        for (int i2 = 0; i2 < numFoes; i2++) {
            updatedForce += coeffs[i2] * eVectors.getEntry(i1, i2) * Math.exp(eVals[i2] * timeStep);
            //                updatedForce+=coeffs[i2]*eVectors.getEntry(i2, i1)*Math.exp(eVals[i2]*timeStep);
            //                updatedForce+=coeffs[i1]*eVectors.getEntry(i2, i1)*Math.exp(eVals[i1]*timeStep);
        }
        forces.get(i1).updateForce(updatedForce);
        currentStep.setForceNumber(updatedForce, i1);
    }
    history.add(currentStep);
    //        eVectors.
    //        this.currentTime++;
    //                Truncator truncator = new Truncator();
    if (true) {
        //            System.out.println("time = " + time);
    }
}

From source file:lanchester.MultiArena2.java

public void step() {
    boolean aboveFloor = true;
    double currentCycle = 0.;
    int numFoes = forces.size();
    if (isMyTurn == null) {
        isMyTurn = new boolean[numFoes][numFoes];
        stanceArray = new int[numFoes][numFoes];
        currentFloor = new double[numFoes][numFoes];
        for (int i1 = 0; i1 < numFoes; i1++) {
            int ind1 = forceMap.get(forces.get(i1));
            for (int i2 = 0; i2 < numFoes; i2++) {
                int ind2 = forceMap.get(forces.get(i2));
                isMyTurn[i1][i2] = true;
                if (i1 == i2) {
                    stanceArray[i1][i2] = AthenaConstants.ALLIED_POSTURE;
                    currentFloor[i1][i2] = 0.;
                } else {
                    stanceArray[i1][i2] = initializeStance(forces.get(i1), forces.get(i2));
                    setFloor(i1, i2);//from  w  w  w  . j av a2  s .  c  o m
                }
            }
        }
    }
    Array2DRowRealMatrix mat = getMat();
    EigenDecomposition eigen = new EigenDecomposition(mat);
    double det = eigen.getDeterminant();
    double[] eVals = eigen.getRealEigenvalues();
    //        for(int i1=0;i1<eVals.length;i1++){
    //            System.out.println("eVals["+i1+"] = "+eVals[i1]);
    //        }
    if (eigen.hasComplexEigenvalues()) {
        System.out.println("Complex eigenvalues");
        for (int i1 = 0; i1 < forces.size(); i1++) {
            AthenaForce f = forces.get(i1);
            System.out.println(f.getName() + " has " + f.getForceSize() + " forces remaining");
        }
    }
    double[] initialNums = getInitialNumbers(forces);
    Array2DRowRealMatrix eVectors = (Array2DRowRealMatrix) eigen.getV();
    LUDecomposition lu = new LUDecomposition(eVectors);
    double det2 = lu.getDeterminant();
    double[] coeffs = new double[numFoes];
    for (int i1 = 0; i1 < numFoes; i1++) {
        Array2DRowRealMatrix tmpMat = (Array2DRowRealMatrix) eVectors.copy();
        tmpMat.setColumn(i1, initialNums);
        LUDecomposition tmpLU = new LUDecomposition(tmpMat);
        double tmpDet = tmpLU.getDeterminant();
        coeffs[i1] = tmpDet / det2;
    }
    aboveFloor = true;
    int cntr = 0;
    int numGone;
    do {
        MultiTimeStep currentStep = new MultiTimeStep(numFoes);
        currentTime += timeStep;
        currentCycle += timeStep;
        currentStep.setTime(currentTime);
        numGone = 0;
        for (int i1 = 0; i1 < numFoes; i1++) {
            double updatedForce = 0.;
            if (forces.get(i1).getForceSize() > lb) {
                for (int i2 = 0; i2 < numFoes; i2++) {
                    //                    updatedForce += coeffs[i2] * eVectors.getEntry(i1, i2) * Math.exp(eVals[i2] * timeStep);
                    updatedForce += coeffs[i2] * eVectors.getEntry(i1, i2) * Math.exp(eVals[i2] * currentCycle);
                    if (updatedForce < 1.) {
                        updatedForce = 0.;
                        numGone++;
                    }
                    //                updatedForce+=coeffs[i2]*eVectors.getEntry(i2, i1)*Math.exp(eVals[i2]*timeStep);
                    //                updatedForce+=coeffs[i1]*eVectors.getEntry(i2, i1)*Math.exp(eVals[i1]*timeStep);
                }
            } else {
                updatedForce = lb / 2.;
                numGone++;
            }
            forces.get(i1).updateForce(updatedForce);
            currentStep.setForceNumber(updatedForce, i1);
        }
        history.add(currentStep);
        aboveFloor = checkAboveFloors();
        cntr++;
    } while (aboveFloor && cntr < 2000 && (numFoes - numGone) > 1);
    for (int i1 = 0; i1 < numFoes; i1++) {
        for (int i2 = 0; i2 < numFoes; i2++) {
            if (i1 == i2) {
                stanceArray[i1][i2] = AthenaConstants.ALLIED_POSTURE;
                currentFloor[i1][i2] = 0.;
            } else {
                stanceArray[i1][i2] = initializeStance(forces.get(i1), forces.get(i2));
                setFloor(i1, i2);
            }
        }
    }

    //        eVectors.
    //        this.currentTime++;
    //                Truncator truncator = new Truncator();
    if (numFoes - numGone == 1) {
        loneSurvivor = true;
        //            System.out.println("time = " + time);
    }
}

From source file:lanchester.MultiArena3.java

public void step() {
    boolean aboveFloor = true;
    double currentCycle = 0.;
    int numFoes = forces.size();
    if (isMyTurn == null) {
        isMyTurn = new boolean[numFoes][numFoes];
        stanceArray = new int[numFoes][numFoes];
        currentFloor = new double[numFoes][numFoes];
        currentCeiling = new double[numFoes][numFoes];
        for (int i1 = 0; i1 < numFoes; i1++) {
            int ind1 = forceMap.get(forces.get(i1));
            for (int i2 = 0; i2 < numFoes; i2++) {
                int ind2 = forceMap.get(forces.get(i2));
                isMyTurn[i1][i2] = true;
                if (i1 == i2) {
                    stanceArray[i1][i2] = AthenaConstants.ALLIED_POSTURE;
                    currentFloor[i1][i2] = 0.;
                    currentCeiling[i1][i2] = 100.;
                } else {
                    stanceArray[i1][i2] = initializeStance(forces.get(i1), forces.get(i2));
                    setFloor(i1, i2);//from  w ww.j  a  v a 2  s .  c o  m
                    setCeiling(i1, i2);
                }
            }
        }
    }
    Array2DRowRealMatrix mat = getMat();
    EigenDecomposition eigen = new EigenDecomposition(mat);
    double det = eigen.getDeterminant();
    double[] eVals = eigen.getRealEigenvalues();
    if (eigen.hasComplexEigenvalues()) {
        System.out.println("Complex eigenvalues");
        for (int i1 = 0; i1 < forces.size(); i1++) {
            MultiForce f = forces.get(i1);
            System.out.println(f.getName() + " has " + f.getNumber() + " forces remaining");
        }
    }
    double[] initialNums = getInitialNumbers(forces);
    Array2DRowRealMatrix eVectors = (Array2DRowRealMatrix) eigen.getV();
    LUDecomposition lu = new LUDecomposition(eVectors);
    double det2 = lu.getDeterminant();
    double[] coeffs = new double[numFoes];
    for (int i1 = 0; i1 < numFoes; i1++) {
        Array2DRowRealMatrix tmpMat = (Array2DRowRealMatrix) eVectors.copy();
        tmpMat.setColumn(i1, initialNums);
        LUDecomposition tmpLU = new LUDecomposition(tmpMat);
        double tmpDet = tmpLU.getDeterminant();
        coeffs[i1] = tmpDet / det2;
    }
    aboveFloor = true;
    boolean belowCeiling = true;
    int cntr = 0;
    int numGone;
    do {
        timeStep = determineTimeStep();
        MultiTimeStep currentStep = new MultiTimeStep(numFoes);
        currentTime += timeStep;
        currentCycle += timeStep;
        currentStep.setTime(currentTime);
        numGone = 0;
        for (int i1 = 0; i1 < numFoes; i1++) {
            double updatedForce = 0.;
            if (stillAlive[i1]) {
                for (int i2 = 0; i2 < numFoes; i2++) {
                    updatedForce += coeffs[i2] * eVectors.getEntry(i1, i2) * Math.exp(eVals[i2] * currentCycle);
                }
                if (updatedForce < 1.) {
                    updatedForce = lb;
                    stillAlive[i1] = false;
                    numGone++;
                }
            } else {
                numGone++;
                updatedForce = lb;
            }
            forces.get(i1).updateForce(updatedForce);
            currentStep.setForceNumber(updatedForce, i1);
        }
        history.add(currentStep);
        aboveFloor = checkAboveFloors();
        belowCeiling = checkBelowCeilings();
        cntr++;
    } while (aboveFloor && belowCeiling && cntr < 2000 && (numFoes - numGone) > 1);
    for (int i1 = 0; i1 < numFoes; i1++) {
        for (int i2 = 0; i2 < numFoes; i2++) {
            if (i1 == i2) {
                stanceArray[i1][i2] = AthenaConstants.ALLIED_POSTURE;
                currentFloor[i1][i2] = 0.;
            } else {
                stanceArray[i1][i2] = initializeStance(forces.get(i1), forces.get(i2));
                setFloor(i1, i2);
            }
        }
    }

    //        eVectors.
    //        this.currentTime++;
    //                Truncator truncator = new Truncator();
    if (numFoes - numGone == 1) {
        loneSurvivor = true;
        //            System.out.println("time = " + time);
    }
}