Changeset 2730 for trunk/src/HuntLocator.cpp

Timestamp:

01/16/13 16:00:28 (11 years ago)

Author:

Takeshi Nakazato

Message:

New Development: No

JIRA Issue: Yes CAS-4770

Ready for Test: Yes

Interface Changes: No

What Interface Changed: Please list interface changes

Test Programs: List test programs

Put in Release Notes: Yes/No?

Module(s): Module Names change impacts.

Description: Describe your changes here...

Rewrite implementations for locator and interpolator.
Documentation (doxygen format) is added to header files.

File:

• trunk/src/HuntLocator.cpp (modified) (3 diffs)

trunk/src/HuntLocator.cpp

@@ -16 +16 @@
 namespace asap {
 
-HuntLocator::HuntLocator(double *v, unsigned int n)
-  : Locator(v, n),
+HuntLocator::HuntLocator()
+  : Locator(),
+    prev_(0)
+{}
+
+HuntLocator::HuntLocator(double *v, unsigned int n, bool copystorage)
+  : Locator(v, n, copystorage),
     prev_(0)
 {}
@@ -28 +33 @@
   if (n_ == 1)
     return 0;
-  bool ascending = (x_[n_-1] >= x_[0]);
-  if (ascending) {
+
+  if (ascending_) {
     if (x <= x_[0])
       return 0;
@@ -44 +49 @@
   unsigned int jl = 0;
   unsigned int ju = n_;
+
+  // hunt phase
   if (prev_ > 0 && prev_ < n_) {
-    // do hunt
     jl = prev_;
-    unsigned int inc = 1;
-    if ((x >= x_[jl]) == ascending) {
-      // hunt up 
-      if (jl >= n_ - 1)
-        return jl;
-      ju = jl + inc;
-      while ((x >= x_[ju]) == ascending) {
-        jl = ju;
-        inc <<= 1;
-        ju = jl + inc;
-        if (ju > n_ - 1) {
-          ju = n_;
-          break;
-        }
-      }
-    }
-    else {
-      // hunt down
-      if (jl == 0) 
-        return jl;
-      ju = jl;
-      jl -= inc;
-      while ((x < x_[jl]) == ascending) {
-        ju = jl;
-        inc <<= 1;
-        if (inc >= ju) {
-          jl = 0;
-          break;
-        }
-        else
-          jl = ju - inc;
-      }
-    }
+    hunt(x, jl, ju);
   }
 
   // final bisection phase
-  unsigned int jm;
-  if (ascending) {
-    while (ju - jl > 1) {
-      jm = (ju + jl) >> 1;
-      if (x > x_[jm])
-        jl = jm;
-      else
-        ju = jm;
-    }
-  }
-  else {
-    while (ju - jl > 1) {
-      jm = (ju + jl) >> 1;
-      if (x < x_[jm])
-        jl = jm;
-      else
-        ju = jm;
-    }
-  }
-  prev_ = jl;
-  return ju;
+  unsigned int j = bisection(x, jl, ju);
+  prev_ = (j > 0) ? j - 1 : 0;
+  return j;
 }