Author: scooter
Date: 2010-10-03 16:27:54 -0700 (Sun, 03 Oct 2010)
New Revision: 22137
Modified:
csplugins/trunk/ucsf/scooter/clusterMaker/src/clusterMaker/algorithms/MCL/RunMCL.java
Log:
Multi-thread the matrix multiplication for MCL
Modified:
csplugins/trunk/ucsf/scooter/clusterMaker/src/clusterMaker/algorithms/MCL/RunMCL.java
===================================================================
---
csplugins/trunk/ucsf/scooter/clusterMaker/src/clusterMaker/algorithms/MCL/RunMCL.java
2010-10-02 01:26:10 UTC (rev 22136)
+++
csplugins/trunk/ucsf/scooter/clusterMaker/src/clusterMaker/algorithms/MCL/RunMCL.java
2010-10-03 23:27:54 UTC (rev 22137)
@@ -1,3 +1,34 @@
+/**
+ * Copyright (c) 2010 The Regents of the University of California.
+ * All rights reserved.
+ *
+ * Redistribution and use in source and binary forms, with or without
+ * modification, are permitted provided that the following conditions
+ * are met:
+ * 1. Redistributions of source code must retain the above copyright
+ * notice, this list of conditions, and the following disclaimer.
+ * 2. Redistributions in binary form must reproduce the above
+ * copyright notice, this list of conditions, and the following
+ * disclaimer in the documentation and/or other materials provided
+ * with the distribution.
+ * 3. Redistributions must acknowledge that this software was
+ * originally developed by the UCSF Computer Graphics Laboratory
+ * under support by the NIH National Center for Research Resources,
+ * grant P41-RR01081.
+ *
+ * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDER "AS IS" AND ANY
+ * EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+ * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
+ * PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS BE LIABLE
+ * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
+ * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT
+ * OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
+ * BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY,
+ * WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE
+ * OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE,
+ * EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+ *
+ */
package clusterMaker.algorithms.MCL;
import java.util.ArrayList;
@@ -4,6 +35,10 @@
import java.util.Arrays;
import java.util.Collection;
import java.util.Comparator;
+import java.util.concurrent.Callable;
+import java.util.concurrent.Executors;
+import java.util.concurrent.ExecutorService;
+import java.util.concurrent.TimeUnit;
import java.util.HashMap;
import java.util.Iterator;
import java.util.List;
@@ -26,6 +61,7 @@
import cern.colt.function.IntIntDoubleFunction;
import cern.colt.matrix.DoubleFactory2D;
+import cern.colt.matrix.DoubleMatrix1D;
import cern.colt.matrix.DoubleMatrix2D;
import cern.colt.matrix.impl.SparseDoubleMatrix2D;
@@ -45,6 +81,7 @@
private DistanceMatrix distanceMatrix = null;
private DoubleMatrix2D matrix = null;
private boolean debug = false;
+ final int NTHREADS = Runtime.getRuntime().availableProcessors();
public RunMCL(DistanceMatrix dMat, double inflationParameter, int
num_iterations,
double clusteringThresh, double maxResidual, CyLogger
logger )
@@ -75,6 +112,8 @@
CyAttributes netAttributes = Cytoscape.getNetworkAttributes();
CyAttributes nodeAttributes = Cytoscape.getNodeAttributes();
+ long startTime = System.currentTimeMillis();
+
// Matrix matrix;
double numClusters;
@@ -96,14 +135,19 @@
{
// Expand
{
+ long t = System.currentTimeMillis();
monitor.setStatus("Iteration: "+(i+1)+"
expanding ");
// debugln("Iteration: "+(i+1)+" expanding ");
// printMatrixInfo(matrix);
- // We really, really want to make sure this is
sparse!
- DoubleMatrix2D newMatrix =
DoubleFactory2D.sparse.make(nodes.size(),nodes.size());
- matrix = matrix.zMult(matrix, newMatrix, 1.0,
0.0, false, false);
+ if (NTHREADS > 1) {
+ matrix = multiplyMatrix(matrix, matrix);
+ } else {
+ DoubleMatrix2D newMatrix =
DoubleFactory2D.sparse.make(matrix.rows(), matrix.columns());
+ matrix = matrix.zMult(matrix,
newMatrix);
+ }
// Normalize
normalize(matrix, clusteringThresh, false);
+ logger.info("Expansion "+(i+1)+" took
"+(System.currentTimeMillis()-t)+"ms");
}
// printMatrix(matrix);
@@ -111,12 +155,14 @@
// Inflate
{
+ long t = System.currentTimeMillis();
monitor.setStatus("Iteration: "+(i+1)+"
inflating");
// debugln("Iteration: "+(i+1)+" inflating");
// printMatrixInfo(matrix);
matrix.forEachNonZero(myPow);
// Normalize
normalize(matrix, clusteringThresh, true);
+ logger.info("Inflation "+(i+1)+" took
"+(System.currentTimeMillis()-t)+"ms");
}
// printMatrix(matrix);
@@ -168,6 +214,8 @@
clusterNumber++;
}
+ logger.info("Total runtime =
"+(System.currentTimeMillis()-startTime)+"ms");
+
Set<NodeCluster>clusters = cMap.keySet();
return new ArrayList(clusters);
}
@@ -278,6 +326,101 @@
if (debug) System.out.print(message);
}
+ private DoubleMatrix2D multiplyMatrix(DoubleMatrix2D A, DoubleMatrix2D
B) {
+ int m = A.rows();
+ int n = A.columns();
+ int p = B.columns();
+
+ // Create views into B
+ final DoubleMatrix1D[] Brows= new DoubleMatrix1D[n];
+ for (int i = n; --i>=0; ) Brows[i] = B.viewRow(i);
+
+ // Create a series of 1D vectors
+ final DoubleMatrix1D[] Crows= new DoubleMatrix1D[n];
+ for (int i = m; --i>=0; ) Crows[i] = B.like1D(m);
+
+ // Create the thread pools
+ final ExecutorService[] threadPools = new
ExecutorService[NTHREADS];
+ for (int pool = 0; pool < threadPools.length; pool++) {
+ threadPools[pool] =
Executors.newFixedThreadPool(1);
+ }
+
+ // final cern.jet.math.PlusMult fun =
cern.jet.math.PlusMult.plusMult(0);
+
+ A.forEachNonZero(
+ new cern.colt.function.IntIntDoubleFunction() {
+ public double apply(int row, int column, double
value) {
+
+ Runnable r = new
ThreadedDotProduct(value, Brows[column], Crows[row]);
+ //r.run();
+ threadPools[row%NTHREADS].submit(r);
+ /*
+ final int frow = row;
+ final int fcolumn = column;
+ final double fvalue = value;
+ threadPool.submit(
+ new Callable <Double>() {
+ public Double call() {
+ final
cern.jet.math.PlusMult fun = cern.jet.math.PlusMult.plusMult(0);
+
fun.multiplicator = fvalue;
+
Crows[frow].assign(Brows[fcolumn], fun);
+ return fvalue;
+ }
+ }
+ );
+ */
+ return value;
+ }
+ }
+ );
+
+ for (int pool = 0; pool < threadPools.length; pool++) {
+ threadPools[pool].shutdown();
+ try {
+ boolean result =
threadPools[pool].awaitTermination(7, TimeUnit.DAYS);
+ } catch (Exception e) {}
+ }
+ // Recreate C
+ return create2DMatrix(Crows);
+ }
+
+ private DoubleMatrix2D create2DMatrix (DoubleMatrix1D[] rows) {
+ int columns = rows[0].size();
+ DoubleMatrix2D C = DoubleFactory2D.sparse.make(rows.length,
columns);
+ for (int row = 0; row < rows.length; row++) {
+ for (int col = 0; col < columns; col++) {
+ double value = rows[row].getQuick(col);
+ if (value != 0.0)
+ C.setQuick(row, col, value);
+ }
+ }
+ return C;
+ }
+
+ private class ThreadedDotProduct implements Runnable {
+ double value;
+ DoubleMatrix1D Bcol;
+ DoubleMatrix1D Crow;
+ // final cern.jet.math.PlusMult fun =
cern.jet.math.PlusMult.plusMult(0);
+
+ ThreadedDotProduct(double value, DoubleMatrix1D Bcol,
+ DoubleMatrix1D Crow) {
+ this.value = value;
+ this.Bcol = Bcol;
+ this.Crow = Crow;
+ }
+
+ public void run() {
+ // fun.multiplicator = value;
+ for (int k = 0; k < Bcol.size(); k++) {
+ if (Bcol.getQuick(k) != 0.0) {
+ Crow.setQuick(k,
Crow.getQuick(k)+Bcol.getQuick(k)*value);
+ }
+ }
+ // Crow.assign(Bcol, fun);
+ }
+ }
+
/**
* The MatrixPow class raises the value of each non-zero cell of the
matrix
* to the power passed in it's constructor.
--
You received this message because you are subscribed to the Google Groups
"cytoscape-cvs" group.
To post to this group, send email to [email protected].
To unsubscribe from this group, send email to
[email protected].
For more options, visit this group at
http://groups.google.com/group/cytoscape-cvs?hl=en.