Fix ghost atom handling

.gitignore

@@ -35,3 +35,8 @@
 /build*/
 /install*/
 /_*/
+
+# project-specific
+/scripts/old/*
+/scripts/*.h
+/scripts/*.cpp

include/dftd_ncoord.h

@@ -200,4 +200,4 @@ extern inline double log_cn_cut(const double cn_max, const double cn);
 
 extern inline double dlog_cn_cut(const double cn_max, const double cn);
 
-}; // namespace dftd4
+} // namespace dftd4

include/dftd_parameters.h

@@ -63,10 +63,10 @@ static const double zeff[119]{
  * These values are actually never used in the code. Only r4r2 is used.
  */
 static const double r2r4[119]{
-  0.0,             //  dummy
-  8.0589,  3.4698, //  H,He
+  0.0,                                                           //  dummy
+  8.0589,  3.4698,                                               //  H,He
   29.0974, 14.8517, 11.8799, 7.8715,  5.5588,  4.7566,  3.8025,
-  3.1036, //  Li-Ne
+  3.1036,                                                        //  Li-Ne
   26.1552, 17.2304, 17.7210, 12.7442, 9.5361,  8.1652,  6.7463,
   5.6004,                                                        //  Na-Ar
   29.2012, 22.3934,                                              //  K,Ca

scripts/README.md

@@ -0,0 +1,3 @@
+# Scripts
+
+Scripts for creating the ORCA source files. This mainly concatenates the corresponding files and adds headers and footers.

scripts/make-qcdftd4.sh

@@ -0,0 +1,317 @@
+#!/bin/bash
+
+INCLUDE="../include"
+SRC="../src"
+
+cat > qcdftd4param.h << 'EOT'
+/*
+ * This is the DFT-D4 equivalent of copyc6, since this is an *empirical*
+ * dispersion correction we can hardly avoid having this big, scary blob
+ * of data lying around somewhere.
+ *
+ * Files like this are usually computer generated so avoid modifying it,
+ * but have the appropriate tool generate it from dftd4.
+ *
+ * Responsible for this mess:
+ * Sebastian Ehlert <[email protected]> (SAW190521)
+ *
+ * Update/Fix for parameters from periodic extension by:
+ * Marvin Friede <[email protected]> (MF161222)
+ *
+ * Extension for Fr, Ra and Actinides by:
+ * Marvin Friede <[email protected]> (MF121223)
+ */
+
+#ifndef QCDFTD4PARAM_H
+#define QCDFTD4PARAM_H
+
+EOT
+
+sed -n '/namespace/,$p' "${INCLUDE}/dftd_parameters.h" >> qcdftd4param.h
+
+# Append the Footer:
+printf "\n#endif // QCDFTD4PARAM_H\n" >> qcdftd4param.h
+
+
+##################################################################
+##################################################################
+
+dosed(){
+	sed -n '/namespace dftd4 {/,/} \/\/ namespace dftd4/p' $1 | \
+	 sed '/namespace dftd4 {/d; /} \/\/ namespace dftd4/d' | \
+	 awk '/./,EOF' > $2
+}
+
+dosed "${INCLUDE}/dftd_cutoff.h" cutoff.txt
+dosed "${INCLUDE}/dftd_model.h" model.txt
+dosed "${INCLUDE}/dftd_dispersion.h" dispersion.txt
+dosed "${INCLUDE}/damping/dftd_rational.h" rational.txt
+dosed "${INCLUDE}/damping/dftd_atm.h" atm.txt
+
+cat > qcdftd4.h << 'EOT'
+/*
+ * D4(EEQ)-ATM implementation
+ *
+ * Sebastian Ehlert <[email protected]> (SAW190521)
+ * Marvin Friede <[email protected]> (MF161222)
+ */
+
+#ifndef __QCDFTD4_H
+#define __QCDFTD4_H
+
+#include "qcinpdat.h" // TGeomInput class
+#include "qcmat1.h"   // TRVector and TRMatrix class
+
+namespace dftd4 {
+
+EOT
+
+######
+
+cat >> qcdftd4.h << 'EOT'
+/* --------------------------------------------------------------------------
+
+  Cutoff (dftd_cutoff.h)
+  https://github.com/dftd4/cpp-d4/blob/main/include/dftd_cutoff.h
+
+  -------------------------------------------------------------------------- */
+
+EOT
+
+cat cutoff.txt >> qcdftd4.h
+rm cutoff.txt
+
+######
+
+cat >> qcdftd4.h << 'EOT'
+/* --------------------------------------------------------------------------
+
+  D4 model (dftd_model.h)
+  https://github.com/dftd4/cpp-d4/blob/main/include/dftd_model.h
+
+  -------------------------------------------------------------------------- */
+
+EOT
+
+cat model.txt >> qcdftd4.h
+rm model.txt
+
+######
+
+cat >> qcdftd4.h << 'EOT'
+/* --------------------------------------------------------------------------
+
+  Dispersion (dftd_dispersion.h)
+  https://github.com/dftd4/cpp-d4/blob/main/include/dftd_dispersion.h
+
+  -------------------------------------------------------------------------- */
+
+EOT
+
+cat dispersion.txt >> qcdftd4.h
+rm dispersion.txt
+
+######
+
+cat >> qcdftd4.h << 'EOT'
+/* --------------------------------------------------------------------------
+
+  Dispersion (damping/dftd_rational.h)
+  https://github.com/dftd4/cpp-d4/blob/main/include/damping/dftd_rational.h
+
+  -------------------------------------------------------------------------- */
+
+EOT
+
+cat rational.txt >> qcdftd4.h
+rm rational.txt
+
+######
+
+cat >> qcdftd4.h << 'EOT'
+/* --------------------------------------------------------------------------
+
+  Dispersion (damping/dftd_atm.h)
+  https://github.com/dftd4/cpp-d4/blob/main/include/damping/dftd_atm.h
+
+  -------------------------------------------------------------------------- */
+
+EOT
+
+cat atm.txt >> qcdftd4.h
+rm atm.txt
+
+######
+
+cat >> qcdftd4.h << 'EOT'
+} // namespace dftd4
+
+/* --------------------------------------------------------------------------------------
+// Calculates the EEQ charges according to the D4 paper
+// https://doi.org/10.1063/1.5090222
+// Bernardo de Souza, 14/09/2023
+ -------------------------------------------------------------------------------------- */
+int CalcEEQCharges(TRMatrix &XYZ, TIVector &ATNO, int NAtoms, int totalcharge, TRVector &q, bool printerror=true);
+
+#endif // __QCDFTD4_H
+
+EOT
+
+
+##################################################################
+##################################################################
+
+
+dosed "${SRC}/dftd_cutoff.cpp" cutoff.txt
+dosed "${SRC}/dftd_model.cpp" model.txt
+dosed "${SRC}/dftd_dispersion.cpp" dispersion.txt
+dosed "${SRC}/damping/dftd_rational.cpp" rational.txt
+dosed "${SRC}/damping/dftd_atm.cpp" atm.txt
+
+cat > qcdftd4.cpp << 'EOT'
+/*
+ * D4(EEQ)-ATM implementation
+ *
+ * Sebastian Ehlert <[email protected]> (SAW190521)
+ * Marvin Friede <[email protected]> (MF161222)
+ */
+#include <cmath>
+#include <limits>
+
+#include "qcinpdat.h" // TGeomInput class
+#include "qcmat2.h"   // BLAS routines
+#include "qcmath.h"   // TVector and TMatrix class
+
+// we cannot avoid one big, scary parameter file...
+#include "qcdftd4param.h"
+
+#include "qceeq.h"
+#include "qcncoord.h"
+
+// always include self
+#include "qcdftd4.h"
+
+/* --------------------------------------------------------------------------------------
+// Calculates the EEQ charges according to the D4 paper
+// https://doi.org/10.1063/1.5090222
+// Bernardo de Souza, 14/09/2023
+ -------------------------------------------------------------------------------------- */
+int CalcEEQCharges(TRMatrix &XYZ, TIVector &ATNO, int NAtoms, int totalcharge, TRVector &q, bool printerror){
+
+  // initialize the charges to zero
+  q.Init();
+
+  // check atomic numbers to guarantee we have all parameters
+  // ghost atoms (ATNO=0) will have no charge
+  for (int i=0;i<NAtoms;i++)
+    if (ATNO(i)>86){
+      if (printerror) printMessage("Atomic number %d detected. EEQ charges can not be calculated.\n",ATNO(i));
+      return 1;
+    }
+
+  TGeomInput molecule;
+  molecule.NAtoms=NAtoms;
+  molecule.CC.CopyMat(XYZ);
+  molecule.ATNO.CopyVec(ATNO);
+  TIVector Index(NAtoms);
+  for (int i=0;i<NAtoms;i++) Index(i)=i;
+  TRMatrix dist;
+  dist.NewMatrix(NAtoms, NAtoms);
+  dftd4::calc_distances(molecule, Index, dist);
+
+  TRMatrix dqdr;
+
+  return dftd4::get_charges(molecule, Index, dist, totalcharge, 25, q, dqdr, false);
+}
+
+namespace dftd4 {
+
+EOT
+
+######
+
+cat >> qcdftd4.cpp << 'EOT'
+/* --------------------------------------------------------------------------
+
+  Cutoff (dftd_cutoff.cpp)
+  https://github.com/dftd4/cpp-d4/blob/main/src/dftd_cutoff.cpp
+
+  -------------------------------------------------------------------------- */
+
+EOT
+
+cat cutoff.txt >> qcdftd4.cpp
+rm cutoff.txt
+
+######
+
+cat >> qcdftd4.cpp << 'EOT'
+/* --------------------------------------------------------------------------
+
+  D4 model (dftd_model.cpp)
+  https://github.com/dftd4/cpp-d4/blob/main/src/dftd_model.cpp
+
+  -------------------------------------------------------------------------- */
+
+EOT
+
+cat model.txt >> qcdftd4.cpp
+rm model.txt
+
+######
+
+cat >> qcdftd4.cpp << 'EOT'
+/* --------------------------------------------------------------------------
+
+  Dispersion (dftd_dispersion.cpp)
+  https://github.com/dftd4/cpp-d4/blob/main/src/dftd_dispersion.cpp
+
+  -------------------------------------------------------------------------- */
+
+EOT
+
+cat dispersion.txt >> qcdftd4.cpp
+rm dispersion.txt
+
+######
+
+cat >> qcdftd4.cpp << 'EOT'
+/* --------------------------------------------------------------------------
+
+  Dispersion (damping/dftd_rational.cpp)
+  https://github.com/dftd4/cpp-d4/blob/main/src/damping/dftd_rational.cpp
+
+  -------------------------------------------------------------------------- */
+
+EOT
+
+cat rational.txt >> qcdftd4.cpp
+rm rational.txt
+
+######
+
+cat >> qcdftd4.cpp << 'EOT'
+/* --------------------------------------------------------------------------
+
+  Dispersion (damping/dftd_atm.cpp)
+  https://github.com/dftd4/cpp-d4/blob/main/src/damping/dftd_atm.cpp
+
+  -------------------------------------------------------------------------- */
+
+EOT
+
+cat atm.txt >> qcdftd4.cpp
+rm atm.txt
+
+######
+
+cat >> qcdftd4.cpp << 'EOT'
+} // namespace dftd4
+EOT
+
+
+##########################################################################
+
+
+sed -i "s/TMolecule/TGeomInput/" *.h
+sed -i "s/TMolecule/TGeomInput/" *.cpp

scripts/make-qcdftd4param.sh

@@ -0,0 +1,34 @@
+#!/bin/bash
+
+INCLUDE="../include"
+
+cat > qcdftd4param.h << 'EOT'
+/*
+ * This is the DFT-D4 equivalent of copyc6, since this is an *empirical*
+ * dispersion correction we can hardly avoid having this big, scary blob
+ * of data lying around somewhere.
+ *
+ * Files like this are usually computer generated so avoid modifying it,
+ * but have the appropriate tool generate it from dftd4.
+ *
+ * Responsible for this mess:
+ * Sebastian Ehlert <[email protected]> (SAW190521)
+ *
+ * Update/Fix for parameters from periodic extension by:
+ * Marvin Friede <[email protected]> (MF161222)
+ *
+ * Extension for Fr, Ra and Actinides by:
+ * Marvin Friede <[email protected]> (MF121223)
+ */
+
+#ifndef QCDFTD4PARAM_H
+#define QCDFTD4PARAM_H
+
+EOT
+
+sed -n '/namespace/,$p' "${INCLUDE}/dftd_parameters.h" >> qcdftd4param.h
+
+printf "\n#endif // QCDFTD4PARAM_H\n" >> qcdftd4param.h
+
+
+sed -i "s/TMolecule/TGeomInput/" qcdftd4param.h