From 0c900e5047d597f52e50917fb011818054129724 Mon Sep 17 00:00:00 2001
From: "Daniel R. Roe" <daniel.r.roe@gmail.com>
Date: Fri, 30 Aug 2024 13:23:42 -0400
Subject: [PATCH] Minor fix for PSF with NUMLP (#1100)

* Make finding PSF tags more robust. Fix issue where PSF read could fail
when NUMLP tag not present

* Hide some debug info

* 6.29.1. Revision bump for PSF NUMLP fix
---
 src/Parm_CharmmPsf.cpp | 49 +++++++++++++++++++++++++++++++++---------
 src/Topology.cpp       | 11 +++++++---
 src/Version.h          |  2 +-
 3 files changed, 48 insertions(+), 14 deletions(-)

diff --git a/src/Parm_CharmmPsf.cpp b/src/Parm_CharmmPsf.cpp
index 1abbc424cd..91a65cc9ba 100644
--- a/src/Parm_CharmmPsf.cpp
+++ b/src/Parm_CharmmPsf.cpp
@@ -48,13 +48,32 @@ int Parm_CharmmPsf::processReadArgs(ArgList& argIn) {
 }
 
 // Parm_CharmmPsf::FindTag()
+/** Scan for the given PSF tag of format <#> <TAG>
+  * \return the <#> value of the tag, or -1 if no tag found.
+  */
 int Parm_CharmmPsf::FindTag(char* tag, const char* target, BufferedLine& infile) {
-  int nval = 0;
-  int tgtsize = strlen( target );
-  while (strncmp(tag,target,tgtsize)!=0) {
-    const char* buffer = infile.Line();
-    if ( buffer == 0 ) return 0;
-    sscanf(buffer,"%i %10s",&nval,tag);
+  int nval = -1;
+  // Check the current line
+  const char* buffer = infile.CurrentLine();
+  //mprintf("DEBUG: Current Line: %s\n", buffer);
+  const char* ptr = strstr(buffer, target);
+  while (ptr == 0) {
+    // Check subsequent lines
+    buffer = infile.Line();
+    if (buffer == 0) return -1;
+    //mprintf("DEBUG: Subsequent line: %s\n", buffer);
+    ptr = strstr(buffer, target);
+  }
+  if (ptr != 0) {
+    //mprintf("DEBUG: Line found: %s\n", buffer);
+    int nscan = sscanf(buffer,"%i %10s", &nval, tag);
+    // Sanity check
+    if (nscan != 2) {
+      mprinterr("Error: Could not get value for PSF tag %s\n", target);
+      return -1;
+    //int tgtsize = strlen( target );
+    //if (strncmp(tag, target, tgtsize) != 0)
+    }
   }
   return nval;
 }
@@ -521,10 +540,20 @@ int Parm_CharmmPsf::ReadParm(FileName const& fname, Topology &parmOut) {
   // Advance to <# lone pairs> <# lone pair hosts> NUMLP NUMLPH
   int numlp = -1;
   int numlph = -1;
-  while (strncmp(tag, "!NUMLP", 6) !=0) {
-    const char* buffer = infile.Line();
-    if ( buffer == 0 ) break;
-    sscanf(buffer,"%i %i %10s", &numlp, &numlph, tag);
+  buffer = infile.CurrentLine();
+  const char* ptr = strstr(buffer, "!NUMLP");
+  while (ptr == 0) {
+    buffer = infile.Line();
+    if (buffer == 0) break;
+    ptr = strstr(buffer, "!NUMLP");
+  }
+  if (ptr != 0) {
+    int nscan = sscanf(buffer, "%i %i %10s", &numlp, &numlph, tag);
+    // Sanity check
+    if (nscan != 3) {
+      mprinterr("Error: Could not scan values for !NUMLP NUMLPH\n");
+      return 1;
+    }
   }
   if (numlp > -1) {
     if (debug_ > 0) mprintf("DEBUG: PSF contains %i lone pairs, %i lone pair hosts.\n", numlp, numlph);
diff --git a/src/Topology.cpp b/src/Topology.cpp
index 678cc0fd01..a7061e3282 100644
--- a/src/Topology.cpp
+++ b/src/Topology.cpp
@@ -826,7 +826,8 @@ int Topology::RemoveBond(int atom1, int atom2)
   * For bonds to H always insert the H second.
   */
 void Topology::AddBond(int atom1, int atom2, int pidxIn) {
-  //mprintf("DEBUG: Enter AddBond(%i, %i, %i)\n", atom1+1, atom2+1, pidxIn);
+  //mprintf("DEBUG: Enter AddBond(%i, %i, %s, %s, %i)\n", atom1+1, atom2+1,
+  //        AtomMaskName(atom1).c_str(), AtomMaskName(atom2).c_str(), pidxIn);
   // Check if atoms are out of range.
   if (WarnOutOfRange(atoms_.size(), atom1, "bond")) return;
   if (WarnOutOfRange(atoms_.size(), atom2, "bond")) return;
@@ -851,9 +852,13 @@ void Topology::AddBond(int atom1, int atom2, int pidxIn) {
           if ( (it->A1() == atom1 && it->A2() == atom2) ||
                (it->A1() == atom2 && it->A2() == atom1) )
           {
-            mprintf("DEBUG: Existing bond found. Existing Idx %i\n", it->Idx());
+            if (debug_ > 0)
+              mprintf("DEBUG: Existing bond found. Existing Idx %i Rk=%f Req=%f\n",
+                      it->Idx(), bondparm_[it->Idx()].Rk(), bondparm_[it->Idx()].Req());
             if (it->Idx() < 0) {
-              mprintf("DEBUG: Adding bond parameter index %i for existing bond.\n", pidxIn);
+              if (debug_ > 0)
+                mprintf("DEBUG: Adding bond parameter index %i Rk=%f Req=%f for existing bond.\n",
+                        pidxIn, bondparm_[pidxIn].Rk(), bondparm_[pidxIn].Req());
               it->SetIdx( pidxIn );
             }
             break;
diff --git a/src/Version.h b/src/Version.h
index a63722cdb3..1c51e948e7 100644
--- a/src/Version.h
+++ b/src/Version.h
@@ -12,7 +12,7 @@
  * Whenever a number that precedes <revision> is incremented, all subsequent
  * numbers should be reset to 0.
  */
-#define CPPTRAJ_INTERNAL_VERSION "V6.29.0"
+#define CPPTRAJ_INTERNAL_VERSION "V6.29.1"
 /// PYTRAJ relies on this
 #define CPPTRAJ_VERSION_STRING CPPTRAJ_INTERNAL_VERSION
 #endif