Angiogenesis

src/arcade/core/util/Graph.java

@@ -31,6 +31,9 @@ public enum Strategy {
     /** Collection of all {@code Edge} objects in a graph. */
     private final Bag allEdges;
 
+    /** Map of {@code Node} objects, for lookup. */
+    private final Map<String, Node> nodes;
+
     /** Map of {@code Node} OUT to bag of {@code Edge} objects. */
     private final Map<Node, Bag> nodeToOutBag;
 
@@ -42,6 +45,7 @@ public Graph() {
         allEdges = new Bag();
         nodeToOutBag = new HashMap<>();
         nodeToInBag = new HashMap<>();
+        nodes = new HashMap<>();
     }
 
     /**
@@ -53,13 +57,32 @@ public void update(Graph graph) {
         allEdges.addAll(graph.allEdges);
         nodeToOutBag.putAll(graph.nodeToOutBag);
         nodeToInBag.putAll(graph.nodeToInBag);
+        updateNodes();
+    }
+
+    /** Updates nodes from graph. */
+    private void updateNodes() {
+        nodes.clear();
+        for (Object obj : allEdges) {
+            Edge edge = (Edge) obj;
+            Node from = edge.getFrom();
+            Node to = edge.getTo();
+
+            if (!nodes.containsKey(from.toString())) {
+                nodes.put(from.toString(), from);
+            }
+            if (!nodes.containsKey(to.toString())) {
+                nodes.put(to.toString(), to);
+            }
+        }
     }
 
     /** Clear edges and nodes from graph. */
     public void clear() {
         allEdges.clear();
         nodeToOutBag.clear();
         nodeToInBag.clear();
+        nodes.clear();
     }
 
     /**
@@ -91,6 +114,28 @@ public boolean contains(Edge edge) {
         return checkEdge(edge);
     }
 
+    /**
+     * Retrieve the node object for the given coordinates. Returns null if no node exists.
+     *
+     * @param x the x coordinate
+     * @param y the y coordinate
+     * @param z the z coordinate
+     * @return the node object
+     */
+    public Node lookup(int x, int y, int z) {
+        return nodes.get("(" + x + "," + y + "," + z + ")");
+    }
+
+    /**
+     * Retrieve the node object for the given coordinates. Returns null if no node exists.
+     *
+     * @param node to lookup an original node
+     * @return the node object
+     */
+    public Node lookup(Node node) {
+        return nodes.get(node.toString());
+    }
+
     /**
      * Gets all edges in the graph.
      *
@@ -255,6 +300,8 @@ public void mergeNodes() {
                     e.setTo(join);
                 }
             }
+
+            nodes.put(join.toString(), join);
         }
     }
 
@@ -278,6 +325,27 @@ public void addEdge(Edge edge) {
         setOutMap(edge.getFrom(), edge);
         setInMap(edge.getTo(), edge);
         setLinks(edge);
+        addNodes(edge);
+    }
+
+    /**
+     * Helper function to adds the nodes from an edge to a graph.
+     *
+     * @param edge the edge to add
+     */
+    private void addNodes(Edge edge) {
+        Node from = edge.getFrom();
+        Node to = edge.getTo();
+        if (!nodes.containsKey(from.toString())) {
+            nodes.put(from.toString(), from);
+        } else {
+            from = nodes.get(from.toString());
+        }
+        if (!nodes.containsKey(to.toString())) {
+            nodes.put(to.toString(), to);
+        } else {
+            to = nodes.get(to.toString());
+        }
     }
 
     /**
@@ -353,6 +421,19 @@ public void removeEdge(Edge edge) {
         unsetOutMap(edge.getFrom(), edge);
         unsetInMap(edge.getTo(), edge);
         unsetLinks(edge);
+        removeNodeIfDetached(edge.getFrom());
+        removeNodeIfDetached(edge.getTo());
+    }
+
+    /**
+     * Remove a node if it is detached from the graph.
+     *
+     * @param node the node to check
+     */
+    private void removeNodeIfDetached(Node node) {
+        if (!nodeToInBag.containsKey(node) && !nodeToOutBag.containsKey(node)) {
+            nodes.remove(node.toString());
+        }
     }
 
     /**
@@ -730,7 +811,8 @@ public static class Edge {
         private final ArrayList<Edge> edgesOut;
 
         /**
-         * Creates an {@code Edge} between two {@link Node} objects.
+         * Creates an {@code Edge} between two {@link Node} objects. Default behavior is to
+         * duplicate nodes.
          *
          * @param from the node the edge is from
          * @param to the node the edge is to
@@ -742,6 +824,21 @@ public Edge(Node from, Node to) {
             edgesOut = new ArrayList<>();
         }
 
+        /**
+         * Creates an {@code Edge} object for graph sites. Used in cases where a new node object is
+         * not needed.
+         *
+         * @param from the node the edge is from
+         * @param to the node the edge is to
+         * @param duplicate {@code true} if nodes should be duplicated, {@code false} otherwise
+         */
+        public Edge(Node from, Node to, boolean duplicate) {
+            this.from = duplicate ? from.duplicate() : from;
+            this.to = duplicate ? to.duplicate() : to;
+            edgesIn = new ArrayList<>();
+            edgesOut = new ArrayList<>();
+        }
+
         /**
          * Gets the node the edge points to based on the calculation strategy (e.g. upstream or
          * downstream).
@@ -822,7 +919,7 @@ public ArrayList<Edge> getEdgesOut() {
          *
          * @return the reversed edge
          */
-        Edge reverse() {
+        public Edge reverse() {
             Node tempTo = to;
             Node tempFrom = from;
             to = tempFrom;

src/arcade/core/util/MatrixArray.java

@@ -0,0 +1,129 @@
+package arcade.core.util;
+
+import java.util.ArrayList;
+import java.util.Collections;
+import arcade.core.util.Matrix.Value;
+
+/**
+ * Container class for array-based sparse matrix representation.
+ *
+ * <p>Class provides a subset of matrix operations needed for solving a system of linear equations
+ * using the successive over-relaxation method in {@link arcade.core.util.Solver}.
+ */
+public class MatrixArray {
+    public MatrixArray(double[][] a) {
+        nRows = a.length;
+        if (nRows == 0) {
+            throw new IllegalArgumentException("MatrixArray ctor: input is empty");
+        }
+        nColumns = a[0].length;
+
+        if (nColumns == 0) {
+            throw new IllegalArgumentException("Matrix array ctor: empty columns");
+        }
+
+        ArrayList<Value> alValues = new ArrayList<>();
+        int nNonZero = 0;
+        for (int row = 0; row < nRows; row++) {
+            if (a[row].length != nColumns) {
+                throw new IllegalArgumentException(
+                        "Matrix array ctor: not all columns are the same length");
+            }
+
+            for (int column = 0; column < nColumns; column++) {
+                if (a[row][column] != 0.) {
+                    alValues.add(new Value(row, column, a[row][column]));
+                    nNonZero++;
+                }
+            }
+        }
+
+        //
+        // No need to sort, since we built them in order.
+        //
+
+        values = new Value[nNonZero];
+        int i = 0;
+        for (Value v : alValues) {
+            values[i] = v;
+            i++;
+        }
+    } // ctor
+
+    public MatrixArray(ArrayList<Value> alValues, int i_nRows, int i_Columns) {
+        nRows = i_nRows;
+        nColumns = i_Columns;
+        values = new Value[alValues.size()];
+
+        Collections.sort(
+                alValues,
+                (v1, v2) -> (v1.i == v2.i ? Integer.compare(v1.j, v2.j) : (v1.i > v2.i ? 1 : -1)));
+
+        int i = 0;
+        for (Value v : alValues) {
+            values[i] = v;
+            i++;
+        }
+    }
+
+    /**
+     * Solves the equation {@code Lx = b} using forward substitution for an array-based sparse
+     * matrix.
+     *
+     * @param b the right-hand side vector
+     * @return the left-hand side vector
+     */
+    public double[] forwardSubstitution(double[] b) {
+
+        int n = b.length;
+        double[] subbed = new double[n];
+        double[] diag = new double[n];
+
+        // Group lower diagonal by row.
+        ArrayList<ArrayList<Value>> rowsL = new ArrayList<ArrayList<Value>>();
+        for (int r = 0; r < n; r++) {
+            rowsL.add(new ArrayList<>());
+        }
+        for (Value v : values) {
+            rowsL.get(v.i).add(v);
+        }
+
+        // Get values along diagonal.
+        for (Value v : values) {
+            if (v.i == v.j) {
+                diag[v.i] = v.v;
+            }
+        }
+
+        // Iterate only through non-zero entries in the lower diagonal matrix.
+        for (int i = 0; i < n; i++) {
+            ArrayList<Value> rowL = rowsL.get(i);
+            double val = 0;
+            for (Value v : rowL) {
+                val += subbed[v.j] * v.v;
+            }
+            val = b[i] - val;
+            subbed[i] = val / diag[i];
+        }
+
+        return subbed;
+    }
+
+    public double[] multiply(double[] b) {
+        if (b.length != nRows) {
+            throw new IllegalArgumentException("MatrixArray.multiply (by a vector): conformation");
+        }
+        double[] multiplied = new double[nRows];
+
+        // Iterate through all entries and multiply.
+        for (Value a : values) {
+            multiplied[a.i] += a.v * b[a.j];
+        }
+
+        return multiplied;
+    }
+
+    int nRows;
+    int nColumns;
+    Value[] values;
+}