better program structure for mess of nested loops

[LMHmedchem]

Hello,

I am developing an algorithm to implement a computational geometry method I am working on. As is typical for me, I haven't fully solved the geometry yet, so I am trying to write a program in a situation where I don't fully know what the program needs to do. I need the program to solve the geometry. Such is research.

The data are data points stored in a vector of row objects. Each point is evaluated in sequence. There is an order to the sequence, but that is not terribly important.

We will call the current point being evaluated "point_B".

Evaluation involves creating and evaluating subsets of point_B plus all possible combinations of three points from the rows that are not point_B. Structurally, point_B is stored in its own object and all of the other points are stored in a vector of objects that has been sorted.

These are combinations were the order is not relevant (1,2,3 is the same as 3,2,1) and no element can occur twice in the same set (1,1,2 is not legal). I need to generate these subsets and then evaluate them for several possible conditions. If any one of the conditions is met, processing stops, the point is labeled, and processing proceeds to the next point.

I have created the possible combinations of 3 points (to go along with point_B) as I always do with a triple loop,

#include <iostream>
#include <vector>

using namespace std;

int main(int argc, char **argv) {

// loop iterators
unsigned int i, j, k;

// point being evaluated
unsigned int point_B = 1;

// vector of row data
vector<unsigned int> points;

// load sample data
points.resize(5);
points[0] = 2;
points[1] = 3;
points[2] = 4;
points[3] = 5;
points[4] = 6;

// triple loop to generate subsets of 3
for(i=0; i<points.size(); i++) {
   for(j=i+1; j<points.size(); j++) {
      for(k=j+1; k<points.size(); k++) {

         // print current subset
         cout << point_B << "," << points[i] << "," << points[j] << "," << points[k] << endl;

      }
   }
}

// main EOF brace
}

The above code generates the expected output,

1,2,3,4
1,2,3,5
1,2,3,6
1,2,4,5
1,2,4,6
1,2,5,6
1,3,4,5
1,3,4,6
1,3,5,6
1,4,5,6

This is pretty straight forward but I thought I should explain the basic structure. Note that in the above code, the points are represented as unsigned int where in the actual code they are objects.

The code gets messy where the sets start to be evaluated. There are at least a dozen different cases that need to be evaluated.

The first is a case where points have been marked to be skipped (not included in the subsets of 3) so this is checked as each level of the subset is created. This just reads a bool in the data point object to see if the point should be processed.

This version of code has the data in primitive objects. Let me know if you want more detail. It's hard to know what folks actually care about so I try to give enough to see what is happening without putting anyone to sleep.

class row_data {
public:
row_data() : check_row(process_point), id(0) { }
   bool process_point;
   unsigned int id;
};

// vector of row data
vector<row_data> points;

// triple loop to generate subsets of 3
for(i=0; i<points.size(); i++) {

   if(points[i].process_point) {

      for(j=i+1; j<points.size(); j++) {

         if(points[j].process_point) {

            for(k=j+1; k<points.size(); k++) {

               if(points[k].process_point) {

                  // print current subset
                  cout << point_B.id << "," << points[i].id <<  "," << points[j].id << "," << points[k].id << endl;

               }   
            }
         }
      }
   }
}

This, again is not too bad, but the code that evaluates the subsets runs 10+ deep in nested brackets, which frankly scares the hell out of me. This is a very linear algorithm, so I am not that worried about creating a rats nest instead of an algorithm, but it seems there must be a more organized and streamlined way of writing this.

In some other applications where I have a complicated series of conditionals I have been able to collect information into a string with a reasonable syntax and then look up the result in a map. I don't think that would work here because this is an algorithm with a long runtime and I don't want to have to evaluate everything for each subset before I know the answer. There are many cases where the first evaluation or so gives me the answer and I would like to break there and skip the rest for that point.

At any rate, can someone give me advice on a better organizational structure for an application that has to consider a great many combinations of factors?

Thanks,

LMHmedchem

[wolle]

This, again is not too bad, but the code that evaluates the subsets runs 10+ deep in nested brackets, which frankly scares the hell out of me.

Rather than a fixed number of nested loops you could use one loop that is called recursively to the wanted depth.

[superbonzo]

you could simply split the combination traversal/filtering/acting logics into different functions/classes.

For example, a purely functional approach could be ( I made it non generic just to ease exposition )

void for_each_combo(Container &points, function<bool(Point const &)> filter, function<void(Point&,Point&,Point&)> action)
{
  for (i = 0; i < points.size(); i++) if (filter(points[i]))
  for (j = i + 1; j < points.size(); j++) if (filter(points[j]))
  for (k = j + 1; k < points.size(); k++) if (filter(points[k]))
    action(points[i], points[j], points[k]);
}

// ...

for_each_combo( points,
  []( Point const& p ){
    // ...
  },
  []( Point& p1, Point& p2, Point& p3 ){
    // ...
  })

in this way, the nested loops are not obfuscated by unrelated statements, filter conditions appear one after the other, the action appears last. If the filter needs to know the nesting level, you can pass it as well ( along with any other 'stack' info ). Of course, if you expect the combo logic to change at compile/runtime, then having a generalized combo traversal with no nested loops ( something in the spirit of std::next_permutation ) might be appropriate. Alternatively, you could also model the traversal/filter/action relationships in a standard OOP way as well.

