Write a C++ parser for sphere/edfb input data

src/sphere/List.h

+#ifndef LIST_OUT_OF_BOUNDS_EXCEPTION
+#define LIST_OUT_OF_BOUNDS_EXCEPTION 1
+#endif
+
+/*! \class List<T>
+ * \brief A class to represent dynamic lists.
+ *
+ * This class helps in the creation and management of dynamic lists of
+ * objects, whose size is not known in advance. List offers the advantage
+ * of saving memory, since only the necessary space will be allocated,
+ * but it has the disadvantage of creating an object which is not contiguous
+ * in memory and, therefore, is very inefficient for subsequent manipulation.
+ *
+ * For this reason, the best use of List objects is to collect all the
+ * desired members and then, once the element number is known, to convert
+ * the List to C array, by calling List.to_array(). This function returns
+ * a contiguous array of type T[SIZE] that can be used for indexed access.
+ */
+template<class T> class List {
+ protected:
+  int size; //!< Size of the List.
+  struct element {T value; element* p_prev;}; //!< List element connector.
+  element *current, //!< Pointer to element affected by last operation.
+    *first, //!< Pointer to the first element in the List.
+    *last; //!< Pointer to the last element in the List.
+
+ public:
+  /*! \fn List(int)
+   * \brief List constructor.
+   *
+   * Use the constructor List<T>(SIZE) to create a new list with a given
+   * size. If the required size is not known in advance, it is recommended
+   * to create a List with SIZE=1 (this is the default behavior) and then
+   * to append the elements dynamically, using List.append(ELEMENT) (where
+   * ELEMENT needs to be a value of type T, corresponding to the class
+   * template specialization). Note that, due to the default behavior, the
+   * following calls are equivalent and they both produce an integer List
+   * with size equal to 1:
+   *
+   * a = List<int>(1);
+   *
+   * b = List<int>();
+   */
+  List(int size = 1) {
+    this->size = size;
+    element* p_prev = NULL;
+    T value = (T)0;
+    current = new element;
+    current->value = value;
+    current->p_prev = p_prev;
+    p_prev = current;
+    first = current;
+    for (int i = 1; i < size; i++) {
+      T value;
+      current = new element;
+      current->value = value;
+      current->p_prev = p_prev;
+      p_prev = current;
+    }
+    last = current;
+  }
+
+  ~List() {
+    current = last;
+    element *old;
+    while (current->p_prev) {
+      old = current;
+      current = old->p_prev;
+      delete old;
+    }
+  }
+  
+  /*! \fn append(T)
+   * \brief Append an element at the end of the list.
+   *
+   * To dynamically create a list whose size is not known in advance,
+   * elements can be appended in an iterative way. Note that element
+   * manipulation is much more effective in a C array than in a List
+   * object. For this reason, after the List has been created, it is
+   * strongly advised to convert it to a C array by calling the function
+   * List.to_array().
+   */
+  void append(T value) {
+    element *p_prev = last;
+    current = new element;
+    current->value = value;
+    current->p_prev = p_prev;
+    last = current;
+    size++;
+  }
+
+  /*! \fn get(int)
+   * \brief Get the element at given index.
+   *
+   * Get the element specified by the index argument. The first element
+   * has index 0 and the last one has index [size - 1].
+   *
+   * \return T: the value of the element at the requested position.
+   * \throw LIST_OUT_OF_BOUNDS_EXCEPTION: raised if the index is out of bounds.
+   */
+  T get(int index) {
+    if (index < 0 || index > size - 1) {
+      throw LIST_OUT_OF_BOUNDS_EXCEPTION;
+    }
+    current = last;
+    for (int i = size - 1; i > index; i--) current = current->p_prev;
+    return current->value;
+  }
+
+  /*! \fn set(int, T)
+   * \brief Set an element by index and value.
+   *
+   * Set the element at the position specified by the index to the value
+   * specified by the value argument.
+   *
+   * \throw LIST_OUT_OF_BOUNDS_EXCEPTION: raised if the index is out of bounds.
+   */
+  void set(int index, T value) {
+    if (index < 0 || index > size - 1) {
+      throw LIST_OUT_OF_BOUNDS_EXCEPTION;
+    }
+    current = last;
+    for (int i = size - 1; i > index; i--) current = current->p_prev;
+    current->value = value;
+  }
+
+  /*! \fn to_array()
+   * \brief Convert the list to a C array.
+   *
+   * The List object is useful to dynamically manage a set of objects
+   * when the number of elements is not known in advance. However, the
+   * resulting object is not contiguosly stored in memory. As a result,
+   * access to specific elements in the middle of the list is not very
+   * effective, because the list needs to be walked every time up to
+   * the desired position. In order to avoid this, List.to_array() makes
+   * a conversion from List to C array, returning a contiguous object,
+   * where indexed access can be used.
+   *
+   * \return T[size]: a C array of type T and size equal to the List size.
+   */
+  T* to_array() {
+    T *result = new T[size];
+    current = last;
+    for (int i = size - 1; i > 0; i--) {
+      result[i] = current->value;
+      current = current->p_prev;
+    }
+    return result;
+  }
+};

