From a6c5371e178bc81d3ae2068cdc90368ad126397e Mon Sep 17 00:00:00 2001
From: Giovanni La Mura <giovanni.lamura@inaf.it>
Date: Thu, 21 Sep 2023 15:04:13 +0200
Subject: [PATCH] Write a C++ parser for sphere/edfb input data
---
src/sphere/List.h | 150 +++++++++++++++++++++++++++++++++++++++++
src/sphere/edfb.cpp | 16 +++++
src/sphere/edfb.h | 160 ++++++++++++++++++++++++++++++++++++++++++++
3 files changed, 326 insertions(+)
create mode 100644 src/sphere/List.h
create mode 100644 src/sphere/edfb.cpp
create mode 100644 src/sphere/edfb.h
diff --git a/src/sphere/List.h b/src/sphere/List.h
new file mode 100644
index 00000000..bed78f2e
--- /dev/null
+++ b/src/sphere/List.h
@@ -0,0 +1,150 @@
+#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;
+ }
+};
diff --git a/src/sphere/edfb.cpp b/src/sphere/edfb.cpp
new file mode 100644
index 00000000..cb0df775
--- /dev/null
+++ b/src/sphere/edfb.cpp
@@ -0,0 +1,16 @@
+/*! \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;
+}
diff --git a/src/sphere/edfb.h b/src/sphere/edfb.h
new file mode 100644
index 00000000..f47e4fc4
--- /dev/null
+++ b/src/sphere/edfb.h
@@ -0,0 +1,160 @@
+/*! \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
--
GitLab