Examins the structure of A for QR decomposition
@param A matrix which is to be decomposed
@return true if the solution is valid or false if the decomposition can't be performed (i.e. requires column pivots)
340054
public boolean process( DMatrixSparseCSC A ) {
        init(A);

        TriangularSolver_DSCC.eliminationTree(A,true,parent,gwork);

        countNonZeroInR(parent);
        countNonZeroInV(parent);

        // if more columns than rows it's possible that Q*R != A. That's because a householder
        // would need to be created that's outside the  m by m Q matrix. In reality it has
        // a partial solution. Column pivot are needed.
        if( m < n ) {
            for (int row = 0; row <m; row++) {
                if( gwork.data[head+row] < 0 ) {
                    return false;
                }
            }
        }
        return true;
    }
1157657
public static boolean decomposeQR_block_col( final int blockLength ,
                                                 final DSubmatrixD1 Y ,
                                                 final double gamma[] )
    {
        int width = Y.col1-Y.col0;
        int height = Y.row1-Y.row0;
        int min = Math.min(width,height);
        for( int i = 0; i < min; i++ ) {
            // compute the householder vector
            if (!computeHouseHolderCol(blockLength, Y, gamma, i))
                return false;

            // apply to rest of the columns in the block
            rank1UpdateMultR_Col(blockLength,Y,i,gamma[Y.col0+i]);
        }

        return true;
    }
1011640
private static boolean[] createInitialR(CompRowMatrix A) {
            boolean[] R = new boolean[A.numRows()];

            int[] rowptr = A.getRowPointers();
            int[] colind = A.getColumnIndices();
            double[] data = A.getData();

            for (int i = 0; i < A.numRows(); ++i) {
                boolean hasOffDiagonal = false;

                for (int j = rowptr[i]; j < rowptr[i + 1]; ++j)
                    if (colind[j] != i && data[j] != 0) {
                        hasOffDiagonal = true;
                        break;
                    }

                R[i] = hasOffDiagonal;
            }

            return R;
        }
217705
public static boolean cholL( DMatrix5x5 A ) {

        A.a11 = Math.sqrt(A.a11);
        A.a12 = 0;
        A.a13 = 0;
        A.a14 = 0;
        A.a15 = 0;
        A.a21 = (A.a21)/A.a11;
        A.a22 = Math.sqrt(A.a22-A.a21*A.a21);
        A.a23 = 0;
        A.a24 = 0;
        A.a25 = 0;
        A.a31 = (A.a31)/A.a11;
        A.a32 = (A.a32-A.a31*A.a21)/A.a22;
        A.a33 = Math.sqrt(A.a33-A.a31*A.a31-A.a32*A.a32);
        A.a34 = 0;
        A.a35 = 0;
        A.a41 = (A.a41)/A.a11;
        A.a42 = (A.a42-A.a41*A.a21)/A.a22;
        A.a43 = (A.a43-A.a41*A.a31-A.a42*A.a32)/A.a33;
        A.a44 = Math.sqrt(A.a44-A.a41*A.a41-A.a42*A.a42-A.a43*A.a43);
        A.a45 = 0;
        A.a51 = (A.a51)/A.a11;
        A.a52 = (A.a52-A.a51*A.a21)/A.a22;
        A.a53 = (A.a53-A.a51*A.a31-A.a52*A.a32)/A.a33;
        A.a54 = (A.a54-A.a51*A.a41-A.a52*A.a42-A.a53*A.a43)/A.a44;
        A.a55 = Math.sqrt(A.a55-A.a51*A.a51-A.a52*A.a52-A.a53*A.a53-A.a54*A.a54);
        return !UtilEjml.isUncountable(A.a55);
    }
