Rename phasorinc, add {C,D}Collection documentation, and main() → main.cpp

.github/workflows/doxygen-gh-pages.yml

@@ -8,7 +8,7 @@ jobs:
   deploy:
     runs-on: ubuntu-latest
     steps:
-      - uses: DenverCoder1/doxygen-github-pages-action@v1.1.0
+      - uses: DenverCoder1/doxygen-github-pages-action@v1.2.0
         with:
           github_token: ${{ secrets.GITHUB_TOKEN }}
           folder: html

README.md

@@ -30,7 +30,7 @@ mkdir build; cd build
 cmake .. \
 # -DMatlab_ROOT_DIR=/usr/local/MATLAB/R2019b/ \
 # -DFFTW_ROOT=/usr/local/fftw3/ \
-# -DCMAKE_INSTALL_PREFIX=$HOME/.local/
+# -DCMAKE_INSTALL_PREFIX=$HOME/.local/ \
 # -DBUILD_TESTING=ON
 make install
 ```
@@ -60,18 +60,48 @@ where lines need to be commented in and the paths modified if cmake cannot
 
 
 ***
+
 ## Usage
 
+Once the executable has been compiled and installed, `tdms` should be in the `PATH`.
+Check that installation worked with
+
+```bash
+$ tdms -h
+```
+
+You can invoke it directly or call it from a MATLAB script.
+We recommend that beginners with MATLAB installed start with the demonstration MATLAB script.
+
 #### To run the demonstration code
 
-Once the executable has been compiled, move into directory _examples/arc_01_,
-launch Matlab and run the Matlab script:
+Move into directory [`examples/arc_01`](./examples/arc_01/),
+launch MATLAB and run the MATLAB script:
 
-_run_pstd_bscan.m_
+[`run_pstd_bscan.m`](./examples/arc_01/run_pstd_bscan.m)
 
 This script will generate the input to the executable, run the executable and
 display sample output.
 
+#### To run standalone
+
+You can also run `tdms` from the command line...
+
+```bash
+$ tdms --help
+Usage:
+tdms [options] infile outfile
+tdms [options] infile gridfile outfile
+Options:
+-h:	Display this help message
+--finite-difference:	Use the finite-difference solver, instead of the pseudo-spectral method.
+-q:	Quiet operation. Silence all logging
+-m:	Minimise output file size by not saving vertex and facet information
+```
+
+The basic workflow is with two arguments, an input file as specified by [`iterate_fdtd_matrix.m`](./tdms/matlab/iteratefdtd_matrix.m), and an output file name to be created.
+
+You can choose two possible solver methods: either pseudo-spectral time-domain (PSTD, the default) or finite-difference (FDTD, with option `--finite-difference`).
 
 #### Parallelism
 
@@ -80,9 +110,10 @@ number of threads can be set with the `OMP_NUM_THREADS` environment variable.
 For example, to use 4 threads, in a bash shell, use:
 
 ```bash
-export OMP_NUM_THREADS=4
+$ export OMP_NUM_THREADS=4
 ```
 
+Before calling the `tdms` executable.
 
 ## Contributing
 

doc/developers.md

@@ -99,7 +99,9 @@ spdlog::debug("Send help");
 ```
 
 ### Compiling on UCL's Myriad cluster
+
 <details>
