add born/matrix support to hybrid pair styles

src/pair.h

@@ -170,10 +170,9 @@ class Pair : protected Pointers {
   }
 
   virtual void born_matrix(int /*i*/, int /*j*/, int /*itype*/, int /*jtype*/, double /*rsq*/,
-               double /*factor_coul*/, double /*factor_lj*/, double& du, double& du2)
+                           double /*factor_coul*/, double /*factor_lj*/, double &du, double &du2)
   {
-    du = 0.0;
+    du = du2 = 0.0;
-    du2 = 0.0;
   }
 
   virtual void settings(int, char **) = 0;

src/pair_hybrid.cpp

@@ -874,6 +874,40 @@ double PairHybrid::single(int i, int j, int itype, int jtype,
   return esum;
 }
 
+/* ----------------------------------------------------------------------
+   call sub-style to compute born matrix interaction
+   error if sub-style does not support born_matrix call
+   since overlay could have multiple sub-styles, sum results explicitly
+------------------------------------------------------------------------- */
+
+void PairHybrid::born_matrix(int i, int j, int itype, int jtype, double rsq,
+                             double factor_coul, double factor_lj,
+                             double &dupair, double &du2pair)
+{
+  if (nmap[itype][jtype] == 0)
+    error->one(FLERR,"Invoked pair born_matrix on pair style none");
+
+  double du, du2;
+  dupair = du2pair = 0.0;
+
+  for (int m = 0; m < nmap[itype][jtype]; m++) {
+    if (rsq < styles[map[itype][jtype][m]]->cutsq[itype][jtype]) {
+      if (styles[map[itype][jtype][m]]->born_matrix_enable == 0)
+        error->one(FLERR,"Pair hybrid sub-style does not support born_matrix call");
+
+      if ((special_lj[map[itype][jtype][m]] != nullptr) ||
+          (special_coul[map[itype][jtype][m]] != nullptr))
+        error->one(FLERR,"Pair hybrid born_matrix calls do not support"
+                   " per sub-style special bond values");
+
+      du = du2 = 0.0;
+      styles[map[itype][jtype][m]]->born_matrix(i,j,itype,jtype,rsq,factor_coul,factor_lj,du,du2);
+      dupair += du;
+      du2pair += du2;
+    }
+  }
+}
+
 /* ----------------------------------------------------------------------
    copy Pair::svector data
 ------------------------------------------------------------------------- */

src/pair_hybrid.h

@@ -49,6 +49,8 @@ class PairHybrid : public Pair {
   void write_restart(FILE *) override;
   void read_restart(FILE *) override;
   double single(int, int, int, int, double, double, double, double &) override;
+  void born_matrix(int, int, int, int, double, double, double, double &, double &) override;
+
   void modify_params(int narg, char **arg) override;
   double memory_usage() override;
 

src/pair_hybrid_scaled.cpp

@@ -413,7 +413,7 @@ double PairHybridScaled::single(int i, int j, int itype, int jtype, double rsq,
           (special_coul[map[itype][jtype][m]] != nullptr))
         error->one(FLERR, "Pair hybrid single() does not support per sub-style special_bond");
 
-      scale = scaleval[map[itype][jtype][m]];
+      double scale = scaleval[map[itype][jtype][m]];
       esum += scale * pstyle->single(i, j, itype, jtype, rsq, factor_coul, factor_lj, fone);
       fforce += scale * fone;
     }
@@ -423,6 +423,57 @@ double PairHybridScaled::single(int i, int j, int itype, int jtype, double rsq,
   return esum;
 }
 
+/* ----------------------------------------------------------------------
+   call sub-style to compute born matrix interaction
+   error if sub-style does not support born_matrix call
+   since overlay could have multiple sub-styles, sum results explicitly
+------------------------------------------------------------------------- */
+
+void PairHybridScaled::born_matrix(int i, int j, int itype, int jtype, double rsq,
+                                   double factor_coul, double factor_lj,
+                                   double &dupair, double &du2pair)
+{
+  if (nmap[itype][jtype] == 0) error->one(FLERR, "Invoked pair born_matrix on pair style none");
+
+  // update scale values from variables where needed
+
+  const int nvars = scalevars.size();
+  if (nvars > 0) {
+    double *vals = new double[nvars];
+    for (int k = 0; k < nvars; ++k) {
+      int m = input->variable->find(scalevars[k].c_str());
+      if (m < 0)
+        error->all(FLERR, "Variable '{}' not found when updating scale factors", scalevars[k]);
+      vals[k] = input->variable->compute_equal(m);
+    }
+    for (int k = 0; k < nstyles; ++k) {
+      if (scaleidx[k] >= 0) scaleval[k] = vals[scaleidx[k]];
+    }
+    delete[] vals;
+  }
+
+  double du, du2, scale;
+  dupair = du2pair = scale = 0.0;
+
+  for (int m = 0; m < nmap[itype][jtype]; m++) {
+    auto pstyle = styles[map[itype][jtype][m]];
+    if (rsq < pstyle->cutsq[itype][jtype]) {
+      if (pstyle->born_matrix_enable == 0)
+        error->one(FLERR, "Pair hybrid sub-style does not support born_matrix call");
+
+      if ((special_lj[map[itype][jtype][m]] != nullptr) ||
+          (special_coul[map[itype][jtype][m]] != nullptr))
+        error->one(FLERR, "Pair hybrid born_matrix() does not support per sub-style special_bond");
+
+      du = du2 = 0.0;
+      scale = scaleval[map[itype][jtype][m]];
+      pstyle->born_matrix(i, j, itype, jtype, rsq, factor_coul, factor_lj, du, du2);
+      dupair += scale*du;
+      du2pair += scale*du2;
+    }
+  }
+}
+
 /* ----------------------------------------------------------------------
    set coeffs for one or more type pairs
 ------------------------------------------------------------------------- */

src/pair_hybrid_scaled.h

@@ -38,6 +38,7 @@ class PairHybridScaled : public PairHybrid {
   void write_restart(FILE *) override;
   void read_restart(FILE *) override;
   double single(int, int, int, int, double, double, double, double &) override;
+  void born_matrix(int, int, int, int, double, double, double, double &, double &) override;
 
   void init_svector() override;
   void copy_svector(int, int) override;