Examples with DoubleMatrix2D cern.colt.matrix.DoubleMatrix2D used on opensource projects

Example 21 with DoubleMatrix2D

use of cern.colt.matrix.DoubleMatrix2D in project tdq-studio-se by Talend.

the class CholeskyDecomposition method XXXsolveBuggy.

/**
 *Solves <tt>A*X = B</tt>; returns <tt>X</tt>.
 *@param  B   A Matrix with as many rows as <tt>A</tt> and any number of columns.
 *@return     <tt>X</tt> so that <tt>L*L'*X = B</tt>.
 *@exception  IllegalArgumentException  if <tt>B.rows() != A.rows()</tt>.
 *@exception  IllegalArgumentException  if <tt>!isSymmetricPositiveDefinite()</tt>.
 */
private DoubleMatrix2D XXXsolveBuggy(DoubleMatrix2D B) {
    cern.jet.math.Functions F = cern.jet.math.Functions.functions;
    if (B.rows() != n) {
        throw new IllegalArgumentException("Matrix row dimensions must agree.");
    }
    if (!isSymmetricPositiveDefinite) {
        throw new IllegalArgumentException("Matrix is not symmetric positive definite.");
    }
    // Copy right hand side.
    DoubleMatrix2D X = B.copy();
    int nx = B.columns();
    // precompute and cache some views to avoid regenerating them time and again
    DoubleMatrix1D[] Xrows = new DoubleMatrix1D[n];
    for (int k = 0; k < n; k++) Xrows[k] = X.viewRow(k);
    // Solve L*Y = B;
    for (int k = 0; k < n; k++) {
        for (int i = k + 1; i < n; i++) {
            // X[i,j] -= X[k,j]*L[i,k]
            Xrows[i].assign(Xrows[k], F.minusMult(L.getQuick(i, k)));
        }
        Xrows[k].assign(F.div(L.getQuick(k, k)));
    }
    // Solve L'*X = Y;
    for (int k = n - 1; k >= 0; k--) {
        Xrows[k].assign(F.div(L.getQuick(k, k)));
        for (int i = 0; i < k; i++) {
            // X[i,j] -= X[k,j]*L[k,i]
            Xrows[i].assign(Xrows[k], F.minusMult(L.getQuick(k, i)));
        }
    }
    return X;
}

Example 22 with DoubleMatrix2D

use of cern.colt.matrix.DoubleMatrix2D in project tdq-studio-se by Talend.

the class Smp method splitStridedNN.

protected DoubleMatrix2D[] splitStridedNN(DoubleMatrix2D A, int threshold, long flops) {
    /*
	determine how to split and parallelize best into blocks
	if more B.columns than tasks --> split B.columns, as follows:
	
			xx|xx|xxx B
			xx|xx|xxx
			xx|xx|xxx
	A
	xxx     xx|xx|xxx C 
	xxx		xx|xx|xxx
	xxx		xx|xx|xxx
	xxx		xx|xx|xxx
	xxx		xx|xx|xxx

	if less B.columns than tasks --> split A.rows, as follows:
	
			xxxxxxx B
			xxxxxxx
			xxxxxxx
	A
	xxx     xxxxxxx C
	xxx     xxxxxxx
	---     -------
	xxx     xxxxxxx
	xxx     xxxxxxx
	---     -------
	xxx     xxxxxxx

	*/
    // long flops = 2L*A.rows()*A.columns()*A.columns();
    // each thread should process at least 30000 flops
    int noOfTasks = (int) Math.min(flops / threshold, this.maxThreads);
    boolean splitHoriz = (A.columns() < noOfTasks);
    // boolean splitHoriz = (A.columns() >= noOfTasks);
    int p = splitHoriz ? A.rows() : A.columns();
    noOfTasks = Math.min(p, noOfTasks);
    if (noOfTasks < 2) {
        // parallelization doesn't pay off (too much start up overhead)
        return null;
    }
    // set up concurrent tasks
    int span = p / noOfTasks;
    final DoubleMatrix2D[] blocks = new DoubleMatrix2D[noOfTasks];
    for (int i = 0; i < noOfTasks; i++) {
        final int offset = i * span;
        // last span may be a bit larger
        if (i == noOfTasks - 1)
            span = p - span * i;
        final DoubleMatrix2D AA, BB, CC;
        if (!splitHoriz) {
            // split B along columns into blocks
            blocks[i] = A.viewPart(0, i, A.rows(), A.columns() - i).viewStrides(1, noOfTasks);
        } else {
            // split A along rows into blocks
            blocks[i] = A.viewPart(i, 0, A.rows() - i, A.columns()).viewStrides(noOfTasks, 1);
        }
    }
    return blocks;
}

Example 23 with DoubleMatrix2D

use of cern.colt.matrix.DoubleMatrix2D in project tdq-studio-se by Talend.

the class QRTest method main.