src/sphere/edfb.cpp

+/*! \brief C++ implementation of EDFB
+ *
+ *  This code aims at replicating the original work-flow in C++.
+ */
+
+#include <cstdio>
+#include <cmath>
+#include <string>
+#include "List.h"
+#include "edfb.h"
+
+int main(int argc, char **argv) {
+  DEDFB input_data = DEDFB(std::string("../../test_data/sphere/DEDFB"));
+  input_data.print();
+  return 0;
+}

src/sphere/edfb.h

+/*! \class DEDFB
+ * \brief A class to represent EDFB configuration data.
+ *
+ * This class replicates the structure of formatted DEDFB input data
+ * files for the EDFB program. The class is built on the SPHERE/DEDFB
+ * file template, but it is meant to grant compliance with the CLUSTER
+ * case too.
+ *
+ * Compatibility between the SPHERE and the CLUSTER cases is planned
+ * through the introduction of dynamic pointers that can fit both cases,
+ * provided that the DEDFB formatting has been properly understood. In
+ * case not, the suggested work-around is to develop specialized classes
+ * for the SPHERE case and the CLUSTER case. In that case, the current
+ * class would already represent a valid implementation for the single
+ * sphere problem.
+ *
+ * It should be noted that the variable names have been left intentionally
+ * equal to the FORTRAN equivalents. In addition they are not documented
+ * yet. Documentation and naming may be revised after the code has been
+ * tested to succesfully reproduce the original work-flow.
+ */
+class DEDFB {
+ protected:
+  int nsph, ies;
+  double exdc, wp, xip;
+  int idfc, nxi, instpc, insn;
+  double *xi_vector;
+  int *iog;
+  int *nshl;
+  double *ros;
+  double *rcf;
+  double *dc0;
+
+public:
+  DEDFB(std::string file_name);
+  void print();
+};
+
+/*! \fn DEDFB(std::string)
+ * \brief Data structure constructor.
+ *
+ * The default procedure to import a configuration file is to build a
+ * configuration object, according to one of the following flavours:
+ *
+ * DEDFB edfb_cfg = DEDFB("FILE_NAME"); // configuration object
+ *
+ * or
+ *
+ * DEDFB *edfb_cfg_ptr = new DEDFB("FILE_NAME"); // pointer version
+ */
+DEDFB::DEDFB(std::string file_name) {
+  FILE *input_file = fopen(file_name.c_str(), "r");
+  int exp_exdc, exp_wp, exp_xip;
+  List<double> xi_list(1);
+  fscanf(input_file, " %d %d", &nsph, &ies);
+  fscanf(
+	 input_file,
+	 " %9lf D%d %9lf D%d %8lf D%d %d %d %d %d",
+	 &exdc, &exp_exdc,
+	 &wp, &exp_wp,
+	 &xip, &exp_xip,
+	 &idfc, &nxi, &instpc, &insn
+	 );
+  exdc *= pow(10.0, exp_exdc);
+  wp *= pow(10.0, exp_wp);
+  xip *= pow(10.0, exp_xip);
+  for (int i_xi = 0; i_xi < nxi; i_xi++) {
+    double xi;
+    int exp_xi;
+    fscanf(input_file, " %9lf D%d", &xi, &exp_xi);
+    xi *= pow(10.0, exp_xi);
+    if (i_xi == 0) {
+      xi_list.set(0, xi);
+    } else {
+      xi_list.append(xi);
+    }
+  }
+  xi_vector = xi_list.to_array();
+  iog = new int[nsph];
+  for (int i_iog = 0; i_iog < nsph; i_iog++) {
+    fscanf(input_file, " %d", (iog + i_iog));
+  }
+  nshl = new int[nsph];
+  ros = new double[nsph];
+  for (int i_nshl = 0; i_nshl < nsph; i_nshl++) {
+    double i_ros;
+    int i_ros_exp;
+    fscanf(input_file, " %d %9lf D%d", (nshl + i_nshl), &i_ros, &i_ros_exp);
+    i_ros *= pow(10.0, i_ros_exp);
+    ros[i_nshl] = i_ros;
+  }
+  rcf = new double[nsph * nshl[0]];
+  for (int i_rcf = 0; i_rcf < nshl[0]; i_rcf++) {
+    for (int j_rcf = 0; j_rcf < nsph; j_rcf++) {
+      double ij_rcf;
+      int ij_rcf_exp;
+      fscanf(input_file, " %8lf D%d", &ij_rcf, &ij_rcf_exp);
+      ij_rcf *= pow(10.0, ij_rcf_exp);
+      rcf[(nsph * i_rcf) + j_rcf] = ij_rcf;
+    }
+  }
+  dc0 = new double[nsph * nshl[0]];
+  for (int i_dc0 = 0; i_dc0 < nsph; i_dc0++) {
+    double ij_dc0;
+    int ij_dc0_exp;
+    fscanf(input_file, " (");
+    for (int j_dc0 = 0; j_dc0 < nshl[0] - 1; j_dc0++) {
+      fscanf(input_file, " %8lf D%d,", &ij_dc0, &ij_dc0_exp);
+      ij_dc0 *= pow(10.0, ij_dc0_exp);
+      dc0[(nshl[0] * i_dc0) + j_dc0] = ij_dc0;
+    }
+    fscanf(input_file, " %8lf D%d)", &ij_dc0, &ij_dc0_exp);
+    ij_dc0 *= pow(10.0, ij_dc0_exp);
+    dc0[(nshl[0] * i_dc0) + nshl[0] - 1] = ij_dc0;
+  }
+  fclose(input_file);
+}
+
+/*! \fn print()
+ * \brief Print the data structure to console.
+ *
+ * This is a simple function to print the data read from a configuration
+ * file to the screen. This operation can be useful for debug and testing.
+ */
+void DEDFB::print() {
+  printf("### CONFIGURATION DATA ###\n");
+  printf("NSPH   = %d\n", nsph);
+  printf("IES    = %d\n", ies);
+  printf("EXDC   = %E\n", exdc);
+  printf("WP     = %E\n", wp);
+  printf("XIP    = %E\n", xip);
+  printf("IDFC   = %d\n", idfc);
+  printf("NXI    = %d\n", nxi);
+  printf("INSTPC = %d\n", instpc);
+  printf("INSN   = %d\n", insn);
+  printf("XIV    = [ %lE", xi_vector[0]);
+  for (int i_xiv = 1; i_xiv < nxi; i_xiv++) printf(", %lE", xi_vector[i_xiv]);
+  printf(" ]\n");
+  printf("IOG    = [ %d", iog[0]);
+  for (int i_iog = 1; i_iog < nsph; i_iog++) printf(", %d", iog[i_iog]);
+  printf(" ]\n");
+  printf("NSHL   = [ %d", nshl[0]);
+  for (int i_nshl = 1; i_nshl < nsph; i_nshl++) printf(", %d", nshl[i_nshl]);
+  printf(" ]\n");
+  printf("ROS    = [ %lE", ros[0]);
+  for (int i_ros = 1; i_ros < nsph; i_ros++) printf(", %lE", ros[i_ros]);
+  printf(" ]\n");
+  printf("RCF    = [ %lE", rcf[0]);
+  for (int i_rcf = 1; i_rcf < nsph * nshl[0]; i_rcf++) printf(", %lE", rcf[i_rcf]);
+  printf(" ]\n");
+  printf("DC0    = [ %lE", dc0[0]);
+  for (int i_dc0 = 1; i_dc0 < nsph * nshl[0]; i_dc0++) printf(", %lE", dc0[i_dc0]);
+  printf(" ]\n");
+  printf("######   END DATA   ######\n");
+  
+}
+
+void inxi(
+  DEDFB *data
+); //!< \brief Initialization process