hi

fastlib/u/garryb/nbr/allnn.cc

@@ -2,6 +2,7 @@
 #include "spbounds.h"
 #include "gnp.h"
 #include "dfs.h"
+#include "nbr_utils.h"
 
 /**
  * An N-Body-Reduce problem.
@@ -28,6 +29,10 @@ class Allnn {
     }
 
    public:
+    void Copy(const Param& other) {
+      dim = other.dim;
+    }
+   
     /**
      * Initialize parameters from a data node (Req NBR).
      */
@@ -63,9 +68,7 @@ class Allnn {
     }
 
    public:
-    void Init(const Param& param,
-        const Vector& q_point, const QPointInfo& q_info,
-        const RNode& r_root) {
+    void Init(const Param& param) {
       distance_sq = DBL_MAX;
       neighbor_i = -1;
     }
@@ -233,9 +236,8 @@ class Allnn {
 int main(int argc, char *argv[]) {
   fx_init(argc, argv);
   
-  DualTreeDepthFirst<Allnn::GNP> dfs;
-  dfs.Init(fx_root);
-  dfs.Begin();
+  nbr_utils::SerialDualTreeMain<Allnn::GNP, DualTreeDepthFirst<Allnn::GNP> >(
+      fx_root, "allnn");
   
   fx_done();
 }

fastlib/u/garryb/nbr/blockdev.cc

@@ -1,16 +1,22 @@
+#include "blockdev.h"
+
+#include <sys/types.h>
+#include <sys/stat.h>
+#include <unistd.h>
+#include <fcntl.h>
 
 void RandomAccessFile::Init(const char *fname, BlockDevice::mode_t mode) {
   int octal_mode;
   
   switch(mode) {
-   BlockDevice::READ:
+   case BlockDevice::READ:
     octal_mode = O_RDONLY;
     break;
-   BlockDevice::MODIFY:
+   case BlockDevice::MODIFY:
     octal_mode = O_RDWR;
     break;
-   BlockDevice::CREATE:
-   BlockDevice::TEMP:
+   case BlockDevice::CREATE:
+   case BlockDevice::TEMP:
     octal_mode = O_RDWR|O_CREAT|O_TRUNC;
     break;
    default: abort();
@@ -22,7 +28,7 @@ void RandomAccessFile::Init(const char *fname, BlockDevice::mode_t mode) {
   }
 }
 
-void RandomAccessFIle::Close() {
+void RandomAccessFile::Close() {
   close(fd_);
 }
 
@@ -42,7 +48,7 @@ void RandomAccessFile::Write(off_t pos, size_t len, const char *buffer) {
       break;
     }
 
-    DEBUG_ASSERT_MSG(written > 0, "error writing %lu bytes", len);
+    DEBUG_ASSERT_MSG(written > 0, "error writing");
 
     buffer += written;
   }
@@ -105,7 +111,7 @@ void DiskBlockDevice::Write(blockid_t blockid,
   file_.Write(off_t(blockid) * n_block_bytes_ + begin, end - begin, data);
 }
 
-blockid_t DiskBlockDevice::AllocBlock() {
+DiskBlockDevice::blockid_t DiskBlockDevice::AllocBlock() {
   blockid_t blockid = n_blocks_;
   n_blocks_ = blockid + 1;
   return blockid;

fastlib/u/garryb/nbr/blockdev.h

@@ -6,6 +6,7 @@
 class BlockActionHandler {
   FORBID_COPY(BlockActionHandler);
  public:
+  BlockActionHandler() {}
   virtual ~BlockActionHandler() {}
   
   virtual void BlockInitFrozen(size_t bytes, char *block) = 0;
@@ -66,7 +67,8 @@ class BlockDevice {
 };
 
 class NullBlockDevice : public BlockDevice {
-  FORBID_COPY(BlankBlockDevice);
+  FORBID_COPY(NullBlockDevice);
+  
  public:
   void Init(blockid_t n_blocks_in, offset_t n_block_bytes_in) {
     n_blocks_ = n_blocks_in;

fastlib/u/garryb/nbr/build.py

@@ -1,15 +1,17 @@
 
+librule(name = "nbr",
+   sources = ["blockdev.cc", "nbr_utils.cc", "cache.cc"],
+   headers = ["blockdev.h", "dfs.h", "nbr_utils.h", "spbounds.h",
+              "cache.h", "gnp.h", "kdtree.h", "spnode.h"],
+   deplibs = ["fastlib:fastlib_int"])
+
 binrule(name = "tkde",
    sources = ["tkde.cc"],
-   headers = ["spnode.h", "spbounds.h", "kdtree.h", "gnp.h", "dfs.h"],
-   deplibs = ["fastlib:fastlib_int"])
-
-
+   deplibs = [":nbr"])
 
 binrule(name = "allnn",
    sources = ["allnn.cc"],
-   headers = ["spnode.h", "spbounds.h", "kdtree.h", "gnp.h", "dfs.h"],
-   deplibs = ["fastlib:fastlib_int"])
+   deplibs = [":nbr"])
 
 
 

fastlib/u/garryb/nbr/cache.cc

@@ -1,7 +1,9 @@
+#include "cache.h"
+
 void SmallCache::Init(BlockDevice *inner_in, BlockActionHandler *handler_in,
     mode_t mode_in) {
   BlockDeviceWrapper::Init(inner_in);
-  metadata_.Init(inner_in->n_blocks);
+  metadata_.Init(inner_in->n_blocks());
   handler_ = handler_in;
   mode_ = mode_in;
 }
@@ -53,9 +55,9 @@ void SmallCache::PerformCacheMiss_(blockid_t blockid) {
   metadata->data = data;
 }
 
-void SmallCache::Writeback_(blockid_block, offset_t begin, offset_t end) {
+void SmallCache::Writeback_(blockid_t blockid, offset_t begin, offset_t end) {
   if (begin != end) {
-    Metadata *metadata = metadata_[blockid];
+    Metadata *metadata = &metadata_[blockid];
     char *data = metadata->data;
     
     if (data) {
@@ -67,7 +69,7 @@ void SmallCache::Writeback_(blockid_block, offset_t begin, offset_t end) {
       }
       
       handler_->BlockRefreeze(n_bytes, buf, buf);
-      inner_->Write(block, begin, end, buf);
+      inner_->Write(blockid, begin, end, buf);
       handler_->BlockThaw(n_bytes, buf);
     }
   }
@@ -81,7 +83,7 @@ void SmallCache::Flush(blockid_t begin_block, offset_t begin_offset,
     if (begin_block == last_block) {
       Writeback_(begin_block, begin_offset, end_offset);
     } else {
-      Writeback_(begin_block, begin_offset, n_block_bytes_ - begin_offest);
+      Writeback_(begin_block, begin_offset, n_block_bytes_ - begin_offset);
       for (blockid_t i = begin_block + 1; i < last_block - 1; i++) {
         Writeback_(i, 0, n_block_bytes_);
       }
@@ -90,7 +92,7 @@ void SmallCache::Flush(blockid_t begin_block, offset_t begin_offset,
   }
 }
 
-void SmallCache::Read(blockid_t block,
+void SmallCache::Read(blockid_t blockid,
     offset_t begin, offset_t end, char *buf) {
   const char *src_buffer = StartRead(blockid) + begin;
   size_t n_bytes = end - begin;
@@ -99,7 +101,7 @@ void SmallCache::Read(blockid_t block,
   StopRead(blockid);
 }
 
-void SmallCache::Write(blockid_t block,
+void SmallCache::Write(blockid_t blockid,
     offset_t begin, offset_t end, const char *buf) {
   char *dest_buffer = StartWrite(blockid) + begin;
   size_t n_bytes = end - begin;

fastlib/u/garryb/nbr/cache.h

@@ -40,13 +40,15 @@ class SmallCache : public BlockDeviceWrapper {
   char *StartWrite(blockid_t blockid);
   void StopRead(blockid_t blockid);
   void StopWrite(blockid_t blockid);
+  void Flush(blockid_t begin_block, offset_t begin_offset,
+    blockid_t last_block, offset_t end_offset);
   void Close();
 
   virtual void Read(blockid_t blockid,
       offset_t begin, offset_t end, char *data);
   virtual void Write(blockid_t blockid,
       offset_t begin, offset_t end, const char *data);
-  
+
   virtual blockid_t AllocBlock() {
     blockid_t blockid = BlockDeviceWrapper::AllocBlock();
     metadata_.Resize(n_blocks());
@@ -55,7 +57,7 @@ class SmallCache : public BlockDeviceWrapper {
 
  private:
   void PerformCacheMiss_(blockid_t blockid);
-  void Writeback_(blockid_block, offset_t begin, offset_t end);
+  void Writeback_(blockid_t blockid, offset_t begin, offset_t end);
 
   Metadata *GetBlock_(blockid_t blockid) {
     Metadata *metadata = &metadata_[blockid];
@@ -79,13 +81,13 @@ class CacheArrayBlockActionHandler : public BlockActionHandler {
 
  public:
   CacheArrayBlockActionHandler() {
-    n_block_elems_ = BIG_BAD_NUMBER;
+    n_elem_bytes_ = BIG_BAD_NUMBER;
     DEBUG_POISON_PTR(default_elem_);
   }
 
   ~CacheArrayBlockActionHandler() {
     mem::Free(default_elem_);
-    n_block_elems_ = BIG_BAD_NUMBER;
+    n_elem_bytes_ = BIG_BAD_NUMBER;
     DEBUG_POISON_PTR(default_elem_);
   }
 
@@ -96,7 +98,7 @@ class CacheArrayBlockActionHandler : public BlockActionHandler {
   }
 
   void BlockInitFrozen(size_t bytes, char *block) {
-    index_t elems = bytes / n_elem_bytes;
+    index_t elems = bytes / n_elem_bytes_;
     for (index_t i = 0; i < elems; i++) {
       mem::CopyBytes(block, default_elem_, n_elem_bytes_);
       block += n_elem_bytes_;
@@ -104,15 +106,16 @@ class CacheArrayBlockActionHandler : public BlockActionHandler {
   }
 
   void BlockRefreeze(size_t bytes, const char *old_location, char *block) {
-    index_t elems = bytes / n_elem_bytes;
+    index_t elems = bytes / n_elem_bytes_;
     for (index_t i = 0; i < elems; i++) {
-      ot::PointerRefreeze(reinterpret_cast<const T*>(old_location), dest);
+      ot::PointerRefreeze(reinterpret_cast<const T*>(old_location), block);
       block += n_elem_bytes_;
+      old_location += n_elem_bytes_;
     }
   }
 
   void BlockThaw(size_t bytes, char *block) {
-    index_t elems = bytes / n_elem_bytes;
+    index_t elems = bytes / n_elem_bytes_;
     for (index_t i = 0; i < elems; i++) {
       ot::PointerThaw<T>(block);
       block += n_elem_bytes_;
@@ -132,7 +135,6 @@ class CacheArray {
 
  public:
   typedef T Element;
-  typedef PointerType<T, AllowWrite>::Type;
 
  private:
   struct Metadata {
@@ -160,18 +162,18 @@ class CacheArray {
   ~CacheArray() {}
 
   /** Reopens another cache array */
-  void Init(CacheArray *other, mode_t mode_in) {
+  void Init(CacheArray *other, BlockDevice::mode_t mode_in) {
     Init(other, mode_in, other->begin_index(), other->end_index());
   }
 
   /** Reopens another cache array */
-  void Init(CacheArray *other, mode_t mode_in, index_t begin_index_in,
-      index_t end_index_in) {
+  void Init(CacheArray *other, BlockDevice::mode_t mode_in,
+      index_t begin_index_in, index_t end_index_in) {
     Init(other->cache_, mode_in, begin_index_in, end_index_in,
         other->n_block_elems_, other->n_elem_bytes_);
   }
 
-  void Init(SmallCache *cache_in, Mode mode_in,
+  void Init(SmallCache *cache_in, BlockDevice::mode_t mode_in,
       index_t begin_index_in, index_t end_index_in,
       index_t n_block_elems_in, size_t n_elem_bytes_in) {
     cache_ = cache_in;
@@ -181,7 +183,7 @@ class CacheArray {
     n_elem_bytes_ = n_elem_bytes_in;
     
     begin_block_ = begin_ / n_block_elems_;
-    last_block_ = (end_ + n_block_elems_ - 1) / n_block_elems_;
+    end_block_ = (end_ + n_block_elems_ - 1) / n_block_elems_;
     
     metadata_.Init(end_block_ - begin_block_);
     
@@ -212,26 +214,26 @@ class CacheArray {
   }
 
   Element *StartWrite(index_t element_id) {
-    DBUG_ASSERT(mode_ != BlockDevice::READ);
+    DEBUG_ASSERT(mode_ != BlockDevice::READ);
     return CheckoutElement_(element_id);
   }
 
-  void StopRead(const Element *ptr, index_t element_id) {
+  void StopRead(index_t element_id) {
     DEBUG_ONLY(BoundsCheck_(element_id));
-    DEBUG_ASSERT(mode_ >= TEMP);
+    DEBUG_ASSERT(mode_ != BlockDevice::READ);
     ReleaseElement_(element_id);
   }
 
-  void StopWrite(Element *ptr, index_t element_id) {
+  void StopWrite(index_t element_id) {
     DEBUG_ONLY(BoundsCheck_(element_id));
-    DEBUG_ASSERT(mode_ >= TEMP);
+    DEBUG_ASSERT(mode_ != BlockDevice::READ);
     ReleaseElement_(element_id);
   }
 
   void Swap(index_t index_a, index_t index_b) {
     DEBUG_ONLY(BoundsCheck_(index_a));
     DEBUG_ONLY(BoundsCheck_(index_b));
-    DEBUG_ASSERT(mode_ >= TEMP);
+    DEBUG_ASSERT(mode_ != BlockDevice::READ);
     char *a = reinterpret_cast<char*>(CheckoutElement_(index_a));
     char *b = reinterpret_cast<char*>(CheckoutElement_(index_b));
     mem::Swap(a, b, n_elem_bytes_);
@@ -242,18 +244,34 @@ class CacheArray {
   void Copy(index_t index_src, index_t index_dest) {
     DEBUG_ONLY(BoundsCheck_(index_src));
     DEBUG_ONLY(BoundsCheck_(index_dest));
-    DEBUG_ASSERT(mode_ >= TEMP);
-    const char *src = reinterpret_cast<char*>(CheckoutElement_(index_a));
-    char *dest = reinterpret_cast<char*>(CheckoutElement_(index_b));
+    DEBUG_ASSERT(mode_ != BlockDevice::READ);
+    const char *src = reinterpret_cast<char*>(CheckoutElement_(index_src));
+    char *dest = reinterpret_cast<char*>(CheckoutElement_(index_dest));
     mem::Copy(dest, src, n_elem_bytes_);
-    ReleaseElement_(index_a);
-    ReleaseElement_(index_b);
+    ReleaseElement_(index_src);
+    ReleaseElement_(index_dest);
   }
 
   /**
    * Flushes all changes.
    */
   void Flush();
+  
+  index_t Alloc() {
+    ++end_index;
+    BlockDevice::blockid_t block = 
+        (end_index + n_block_elems_ - 1) / n_block_elems_;
+    if (block != end_block_) {
+      end_block_ = block;
+      metadata_.Resize(end_block_ - begin_block_);
+      // Okay, notify the lower layers we're allocating.
+      index_t block_allocated = cache_->AllocBlock();
+      DEBUG_ASSERT_MSG(block_allocated == end_block_,
+          "Distributed data structure creation "
+          "is not yet supported by CacheArray.");
+    }
+    return end_index - 1;
+  }
 
  private:
   void BoundsCheck_(index_t element_id) {
@@ -275,7 +293,7 @@ class CacheArray {
     DEBUG_ONLY(metadata_->lock_count++);
     Element *ptr = reinterpret_cast<Element*>(metadata_->data + offset);
 
-    if (unlikely(!metadata_->data)) {
+    if (unlikely(metadata_->data == NULL)) {
       return HandleCacheMiss_(element_id);
     } else {
       return ptr;
@@ -290,12 +308,13 @@ class CacheArray {
 
 template<typename T>
 void CacheArray<T>::Flush() {
-  for (blockid_t block = begin_block_; block < end_block_; block++) {
-    if (metadata_[block - begin_block_].data) {
+  for (BlockDevice::blockid_t blockid = begin_block_;
+      blockid < end_block_; blockid++) {
+    if (metadata_[blockid - begin_block_].data) {
       if (mode_ != BlockDevice::READ) {
-        cache_->StopWrite(block);
+        cache_->StopWrite(blockid);
       } else {
-        cache_->StopRead(block);
+        cache_->StopRead(blockid);
       }
     }
   }
@@ -322,30 +341,6 @@ typename CacheArray<T>::Element* CacheArray<T>::HandleCacheMiss_(
   return reinterpret_cast<Element*>(metadata->data + offset);
 }
 
-template<typename T>
-class CacheAutoPtrConst {
- private:
-  const T* ptr_;
-  index_t i;
-  
- public:
-  CacheAutoPtrConst(const CacheArray<T>& array, ) {
-    ptr_ = ;
-  }
-  ~CacheAutoPtrConst() {
-    arrayi;
-  }
-  
-  operator const T* () {
-    return ptr_;
-  }
-  
-  const T* ptr() const {
-    return ptr_;
-  }
-};
-
-
 template<typename T>
 class TempCacheArray : public CacheArray<T> {
  private:
@@ -358,20 +353,19 @@ class TempCacheArray : public CacheArray<T> {
 
   /** Creates a blank, temporary cached array */
   void Init(const T& default_obj,
-      index_t n_elems,
-      unsigned int n_block_elems) {
-    CacheArrayBlockActionHandler<QMutableInfo> *handler
-        = new CacheArrayBlockActionHandler<QMutableInfo>;
+      index_t n_elems_in,
+      unsigned int n_block_elems_in) {
+    CacheArrayBlockActionHandler<T> *handler =
+        new CacheArrayBlockActionHandler<T>;
     handler->Init(default_obj);
 
-    null_device_.Init((n_elems + n_block_elems + 1) / n_block_elems,
-        n_block_elems);
+    null_device_.Init((n_elems_in + n_block_elems_in + 1) / n_block_elems_in,
+        n_block_elems_in * handler->n_elem_bytes());
     underlying_cache_.Init(&null_device_, handler, BlockDevice::TEMP);
     
     CacheArray<T>::Init(&underlying_cache_, BlockDevice::TEMP,
-        0, q_nodes_.end_index(), q_nodes_.n_block_elems(),
-        mutables_handler->n_elem_bytes());
+        0, n_elems_in, n_block_elems_in, handler->n_elem_bytes());
   }
-}
+};
 
 #endif

fastlib/u/garryb/nbr/dfs.h

@@ -38,8 +38,21 @@ class DualTreeDepthFirst {
   uint64 n_recurse_;
   
  public:
-  void Init(datanode *datanode);
+  void Init(
+      struct datanode *datanode_in,
+      const typename GNP::Param& param_in,
+      CacheArray<typename GNP::Point> *q_points,
+      CacheArray<typename GNP::QPointInfo> *q_point_infos,
+      CacheArray<typename GNP::QNode> *q_nodes,
+      CacheArray<typename GNP::Point> *r_points,
+      CacheArray<typename GNP::RPointInfo> *r_point_infos,
+      CacheArray<typename GNP::RNode> *r_nodes,
+      CacheArray<typename GNP::QResult> *q_results);
   void Begin();
+  
+  const typename GNP::GlobalResult& global_result() const {
+    return global_result_;
+  }
 
  private:
   void Pair_(index_t q_node_i, index_t r_node_i,
@@ -56,6 +69,7 @@ class DualTreeDepthFirst {
 template<typename GNP>
 void DualTreeDepthFirst<GNP>::Init(
     struct datanode *datanode_in,
+    const typename GNP::Param& param_in,
     CacheArray<typename GNP::Point> *q_points,
     CacheArray<typename GNP::QPointInfo> *q_point_infos,
     CacheArray<typename GNP::QNode> *q_nodes,
@@ -63,7 +77,7 @@ void DualTreeDepthFirst<GNP>::Init(
     CacheArray<typename GNP::RPointInfo> *r_point_infos,
     CacheArray<typename GNP::RNode> *r_nodes,
     CacheArray<typename GNP::QResult> *q_results) {
-  param_ = WALDO;
+  param_.Copy(param_in);
   
   q_points_.Init(q_points, BlockDevice::READ);
   q_point_infos_.Init(q_point_infos, BlockDevice::READ);
@@ -78,75 +92,10 @@ void DualTreeDepthFirst<GNP>::Init(
   default_mutable.postponed.Init(param_);
   q_mutables_.Init(default_mutable, q_nodes_.n_block_elems());
   
-  datanode_ = datanode_in;
-  do_naive_ = fx_param_bool(datanode_, "do_naive", false);
-}
-
-template<typename GNP, typename Solver>
-void DualTreeMain(datanode *datanode) {
-  datanode_ = datanode;
-
-  fx_timer_start(datanode_, "q_matrix");
-  cache_utils::Load(fx_param_str(datanode, "q"), );
-  fx_timer_stop(datanode_, "q_matrix");
-
-  TempCacheArray<typename GNP::Point> q_points;
-  TempCacheArray<typename GNP::QPointInfo> q_point_infos;
-  TempCacheArray<typename GNP::QNode> q_nodes;
-  TempCacheArray<typename GNP::Point> r_points;
-  TempCacheArray<typename GNP::RPointInfo> r_point_infos;
-  TempCacheArray<typename GNP::RNode> r_nodes;
-  TempCacheArray<typename GNP::QResult> q_results;
-  
-  Solver solver;
-  solver.Init(&q_points, &q_point_infos, &q_nodes,
-      &r_points, &r_point_infos, &r_nodes, &q_result);
-  
-  fx_timer_start(datanode_, "q_matrix");
-  Matrix q_matrix;
-  ASSERT_PASS(data::Load(fx_param_str_req(datanode, "q"), &q_matrix));
-  fx_timer_stop(datanode_, "q_matrix");
-  
-  fx_timer_start(datanode_, "r_matrix");
-  Matrix r_matrix;
-  ASSERT_PASS(data::Load(fx_param_str_req(datanode, "r"), &r_matrix));
-  fx_timer_stop(datanode_, "r_matrix");
-  
-  param_.Init(fx_submodule(datanode, "param", "param"),
-      q_matrix, r_matrix);
-  
-  ArrayList<typename GNP::QPointInfo> q_point_info;
-  q_point_info.Init(q_matrix.n_cols());
-  // TODO: Read info?
-  
-  ArrayList<typename GNP::RPointInfo> r_point_info;
-  r_point_info.Init(r_matrix.n_cols());
-  // TODO: Read info?
-  
-  fx_timer_start(datanode_, "q_tree");
-  q_tree_.Init(&param_, q_matrix, q_point_info);
-  q_tree_.Build();
-  fx_timer_stop(datanode_, "q_tree");
-  
-  fx_timer_start(datanode_, "r_tree");
-  r_tree_.Init(&param_, r_matrix, r_point_info);
-  r_tree_.Build();
-  fx_timer_stop(datanode_, "r_tree");
-  
-  q_results_.Init(q_tree_.points().size());
-  for (index_t i = 0; i < q_tree_.points().size(); i++) {
-    q_results_[i].Init(param_, q_tree_.points()[i], q_tree_.point_info()[i],
-        r_tree_.nodes()[0]);
-  }
-  
-  q_mutables_.Init(q_tree_.nodes().size());
-  for (index_t i = 0; i < q_tree_.nodes().size(); i++) {
-    q_mutables_[i].mass_result.Init(param_);
-    q_mutables_[i].postponed.Init(param_);
-  }
-  
   global_result_.Init(param_);
   
+  datanode_ = datanode_in;
+  do_naive_ = fx_param_bool(datanode_, "do_naive", false);
 }
 
 template<typename GNP>
@@ -178,13 +127,13 @@ void DualTreeDepthFirst<GNP>::Begin() {
       PushDown_(0);
     }
   }
-  
+
   q_nodes_.StopRead(0);
   q_mutables_.StopWrite(0);
   r_nodes_.StopRead(0);
-  
+
   fx_timer_stop(datanode_, "execute");
-  
+
   DEBUG_ONLY(fx_format_result(datanode_, "naive_ratio", "%f",
       1.0 * n_naive_ / q_root->count() / r_root->count()));
   DEBUG_ONLY(fx_format_result(datanode_, "naive_per_query", "%f",
@@ -197,10 +146,13 @@ void DualTreeDepthFirst<GNP>::Begin() {
       1.0 * n_recurse_ / q_root->count() / r_root->count()));
   DEBUG_ONLY(fx_format_result(datanode_, "recurse_per_query", "%f",
       1.0 * n_recurse_ / q_root->count()));
-  
+
   if (fx_param_bool(datanode_, "print", 0)) {
     ot::Print(q_results_);
   }
+
+  q_mutables_.Flush();
+  q_results_.Flush();
 }
 
 template<typename GNP>
@@ -312,9 +264,9 @@ void DualTreeDepthFirst<GNP>::Pair_(index_t q_node_i, index_t r_node_i,
       const typename GNP::RNode *r_child2 = r_nodes_.StartRead(r_child2_i);
       
       double r_child1_h = GNP::Algorithm::Heuristic(
-          param_, *q_node, *rnode_(r_child1_i));
+          param_, *q_node, *r_child1);
       double r_child2_h = GNP::Algorithm::Heuristic(
-          param_, *q_node, *rnode_(r_child2_i));
+          param_, *q_node, *r_child2);
 
       if (unlikely(r_child2_h < r_child1_h)) {
         const typename GNP::RNode *r_child_t = r_child1;
@@ -367,12 +319,12 @@ void DualTreeDepthFirst<GNP>::BaseCase_(
     const typename GNP::QMassResult& exclusive_unvisited,
     QMutableInfo *q_node_mut) {
   typename GNP::PairVisitor visitor;
-  const typename GNP::Point *r_points[r_node->count()];
-  const typename GNP::RPointInfo *r_infos[r_node->count()];
+  const typename GNP::Point *r_local_points[r_node->count()];
+  const typename GNP::RPointInfo *r_local_infos[r_node->count()];
 
   for (index_t r_i_rel = 0; r_i_rel < r_node->count(); ++r_i_rel) {
-    r_points[r_i_rel] = r_points_.StartRead(r_i_rel + r_node->begin());
-    r_infos[r_i_rel] = r_point_infos_.StartRead(r_i_rel + r_node->begin());
+    r_local_points[r_i_rel] = r_points_.StartRead(r_i_rel + r_node->begin());
+    r_local_infos[r_i_rel] = r_point_infos_.StartRead(r_i_rel + r_node->begin());
   }
 
   DEBUG_ONLY(n_pre_naive_ += q_node->count() * r_node->count());
@@ -380,8 +332,8 @@ void DualTreeDepthFirst<GNP>::BaseCase_(
   visitor.Init(param_);
 
   q_node_mut->mass_result.StartReaccumulate(param_, *q_node);
-  DEBUG_ASSERT_MSG(q_node->count() != 0, "index %d, count = %d, begin = %d",
-      int(q_node - qnode_(0)), q_node->begin(), q_node->count());
+  //DEBUG_ASSERT_MSG(q_node->count() != 0, "index %d, count = %d, begin = %d",
+  //    int(q_node - qnode_(0)), q_node->begin(), q_node->count());
 
   for (index_t q_i = q_node->begin(); q_i < q_node->end(); ++q_i) {
     const typename GNP::Point *q_point = q_points_.StartRead(q_i);
@@ -396,7 +348,7 @@ void DualTreeDepthFirst<GNP>::BaseCase_(
 
       for (index_t r_i_rel = 0; r_i_rel < r_count; ++r_i_rel) {
         visitor.VisitPair(param_, *q_point, *q_info, q_i,
-            *r_point, *r_info, r_i_rel + r_node->begin());
+            r_local_points[i], r_local_infos[i], r_i_rel + r_node->begin());
       }
 
       visitor.FinishVisitingQueryPoint(param_, *q_point, *q_info, *r_node,
@@ -413,8 +365,8 @@ void DualTreeDepthFirst<GNP>::BaseCase_(
   }
 
   for (index_t r_i_rel = 0; r_i_rel < r_node->count(); ++r_i_rel) {
-    r_points_.StopRead(r_i_rel + r_node->begin());
-    r_point_infos_.StopRead(r_i_rel + r_node->begin());
+    r_local_points.StopRead(r_i_rel + r_node->begin());
+    r_local_point_infos.StopRead(r_i_rel + r_node->begin());
   }
 
   q_node_mut->mass_result.FinishReaccumulate(param_, *q_node);

fastlib/u/garryb/nbr/kdtree.h

@@ -13,6 +13,7 @@
 
 #include "spnode.h"
 #include "spbounds.h"
+#include "cache.h"
 
 #include "base/common.h"
 #include "col/arraylist.h"
@@ -28,29 +29,27 @@ class KdTreeMidpointBuilder {
   typedef TNode Node;
   typedef typename TNode::Bound Bound;
   typedef TParam Param;
-  
+
  private:
   const Param* param_;
   CacheArray<Vector> points_;
   CacheArray<PointInfo> point_infos_;
   CacheArray<Node> nodes_;
-  CacheArrayAllocator allocator_;
   index_t leaf_size_;
   index_t dim_;
-  
+
  public:
   void InitBuild(
       struct datanode *module,
       const Param* param_in_,
-      CacheArray<Vector> *points_in,
-      CacheArray<PointInfo> *point_infos_in,
+      CacheArray<Vector> *points_inout,
+      CacheArray<PointInfo> *point_infos_inout,
       CacheArray<Node> *nodes_out) {
     param_ = param_in_;
 
-    points_.Init(points_in, BlockDevice::MODIFY);
-    point_infos_.Init(nodes_, BlockDevice::MODIFY);
-    nodes_.Init(nodes_out_, BlockDevice::CREATE);
-    allocator_.Init(&nodes_);
+    points_.Init(points_inout, BlockDevice::MODIFY);
+    point_infos_.Init(point_infos_inout, BlockDevice::MODIFY);
+    nodes_.Init(nodes_out, BlockDevice::CREATE);
 
     const Vector *first_point = points_.StartRead(points_.begin_index());
     dim_ = first_point->length();
@@ -111,7 +110,7 @@ index_t KdTreeMidpointBuilder<TPointInfo, TNode, TParam>::Partition_(
       left++;
     }
 
-    while (1)
+    while (1) {
       right_v = points_.StartWrite(right);
       if (right_v->get(split_dim) < splitvalue || unlikely(left > right)) {
         break;
@@ -125,7 +124,7 @@ index_t KdTreeMidpointBuilder<TPointInfo, TNode, TParam>::Partition_(
       break;
     }
 
-    *left_v.SwapValues(right_v);
+    left_v->SwapValues(right_v);
     // TODO: If point info has pointers this will incur bad cache performance
     // In the future we rely on OT frozen storage
     point_infos_.Swap(left, right);
@@ -171,8 +170,8 @@ void KdTreeMidpointBuilder<TPointInfo, TNode, TParam>::KdTreeMidpointBuilder::Bu
     double split_val = node->bound().get(split_dim).mid();
 
     if (max_width != 0) {
-      index_t left_i = allocator_.Alloc();
-      index_t right_i = allocator_.Alloc();
+      index_t left_i = nodes_.Alloc();
+      index_t right_i = nodes_.Alloc();
       Node *left = nodes_.StartWrite(left_i);
       Node *right = nodes_.StartWrite(right_i);
 
@@ -219,7 +218,11 @@ void KdTreeMidpointBuilder<TPointInfo, TNode, TParam>::KdTreeMidpointBuilder::Bu
     node->set_leaf();
   
     for (index_t i = node->begin(); i < node->end(); i++) {
-      node->stat().Accumulate(*param_, points_[i], point_infos_[i]);
+      const Vector *point = points_.StartRead(i);
+      const PointInfo *point_info = point_infos_.StopRead(i);
+      node->stat().Accumulate(*param_, point, point_info);
+      points_.StopRead(i);
+      point_infos_.StopRead(i);
     }
     node->stat().Postprocess(*param_, node->bound(), node->count());
   }
@@ -228,16 +231,19 @@ void KdTreeMidpointBuilder<TPointInfo, TNode, TParam>::KdTreeMidpointBuilder::Bu
 }
 
 template<typename TPointInfo, typename TNode, typename TParam>
-void KdTreeMidpointBuilder<TPointInfo, TNode, TParam>::Build() {
-  index_t node_i = 0;
+void KdTreeMidpointBuilder<TPointInfo, TNode, TParam>::Build_() {
+  index_t node_i = nodes_.Alloc();
+  Node *node = nodes_.StartWrite(node_i);
+  
+  DEBUG_SAME_INT(node_i, 0);
 
   nodes_.AddBack();
-  node->Init(0, points_.size());
+  node->Init(points_.begin_index(), points_.end_index());
   node->bound().Init(dim_);
   
-  FindBoundingBox_(node_i, points_.size(), &node->bound());
+  FindBoundingBox_(node->begin(), node->end(), &node->bound());
 
-  Build_(0);
+  Build_(node_i);
 }
 
 #endif

fastlib/u/garryb/nbr/nbr_utils.cc

@@ -1,13 +1,14 @@
 #include "nbr_utils.h"
 
-success_t nbr_utils::Load(const char *fname, TempCacheArray<Vector> *cache_out) {
+success_t nbr_utils::Load(const char *fname,
+    TempCacheArray<Vector> *cache_out, index_t vectors_per_block) {
   Matrix matrix;
   Vector first_row;
   success_t success = data::Load(fname, &matrix);
   
   matrix.MakeColumnVector(0, &first_row);
   cache_out->Init(first_row, matrix.n_cols(),
-      max(matrix.n_cols() + 1, 256));
+      max(1, min(matrix.n_cols(), vectors_per_block)));
   
   return success;
 }

fastlib/u/garryb/nbr/nbr_utils.h

@@ -1,19 +1,64 @@
 #ifndef NBR_UTILS_H
 #define NBR_UTILS_H
 
+#include "kdtree.h"
+
 namespace nbr_utils {
-  success_t Load(const char *fname, TempCacheArray<Vector> *cache_out);
+  success_t Load(const char *fname, TempCacheArray<Vector> *cache_out,
+      index_t vectors_per_block);
 
   template<typename PointInfo, typename Node, typename Param>
   success_t LoadKdTree(struct datanode *datanode,
       const Param& param,
-      CacheArray<PointInfo> *point_infos_to_rearrange,
+      TempCacheArray<PointInfo> *point_infos_out,
       TempCacheArray<Vector> *points_out,
       TempCacheArray<Node> *nodes_out) {
-    nbr_utils::Load(fx_param_str_req(datanode, ""), points_out);
-    KdTreeMidpointBuilder<TPointInfo, TNode, TParam> builder;
+    index_t vectors_per_block = fx_param_int(
+        datanode, "vectors_per_block", 256);
+    
+    fx_timer_start(datanode, "read");
+    nbr_utils::Load(fx_param_str_req(datanode, ""), points_out,
+        vectors_per_block);
+    fx_timer_stop(datanode, "read");
+    
+    PointInfo blank_info;
+    point_infos_out->Init(blank_info,
+        points_out->end_index(), vectors_per_block);
+    
+    fx_timer_start(datanode, "tree");
+    KdTreeMidpointBuilder<PointInfo, Node, Param> builder;
     builder.InitBuild(datanode, &param, points_out,
-      point_infos_to_rearrange, nodes_out);
+      point_infos_out, nodes_out);
+    fx_timer_stop(datanode, "tree");
+  }
+  
+  template<typename GNP, typename Solver>
+  void SerialDualTreeMain(datanode *datanode, const char *gnp_name) {
+    typename GNP::Param param;
+    param.Init(fx_submodule(datanode, gnp_name, gnp_name));
+
+    TempCacheArray<typename GNP::Point> q_points;
+    TempCacheArray<typename GNP::QPointInfo> q_point_infos;
+    TempCacheArray<typename GNP::QNode> q_nodes;
+    TempCacheArray<typename GNP::Point> r_points;
+    TempCacheArray<typename GNP::RPointInfo> r_point_infos;
+    TempCacheArray<typename GNP::RNode> r_nodes;
+    TempCacheArray<typename GNP::QResult> q_results;
+
+    nbr_utils::LoadKdTree(fx_param_str_req(datanode, "q"),
+        param, &q_point_infos, &q_points, &q_nodes);
+    nbr_utils::LoadKdTree(fx_param_str_req(datanode, "r"),
+        param, &r_point_infos, &r_points, &r_nodes);
+
+    typename GNP::QResult default_result;
+    default_result.Init(param);
+    q_results.Init();
+
+    Solver solver;
+    solver.Init(fx_submodule(datanode, "solver", "solver"), param,
+        &q_points, &q_point_infos, &q_nodes,
+        &r_points, &r_point_infos, &r_nodes, &q_results);
+    solver.Begin();
   }
 };
 

fastlib/u/garryb/nbr/paper/paper.tex

@@ -173,6 +173,8 @@ Although we briefly mention strategies for massive parallelization, this paper f
       - their algorithms assume only *intrinsic* prunes
     - we also can't predict load balancing
     - since I'm actually *generalizing* their work, I have to emphasize and back up this claim
+- "Oh yeah"
+  - Define "map" operator and the rest of our notation
 \end{verbatim}
 
 \section{Dual-Tree Algorithm Framework}

fastlib/u/garryb/nbr/tkde.cc

@@ -72,10 +72,17 @@ class Tkde {
 
     OT_DEF(Param) {
       OT_MY_OBJECT(kernel);
+      OT_MY_OBJECT(thresh);
       OT_MY_OBJECT(dim);
     }
 
    public:
+    void Copy(const Param& other) {
+      kernel.Copy(other.kernel);
+      thresh = other.thresh;
+      dim = other.dim;
+    }
+    
     /**
      * Initialize parameters from a data node (Req NBR).
      */
@@ -356,9 +363,7 @@ class Tkde {
     }
 
    public:
-    void Init(const Param& param,
-        const Vector& q_point, const QPointInfo& q_info,
-        const RNode& r_root) {
+    void Init(const Param& param) {
       density = 0;
       label = LAB_EITHER;
       DEBUG_ONLY(n_r = 0);