[LMHmedchem]

I have implemented some of the suggestions and now have more reasonable looking code. Mostly I just moved some of the filtering to functions.

The main function loops over all of the data passing each point to a process function along with an object to contain the results of processing.

in main.cpp

class row_data {
public:
   // initialize class members
   row_data()
      : point_B(false),
        c_set(false),
        i_hull(false),
        inconclusive(false),
        special_case(false),
        id(0),
        c_distance(0.0),
        v_magnitude(0.0),
        row_name("") { }

   bool point_B, c_set, i_hull, inconclusive;
   unsigned int id;
   double c_distance, v_magnitude;
   string row_name;
   vector<string> input_string_data;
   vector<int> int_data;
   vector<double> double_data;
   vector<double> point_coordinates;
};

// returned by evaluation functions with results of evaluation
class result_object {
public:
   // initialize class members
   result_object()
      : c_set(false),
        i_hull(false),
        inconclusive(false),
        special_case(false) { }

   bool c_set, i_hull, inconclusive, special_case;
   vector<unsigned int>special_case_set;
};


// loop iterator
unsigned int B;

// container for input data
vector<row_data> data_rows;

// loop over each row of data
for(B=0; B<data_rows.size(); B++) {

   // create return value
   result_object current_result;

   // set bool to skip the current row in processing
   data_rows[B].point_B = true;

   // call function to process each point
   process_point(current_result, data_rows[B], data_rows);

   // evaluate current_result object

}

The process point function creates all possible combinations of the rows in the data set excluding the current point_B. This involves loading all of the data from the points into an object and then passing that object to an evaluation function. The results of the evaluation function are stored in the result object.

in process_point.cpp

void process_point(result_object &current_result, row_data& point_B, vector<row_data>& rows_data){
   class plane_set {
   public:
      // initialize class members
      plane_set()
         : ct_GT_90(0),
           p_intersection_sum(0.0) { }
    
      unsigned int ct_GT_90;
      double _intersection_sum;
      vector<double> angle_BC_values;
      vector_data sum_vector, vector_BA, vector_BD, vector_BE;
      row_data point_B;
   };
    
   // code here creates a vector from point B to each other point
   // these vectors are stored in a vector of vector_data obejcts
   vector<vector_data> vecs_from_pt_B;
    
    // a sum vector is created adding all of the vectors in vecs_from_pt_B[]
   vector<double> sum_vector;
    
   // create plane set here and pass to evaluate_plane_set()
   plane_set current_plane_set;
   // add current point_B
   current_plane_set.point_B = point_B;
   // add current sum vector
   current_plane_set.sum_vector = sum_vector;

   // loop iterators
   unsigned int i, j, k;

   // loop on all combinations of vector BA, BD, and BE
   for(i=0; i<vecs_from_pt_B.size(); i++) {

      // add current vector_BA
      current_plane_set.vector_BA = vecs_from_pt_B[i];

      for(j=i+1; j<vecs_from_pt_B.size(); j++) {

         // add current vector_BD
         current_plane_set.vector_BD = vecs_from_pt_B[j];

         // loop on all vectors other than current i/j plane set
         for(k=j+1; k<vecs_from_pt_B.size(); k++) {

            // add current vector_BE
            current_plane_set.vector_BE = vecs_from_pt_B[k];

            // reinitialize values calculated in evaluate function
            current_plane_set.ct_GT_90 = 0;
            current_plane_set.p_intersection_sum = 0.0;
            // clear vector of BC angles
            current_plane_set.angle_BC_values.clear();

            // call function to evaluate plane set
            evaluate_plane_set(current_result, current_plane_set);

            // if check determines current_result.i_hull == true, return, point is finished
            if(current_result.i_hull) { return; }
            // if check determines current_result.c_set == true, return, point is finished
            else if(current_result.c_set) { return; }

            // if point has not been found to be one of the above, keep looping on subsets

         }
      }
   }

   // if all combinations have been checked without returning i_hull==true, point is c_set
   current_result.c_set = true;

return;
}

This is more reasonable than before and performance is also reasonable. I am doing more copping than I would like with statements like,

current_plane_set.vector_BA = vecs_from_pt_B[i];

where the vector object that resides in vecs_from_pt_B[] is being copped to the plane set object. It would be nicer I think to do this with a pointer or by reference. Is there a way I could make current_plane_set.vector_BA &vecs_from_pt_B[i] or *vecs_from_pt_B[i] instead of it being a copy?

I understand the thrust of this suggestion, but I don't recognize the syntax.

void for_each_combo(Container &points, function<bool(Point const &)> filter, function<void(Point&,Point&,Point&)> action) {
  for (i = 0; i < points.size(); i++) if (filter(points[i]))
  for (j = i + 1; j < points.size(); j++) if (filter(points[j]))
  for (k = j + 1; k < points.size(); k++) if (filter(points[k]))
    action(points[i], points[j], points[k]);
}

// ...

for_each_combo( points,
  []( Point const& p ){
    // ...
  },
  []( Point& p1, Point& p2, Point& p3 ){
    // ...
  })

I don't see where the filter and action functions get called. I also don't think I have ever used a function that included another function in it's definition. Is there a reference where I can do a bit of reading on what you are suggesting here?

LMHmedchem