public static void main(String[] args) {
    // For COLT
    DoubleMatrix2D xmatrix, ymatrix, zmatrix;
    DoubleFactory2D myfactory;
    myfactory = DoubleFactory2D.dense;
    xmatrix = myfactory.make(8, 2);
    ymatrix = myfactory.make(8, 1);
    xmatrix.set(0, 0, 1);
    xmatrix.set(1, 0, 1);
    xmatrix.set(2, 0, 1);
    xmatrix.set(3, 0, 1);
    xmatrix.set(4, 0, 1);
    xmatrix.set(5, 0, 1);
    xmatrix.set(6, 0, 1);
    xmatrix.set(7, 0, 1);
    xmatrix.set(0, 1, 80);
    xmatrix.set(1, 1, 220);
    xmatrix.set(2, 1, 140);
    xmatrix.set(3, 1, 120);
    xmatrix.set(4, 1, 180);
    xmatrix.set(5, 1, 100);
    xmatrix.set(6, 1, 200);
    xmatrix.set(7, 1, 160);
    ymatrix.set(0, 0, 0.6);
    ymatrix.set(1, 0, 6.70);
    ymatrix.set(2, 0, 5.30);
    ymatrix.set(3, 0, 4.00);
    ymatrix.set(4, 0, 6.55);
    ymatrix.set(5, 0, 2.15);
    ymatrix.set(6, 0, 6.60);
    ymatrix.set(7, 0, 5.75);
    Algebra myAlgebra = new Algebra();
    zmatrix = myAlgebra.solve(xmatrix, ymatrix);
    System.err.println(xmatrix);
    System.err.println(ymatrix);
    System.err.println(zmatrix);
/*
       // For JAMA
       Matrix amatrix,bmatrix,cmatrix;
       amatrix = new Matrix(8,2);
       bmatrix = new Matrix(8,1);

       amatrix.set(0,0,1);
       amatrix.set(1,0,1);
       amatrix.set(2,0,1);
       amatrix.set(3,0,1);
       amatrix.set(4,0,1);
       amatrix.set(5,0,1);
       amatrix.set(6,0,1);
       amatrix.set(7,0,1);

       amatrix.set(0,1,80);
       amatrix.set(1,1,220);
       amatrix.set(2,1,140);
       amatrix.set(3,1,120);
       amatrix.set(4,1,180);
       amatrix.set(5,1,100);
       amatrix.set(6,1,200);
       amatrix.set(7,1,160);

       bmatrix.set(0,0,0.6);
       bmatrix.set(1,0,6.70);
       bmatrix.set(2,0,5.30);
       bmatrix.set(3,0,4.00);
       bmatrix.set(4,0,6.55);
       bmatrix.set(5,0,2.15);
       bmatrix.set(6,0,6.60);
       bmatrix.set(7,0,5.75);

       cmatrix = amatrix.solve(bmatrix);
       amatrix.print(8,5);
       bmatrix.print(8,5);
       cmatrix.print(8,5);
		*/
}

Example 24 with DoubleMatrix2D

use of cern.colt.matrix.DoubleMatrix2D in project tdq-studio-se by Talend.

the class Benchmark method benchmark.

/**
 * Runs a bench on matrices holding double elements.
 */
