From 3a53ae3660f2aef8c6ed0a8131d139a20267b31a Mon Sep 17 00:00:00 2001
From: Ray Speth <speth@mit.edu>
Date: Tue, 5 Mar 2013 17:02:56 +0000
Subject: [PATCH] [1D] Fixed infinite loops in grid refinement / pruning

Improved determination of the conditions under which points can be removed
without violating the 'ratio' constraint. Previously, required points were being
removed from the grid.
---
 src/oneD/refine.cpp | 22 ++++++++++++++++++++--
 1 file changed, 20 insertions(+), 2 deletions(-)

diff --git a/src/oneD/refine.cpp b/src/oneD/refine.cpp
index b0d6fcbd8..48044c1ac 100644
--- a/src/oneD/refine.cpp
+++ b/src/oneD/refine.cpp
@@ -134,19 +134,37 @@ int Refiner::analyze(size_t n, const doublereal* z,
         }
     }
 
+    // Refine based on properties of the grid itself
     for (size_t j = 1; j < n-1; j++) {
+        // Add a new point if the ratio with left interval is too large
         if (dz[j] > m_ratio*dz[j-1]) {
             m_loc[j] = 1;
             m_c["point "+int2str(j)] = 1;
+            m_keep[j-1] = 1;
+            m_keep[j] = 1;
+            m_keep[j+1] = 1;
+            m_keep[j+2] = 1;
         }
+
+        // Add a point if the ratio with right interval is too large
         if (dz[j] < dz[j-1]/m_ratio) {
             m_loc[j-1] = 1;
             m_c["point "+int2str(j-1)] = 1;
+            m_keep[j-2] = 1;
+            m_keep[j-1] = 1;
+            m_keep[j] = 1;
+            m_keep[j+1] = 1;
         }
-        if (j > 1 && z[j+1]-z[j] > m_ratio * dz[j-2]) {
+
+        // Keep the point if removing would make the ratio with the left
+        // interval too large.
+        if (j > 1 && z[j+1]-z[j-1] > m_ratio * dz[j-2]) {
             m_keep[j] = 1;
         }
-        if (j < n-2 && z[j+1]-z[j] > m_ratio * dz[j+1]) {
+
+        // Keep the point if removing would make the ratio with the right
+        // interval too large.
+        if (j < n-2 && z[j+1]-z[j-1] > m_ratio * dz[j+1]) {
             m_keep[j] = 1;
         }