Minor corrections (#23)

Refactoring and restructuring.

.gitignore

@@ -32,5 +32,6 @@
 *.app
 
 # standard build directories
-build
-_build
+/build*/
+/install*/
+/_*/

src/program_dftd.cpp → app/main.cpp

@@ -106,13 +106,13 @@ void help() {
            be verbose
    -h, --help
            print this message
-
+   -g, --grad
+           calculate dispersion gradient
 )";
 }
 
 int main(int argc, char **argv) {
   std::string func;
-  bool lverbose{false};
   bool lmbd{true}, lgrad{false};
   dftd4::dparam par; // damping parameter for DFT-D4 calculation
   dftd4::TMolecule mol;
@@ -137,12 +137,21 @@ int main(int argc, char **argv) {
   }
 
   // get other flags
-  if (args.getflag("-v") || args.getflag("--verbose")) { lverbose = true; }
+  // bool lverbose{false};
+  // if (args.getflag("-v") || args.getflag("--verbose")) { lverbose = true; }
   if (args.getflag("-g") || args.getflag("--grad")) { lgrad = true; }
+
   if (args.getflag("--func")) {
     func = args.getopt("--func");
-    dftd4::d4par(func, par, lmbd);
+  } else if (args.getflag("-f")) {
+    func = args.getopt("-f");
+  } else {
+    printf("WARNING: No functional given (via '--func'). Defaulting to PBE.\n\n"
+    );
+    func = "pbe";
   }
+  dftd4::d4par(func, par, lmbd);
+
   // last argument is assumed to filename since
   // dftd4 [options] <file>
   std::string fname{argv[argc - 1]};
@@ -154,15 +163,6 @@ int main(int argc, char **argv) {
   dftd4::TCutoff cutoff;
   dftd4::TD4Model d4;
 
-  // masking (nothing excluded)
-  dftd4::TVector<int> realIdx;
-  realIdx.NewVec(mol.NAtoms);
-  int nat = 0;
-  for (int i = 0; i != mol.NAtoms; i++) {
-    realIdx(i) = nat;
-    nat++;
-  }
-
   // analytical gradient
   double *d4grad;
   if (lgrad) {

src/dftd_readxyz.cpp → app/readxyz.cpp

@@ -15,7 +15,6 @@
  * You should have received a copy of the GNU Lesser General Public License
  * along with cpp-d4.  If not, see <https://www.gnu.org/licenses/>.
  */
-
 #include "dftd_readxyz.h"
 
 #include <algorithm>
@@ -66,19 +65,6 @@ void read_xyzfile(const std::string &name, dftd4::TMolecule &mol) {
 
 int element(const std::string &sym) {
   char elem[3]{"  "};
-  char pse[118][3]{
-    "h ", "he", "li", "be", "b ", "c ", "n ", "o ", "f ", "ne", "na", "mg",
-    "al", "si", "p ", "s ", "cl", "ar", "k ", "ca", "sc", "ti", "v ", "cr",
-    "mn", "fe", "co", "ni", "cu", "zn", "ga", "ge", "as", "se", "br", "kr",
-    "rb", "sr", "y ", "zr", "nb", "mo", "tc", "ru", "rh", "pd", "ag", "cd",
-    "in", "sn", "sb", "te", "i ", "xe", "cs", "ba", "la", "ce", "pr", "nd",
-    "pm", "sm", "eu", "gd", "tb", "dy", "ho", "er", "tm", "yb", "lu", "hf",
-    "ta", "w ", "re", "os", "ir", "pt", "au", "hg", "tl", "pb", "bi", "po",
-    "at", "rn", "fr", "ra", "ac", "th", "pa", "u ", "np", "pu", "am", "cm",
-    "bk", "cf", "es", "fm", "md", "no", "lr", "rf", "db", "sg", "bh", "hs",
-    "mt", "ds", "rg", "cn", "nh", "fl", "mc", "lv", "ts", "og",
-  };
-
   std::transform(sym.begin(), sym.end(), elem, ::tolower);
 
   for (int i = 0; i != 118; i++) {

include/dftd_cblas.h

@@ -26,7 +26,7 @@ namespace dftd4 {
 
 /**
  * @brief General matrix vector multiplication (`C = alpha * A * V + C`).
- * 
+ *
  * @param C Result vector C. Modified in-place.
  * @param A Matrix A.
  * @param V Vector V.
@@ -80,10 +80,9 @@ inline int BLAS_Add_Mat_x_Vec(
   return EXIT_FAILURE;
 };
 
-
 /**
  * @brief General matrix-matrix multiplication (`C = alpha * A * B + C`).
- * 
+ *
  * @param C Result matrix C. Modified in-place.
  * @param A Matrix A.
  * @param B Matrix B.
@@ -205,7 +204,7 @@ inline int BLAS_Add_Mat_x_Mat(
 
 /**
  * @brief Compute inverse of a matrix using LU decomposition.
- * 
+ *
  * @param a Matrix a.
  * @return Exit code.
  */

include/dftd_cutoff.h

@@ -33,18 +33,18 @@ static const double disp3_default = 40.0;
 
 // Collection of real space cutoffs.
 class TCutoff {
-    public:
-  double disp2;
-  double disp3;
-  double cn;
-  double cn_eeq;
-
-  explicit TCutoff(
-    double cut_disp2 = disp2_default,
-    double cut_disp3 = disp3_default,
-    double cut_cn = cn_default,
-    double cut_cn_eeq = cn_eeq_default
-  );
+  public:
+    double disp2;
+    double disp3;
+    double cn;
+    double cn_eeq;
+
+    explicit TCutoff(
+      double cut_disp2 = disp2_default,
+      double cut_disp3 = disp3_default,
+      double cut_cn = cn_default,
+      double cut_cn_eeq = cn_eeq_default
+    );
 };
 
 } // namespace dftd4

include/dftd_geometry.h

@@ -23,25 +23,25 @@
 namespace dftd4 {
 // Input of the molecular geometry
 class TMolecule {
-    public:
-  int NAtoms;
-  TMatrix<double> xyz; // Cartesian Coordinates: (NAtoms x 3)-matrix
-  TVector<int> at;     // atomic numbers
+  public:
+    int NAtoms;
+    TMatrix<double> xyz; // Cartesian Coordinates: (NAtoms x 3)-matrix
+    TVector<int> at;     // atomic numbers
 
-  TMolecule() { NAtoms = 0; }
-  ~TMolecule() { FreeMemory(); }
+    TMolecule() { NAtoms = 0; }
+    ~TMolecule() { FreeMemory(); }
 
-  void GetMemory(int NumAt_) {
-    FreeMemory();
-    NAtoms = NumAt_;
-    xyz.New(NAtoms, 3);
-    at.New(NAtoms);
-  }
+    void GetMemory(int NumAt_) {
+      FreeMemory();
+      NAtoms = NumAt_;
+      xyz.New(NAtoms, 3);
+      at.New(NAtoms);
+    }
 
-  void FreeMemory(void) {
-    xyz.Delete();
-    at.Delete();
-  }
+    void FreeMemory(void) {
+      xyz.Delete();
+      at.Delete();
+    }
 };
 
 } // namespace dftd4

include/dftd_readxyz.h

@@ -21,6 +21,19 @@
 
 #include "dftd_geometry.h"
 
+static const char pse[118][3]{
+  "h ", "he", "li", "be", "b ", "c ", "n ", "o ", "f ", "ne", "na", "mg",
+  "al", "si", "p ", "s ", "cl", "ar", "k ", "ca", "sc", "ti", "v ", "cr",
+  "mn", "fe", "co", "ni", "cu", "zn", "ga", "ge", "as", "se", "br", "kr",
+  "rb", "sr", "y ", "zr", "nb", "mo", "tc", "ru", "rh", "pd", "ag", "cd",
+  "in", "sn", "sb", "te", "i ", "xe", "cs", "ba", "la", "ce", "pr", "nd",
+  "pm", "sm", "eu", "gd", "tb", "dy", "ho", "er", "tm", "yb", "lu", "hf",
+  "ta", "w ", "re", "os", "ir", "pt", "au", "hg", "tl", "pb", "bi", "po",
+  "at", "rn", "fr", "ra", "ac", "th", "pa", "u ", "np", "pu", "am", "cm",
+  "bk", "cf", "es", "fm", "md", "no", "lr", "rf", "db", "sg", "bh", "hs",
+  "mt", "ds", "rg", "cn", "nh", "fl", "mc", "lv", "ts", "og",
+};
+
 extern void read_xyzfile(const std::string &, dftd4::TMolecule &);
 
 extern int element(const std::string &);

meson.build

@@ -20,7 +20,10 @@ project(
   'cpp',
   version: '2.0',
   license: 'LGPL3-or-later',
-  meson_version: '>=0.45',
+  meson_version: '>=0.46',
+  default_options: [
+    'default_library=both',
+  ],
 )
 
 cxx = meson.get_compiler('cpp')
@@ -41,8 +44,8 @@ srcs = files(
   'src/dftd_eeq.cpp',
   'src/dftd_model.cpp',
   'src/dftd_ncoord.cpp',
-  'src/damping/atm.cpp',
-  'src/damping/rational.cpp',
+  'src/damping/dftd_atm.cpp',
+  'src/damping/dftd_rational.cpp',
 )
 
 cpp_d4_inc = include_directories('include')
@@ -66,8 +69,8 @@ install_headers(
   'include/dftd_model.h',
   'include/dftd_ncoord.h',
   'include/dftd_parameters.h',
-  'include/damping/atm.h',
-  'include/damping/rational.h',
+  'include/damping/dftd_atm.h',
+  'include/damping/dftd_rational.h',
 )
 
 cpp_d4_dep = declare_dependency(
@@ -79,8 +82,8 @@ cpp_d4_dep = declare_dependency(
 executable(
   meson.project_name(),
   sources: files(
-    'src/program_dftd.cpp',
-    'src/dftd_readxyz.cpp',
+    'app/main.cpp',
+    'app/readxyz.cpp',
   ),
   dependencies: cpp_d4_dep,
   install: true,

src/damping/atm.cpp → src/damping/dftd_atm.cpp

@@ -30,7 +30,7 @@
 #include "dftd_matrix.h"
 #include "dftd_ncoord.h"
 #include "dftd_parameters.h"
-#include "damping/atm.h"
+#include "damping/dftd_atm.h"
 
 namespace dftd4 {
 

src/damping/rational.cpp → src/damping/dftd_rational.cpp

@@ -25,8 +25,8 @@
 #include "dftd_matrix.h"
 #include "dftd_ncoord.h"
 #include "dftd_parameters.h"
-#include "damping/atm.h"
-#include "damping/rational.h"
+#include "damping/dftd_atm.h"
+#include "damping/dftd_rational.h"
 
 
 namespace dftd4 {

src/dftd_dispersion.cpp

@@ -21,8 +21,8 @@
  */
 #include <cmath>
 
-#include "damping/atm.h"
-#include "damping/rational.h"
+#include "damping/dftd_atm.h"
+#include "damping/dftd_rational.h"
 #include "dftd_cblas.h"
 #include "dftd_dispersion.h"
 #include "dftd_eeq.h"

src/dftd_model.cpp

@@ -258,21 +258,21 @@ int TD4Model::get_atomic_c6(
 
 int TD4Model::set_refalpha_eeq(const TMolecule &mol, TMatrix<double> &alpha)
   const {
-  int iat{0}, is{0};
+  int izp{0}, is{0};
   double iz{0.0}, aiw{0.0};
 
-  for (int i = 0; i != mol.NAtoms; i++) {
-    iat = mol.at(i);
-    for (int ir = 0; ir != refn[iat]; ir++) {
-      is = refsys[iat][ir];
+  for (int iat = 0; iat != mol.NAtoms; iat++) {
+    izp = mol.at(iat);
+    for (int ir = 0; ir != refn[izp]; ir++) {
+      is = refsys[izp][ir];
       iz = zeff[is];
       for (int k = 0; k != 23; k++) {
         aiw = secscale[is] * secalpha[is][k] *
-              zeta(ga, gam[is] * gc, iz, refsq[iat][ir] + iz);
-        alpha(i, 23 * ir + k) = std::max(
+              zeta(ga, gam[is] * gc, iz, refsq[izp][ir] + iz);
+        alpha(iat, 23 * ir + k) = std::max(
           0.0,
-          refascale[iat][ir] *
-            (refalpha[iat][23 * ir + k] - refscount[iat][ir] * aiw)
+          refascale[izp][ir] *
+            (refalpha[izp][23 * ir + k] - refscount[izp][ir] * aiw)
         );
       }
     }

src/dftd_ncoord.cpp

@@ -23,12 +23,11 @@
  * This module works on a distance matrix to avoid recalculating
  * the distances every time.
  */
-#include "dftd_ncoord.h"
-
 #include <cmath>
 #include <iostream>
 
 #include "dftd_geometry.h"
+#include "dftd_ncoord.h"
 #include "dftd_matrix.h"
 
 namespace dftd4 {