860019
public void decompose( DMatrixRMaj A ) {

        this.QR = A.copy();

        int N = Math.min(A.numCols,A.numRows);

        gammas = new double[ A.numCols ];

        DMatrixRMaj A_small = new DMatrixRMaj(A.numRows,A.numCols);
        DMatrixRMaj A_mod = new DMatrixRMaj(A.numRows,A.numCols);
        DMatrixRMaj v = new DMatrixRMaj(A.numRows,1);
        DMatrixRMaj Q_k = new DMatrixRMaj(A.numRows,A.numRows);

        for( int i = 0; i < N; i++ ) {
            // reshape temporary variables
            A_small.reshape(QR.numRows-i,QR.numCols-i);
            v.reshape(A_small.numRows,1);
            Q_k.reshape(v.getNumElements(),v.getNumElements());

            // use extract matrix to get the column that is to be zeroed
            CommonOps_DDRM.extract(QR,i,QR.numRows,i,i+1,v,0,0);

            double max = CommonOps_DDRM.elementMaxAbs(v);

            if( max > 0 && v.getNumElements() > 1 ) {
                // normalize to reduce overflow issues
                CommonOps_DDRM.divide(v,max);

                // compute the magnitude of the vector
                double tau = NormOps_DDRM.normF(v);

                if( v.get(0) < 0 )
                    tau *= -1.0;

                double u_0 = v.get(0) + tau;
                double gamma = u_0/tau;

                CommonOps_DDRM.divide(v,u_0);
                v.set(0,1.0);

                // extract the submatrix of A which is being operated on
                CommonOps_DDRM.extract(QR,i,QR.numRows,i,QR.numCols,A_small,0,0);

                // A = (I - &gamma;*u*u<sup>T</sup>)A
                CommonOps_DDRM.setIdentity(Q_k);
                CommonOps_DDRM.multAddTransB(-gamma,v,v,Q_k);
                CommonOps_DDRM.mult(Q_k,A_small,A_mod);

                // save the results
                CommonOps_DDRM.insert(A_mod, QR, i,i);
                CommonOps_DDRM.insert(v, QR, i,i);
                QR.unsafe_set(i,i,-tau*max);

                // save gamma for recomputing Q later on
                gammas[i] = gamma;
            }
        }
    }
359650
public static void abs(DMatrixD1 a ) {
        final int length = a.getNumElements();

        for ( int i = 0; i < length; i++ ) {
            a.data[i] = Math.abs(a.data[i]);
        }
    }
819024
public static void standard(DMatrix1Row A, DMatrix1Row A_tran)
    {
        int index = 0;
        for( int i = 0; i < A_tran.numRows; i++ ) {
            int index2 = i;

            int end = index + A_tran.numCols;
            while( index < end ) {
                A_tran.data[index++ ] = A.data[ index2 ];
                index2 += A.numCols;
            }
        }
    }
779188
@Override
    public DMatrixRMaj getQ(DMatrixRMaj Q , boolean compact ) {
        if( compact ) {
            if( Q == null ) {
                Q = CommonOps_DDRM.identity(numRows,minLength);
            } else {
                if( Q.numRows != numRows || Q.numCols != minLength ) {
                    throw new IllegalArgumentException("Unexpected matrix dimension.");
                } else {
                    CommonOps_DDRM.setIdentity(Q);
                }
            }
        } else {
            if( Q == null ) {
                Q = CommonOps_DDRM.identity(numRows);
            } else {
                if( Q.numRows != numRows || Q.numCols != numRows ) {
                    throw new IllegalArgumentException("Unexpected matrix dimension.");
                } else {
                    CommonOps_DDRM.setIdentity(Q);
                }
            }
        }

        for( int j = minLength-1; j >= 0; j-- ) {
            u[j] = 1;
            for( int i = j+1; i < numRows; i++ ) {
                u[i] = QR.get(i,j);
            }
            QrHelperFunctions_DDRM.rank1UpdateMultR(Q, u, gammas[j], j, j, numRows, v);
        }

        return Q;
    }
548487
@Override
    public boolean setA(ZMatrixRMaj A) {
        if( A.numRows > maxRows || A.numCols > maxCols ) {
            setMaxSize(A.numRows,A.numCols);
        }

        _setA(A);
        if( !decomposer.decompose(A) )
            return false;

        Q.reshape(numRows,numRows);
        R.reshape(numRows,numCols);
        decomposer.getQ(Q,false);
        decomposer.getR(R,false);
        CommonOps_ZDRM.transposeConjugate(Q,Qt);

        return true;
    }