+
   > **Warning**
   > These instructions are a bit experimental. Please use with care (and report anything that's wrong here)!
 
@@ -134,13 +136,13 @@ spdlog::debug("Send help");
   CApath: none
   ```
   it's because the MATLAB module is interfering with the SSL certificates (and we clone over https by default). This issue is known and reported. As a workaround, we've added the build option `-DGIT_SSH=ON` to switch to `git clone` over ssh instead.
+
 </details>
 
 
 ## Where's the main?
 
-The C++ `main` function is in openandorder.cpp <!-- words with a dot in them are assumed to be files so this will hyperlink to openandorder.cpp iff *that* file is also documented. --> however this only really does file I/O and setup.
-The main FDTD algorithm code is in iterator.cpp <!-- won't be linked as an undocumented file doesn't exist for doxygen... this is fine, we can link to the real file in github.--> and classes included therein.
+The C++ `main` function is in main.cpp however the main algorithm code is execute_simulation() in iterator.cpp
 
 ## Testing
 

tdms/CMakeLists.txt

@@ -83,32 +83,10 @@ endif()
 
 
 # TDMS target -----------------------------------------------------------------
-set(SOURCES
-    src/fields/base.cpp
-    src/fields/electric.cpp
-    src/fields/magnetic.cpp
-    src/fields/split.cpp
-    src/fields/td_field_exporter_2d.cpp
-    src/argument_parser.cpp
-    src/array_init.cpp
-    src/arrays.cpp
-    src/dimensions.cpp
-    src/fdtd_grid_initialiser.cpp
-    src/grid_labels.cpp
-    src/interface.cpp
-    src/interpolate.cpp
-    src/iterator.cpp
-    src/interpolation_methods.cpp
-    src/matrix_collection.cpp
-    src/mesh_base.cpp
-    src/numerical_derivative.cpp
-    src/shapes.cpp
-    src/simulation_parameters.cpp
-    src/source.cpp
-    src/timer.cpp
-    src/utils.cpp
-    matlabio/matlabio.cpp
-    )
+file(GLOB_RECURSE SOURCES "src/*.cpp")
+
+# TODO: can delete this line when matlabio is removed.
+set(SOURCES ${SOURCES} matlabio/matlabio.cpp)
 
 if (BUILD_TESTING)
     test_target()

tdms/cmake/targets.cmake

@@ -1,7 +1,7 @@
 function(release_target)
     add_executable(tdms)
 
-    target_sources(tdms PRIVATE ${SOURCES} src/openandorder.cpp)
+    target_sources(tdms PRIVATE ${SOURCES} src/main.cpp)
 
     target_link_libraries(tdms
             FFTW::Double
@@ -36,7 +36,7 @@ function(test_target)
 
     add_library(tdms_lib SHARED)
     target_sources(tdms_lib PUBLIC ${SOURCES})
-    add_executable(tdms "src/openandorder.cpp")
+    add_executable(tdms "src/main.cpp")
 
     target_link_libraries(tdms_lib LINK_PUBLIC
             FFTW::Double

tdms/include/arrays.h

@@ -91,14 +91,26 @@ class MaterialCollection{
   static void init_xyz_vectors(const mxArray *ptr, XYZVectors &arrays, const std::string &prefix);
 };
 
+/**
+ * @brief A class to encapsulate collection of algebraic terms in the
+ * discretized forms of Maxwells equations for E fields. The symbol chosen in
+ * the original reference is \f$C\f$.
+ *
+ * @details Algebraic terms \f$C_{a,b,c}\f$ defined in Section 4.2 of Munro, P,.
+ * "Application of numerical methods to high numerical aperture imaging", 2006,
+ * PhD thesis, Imperial College London.
+ *
+ * The definitions are equations 4.13, 4.14 (pp 82-3). Part of Maxwell's E-field
+ * equations in equations 4.7-9.
+ */
 class CCollectionBase {
 public:
   XYZVectors a;
   XYZVectors b;
   XYZVectors c;
 };
 
-// TODO: docstring
+/*! @copydoc CCollectionBase */
 class CCollection : public CCollectionBase {
 private:
   void init_xyz_vectors(const mxArray *ptr, XYZVectors &arrays, const std::string &prefix);
@@ -110,19 +122,31 @@ class CCollection : public CCollectionBase {
   explicit CCollection(const mxArray *ptr);
 };
 
-// TODO: docstring
+/*! @copydoc CCollectionBase */
 class CMaterial : public CCollectionBase, MaterialCollection {
 public:
   explicit CMaterial(const mxArray *ptr);
 };
 
+/**
+ * @brief A class to encapsulate collection of algebraic terms in the
+ * discretized forms of Maxwells equations for H fields. The symbol chosen in
+ * the original reference is \f$D\f$.
+ *
+ * @details Algebraic terms \f$D_{a,b}\f$ defined in Section 4.2 of Munro, P,.
+ * "Application of numerical methods to high numerical aperture imaging", 2006,
+ * PhD thesis, Imperial College London.
+ *
+ * The definitions are equations 4.15, 4.16 (pp 82-3). Part of Maxwell's H-field
+ * equations in equations 4.10-12.
+ */
 class DCollectionBase {
 public:
   XYZVectors a;
   XYZVectors b;
 };
 
-// TODO: docstring
+/*! @copydoc DCollectionBase */
 class DCollection: public DCollectionBase{
 private:
   static void init_xyz_vectors(const mxArray *ptr, XYZVectors &arrays, const std::string &prefix);
@@ -131,7 +155,7 @@ class DCollection: public DCollectionBase{
   explicit DCollection(const mxArray *ptr);
 };
 
-// TODO: docstring
+/*! @copydoc DCollectionBase */
 class DMaterial : public DCollectionBase, MaterialCollection {
 public:
   explicit DMaterial(const mxArray *ptr);

tdms/include/mesh_base.h

@@ -94,6 +94,6 @@ void triangulateCuboidSkip(int I0, int I1, int J0, int J1, int K0, int K1, mxArr
 void conciseTriangulateCuboid(int I0, int I1, int J0, int J1, int K0, int K1, mxArray **vertices,
                               mxArray **facets);
 void conciseTriangulateCuboidSkip(int I0, int I1, int J0, int J1, int K0, int K1,
-                                  PhasorInc &phasorinc, mxArray **vertices, mxArray **facets);
+                                  SurfaceSpacingStride &spacing_stride, mxArray **vertices, mxArray **facets);
 
 void conciseCreateBoundary(int I0, int I1,int K0, int K1, mxArray **vertices, mxArray ** facets);

tdms/include/simulation_parameters.h

@@ -34,9 +34,9 @@ enum RunMode{
 };
 
 enum Dimension{
-  THREE,  // Full dimensionality - compute all H and E components
-  TE,     // Transverse electric - only compute Ex, Ey, and Hz components
-  TM      // Transverse magnetic - only compute Hx, Hy, and Ez components
+  THREE,  //< Full dimensionality - compute all H and E components
+  TE,     //< Transverse electric - only compute Ex, Ey, and Hz components
+  TM      //< Transverse magnetic - only compute Hx, Hy, and Ez components
 };
 
 enum InterpolationMethod{
@@ -50,7 +50,7 @@ enum InterpolationMethod{
  * to extract the phasors on. the surface i.e. x = 2 means extract from
  * every 2nd Yee cell.
  */
-struct PhasorInc{
+struct SurfaceSpacingStride{
   int x = 1;
   int y = 1;
   int z = 1;
@@ -74,44 +74,44 @@ class SimulationParameters{
 public:
     SimulationParameters();
 
-    double       omega_an      = 0.0;       // Angular ω
-    unsigned int Nt            = 0;         // Number of simulation steps
-    unsigned int Np            = 0;         // Number of times to extract the phasors
-    unsigned int Npe           = 0;         // Extract the phasors every this number of iterations
-    unsigned int start_tind    = 0;         // Starting iteration number for the time steps
-    double       dt            = 0.0;       // Time step
-    bool         has_tdfdir    = false;     // Is the tdfdir (time domain field directory) defined?
-    bool         is_multilayer = false;     // Is this simulation of a multilayer?
-    bool         is_disp_ml    = false;     // Is the material dispersive?
-    double       to_l          = 0.0;       // Time delay of pulse
-    double       hwhm          = 0.0;       // Half width at half max of pulse
-    PerfectlyMatchedLayer pml;              // Perfectly matched layer struct with size attributes
-    bool         exphasorsvolume = false;   // Should phasors be extracted in the whole volume?
-    bool         exphasorssurface = false;  // Should phasors be extracted on a surface?
-    bool         intphasorssurface = false; // Should phasors be extracted/interpolated?
-    RunMode      run_mode     = complete;   // Run mode
-    SourceMode   source_mode  = pulsed;     // Source mode
-    Dimension    dimension    = THREE;      // Dimensions to calculate in
-    bool         is_structure = false;      // Has a grating structure been defined?
-    bool         exdetintegral = false;     // TODO: detector sensitivity evaluation ?
-    int          k_det_obs    = 0;          // TODO: no idea what this is?!
-    double       z_obs        = 0.0;        // Value of the input grid_labels_z at k_det_obs
+    double       omega_an      = 0.0;       //< Angular ω
+    unsigned int Nt            = 0;         //< Number of simulation steps
+    unsigned int Np            = 0;         //< Number of times to extract the phasors
+    unsigned int Npe           = 0;         //< Extract the phasors every this number of iterations
+    unsigned int start_tind    = 0;         //< Starting iteration number for the time steps
+    double       dt            = 0.0;       //< Time step
+    bool         has_tdfdir    = false;     //< Is the tdfdir (time domain field directory) defined?
+    bool         is_multilayer = false;     //< Is this simulation of a multilayer?
+    bool         is_disp_ml    = false;     //< Is the material dispersive?
+    double       to_l          = 0.0;       //< Time delay of pulse
+    double       hwhm          = 0.0;       //< Half width at half max of pulse
+    PerfectlyMatchedLayer pml;              //< Perfectly matched layer struct with size attributes
+    bool         exphasorsvolume = false;   //< Should phasors be extracted in the whole volume?
+    bool         exphasorssurface = false;  //< Should phasors be extracted on a surface?
+    bool         intphasorssurface = false; //< Should phasors be extracted/interpolated?
+    RunMode      run_mode     = complete;   //< Run mode
+    SourceMode   source_mode  = pulsed;     //< Source mode
+    Dimension    dimension    = THREE;      //< Dimensions to calculate in
+    bool         is_structure = false;      //< Has a grating structure been defined?
+    bool         exdetintegral = false;     //< TODO: detector sensitivity evaluation ?
+    int          k_det_obs    = 0;          //< TODO: no idea what this is?!
+    double       z_obs        = 0.0;        //< Value of the input grid_labels_z at k_det_obs
     bool         air_interface_present = false;
-    double       air_interface = 0.0;       // TODO: what is this?!
-    bool         interp_mat_props = false;  // Should the material properties be interpolated?
-    InterpolationMethod interp_method = cubic; // Type of surface field interpolation to do
-    bool         exi_present = false;       // Is the time dependent x incident field present?
-    bool         eyi_present = false;       // Is the time dependent y incident field present?
-    PhasorInc    phasorinc;                 // Surface stride for extracting phasors
-    YeeCellDimensions delta;                // Yee cell dimensions (dx, dy, dz)
+    double       air_interface = 0.0;       //< TODO: what is this?!
+    bool         interp_mat_props = false;  //< Should the material properties be interpolated?
+    InterpolationMethod interp_method = cubic; //< Type of surface field interpolation to do
+    bool         exi_present = false;       //< Is the time dependent x incident field present?
+    bool         eyi_present = false;       //< Is the time dependent y incident field present?
+    SurfaceSpacingStride spacing_stride;    //< Surface stride for extracting phasors (in matlab this is called 'phasorinc')
+    YeeCellDimensions delta;                //< Yee cell dimensions (dx, dy, dz)
 
     void set_run_mode(std::string mode_string);
 
     void set_source_mode(std::string mode_string);
 
     void set_dimension(std::string mode_string);
 
-    void set_phasorinc(const double* vector);
+    void set_spacing_stride(const double* vector);
 
     void set_Np(FrequencyExtractVector &f_ex_vec);
 };

tdms/src/argument_parser.cpp

@@ -47,8 +47,8 @@ ArgumentNamespace ArgumentParser::parse_args(int n_args, char *arg_ptrs[]) {
 
 void ArgumentParser::print_help_message(){
   fprintf(stdout,"Usage:\n"
-                 "openandorder [options] infile outfile\n"
-                 "openandorder [options] infile gridfile outfile\n"
+                 "tdms [options] infile outfile\n"
+                 "tdms [options] infile gridfile outfile\n"
                  "Options:\n"
                  "-h:\tDisplay this help message\n"
                  "--finite-difference:\tUse the finite-difference solver, instead of the pseudo-spectral method.\n"

tdms/src/iterator.cpp

@@ -412,7 +412,7 @@ void execute_simulation(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs
   input_counter++;
   //fprintf(stderr,"Got   phasorsurface\n");
 
-  params.set_phasorinc(mxGetPr(prhs[input_counter++]));
+  params.set_spacing_stride(mxGetPr(prhs[input_counter++]));
   params.set_dimension(string_in(prhs[input_counter++], "dimension"));
 
   /*Get conductive_aux */
@@ -582,7 +582,7 @@ void execute_simulation(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs
                             &mx_surface_facets);
     else
       conciseTriangulateCuboidSkip(cuboid[0], cuboid[1], cuboid[2], cuboid[3], cuboid[4], cuboid[5],
-                                   params.phasorinc, &mx_surface_vertices,
+                                   params.spacing_stride, &mx_surface_vertices,
                                    &mx_surface_facets);
     //fprintf(stderr,"Qos 00a:\n");
     //we don't need the facets so destroy the matrix now to save memory
@@ -4636,7 +4636,7 @@ void execute_simulation(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs
                             &mx_surface_facets);
     else
       conciseTriangulateCuboidSkip(cuboid[0], cuboid[1], cuboid[2], cuboid[3], cuboid[4], cuboid[5],
-                                   params.phasorinc, &dummy_vertex_list,
+                                   params.spacing_stride, &dummy_vertex_list,
                                    &mx_surface_facets);
     mxDestroyArray(dummy_vertex_list);
     mxArray *vertex_list;