make fix qeq/shielded fully compatible with pair style reax/c as a replacement for fix qeq/reax

src/QEQ/fix_qeq.cpp

@@ -100,6 +100,7 @@ FixQEq::FixQEq(LAMMPS *lmp, int narg, char **arg) :
   q1 = NULL;
   q2 = NULL;
   streitz_flag = 0;
+  reax_flag = 0;
   qv = NULL;
 
   comm_forward = comm_reverse = 1;
@@ -117,6 +118,8 @@ FixQEq::FixQEq(LAMMPS *lmp, int narg, char **arg) :
 
   if (strcmp(arg[7],"coul/streitz") == 0) {
     streitz_flag = 1;
+  } else if (strcmp(arg[7],"reax/c") == 0) {
+    reax_flag = 1;
   } else {
     read_file(arg[7]);
   }
@@ -138,7 +141,7 @@ FixQEq::~FixQEq()
 
   memory->destroy(shld);
 
-  if (!streitz_flag) {
+  if (!streitz_flag && !reax_flag) {
     memory->destroy(chi);
     memory->destroy(eta);
     memory->destroy(gamma);

src/QEQ/fix_qeq.h

@@ -70,7 +70,7 @@ class FixQEq : public Fix {
   double *chi,*eta,*gamma,*zeta,*zcore;  // qeq parameters
   double *chizj;
   double **shld;
-  int streitz_flag;
+  int streitz_flag, reax_flag;
 
   bigint ngroup;
 

src/QEQ/fix_qeq_shielded.cpp

@@ -29,6 +29,7 @@
 #include "update.h"
 #include "force.h"
 #include "group.h"
+#include "pair.h"
 #include "kspace.h"
 #include "respa.h"
 #include "memory.h"
@@ -39,7 +40,9 @@ using namespace LAMMPS_NS;
 /* ---------------------------------------------------------------------- */
 
 FixQEqShielded::FixQEqShielded(LAMMPS *lmp, int narg, char **arg) :
-  FixQEq(lmp, narg, arg) {}
+  FixQEq(lmp, narg, arg) {
+  if (reax_flag) extract_reax();
+}
 
 /* ---------------------------------------------------------------------- */
 
@@ -73,6 +76,22 @@ void FixQEqShielded::init()
 
 }
 
+/* ---------------------------------------------------------------------- */
+
+void FixQEqShielded::extract_reax()
+{
+  Pair *pair = force->pair_match("reax/c",1);
+  if (pair == NULL) error->all(FLERR,"No pair reax/c for fix qeq/shielded");
+  int tmp;
+  chi = (double *) pair->extract("chi",tmp);
+  eta = (double *) pair->extract("eta",tmp);
+  gamma = (double *) pair->extract("gamma",tmp);
+  if (chi == NULL || eta == NULL || gamma == NULL)
+    error->all(FLERR,
+        "Fix qeq/slater could not extract params from pair reax/c");
+}
+
+
 /* ---------------------------------------------------------------------- */
 
 void FixQEqShielded::init_shielding()

src/QEQ/fix_qeq_shielded.h

@@ -32,6 +32,7 @@ class FixQEqShielded : public FixQEq {
   void pre_force(int);
 
  private:
+  void extract_reax();
   void init_shielding();
   void init_matvec();
   void compute_H();

src/USER-REAXC/pair_reaxc.cpp

@@ -361,7 +361,8 @@ void PairReaxC::init_style( )
 
   int iqeq;
   for (iqeq = 0; iqeq < modify->nfix; iqeq++)
-    if (strstr(modify->fix[iqeq]->style,"qeq/reax")) break;
+    if (strstr(modify->fix[iqeq]->style,"qeq/reax")
+       || strstr(modify->fix[iqeq]->style,"qeq/shielded")) break;
   if (iqeq == modify->nfix && qeqflag == 1)
     error->all(FLERR,"Pair reax/c requires use of fix qeq/reax");