227785
public static boolean cholU( DMatrix5x5 A ) {

        A.a11 = Math.sqrt(A.a11);
        A.a21 = 0;
        A.a31 = 0;
        A.a41 = 0;
        A.a51 = 0;
        A.a12 = (A.a12)/A.a11;
        A.a22 = Math.sqrt(A.a22-A.a12*A.a12);
        A.a32 = 0;
        A.a42 = 0;
        A.a52 = 0;
        A.a13 = (A.a13)/A.a11;
        A.a23 = (A.a23-A.a12*A.a13)/A.a22;
        A.a33 = Math.sqrt(A.a33-A.a13*A.a13-A.a23*A.a23);
        A.a43 = 0;
        A.a53 = 0;
        A.a14 = (A.a14)/A.a11;
        A.a24 = (A.a24-A.a12*A.a14)/A.a22;
        A.a34 = (A.a34-A.a13*A.a14-A.a23*A.a24)/A.a33;
        A.a44 = Math.sqrt(A.a44-A.a14*A.a14-A.a24*A.a24-A.a34*A.a34);
        A.a54 = 0;
        A.a15 = (A.a15)/A.a11;
        A.a25 = (A.a25-A.a12*A.a15)/A.a22;
        A.a35 = (A.a35-A.a13*A.a15-A.a23*A.a25)/A.a33;
        A.a45 = (A.a45-A.a14*A.a15-A.a24*A.a25-A.a34*A.a35)/A.a44;
        A.a55 = Math.sqrt(A.a55-A.a15*A.a15-A.a25*A.a25-A.a35*A.a35-A.a45*A.a45);
        return !UtilEjml.isUncountable(A.a55);
    }
661185
public static boolean cholL( DMatrix2x2 A ) {

        A.a11 = Math.sqrt(A.a11);
        A.a12 = 0;
        A.a21 = (A.a21)/A.a11;
        A.a22 = Math.sqrt(A.a22-A.a21*A.a21);
        return !UtilEjml.isUncountable(A.a22);
    }
993831
protected boolean computeH(DMatrixRMaj A, DMatrixRMaj H) {

		if( !solverNullspace.process(A.copy(),1,H) )
			return true;

		H.numRows = 3;
		H.numCols = 3;

		return false;
	}
480943
@Override
    public boolean setA(DMatrixRBlock A) {
        // Extract a lower triangular solution
        if( !decomposer.decompose(A) )
            return false;

        blockLength = A.blockLength;

        return true;
    }
784154
@Override
	public DoubleMatrix1D[] solve() throws Exception {

		final CholeskyFactorization HFact = new CholeskyFactorization(H, new Matrix1NornRescaler());
		boolean isHReducible = true;
		try{
			HFact.factorize();
		}catch(Exception e){
			isHReducible = false;
		}


        final DoubleMatrix1D v;// dim equals cols of A
        final DoubleMatrix1D w;// dim equals rank of A

		if (isHReducible) {
			// Solving KKT system via elimination
			DoubleMatrix1D HInvg = HFact.solve(g);
			if (A != null) {
				DoubleMatrix2D HInvAT = HFact.solve(AT);
				DoubleMatrix2D MenoSLower = ColtUtils.subdiagonalMultiply(A, HInvAT);
				DoubleMatrix1D AHInvg = ALG.mult(A, HInvg);
				
				CholeskyFactorization MSFact = new CholeskyFactorization(MenoSLower, new Matrix1NornRescaler());
				MSFact.factorize();
				if(h == null){
					w = MSFact.solve(ColtUtils.scalarMult(AHInvg, -1));
				}else{
					w = MSFact.solve(ColtUtils.add(h, AHInvg, -1));
				}
				
				v = HInvg.assign(ALG.mult(HInvAT, w), Functions.plus).assign(Mult.mult(-1));
			} else {
				w = null;
				v = HInvg.assign(Mult.mult(-1));
			}
		} else {
			// H is singular
			// Solving the full KKT system
			if(A!=null){
				DoubleMatrix1D[] fullSol =  this.solveAugmentedKKT();
				v = fullSol[0];
				w = fullSol[1];
			}else{
				//@TODO: try with rescaled H
				throw new Exception("KKT solution failed");
			}
		}

		// solution checking
		if (this.checkKKTSolutionAccuracy && !this.checkKKTSolutionAccuracy(v, w)) {
			log.error("KKT solution failed");
			throw new Exception("KKT solution failed");
		}

		DoubleMatrix1D[] ret = new DoubleMatrix1D[2];
		ret[0] = v;
		ret[1] = w;
		return ret;
	}
