octave-nkf: src/DLD-FUNCTIONS/symrcm.cc comparison

comparison src/DLD-FUNCTIONS/symrcm.cc @ 11586:12df7854fa7c

strip trailing whitespace from source files

author	John W. Eaton <jwe@octave.org>
date	Thu, 20 Jan 2011 17:24:59 -0500
parents	01f703952eff
children	72c96de7a403

comparison

equal deleted inserted replaced

-:1473d0cf86d2
+:12df7854fa7c
 // Queues Q have a fixed maximum size N (rows,cols of the matrix) and are
 // stored in an array. qh and qt point to queue head and tail.
 // Enqueue operation (adds a node "o" at the tail)
 inline static void
 Q_enq (CMK_Node *Q, octave_idx_type N, octave_idx_type& qt, const CMK_Node& o)
 {
 Q[qt] = o;
 qt = (qt + 1) % (N + 1);
 }
 // Dequeue operation (removes a node from the head)
 inline static CMK_Node
 Q_deq (CMK_Node * Q, octave_idx_type N, octave_idx_type& qh)
 {
 CMK_Node r = Q[qh];
 qh = (qh + 1) % (N + 1);
 return r;
 #define PARENT(i)       (((i) - 1) >> 1)        // = floor(((i)-1)/2)
 // Builds a min-heap (the root contains the smallest element). A is an array
 // with the graph's nodes, i is a starting position, size is the length of A.
 static void
 H_heapify_min (CMK_Node *A, octave_idx_type i, octave_idx_type size)
 {
 octave_idx_type j = i;
 for (;;)
 {
 CMK_Node tmp = A[j];
 A[j] = A[smallest];
 A[smallest] = tmp;
 j = smallest;
 }
 else
 break;
 }
 }
 // Heap operation insert. Running time is O(log(n))
 static void
 H_insert (CMK_Node *H, octave_idx_type& h, const CMK_Node& o)
 {
 octave_idx_type i = h++;
 H[i] = o;
 if (i == 0)
 return;
 do
 {
 octave_idx_type p = PARENT(i);
 if (H[i].deg < H[p].deg)
 H[i] = H[p];
 H[p] = tmp;
 i = p;
 }
 else
 break;
 }
 while (i > 0);
 }
 // Heap operation remove-min. Removes the smalles element in O(1) and
 // reorganizes the heap optionally in O(log(n))
 inline static CMK_Node
 H_remove_min (CMK_Node *H, octave_idx_type& h, int reorg/*=1*/)
 {
 CMK_Node r = H[0];
 H[0] = H[--h];
 if (reorg)
 H_heapify_min(H, 0, h);
 return r;
 }
 // Predicate (heap empty)
 #define H_empty(H, h)   ((h) == 0)
 // Helper function for the Cuthill-McKee algorithm. Tries to determine a
 // pseudo-peripheral node of the graph as starting node.
 static octave_idx_type
 find_starting_node (octave_idx_type N, const octave_idx_type *ridx,
 const octave_idx_type *cidx, const octave_idx_type *ridx2,
 const octave_idx_type *cidx2, octave_idx_type *D,
 octave_idx_type start)
 {
 CMK_Node w;
 OCTAVE_LOCAL_BUFFER (CMK_Node, Q, N+1);
 }
 return x.id;
 }
 // Calculates the node's degrees. This means counting the non-zero elements
 // in the symmetric matrix' rows. This works for non-symmetric matrices
 // as well.
 static octave_idx_type
 calc_degrees (octave_idx_type N, const octave_idx_type *ridx,
 const octave_idx_type *cidx, octave_idx_type *D)
 {
 octave_idx_type max_deg = 0;
 for (octave_idx_type i = 0; i < N; i++)
 D[i] = 0;
 for (octave_idx_type j = 0; j < N; j++)
 {
 for (octave_idx_type i = cidx[j]; i < cidx[j+1]; i++)
 {
 OCTAVE_QUIT;
 octave_idx_type k = ridx[i];
 // there is a non-zero element (k,j)
 D[k]++;
 if (D[k] > max_deg)
 max_deg = D[k];
 // if there is no element (j,k) there is one in
 // the symmetric matrix:
 if (k != j)
 {
 if (! found)
 {
 // A(j,k) == 0
 D[j]++;
 if (D[j] > max_deg)
 max_deg = D[j];
 }
 }
 }
 }
 }
 // Transpose of the structure of a square sparse matrix
 static void
 transpose (octave_idx_type N, const octave_idx_type *ridx,
 const octave_idx_type *cidx, octave_idx_type *ridx2,
 octave_idx_type *cidx2)
 {
 octave_idx_type nz = cidx[N];
 OCTAVE_LOCAL_BUFFER (octave_idx_type, w, N + 1);
 return retval;
 }
 octave_value arg = args(0);
 // the parameter of the matrix is converted into a sparse matrix
 //(if necessary)
 octave_idx_type *cidx;
 octave_idx_type *ridx;
 SparseMatrix Ar;
 SparseComplexMatrix Ac;
 // if none of the nodes has a degree > 0 (a matrix of zeros)
 // the return value corresponds to the identity permutation
 if (max_deg == 0)
 {
 for (octave_idx_type i = 0; i < N; i++)
 P(i) = i;
 return octave_value (P);
 }
 // a heap for the a node's neighbors. The number of neighbors is
 // bandwidths.
 octave_idx_type B = 0;
 // mark all nodes as unvisited; with the exception of the nodes
 // that have degree==0 and build a CC of the graph.
 boolNDArray btmp (dim_vector (1, N), false);
 bool *visit = btmp.fortran_vec ();
 do
 {
 // locate an unvisited starting node of the graph
 octave_idx_type i;
 for (i = 0; i < N; i++)
 if (! visit[i])
 break;
 // locate a probably better starting node
 v.id = find_starting_node (N, ridx, cidx, ridx2, cidx2, D, i);
 // mark the node as visited and enqueue it (a starting node
 // for the BFS). Since the node will be a root of a spanning
 // tree, its dist is 0.
 v.deg = D[v.id];
 v.dist = 0;
 octave_idx_type Bsub = 0;
 // min. dist. to the root is 0
 octave_idx_type level = 0;
 // the root is the first/only node on level 0
 octave_idx_type level_N = 1;
 while (! Q_empty (Q, N, qh, qt))
 {
 v = Q_deq (Q, N, qh);
 i = v.id;
 c++;
 // for computing the inverse permutation P where
 // A(inv(P),inv(P)) or P'*A*P is banded
 //         P(i) = c;
 // for computing permutation P where
 // A(P(i),P(j)) or P*A*P' is banded
 P(c) = i;
 // put all unvisited neighbors j of node i on the heap
 {
 OCTAVE_QUIT;
 // locate a neighbor of i with minimal degree in O(log(N))
 v = H_remove_min(S, s, 1);
 // entered the BFS a new level?
 if (v.dist > level)
 {
 // adjustment of bandwith:
 // "[...] the minimum bandwidth that
 // can be obtained [...] is the
 //  maximum number of nodes per level"
 if (Bsub < level_N)
 Bsub = level_N;
 level = v.dist;
 // v is the first node on the new level
 level_N = 1;
 }
 else
 {
 // there is no new level but another node on
 // this level:
 level_N++;
 }
 // enqueue v in O(1)
 Q_enq (Q, N, qt, v);
 }
 // synchronize the bandwidth with level_N once again:
 if (Bsub < level_N)
 Bsub = level_N;
 }
 // finish of BFS. If there are still unvisited nodes in the graph

Mercurial > octave-nkf

comparison src/DLD-FUNCTIONS/symrcm.cc @ 11586:12df7854fa7c