X-Git-Url: http://git.nguyen.vg/gitweb/?p=SXSI%2Flibbp.git;a=blobdiff_plain;f=bp.c;fp=bp.c;h=b85508bd8d5c38cb2e2e6cb4114b90f436442282;hp=0000000000000000000000000000000000000000;hb=b94f8d72735df125b191bf5a49cba0c037278787;hpb=657b1bc3dc283b45b9cceab5f1825a06496146cd

diff --git a/bp.c b/bp.c
new file mode 100644
index 0000000..b85508b
--- /dev/null
+++ b/bp.c
@@ -0,0 +1,1161 @@
+#include "bp.h"
+
+#define max(a,b) \
+  ({ __typeof__ (a) _a = (a); \
+  __typeof__ (b) _b = (b); \
+  _a > _b ? _a : _b; })
+
+
+#define min(a,b) \
+  ({ __typeof__ (a) _a = (a); \
+  __typeof__ (b) _b = (b); \
+   _a <= _b ? _a : _b; })
+
+//#define CHECK
+#define RANDOM
+
+int msize=0;
+#define mymalloc(p,n,f) {				\
+  p = (__typeof__(p)) malloc((n)*sizeof(*p));	        \
+if ((p)==NULL) {printf("not enough memory (%d bytes) in line %d\n",msize,__LINE__); \
+  exit(1);};								\
+msize += (f)*(n)*sizeof(*p);						\
+;}
+
+int postorder_select_bsearch(bp *b,int s);
+
+int naive_depth(bp *b, int s)
+{
+  int i,d;
+  if (s < 0) return 0;
+  d = 0;
+  for (i=0; i<=s; i++) {
+    if (getbit(b->B,i)==OP) {
+      d++;
+    } else {
+      d--;
+    }
+  }
+  return d;
+}
+
+void printbp(bp *b, int s, int t)
+{
+  int i,c,d;
+  d = 0;
+  for (i=s; i<=t; i++) {
+    if (getbit(b->B,i)==OP) {
+      c = '(';
+      d++;
+    } else {
+      c = ')';
+      d--;
+    }
+    putchar(c);
+  }
+}
+
+int *matchtbl,*parenttbl;
+void make_naivetbl(pb *B,int n)
+{
+  int i,v;
+  int *stack,s;
+
+  mymalloc(matchtbl,n,0);
+  mymalloc(parenttbl,n,0);
+  mymalloc(stack,n,0);
+
+  for (i=0; i<n; i++) matchtbl[i] = parenttbl[i] = -1;
+
+  s = 0;
+  v = 0;
+  stack[0] = -1;
+  for (i=0; i<n; i++) {
+    if (getbit(B,i)==OP) {
+      v++;
+      if (v > 0) {
+	parenttbl[i] = stack[v-1];
+	stack[v] = i;
+      }
+    } else {
+      if (v > 0) {
+	matchtbl[stack[v]] = i;  // close
+	matchtbl[i] = stack[v];   // open
+      }
+      v--;
+    }
+  }
+
+  free(stack);
+}
+
+int popCount[1<<ETW];
+int fwdtbl[(2*ETW+1)*(1<<ETW)];
+int bwdtbl[(2*ETW+1)*(1<<ETW)];
+#if 0
+int mintbl_li[1<<ETW], mintbl_lv[1<<ETW];
+int mintbl_ri[1<<ETW], mintbl_rv[1<<ETW];
+int maxtbl_li[1<<ETW], maxtbl_lv[1<<ETW];
+int maxtbl_ri[1<<ETW], maxtbl_rv[1<<ETW];
+#endif
+int minmaxtbl_i[4][1<<ETW], minmaxtbl_v[4][1<<ETW];
+int degtbl[1<<ETW];
+int degtbl2[(2*ETW+1)*(1<<ETW)];
+int childtbl[(ETW)*(1<<ETW)];
+int childtbl2[2*ETW+1][ETW][(1<<ETW)];
+int depthtbl[(2*ETW+1)*(1<<ETW)];
+int inited = 0;
+void make_matchtbl(void)
+{
+  int i,j,x,r;
+  int m,M;
+  pb buf[1];
+  int deg;
+  if (inited)
+    return;
+  inited = 1;
+  for (x = 0; x < (1<<ETW); x++) {
+    setbits(buf,0,ETW,x);
+    for (r=-ETW; r<=ETW; r++) fwdtbl[((r+ETW)<<ETW)+x] = ETW;
+    for (r=-ETW; r<=ETW; r++) bwdtbl[((r+ETW)<<ETW)+x] = ETW;
+    for (r=-ETW; r<=ETW; r++) degtbl2[((r+ETW)<<ETW)+x] = 0;
+    for (r=-ETW; r<=ETW; r++) depthtbl[((r+ETW)<<ETW)+x] = 0;
+
+    r = 0;
+    for (i=0; i<ETW; i++) {
+      if (getbit(buf,i)==OP) {
+	r++;
+      } else {
+	r--;
+      }
+      if (fwdtbl[((r+ETW)<<ETW)+x] == ETW) fwdtbl[((r+ETW)<<ETW)+x] = i;
+    }
+
+    r = 0;
+    for (i=ETW-1; i>=0; i--) {
+      if (getbit(buf,i)==OP) {
+	r--;
+      } else {
+	r++;
+      }
+      if (bwdtbl[((r+ETW)<<ETW)+x] == ETW) bwdtbl[((r+ETW)<<ETW)+x] = ETW-1-i;
+    }
+
+    r = 0;
+    for (i=0; i<ETW; i++) {
+      if (getbit(buf,i)==OP) {
+	r++;
+      } else {
+	r--;
+      }
+      depthtbl[((r+ETW)<<ETW)+x] += (1<<(ETW-1));
+    }
+
+    r = 0;
+    for (i=0; i<ETW; i++) {
+      if (getbit(buf,i)==OP) r++;
+    }
+    popCount[x] = r;
+
+    r = 0;  m = 0;  M = 0;
+    m = ETW+1;  M = -ETW-1;
+    //maxtbl_lv[x] = -ETW-1;
+    //mintbl_lv[x] = ETW+1;
+    minmaxtbl_v[OPT_MAX | OPT_LEFT][x] = -ETW-1;
+    minmaxtbl_v[OPT_MIN | OPT_LEFT][x] = ETW+1;
+    deg = 0;
+    for (i=0; i<ETW; i++) {
+      if (getbit(buf,i)==OP) {
+	r++;
+	if (r > M) {
+	  M = r;  
+	  //maxtbl_li[x] = i;  maxtbl_lv[x] = r;
+	  minmaxtbl_i[OPT_MAX | OPT_LEFT][x] = i;
+	  minmaxtbl_v[OPT_MAX | OPT_LEFT][x] = r;
+	}
+      } else {
+	r--;
+	if (r == m) {
+	  deg++;
+	  childtbl[((deg-1)<<ETW) + x] = i;
+	}
+	if (r < m) {
+	  m = r;
+	  //mintbl_li[x] = i;  mintbl_lv[x] = r;
+	  minmaxtbl_i[OPT_MIN | OPT_LEFT][x] = i;
+	  minmaxtbl_v[OPT_MIN | OPT_LEFT][x] = r;
+	  deg = 1;
+	  childtbl[((deg-1)<<ETW) + x] = i;
+	}
+      }
+      if (r <= m) degtbl2[((r+ETW)<<ETW)+x]++;
+    }
+    degtbl[x] = deg;
+
+    r = 0;  m = 0;  M = 0;
+    //maxtbl_rv[x] = -ETW-1;
+    //mintbl_rv[x] = ETW+1;
+    minmaxtbl_v[OPT_MAX | OPT_RIGHT][x] = -ETW-1;
+    minmaxtbl_v[OPT_MIN | OPT_RIGHT][x] = ETW+1;
+    for (i=0; i<ETW; i++) {
+      if (getbit(buf,i)==OP) {
+	r++;
+	if (r >= M) {
+	  M = r;  
+	  //maxtbl_ri[x] = i;  maxtbl_rv[x] = r;
+	  minmaxtbl_i[OPT_MAX | OPT_RIGHT][x] = i;
+	  minmaxtbl_v[OPT_MAX | OPT_RIGHT][x] = r;
+	}
+      } else {
+	r--;
+	if (r <= m) {
+	  m = r;
+	  //mintbl_ri[x] = i;  mintbl_rv[x] = r;
+	  minmaxtbl_i[OPT_MIN | OPT_RIGHT][x] = i;
+	  minmaxtbl_v[OPT_MIN | OPT_RIGHT][x] = r;
+	}
+      }
+    }
+
+    for (i = 0; i < ETW; i++) {
+      for (j = -ETW; j <= ETW; j++) {
+	childtbl2[j+ETW][i][x] = -1;
+      }
+    }
+
+    for (j=-ETW; j<=ETW; j++) {
+      int ith;
+      ith = 0;
+      r = 0;
+      for (i = 0; i < ETW; i++) {
+	if (getbit(buf,i)==OP) {
+	  r++;
+	} else {
+	  r--;
+	  if (r < j) break;
+	  if (r == j) {
+	    ith++;
+	    childtbl2[j+ETW][ith-1][x] = i;
+	  }
+	}
+      }
+    }
+  }
+
+}
+
+
+int bp_construct(bp *b,int n, pb *B, int opt)
+{
+  int i,j,d;
+  int m,M,ds;
+  int ns,nm;
+  byte *sm, *sM;
+  byte *sd;
+  int *mm, *mM;
+  int *md;
+  int m_ofs;
+  int r; // # of minimum values
+
+#if 0
+  if (SB % D != 0) {
+    printf("warning: SB=%d should be a multiple of D=%d\n",SB,D);
+    // not necessarily?
+  }
+  if (SB % RRR != 0) {
+    printf("warning: SB=%d should be a multiple of RRR=%d\n",SB,RRR);
+  }
+#endif
+
+  b->B = B;
+  b->n = n;
+  b->opt = opt;
+  b->idx_size = 0;
+  b->sm = NULL;
+  b->sM = NULL;
+  b->sd = NULL;
+  b->mm = NULL;
+  b->mM = NULL;
+  b->md = NULL;
+  b->da_leaf = NULL;
+  b->da_inorder = NULL;
+  b->da_postorder = NULL;
+  b->da_dfuds_leaf = NULL;
+  mymalloc(b->da,1,0);
+  darray_construct(b->da,n,B, opt & OPT_FAST_PREORDER_SELECT);
+  b->idx_size += b->da->idx_size;
+  //Kim: comment this and the following, they polute the printing of the xpath library
+  //printf("preorder rank/select table: %d bytes (%1.2f bpc)\n",b->da->idx_size,(double)b->da->idx_size*8/n);
+
+  make_matchtbl();
+
+  ns = (n+SB-1)/SB;
+  mymalloc(sm, ns, 0);  b->idx_size += ns * sizeof(*sm);
+  mymalloc(sM, ns, 0);  b->idx_size += ns * sizeof(*sM);
+  b->sm = sm;
+  b->sM = sM;
+  if (opt & OPT_DEGREE) {
+    mymalloc(sd, ns, 0);  b->idx_size += ns * sizeof(*sd);
+    b->sd = sd;
+    //printf("SB degree table: %d bytes (%1.2f bpc)\n",ns * sizeof(*sd), (double)ns * sizeof(*sd) * 8/n);
+  }
+  //printf("SB table: %d bytes (%1.2f bpc)\n",ns * sizeof(*sm) * 2, (double)ns * sizeof(*sm)*2 * 8/n);
+
+  for (i=0; i<n; i++) {
+    if (i % SB == 0) {
+      ds = depth(b,i);
+      m = M = ds;
+      r = 1;
+    } else {
+      d = depth(b,i);
+      if (d == m) r++;
+      if (d < m) {
+	m = d;
+	r = 1;
+      }
+      if (d > M) M = d;
+    }
+    if (i % SB == SB-1 || i==n-1) {
+      ds = depth(b,(i/SB)*SB-1);
+      if (m - ds + SB < 0 || m - ds + SB > 255) {
+	printf("error m=%d ds=%d\n",m,ds);
+      }
+      if (M - ds + 1 < 0 || M - ds + 1 > 255) {
+	printf("error M=%d ds=%d\n",M,ds);
+      }
+      sm[i/SB] = m - ds + SB;
+      sM[i/SB] = M - ds + 1;
+      if (opt & OPT_DEGREE) sd[i/SB] = r;
+    }
+  }
+
+#if 0
+  printf("sd: ");
+  for (i=0;i<n/SB;i++) printf("%d ",sd[i]);
+  printf("\n");
+#endif
+
+
+  nm = (n+MB-1)/MB;
+  m_ofs = 1 << blog(nm-1);
+  b->m_ofs = m_ofs;
+
+  mymalloc(mm, nm + m_ofs, 0);  b->idx_size += (nm+m_ofs) * sizeof(*mm);
+  mymalloc(mM, nm + m_ofs, 0);  b->idx_size += (nm+m_ofs) * sizeof(*mM);
+  b->mm = mm;
+  b->mM = mM;
+  if (opt & OPT_DEGREE) {
+    mymalloc(md, nm + m_ofs, 0);  b->idx_size += (nm+m_ofs) * sizeof(*md);
+    b->md = md;
+    //printf("MB degree table: %d bytes (%1.2f bpc)\n",(nm+m_ofs) * sizeof(*md), (double)(nm+m_ofs) * sizeof(*md) * 8/n);
+  }
+  //printf("MB table: %d bytes (%1.2f bpc)\n",(nm+m_ofs) * sizeof(*mm) * 2, (double)(nm+m_ofs) * sizeof(*mm)*2 * 8/n);
+
+  for (i=0; i<n; i++) {
+    d = depth(b,i);
+    if (i % MB == 0) {
+      m = M = d;
+      r = 1;
+    } else {
+      if (d == m) r++;
+      if (d < m) {
+	m = d;
+	r = 1;
+      }
+      if (d > M) M = d;
+    }
+    if (i % MB == MB-1 || i==n-1) {
+      mm[m_ofs+ i/MB] = m;
+      mM[m_ofs+ i/MB] = M;
+      if (opt & OPT_DEGREE) md[m_ofs+ i/MB] = r;
+    }
+  }
+  
+  for (j=m_ofs-1; j > 0; j--) {
+    m = 0;
+    if (j*2 < nm + m_ofs) m = mm[j*2];
+    if (j*2+1 < nm + m_ofs) m = min(m,mm[j*2+1]);
+    M = 0;
+    if (j*2 < nm + m_ofs) M = mM[j*2];
+    if (j*2+1 < nm + m_ofs) M = max(M,mM[j*2+1]);
+    mm[j] = m;  mM[j] = M;
+    if (opt & OPT_DEGREE) {
+      d = 0;
+      if (j*2 < nm + m_ofs) d = md[j*2];
+      if (j*2+1 < nm + m_ofs) {
+	if (mm[j*2] == mm[j*2+1]) d += md[j*2+1];
+	if (mm[j*2] > mm[j*2+1]) d = md[j*2+1];
+      }
+      md[j] = d;
+    }
+  }
+  mm[0] = -1;
+  mM[0] = mM[1];
+  if (opt & OPT_DEGREE) {
+    md[0] = -1;
+  }
+
+
+#if 0
+  printf("md: ");
+  for (i=0;i<m_ofs + n/MB;i++) printf("%d ",md[i]);
+  printf("\n");
+#endif
+
+  if (opt & OPT_LEAF) {
+    mymalloc(b->da_leaf,1,0);
+    darray_pat_construct(b->da_leaf, n, B, 2, 0x2, opt & OPT_FAST_LEAF_SELECT);
+    //printf("leaf rank/select table: %d bytes (%1.2f bpc)\n",b->da_leaf->idx_size,(double)b->da_leaf->idx_size*8/n);
+    b->idx_size += b->da_leaf->idx_size;  
+  } else {
+    b->da_leaf = NULL;
+  }
+
+  if (opt & OPT_INORDER) {
+    mymalloc(b->da_inorder,1,0);
+    darray_pat_construct(b->da_inorder, n, B, 2, 0x1, opt & OPT_FAST_INORDER_SELECT);
+    //printf("inorder rank/select table: %d bytes (%1.2f bpc)\n",b->da_inorder->idx_size,(double)b->da_inorder->idx_size*8/n);
+    b->idx_size += b->da_inorder->idx_size;
+  } else {
+    b->da_inorder = NULL;
+  }
+
+  if (opt & OPT_FAST_POSTORDER_SELECT) {
+    mymalloc(b->da_postorder,1,0);
+    darray_pat_construct(b->da_postorder, n, B, 1, 0x0, (opt & OPT_FAST_POSTORDER_SELECT) | OPT_NO_RANK);
+    //printf("postorder rank/select table: %d bytes (%1.2f bpc)\n",b->da_postorder->idx_size,(double)b->da_postorder->idx_size*8/n);
+    b->idx_size += b->da_postorder->idx_size;
+  } else {
+    b->da_postorder = NULL;
+  }
+
+  if (opt & OPT_DFUDS_LEAF) {
+    mymalloc(b->da_dfuds_leaf,1,0);
+    darray_pat_construct(b->da_dfuds_leaf, n, B, 2, 0x0, opt & OPT_FAST_DFUDS_LEAF_SELECT);
+    //printf("dfuds leaf rank/select table: %d bytes (%1.2f bpc)\n",b->da_dfuds_leaf->idx_size,(double)b->da_dfuds_leaf->idx_size*8/n);
+    b->idx_size += b->da_dfuds_leaf->idx_size;
+  } else {
+    b->da_dfuds_leaf = NULL;
+  }
+
+  return 0;
+}
+
+// destroyTree: frees the memory of tree.
+void destroyTree(bp *b) {
+   if (!b) return; // nothing to free
+
+   destroyDarray(b->da); // destroys da data structure
+   if (b->da) free(b->da);     
+   
+   if (b->sm) free(b->sm);
+   if (b->sM) free(b->sM);
+   if (b->sd) free(b->sd);
+   if (b->mm) free(b->mm);
+   if (b->mM) free(b->mM);
+   if (b->md) free(b->md);
+   
+   destroyDarray(b->da_leaf);
+   if (b->da_leaf) free(b->da_leaf);
+   
+   destroyDarray(b->da_inorder);
+   if (b->da_inorder) free(b->da_inorder);
+   
+   destroyDarray(b->da_postorder);
+   if (b->da_postorder) free(b->da_postorder);
+   
+   destroyDarray(b->da_dfuds_leaf);
+   if (b->da_dfuds_leaf) free(b->da_dfuds_leaf);
+}
+
+
+// saveTree: saves parentheses data structure to file
+// By Diego Arroyuelo
+void saveTree(bp *b, FILE *fp) {
+
+   if (fwrite(&(b->n), sizeof(int), 1, fp) != 1) {
+      printf("Error: cannot save number of parentheses to file\n");
+      exit(1);
+   }
+
+   if (fwrite(b->B, sizeof(pb), (b->n+D-1)/D, fp) != ((b->n+D-1)/D)) {
+      printf("Error: cannot save parentheses sequence to file\n");
+      exit(1);
+   }
+
+   if (fwrite(&(b->opt), sizeof(int), 1, fp) != 1) {
+      printf("Error: cannot save opt in parentheses to file\n");
+      exit(1);
+   }
+}
+
+// loadTree: load parentheses data structure from file
+// By Diego Arroyuelo
+void loadTree(bp *b, FILE *fp) {
+
+   pb *B;
+   int n, opt;
+
+   if (fread(&n, sizeof(int), 1, fp) != 1) {
+      printf("Error: cannot read number of parentheses from file\n");
+      exit(1);
+   }
+
+   mymalloc(B,(n+D-1)/D,0);
+   
+   if (fread(B, sizeof(pb), (n+D-1)/D, fp) != ((n+D-1)/D)) {
+      printf("Error: cannot read parentheses sequence from file\n");
+      exit(1);
+   }
+
+   if (fread(&opt, sizeof(int), 1, fp) != 1) {
+      printf("Error: cannot read opt in parentheses from file\n");
+      exit(1);
+   }
+   
+   bp_construct(b, n, B, opt); 
+   
+}
+
+
+
+int naive_fwd_excess(bp *b,int s, int rel)
+{
+  int i,v,n;
+  pb *B;
+  n = b->n;  B = b->B;
+  v = 0;
+  for (i=s+1; i<n; i++) {
+    if (getbit(B,i)==OP) {
+      v++;
+    } else {
+      v--;
+    }
+    if (v == rel) return i;
+  }
+  return -1;
+}
+
+int naive_bwd_excess(bp *b,int s, int rel)
+{
+  int i,v;
+  pb *B;
+  B = b->B;
+  v = 0;
+  for (i=s; i>=0; i--) {
+    if (getbit(B,i)==OP) {
+      v--;
+    } else {
+      v++;
+    }
+    if (v == rel) return i-1;
+  }
+  return -2;
+}
+
+int naive_search_SB_l(bp *b, int i, int rel)
+{
+  int il,v;
+
+  il = (i / SB) * SB;
+  for (; i>=il; i--) {
+    if (getbit(b->B,i)==OP) {
+      rel++;
+    } else {
+      rel--;
+    }
+    if (rel == 0) return i-1;
+  }
+  if (i < 0) return -2;
+  return -3;
+}
+
+int naive_rmq(bp *b, int s, int t,int opt)
+{
+  int d,i,dm,im;
+
+  if (opt & OPT_RIGHT) {
+    d = dm = depth(b,t);  im = t;
+    i = t-1;
+    while (i >= s) {
+      if (getbit(b->B,i+1)==CP) {
+	d++;
+	if (opt & OPT_MAX) {
+	  if (d > dm) {
+	    dm = d;  im = i;
+	  }
+	}
+      } else {
+	d--;
+	if (!(opt & OPT_MAX)) {
+	  if (d < dm) {
+	    dm = d;  im = i;
+	  }
+	}
+      }
+      i--;
+    }
+  } else {
+    d = dm = depth(b,s);  im = s;
+    i = s+1;
+    while (i <= t) {
+      if (getbit(b->B,i)==OP) {
+	d++;
+	if (opt & OPT_MAX) {
+	  if (d > dm) {
+	    dm = d;  im = i;
+	  }
+	}
+      } else {
+	d--;
+	if (!(opt & OPT_MAX)) {
+	  if (d < dm) {
+	    dm = d;  im = i;
+	  }
+	}
+      }
+      i++;
+    }
+  }
+  return im;
+}
+
+int root_node(bp *b)
+{
+  return 0;
+}
+
+
+int rank_open(bp *b, int s)
+{
+  return darray_rank(b->da,s);
+}
+
+int rank_close(bp *b, int s)
+{
+  return s+1 - darray_rank(b->da,s);
+}
+
+int select_open(bp *b, int s)
+{
+  if (b->opt & OPT_FAST_PREORDER_SELECT) {
+    return darray_select(b->da,s,1);
+  } else {
+    return darray_select_bsearch(b->da,s,getpat_preorder);
+  }
+}
+
+int select_close(bp *b, int s)
+{
+  if (b->opt & OPT_FAST_POSTORDER_SELECT) {
+    return darray_pat_select(b->da_postorder,s,getpat_postorder);
+  } else {
+    return postorder_select_bsearch(b,s);
+  }
+}
+
+///////////////////////////////////////////
+//  find_close(bp *b,int s)
+//    returns the matching close parenthesis of s
+///////////////////////////////////////////
+int find_close(bp *b,int s)
+{
+  return fwd_excess(b,s,-1);
+}
+
+///////////////////////////////////////////
+//  find_open(bp *b,int s)
+//    returns the matching open parenthesis of s
+///////////////////////////////////////////
+int find_open(bp *b,int s)
+{
+  int r;
+  r = bwd_excess(b,s,0);
+  if (r >= -1) return r+1;
+  return -1;
+}
+
+///////////////////////////////////////////
+//  parent(bp *b,int s)
+//    returns the parent of s
+//            -1 if s is the root
+///////////////////////////////////////////
+int parent(bp *b,int s)
+{
+  int r;
+  r = bwd_excess(b,s,-2);
+  if (r >= -1) return r+1;
+  return -1;
+}
+
+int enclose(bp *b,int s)
+{
+  return parent(b,s);
+}
+
+///////////////////////////////////////////
+//  level_ancestor(bp *b,int s,int d)
+//    returns the ancestor of s with relative depth d (d < 0)
+//            -1 if no such node
+///////////////////////////////////////////
+int level_ancestor(bp *b,int s,int d)
+{
+  int r;
+  r = bwd_excess(b,s,d-1);
+  if (r >= -1) return r+1;
+  return -1;
+}
+
+///////////////////////////////////////////
+//  lca(bp *b, int s, int t)
+//    returns the lowest common ancestor of s and t
+///////////////////////////////////////////
+int lca(bp *b, int s, int t)
+{
+  return parent(b,rmq(b,s,t,0)+1);
+}
+
+
+///////////////////////////////////////////
+//  preorder_rank(bp *b,int s)
+//    returns the preorder (>= 1) of node s (s >= 0)
+///////////////////////////////////////////
+int preorder_rank(bp *b,int s)
+{
+  return darray_rank(b->da,s);
+}
+
+///////////////////////////////////////////
+//  preorder_select(bp *b,int s)
+//    returns the node with preorder s (s >= 1)
+//            -1 if no such node
+///////////////////////////////////////////
+int preorder_select(bp *b,int s)
+{
+  //  no error handling
+  if (b->opt & OPT_FAST_PREORDER_SELECT) {
+    return darray_select(b->da,s,1);
+  } else {
+    return darray_select_bsearch(b->da,s,getpat_preorder);
+  }
+}
+
+///////////////////////////////////////////
+//  postorder_rank(bp *b,int s)
+//    returns the postorder (>= 1) of node s (s >= 0)
+//            -1 if s-th bit is not OP
+///////////////////////////////////////////
+int postorder_rank(bp *b,int s)
+{
+  int t;
+  if (inspect(b,s) == CP) return -1;
+  t = find_close(b,s);
+  //  return t+1 - darray_rank(b->da,t);
+  return rank_close(b,t);
+}
+
+int postorder_select_bsearch(bp *b,int s)
+{
+  int l,r,m;
+
+  if (s == 0) return -1;
+
+  if (s > b->da->n - b->da->m) {
+    return -1;
+  }
+  l = 0;  r = b->da->n - 1;
+
+  while (l < r) {
+    m = (l+r)/2;
+    //printf("m=%d rank=%d s=%d\n",m,m+1 - darray_rank(b->da,m),s);
+    if (m+1 - darray_rank(b->da,m) >= s) {
+      r = m;
+    } else {
+      l = m+1;
+    }
+  }
+  return l;
+}
+
+///////////////////////////////////////////
+//  postorder_select(bp *b,int s)
+//    returns the position of CP of the node with postorder s (>= 1)
+///////////////////////////////////////////
+int postorder_select(bp *b,int s)
+{
+#if 0
+  if (b->opt & OPT_FAST_POSTORDER_SELECT) {
+    return darray_pat_select(b->da_postorder,s,getpat_postorder);
+  } else {
+    return postorder_select_bsearch(b->da,s);
+  }
+#else
+  return select_close(b,s);
+#endif
+}
+
+///////////////////////////////////////////
+//  leaf_rank(bp *b,int s)
+//    returns the number of leaves to the left of s
+///////////////////////////////////////////
+int leaf_rank(bp *b,int s)
+{
+  if ((b->opt & OPT_LEAF) == 0) {
+    printf("leaf_rank: error!!!   not supported\n");
+    return -1;
+  }
+  if (s >= b->n-1) {
+    s = b->n-2;
+  }
+  return darray_pat_rank(b->da_leaf,s,getpat_leaf);
+}
+
+///////////////////////////////////////////
+//  leaf_select(bp *b,int s)
+//    returns the position of s-th leaf
+///////////////////////////////////////////
+int leaf_select(bp *b,int s)
+{
+  if ((b->opt & OPT_LEAF) == 0) {
+    printf("leaf_select: error!!!   not supported\n");
+    return -1;
+  }
+  if (s > b->da_leaf->m) return -1;
+  if (b->opt & OPT_FAST_LEAF_SELECT) {
+    return darray_pat_select(b->da_leaf,s,getpat_leaf);
+  } else {
+    return darray_select_bsearch(b->da_leaf,s,getpat_leaf);
+  }
+}
+
+
+///////////////////////////////////////////
+//  inorder_rank(bp *b,int s)
+//    returns the number of ")("  (s >= 0)
+///////////////////////////////////////////
+int inorder_rank(bp *b,int s)
+{
+  if ((b->opt & OPT_INORDER) == 0) {
+    printf("inorder_rank: error!!!   not supported\n");
+    return -1;
+  }
+  if (s >= b->n-1) {
+    s = b->n-2;
+  }
+  return darray_pat_rank(b->da_inorder,s,getpat_inorder);
+}
+
+///////////////////////////////////////////
+//  inorder_select(bp *b,int s)
+//    returns the s-th position of ")("  (s >= 1)
+///////////////////////////////////////////
+int inorder_select(bp *b,int s)
+{
+  if ((b->opt & OPT_INORDER) == 0) {
+    printf("inorder_select: error!!!   not supported\n");
+    return -1;
+  }
+  if (b->opt & OPT_FAST_INORDER_SELECT) {
+    return darray_pat_select(b->da_inorder,s,getpat_inorder);
+  } else {
+    return darray_select_bsearch(b->da_inorder,s,getpat_inorder);
+  }
+}
+
+///////////////////////////////////////////
+//  leftmost_leaf(bp *b, int s)
+///////////////////////////////////////////
+int leftmost_leaf(bp *b, int s)
+{
+  if ((b->opt & OPT_LEAF) == 0) {
+    printf("leftmost_leaf: error!!!   not supported\n");
+    return -1;
+  }
+  return leaf_select(b,leaf_rank(b,s)+1);
+}
+
+///////////////////////////////////////////
+//  rightmost_leaf(bp *b, int s)
+///////////////////////////////////////////
+int rightmost_leaf(bp *b, int s)
+{
+  int t;
+  if ((b->opt & OPT_LEAF) == 0) {
+    printf("leftmost_leaf: error!!!   not supported\n");
+    return -1;
+  }
+  t = find_close(b,s);
+  return leaf_select(b,leaf_rank(b,t));
+}
+
+
+
+///////////////////////////////////////////
+//  inspect(bp *b, int s)
+//    returns OP (==1) or CP (==0) at s-th bit (0 <= s < n)
+///////////////////////////////////////////
+int inspect(bp *b, int s)
+{
+  if (s < 0 || s >= b->n) {
+    printf("inspect: error s=%d is out of [%d,%d]\n",s,0,b->n-1);
+  }
+  return getbit(b->B,s);
+}
+
+int isleaf(bp *b, int s)
+{
+  if (inspect(b,s) != OP) {
+    printf("isleaf: error!!! B[%d] = OP\n",s);
+  }
+  if (inspect(b,s+1) == CP) return 1;
+  else return 0;
+}
+
+
+///////////////////////////////////////////
+//  subtree_size(bp *b, int s)
+//    returns the number of nodes in the subtree of s
+///////////////////////////////////////////
+int subtree_size(bp *b, int s)
+{
+  return (find_close(b,s) - s + 1) / 2;
+}
+
+///////////////////////////////////////////
+//  first_child(bp *b, int s)
+//    returns the first child
+//            -1 if s is a leaf
+///////////////////////////////////////////
+int first_child(bp *b, int s)
+{
+  if (inspect(b,s+1) == CP) return -1;
+  return s+1;
+}
+
+///////////////////////////////////////////
+//  next_sibling(bp *b,int s)
+//    returns the next sibling of parent(s)
+//            -1 if s is the last child
+//////////////////////////////////////////
+int next_sibling(bp *b, int s)
+{
+  int t;
+  t = find_close(b,s)+1;
+  if (t >= b->n) {
+    printf("next_sibling: error s=%d t=%d\n",s,t);
+  }
+  if (inspect(b,t) == CP) return -1;
+  return t;
+}
+
+///////////////////////////////////////////
+//  prev_sibling(bp *b,int s)
+//    returns the previous sibling of parent(s)
+//            -1 if s is the first child
+//////////////////////////////////////////
+int prev_sibling(bp *b, int s)
+{
+  int t;
+  if (s < 0) {
+    printf("prev_sibling: error s=%d\n",s);
+  }
+  if (s == 0) return -1;
+  if (inspect(b,s-1) == OP) return -1;
+  t = find_open(b,s-1);
+  return t;
+}
+
+///////////////////////////////////////////
+//  deepest_node(bp *b,int s)
+//    returns the first node with the largest depth in the subtree of s
+///////////////////////////////////////////
+int deepest_node(bp *b,int s)
+{
+  int t,m;
+  t = find_close(b,s);
+  m = rmq(b,s,t, OPT_MAX);
+  return m;
+}
+
+///////////////////////////////////////////
+//  subtree_height(bp *b,int s)
+//    returns the height of the subtree of s
+//            0 if s is a leaf
+///////////////////////////////////////////
+int subtree_height(bp *b,int s)
+{
+  int t;
+  t = deepest_node(b,s);
+  return depth(b,t) - depth(b,s);
+}
+
+int naive_degree(bp *b, int s)
+{
+  int t,d;
+  d = 0;
+  t = first_child(b,s);
+  while (t >= 0) {
+    d++;
+    t = next_sibling(b,t);
+  }
+  return d;
+}
+
+///////////////////////////////////////////
+//  degree(bp *b, int s)
+//    returns the number of children of s
+//            0 if s is a leaf
+///////////////////////////////////////////
+int degree(bp *b, int s)
+{
+  if (b->opt & OPT_DEGREE) {
+    return fast_degree(b,s,b->n,0);
+  } else {
+    return naive_degree(b,s);
+  }
+}
+
+int naive_child(bp *b, int s, int d)
+{
+  int t,i;
+  t = first_child(b,s);
+  for (i = 1; i < d; i++) {
+    if (t == -1) break;
+    t = next_sibling(b,t);
+  }
+  return t;
+}
+
+///////////////////////////////////////////
+//  child(bp *b, int s, int d)
+//    returns the d-th child of s  (1 <= d <= degree(s))
+//            -1 if no such node
+///////////////////////////////////////////
+int child(bp *b, int s, int d)
+{
+  int r;
+  if (b->opt & OPT_DEGREE) {
+    //return find_open(b,fast_degree(b,s,b->n,d));
+    if (d==1) return first_child(b,s);
+    r = fast_degree(b,s,b->n,d-1)+1;
+    if (inspect(b,r) == CP) return -1;
+    return r;
+  } else {
+    return naive_child(b,s,d);
+  }
+    
+}
+
+
+int naive_child_rank(bp *b, int t)
+{
+  int v,d;
+  d = 0;
+  while (t != -1) {
+    d++;
+    t = prev_sibling(b,t);
+  }
+  return d;
+}
+
+///////////////////////////////////////////
+//  child_rank(bp *b, int t)
+//    returns d if t is the d-th child of the parent of t (d >= 1)
+//            1 if t is the root
+///////////////////////////////////////////
+int child_rank(bp *b, int t)
+{
+  int r;
+  if (t == root_node(b)) return 1;
+  if (b->opt & OPT_DEGREE) {
+    r = parent(b,t);
+    return fast_degree(b,r,t,0)+1;
+  } else {
+    return naive_child_rank(b,t);
+  }
+}
+
+
+
+///////////////////////////////////////////
+//  is_ancestor(bp *b, int s, int t)
+//     returns 1 if s is an ancestor of t
+//             0 otherwise
+///////////////////////////////////////////
+int is_ancestor(bp *b, int s, int t)
+{
+  int v;
+  v = find_close(b,s);
+  if (s <= t && t <= v) return 1;
+  return 0;
+}
+
+///////////////////////////////////////////
+//  distance(bp *b, int s, int t)
+//    returns the length of the shortest path from s to t in the tree
+///////////////////////////////////////////
+int distance(bp *b, int s, int t)
+{
+  int v,d;
+  v = lca(b,s,t);
+  d = depth(b,v);
+  return (depth(b,s) - d) + (depth(b,t) - d);
+}
+
+///////////////////////////////////////////
+//  level_next(bp *b, int d)
+///////////////////////////////////////////
+int level_next(bp *b,int s)
+{
+  int t;
+  t = fwd_excess(b,s,0);
+  return t;
+}
+
+///////////////////////////////////////////
+//  level_prev(bp *b, int d)
+///////////////////////////////////////////
+int level_prev(bp *b,int s)
+{
+  int t;
+  t = bwd_excess(b,s,0);
+  return t;
+}
+
+///////////////////////////////////////////
+//  level_leftmost(bp *b, int d)
+///////////////////////////////////////////
+int level_leftmost(bp *b, int d)
+{
+  int t;
+  if (d < 1) return -1;
+  if (d == 1) return 0;
+  t = fwd_excess(b,0,d);
+  return t;
+}
+
+///////////////////////////////////////////
+//  level_rigthmost(bp *b, int d)
+///////////////////////////////////////////
+int level_rigthmost(bp *b, int d)
+{
+  int t;
+  if (d < 1) return -1;
+  if (d == 1) return 0;
+  t = bwd_excess(b,0,d-1);
+  return find_open(b,t);
+}
+
+///////////////////////////////////////////
+//  leaf_size(bp *b, int s)
+///////////////////////////////////////////
+int leaf_size(bp *b, int s)
+{
+  int t;
+  if ((b->opt & OPT_LEAF) == 0) {
+    printf("leaf_size: error!!!   not supported\n");
+    return -1;
+  }
+  t = find_close(b,s);
+  return leaf_rank(b,t) - leaf_rank(b,s);
+}