780464
@Override
    public boolean setA(DMatrixRBlock A) {
        if( A.numRows < A.numCols )
            throw new IllegalArgumentException("Number of rows must be more than or equal to the number of columns.  " +
                    "Can't solve an underdetermined system.");

        if( !decomposer.decompose(A))
            return false;

        this.QR = decomposer.getQR();

        return true;
    }
324694
public Matrix solve(Matrix rhs) {
        if (M != N || rhs.M != N || rhs.N != 1) {
            throw new RuntimeException("Illegal matrix dimensions.");
        }

        // create copies of the data
        Matrix A = new Matrix(this);
        Matrix b = new Matrix(rhs);

        // Gaussian elimination with partial pivoting
        for (int i = 0; i < N; i++) {

            // find pivot row and swap
            int max = i;
            for (int j = i + 1; j < N; j++) {
                if (Math.abs(A.data[j][i]) > Math.abs(A.data[max][i])) {
                    max = j;
                }
            }
            A.swap(i, max);
            b.swap(i, max);

            // singular
            if (A.data[i][i] == 0.0) {
                throw new RuntimeException("Matrix is singular.");
            }

            // pivot within b
            for (int j = i + 1; j < N; j++) {
                b.data[j][0] -= b.data[i][0] * A.data[j][i] / A.data[i][i];
            }

            // pivot within A
            for (int j = i + 1; j < N; j++) {
                double m = A.data[j][i] / A.data[i][i];
                for (int k = i + 1; k < N; k++) {
                    A.data[j][k] -= A.data[i][k] * m;
                }
                A.data[j][i] = 0.0;
            }
        }

        // back substitution
        Matrix x = new Matrix(N, 1);
        for (int j = N - 1; j >= 0; j--) {
            double t = 0.0;
            for (int k = j + 1; k < N; k++) {
                t += A.data[j][k] * x.data[k][0];
            }
            x.data[j][0] = (b.data[j][0] - t) / A.data[j][j];
        }
        return x;

    }
726928
public static boolean cholL( DMatrix6x6 A ) {

        A.a11 = Math.sqrt(A.a11);
        A.a12 = 0;
        A.a13 = 0;
        A.a14 = 0;
        A.a15 = 0;
        A.a16 = 0;
        A.a21 = (A.a21)/A.a11;
        A.a22 = Math.sqrt(A.a22-A.a21*A.a21);
        A.a23 = 0;
        A.a24 = 0;
        A.a25 = 0;
        A.a26 = 0;
        A.a31 = (A.a31)/A.a11;
        A.a32 = (A.a32-A.a31*A.a21)/A.a22;
        A.a33 = Math.sqrt(A.a33-A.a31*A.a31-A.a32*A.a32);
        A.a34 = 0;
        A.a35 = 0;
        A.a36 = 0;
        A.a41 = (A.a41)/A.a11;
        A.a42 = (A.a42-A.a41*A.a21)/A.a22;
        A.a43 = (A.a43-A.a41*A.a31-A.a42*A.a32)/A.a33;
        A.a44 = Math.sqrt(A.a44-A.a41*A.a41-A.a42*A.a42-A.a43*A.a43);
        A.a45 = 0;
        A.a46 = 0;
        A.a51 = (A.a51)/A.a11;
        A.a52 = (A.a52-A.a51*A.a21)/A.a22;
        A.a53 = (A.a53-A.a51*A.a31-A.a52*A.a32)/A.a33;
        A.a54 = (A.a54-A.a51*A.a41-A.a52*A.a42-A.a53*A.a43)/A.a44;
        A.a55 = Math.sqrt(A.a55-A.a51*A.a51-A.a52*A.a52-A.a53*A.a53-A.a54*A.a54);
        A.a56 = 0;
        A.a61 = (A.a61)/A.a11;
        A.a62 = (A.a62-A.a61*A.a21)/A.a22;
        A.a63 = (A.a63-A.a61*A.a31-A.a62*A.a32)/A.a33;
        A.a64 = (A.a64-A.a61*A.a41-A.a62*A.a42-A.a63*A.a43)/A.a44;
        A.a65 = (A.a65-A.a61*A.a51-A.a62*A.a52-A.a63*A.a53-A.a64*A.a54)/A.a55;
        A.a66 = Math.sqrt(A.a66-A.a61*A.a61-A.a62*A.a62-A.a63*A.a63-A.a64*A.a64-A.a65*A.a65);
        return !UtilEjml.isUncountable(A.a66);
    }