simple-squiggle

A restricted subset of Squiggle
Log | Files | Refs | README

csDfs.js (2434B)


      1 import { csMarked } from './csMarked.js';
      2 import { csMark } from './csMark.js';
      3 import { csUnflip } from './csUnflip.js';
      4 /**
      5  * Depth-first search computes the nonzero pattern xi of the directed graph G (Matrix) starting
      6  * at nodes in B (see csReach()).
      7  *
      8  * @param {Number}  j               The starting node for the DFS algorithm
      9  * @param {Matrix}  g               The G matrix to search, ptr array modified, then restored
     10  * @param {Number}  top             Start index in stack xi[top..n-1]
     11  * @param {Number}  k               The kth column in B
     12  * @param {Array}   xi              The nonzero pattern xi[top] .. xi[n - 1], an array of size = 2 * n
     13  *                                  The first n entries is the nonzero pattern, the last n entries is the stack
     14  * @param {Array}   pinv            The inverse row permutation vector, must be null for L * x = b
     15  *
     16  * @return {Number}                 New value of top
     17  *
     18  * Reference: http://faculty.cse.tamu.edu/davis/publications.html
     19  */
     20 
     21 export function csDfs(j, g, top, xi, pinv) {
     22   // g arrays
     23   var index = g._index;
     24   var ptr = g._ptr;
     25   var size = g._size; // columns
     26 
     27   var n = size[1]; // vars
     28 
     29   var i, p, p2; // initialize head
     30 
     31   var head = 0; // initialize the recursion stack
     32 
     33   xi[0] = j; // loop
     34 
     35   while (head >= 0) {
     36     // get j from the top of the recursion stack
     37     j = xi[head]; // apply permutation vector
     38 
     39     var jnew = pinv ? pinv[j] : j; // check node j is marked
     40 
     41     if (!csMarked(ptr, j)) {
     42       // mark node j as visited
     43       csMark(ptr, j); // update stack (last n entries in xi)
     44 
     45       xi[n + head] = jnew < 0 ? 0 : csUnflip(ptr[jnew]);
     46     } // node j done if no unvisited neighbors
     47 
     48 
     49     var done = 1; // examine all neighbors of j, stack (last n entries in xi)
     50 
     51     for (p = xi[n + head], p2 = jnew < 0 ? 0 : csUnflip(ptr[jnew + 1]); p < p2; p++) {
     52       // consider neighbor node i
     53       i = index[p]; // check we have visited node i, skip it
     54 
     55       if (csMarked(ptr, i)) {
     56         continue;
     57       } // pause depth-first search of node j, update stack (last n entries in xi)
     58 
     59 
     60       xi[n + head] = p; // start dfs at node i
     61 
     62       xi[++head] = i; // node j is not done
     63 
     64       done = 0; // break, to start dfs(i)
     65 
     66       break;
     67     } // check depth-first search at node j is done
     68 
     69 
     70     if (done) {
     71       // remove j from the recursion stack
     72       head--; // and place in the output stack
     73 
     74       xi[--top] = j;
     75     }
     76   }
     77 
     78   return top;
     79 }