public static void benchmark(int runs, int size, String kind, boolean print, int initialCapacity, double minLoadFactor, double maxLoadFactor, double percentNonZero) {
    // certain loops need to be constructed so that the jitter can't optimize them away and we get fantastic numbers.
    // this involves primarly read-loops
    cern.colt.Timer timer1 = new cern.colt.Timer();
    cern.colt.Timer timer2 = new cern.colt.Timer();
    cern.colt.Timer timer3 = new cern.colt.Timer();
    cern.colt.Timer timer4 = new cern.colt.Timer();
    cern.colt.Timer timer5 = new cern.colt.Timer();
    cern.colt.Timer timer6 = new cern.colt.Timer();
    DoubleMatrix2D matrix = null;
    if (kind.equals("sparse"))
        matrix = new SparseDoubleMatrix2D(size, size, initialCapacity, minLoadFactor, maxLoadFactor);
    else if (kind.equals("dense"))
        matrix = cern.colt.matrix.DoubleFactory2D.dense.make(size, size);
    else
        // else if (kind.equals("denseArray")) matrix = new DoubleArrayMatrix2D(size,size);
        throw new RuntimeException("unknown kind");
    System.out.println("\nNow initializing...");
    // Matrix AJ = new Matrix(columnwise,3);
    // Basic.random(matrix, new cern.jet.random.Uniform(new cern.jet.random.engine.MersenneTwister()));
    double value = 2;
    DoubleMatrix2D tmp = DoubleFactory2D.dense.sample(matrix.rows(), matrix.columns(), value, percentNonZero);
    matrix.assign(tmp);
    tmp = null;
    /*
	long NN = matrix.size();
	int nn = (int) (NN*percentNonZero);
	long[] nonZeroIndexes = new long[nn];
	cern.jet.random.sampling.RandomSampler sampler = new cern.jet.random.sampling.RandomSampler(nn,NN,0,new cern.jet.random.engine.MersenneTwister());
	sampler.nextBlock(nn,nonZeroIndexes,0);
	for (int i=nn; --i >=0; ) {
		int row = (int) (nonZeroIndexes[i]/size);
		int column = (int) (nonZeroIndexes[i]%size);
		matrix.set(row,column, value);
	}
	*/
    /*
	timer1.start();
	for (int i=0; i<runs; i++) {
		LUDecomposition LU = new LUDecomposition(matrix);
	}
	timer1.stop();
	timer1.display();

	{
		Jama.Matrix jmatrix = new Jama.Matrix(matrix.toArray());
		timer2.start();
		for (int i=0; i<runs; i++) {
			Jama.LUDecomposition LU = new Jama.LUDecomposition(jmatrix);
		}
		timer2.stop();
		timer2.display();
	}
	*/
    System.out.println("\ntesting...");
    if (print)
        System.out.println(matrix);
    DoubleMatrix2D dense = DoubleFactory2D.dense.make(size, size);
    dense.assign(matrix);
    if (!dense.equals(matrix))
        throw new InternalError();
    DoubleMatrix2D ADense = dense.copy();
    DoubleMatrix2D BDense = dense.copy();
    DoubleMatrix2D CDense = dense.copy();
    ADense.zMult(BDense, CDense);
    System.out.println("\nNext testing...");
    /*
	{
		timer6.start();
		double a = cubicLoop(runs,size);
		timer6.stop();
		timer6.display();
		System.out.println(a);
	}
	*/
    {
        DoubleMatrix2D A = matrix.copy();
        DoubleMatrix2D B = matrix.copy();
        // DoubleMatrix2D C = Basic.product(A,B);
        DoubleMatrix2D C = matrix.copy();
        A.zMult(B, C);
        if (!(C.equals(CDense)))
            throw new InternalError();
        C.assign(matrix);
        System.out.println("\nNow benchmarking...");
        timer3.start();
        for (int i = 0; i < runs; i++) {
            A.zMult(B, C);
        }
        timer3.stop();
        timer3.display();
        int m = A.rows();
        int n = A.columns();
        int p = B.rows();
        int reps = runs;
        double mflops = 1.0e-3 * (2.0 * m * n * p * reps) / timer3.millis();
        System.out.println("mflops: " + mflops);
    }
    if (print)
        System.out.println(matrix);
    System.out.println("bye bye.");
}

Example 25 with DoubleMatrix2D

use of cern.colt.matrix.DoubleMatrix2D in project tdq-studio-se by Talend.

the class DenseDoubleMatrix2D method assign.

/**
 * Replaces all cell values of the receiver with the values of another matrix.
 * Both matrices must have the same number of rows and columns.
 * If both matrices share the same cells (as is the case if they are views derived from the same matrix) and intersect in an ambiguous way, then replaces <i>as if</i> using an intermediate auxiliary deep copy of <tt>other</tt>.
 *
 * @param     source   the source matrix to copy from (may be identical to the receiver).
 * @return <tt>this</tt> (for convenience only).
 * @throws	IllegalArgumentException if <tt>columns() != source.columns() || rows() != source.rows()</tt>
 */
public DoubleMatrix2D assign(DoubleMatrix2D source) {
    // overriden for performance only
    if (!(source instanceof DenseDoubleMatrix2D)) {
        return super.assign(source);
    }
    DenseDoubleMatrix2D other = (DenseDoubleMatrix2D) source;
    // nothing to do
    if (other == this)
        return this;
    checkShape(other);
    if (this.isNoView && other.isNoView) {
        // quickest
        System.arraycopy(other.elements, 0, this.elements, 0, this.elements.length);
        return this;
    }
    if (haveSharedCells(other)) {
        DoubleMatrix2D c = other.copy();
        if (!(c instanceof DenseDoubleMatrix2D)) {
            // should not happen
            return super.assign(other);
        }
        other = (DenseDoubleMatrix2D) c;
    }
    final double[] elems = this.elements;
    final double[] otherElems = other.elements;
    if (elems == null || otherElems == null)
        throw new InternalError();
    int cs = this.columnStride;
    int ocs = other.columnStride;
    int rs = this.rowStride;
    int ors = other.rowStride;
    int otherIndex = other.index(0, 0);
    int index = index(0, 0);
    for (int row = rows; --row >= 0; ) {
        for (int i = index, j = otherIndex, column = columns; --column >= 0; ) {
            elems[i] = otherElems[j];
            i += cs;
            j += ocs;
        }
        index += rs;
        otherIndex += ors;
    }
    return this;
}