diff --git a/doc/src/fix_adapt.txt b/doc/src/fix_adapt.txt
index a35357a7ec811236d9d5fe5cb17865423bf7b700..d7c32bef3d0d360fb52cc098684ff8dc9ba3c08b 100644
--- a/doc/src/fix_adapt.txt
+++ b/doc/src/fix_adapt.txt
@@ -22,6 +22,11 @@ attribute = {pair} or {kspace} or {atom} :l
pparam = parameter to adapt over time
I,J = type pair(s) to set parameter for
v_name = variable with name that calculates value of pparam
+ {bond} args = bstyle bparam I v_name
+ bstyle = bond style name, e.g. harmonic
+ bparam = parameter to adapt over time
+ I = type bond to set parameter for
+ v_name = variable with name that calculates value of bparam
{kspace} arg = v_name
v_name = variable with name that calculates scale factor on K-space terms
{atom} args = aparam v_name
@@ -42,7 +47,10 @@ keyword = {scale} or {reset} :l
fix 1 all adapt 1 pair soft a 1 1 v_prefactor
fix 1 all adapt 1 pair soft a 2* 3 v_prefactor
fix 1 all adapt 1 pair lj/cut epsilon * * v_scale1 coul/cut scale 3 3 v_scale2 scale yes reset yes
-fix 1 all adapt 10 atom diameter v_size :pre
+fix 1 all adapt 10 atom diameter v_size
+
+variable ramp_up equal "ramp(0.01,0.5)"
+fix stretch all adapt 1 bond harmonic r0 1 v_ramp_up :pre
[Description:]
@@ -192,6 +200,19 @@ fix 1 all adapt 1 pair soft a * * v_prefactor :pre
:line
+The {bond} keyword uses the specified variable to change the value of
+a bond coefficient over time, very similar to how the {pair} keyword
+operates. The only difference is that now a bond coefficient for a
+given bond type is adapted.
+
+Currently {bond} does not support bond_style hybrid nor bond_style
+hybrid/overlay as bond styles. The only bonds that currently are
+working with fix_adapt are
+
+"harmonic"_bond_harmonic.html: k,r0: type bonds :tb(c=3,s=:)
+
+:line
+
The {kspace} keyword used the specified variable as a scale factor on
the energy, forces, virial calculated by whatever K-Space solver is
defined by the "kspace_style"_kspace_style.html command. If the
diff --git a/src/MOLECULE/bond_harmonic.cpp b/src/MOLECULE/bond_harmonic.cpp
index f164a51de4309d57cff2d06a90ffa6b67c203639..0763d7d3e24667be3c3d81adbd236c16c3a54c7b 100644
--- a/src/MOLECULE/bond_harmonic.cpp
+++ b/src/MOLECULE/bond_harmonic.cpp
@@ -13,6 +13,7 @@
#include <math.h>
#include <stdlib.h>
+#include <string.h>
#include "bond_harmonic.h"
#include "atom.h"
#include "neighbor.h"
@@ -26,7 +27,10 @@ using namespace LAMMPS_NS;
/* ---------------------------------------------------------------------- */
-BondHarmonic::BondHarmonic(LAMMPS *lmp) : Bond(lmp) {}
+BondHarmonic::BondHarmonic(LAMMPS *lmp) : Bond(lmp)
+{
+ reinitflag = 1;
+}
/* ---------------------------------------------------------------------- */
@@ -196,3 +200,16 @@ double BondHarmonic::single(int type, double rsq, int i, int j,
if (r > 0.0) fforce = -2.0*rk/r;
return rk*dr;
}
+
+/* ----------------------------------------------------------------------
+ Return ptr to internal members upon request.
+------------------------------------------------------------------------ */
+void *BondHarmonic::extract( char *str, int &dim )
+{
+ dim = 1;
+ if( strcmp(str,"kappa")==0) return (void*) k;
+ if( strcmp(str,"r0")==0) return (void*) r0;
+ return NULL;
+}
+
+
diff --git a/src/MOLECULE/bond_harmonic.h b/src/MOLECULE/bond_harmonic.h
index 7c7125b04ce233b5543b552ceb65b89ffea146f9..a0fd24577a8988e18fcb7a676a0041c353c52191 100644
--- a/src/MOLECULE/bond_harmonic.h
+++ b/src/MOLECULE/bond_harmonic.h
@@ -36,6 +36,7 @@ class BondHarmonic : public Bond {
void read_restart(FILE *);
void write_data(FILE *);
double single(int, double, int, int, double &);
+ virtual void *extract(char *, int &);
protected:
double *k,*r0;
diff --git a/src/bond.cpp b/src/bond.cpp
index 5a33f107cf44ade599dff8b5a9a9b81f00d04386..825ff1b199a48b854f5b7d716eac80c77228345e 100644
--- a/src/bond.cpp
+++ b/src/bond.cpp
@@ -292,3 +292,14 @@ double Bond::memory_usage()
bytes += comm->nthreads*maxvatom*6 * sizeof(double);
return bytes;
}
+
+/* -----------------------------------------------------------------------
+ Reset all type-based bond params via init.
+-------------------------------------------------------------------------- */
+void Bond::reinit()
+{
+ if (!reinitflag)
+ error->all(FLERR,"Fix adapt interface to this bond style not supported");
+
+ init();
+}
diff --git a/src/bond.h b/src/bond.h
index 41604387a3f8e056a92e7b08a25007a062ca3021..29de7ad7d2e5c907632f63c987af15e3ec75328d 100644
--- a/src/bond.h
+++ b/src/bond.h
@@ -30,6 +30,8 @@ class Bond : protected Pointers {
double virial[6]; // accumulated virial
double *eatom,**vatom; // accumulated per-atom energy/virial
+ int reinitflag; // 1 if compatible with fix adapt and alike
+
// KOKKOS host/device flag and data masks
ExecutionSpace execution_space;
@@ -49,6 +51,8 @@ class Bond : protected Pointers {
virtual void write_data(FILE *) {}
virtual double single(int, double, int, int, double &) = 0;
virtual double memory_usage();
+ virtual void *extract(char *, int &) {return NULL;}
+ virtual void reinit();
void write_file(int, char**);
diff --git a/src/fix_adapt.cpp b/src/fix_adapt.cpp
index 4c7eb5f218e4a7d0a1e83996ca10cbe6d74c8b07..07559c77198ec79148e07250cae2bd1410351db5 100644
--- a/src/fix_adapt.cpp
+++ b/src/fix_adapt.cpp
@@ -16,6 +16,7 @@
#include <stdlib.h>
#include "fix_adapt.h"
#include "atom.h"
+#include "bond.h"
#include "update.h"
#include "group.h"
#include "modify.h"
@@ -35,7 +36,7 @@ using namespace LAMMPS_NS;
using namespace FixConst;
using namespace MathConst;
-enum{PAIR,KSPACE,ATOM};
+enum{PAIR,KSPACE,ATOM,BOND};
enum{DIAMETER,CHARGE};
/* ---------------------------------------------------------------------- */
@@ -68,6 +69,10 @@ nadapt(0), id_fix_diam(NULL), id_fix_chg(NULL), adapt(NULL)
if (iarg+3 > narg) error->all(FLERR,"Illegal fix adapt command");
nadapt++;
iarg += 3;
+ } else if (strcmp(arg[iarg],"bond") == 0 ){
+ if (iarg+5 > narg) error->all(FLERR,"Illegal fix adapt command");
+ nadapt++;
+ iarg += 5;
} else break;
}
@@ -103,6 +108,25 @@ nadapt(0), id_fix_diam(NULL), id_fix_chg(NULL), adapt(NULL)
} else error->all(FLERR,"Illegal fix adapt command");
nadapt++;
iarg += 6;
+ } else if (strcmp(arg[iarg],"bond") == 0 ){
+ if (iarg+5 > narg) error->all(FLERR, "Illegal fix adapt command");
+ adapt[nadapt].which = BOND;
+ int n = strlen(arg[iarg+1]) + 1;
+ adapt[nadapt].bstyle = new char[n];
+ strcpy(adapt[nadapt].bstyle,arg[iarg+1]);
+ n = strlen(arg[iarg+2]) + 1;
+ adapt[nadapt].bparam = new char[n];
+ adapt[nadapt].bond = NULL;
+ strcpy(adapt[nadapt].bparam,arg[iarg+2]);
+ force->bounds(FLERR,arg[iarg+3],atom->ntypes,
+ adapt[nadapt].ilo,adapt[nadapt].ihi);
+ if (strstr(arg[iarg+4],"v_") == arg[iarg+4]) {
+ n = strlen(&arg[iarg+4][2]) + 1;
+ adapt[nadapt].var = new char[n];
+ strcpy(adapt[nadapt].var,&arg[iarg+4][2]);
+ } else error->all(FLERR,"Illegal fix adapt command");
+ nadapt++;
+ iarg += 5;
} else if (strcmp(arg[iarg],"kspace") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal fix adapt command");
adapt[nadapt].which = KSPACE;
@@ -160,6 +184,13 @@ nadapt(0), id_fix_diam(NULL), id_fix_chg(NULL), adapt(NULL)
for (int m = 0; m < nadapt; m++)
if (adapt[m].which == PAIR)
memory->create(adapt[m].array_orig,n+1,n+1,"adapt:array_orig");
+
+ // allocate bond style arrays:
+ n = atom->nbondtypes;
+ for (int m = 0; m < nadapt; ++m)
+ if (adapt[m].which == BOND)
+ // For now just use same storage and fake it to be one-dimensional:
+ memory->create(adapt[m].vector_orig,n+1,"adapt:vector_orig");
}
/* ---------------------------------------------------------------------- */
@@ -172,6 +203,10 @@ FixAdapt::~FixAdapt()
delete [] adapt[m].pstyle;
delete [] adapt[m].pparam;
memory->destroy(adapt[m].array_orig);
+ } else if (adapt[m].which == BOND) {
+ delete [] adapt[m].bstyle;
+ delete [] adapt[m].bparam;
+ memory->destroy(adapt[m].vector_orig);
}
}
delete [] adapt;
@@ -282,6 +317,7 @@ void FixAdapt::init()
// setup and error checks
anypair = 0;
+ anybond = 0;
for (int m = 0; m < nadapt; m++) {
Adapt *ad = &adapt[m];
@@ -350,7 +386,44 @@ void FixAdapt::init()
}
delete [] pstyle;
+ } else if (ad->which == BOND){
+ ad->bond = NULL;
+ anybond = 1;
+ // Use same routines from pair to strip any suffices:
+
+ int n = strlen(ad->bstyle) + 1;
+ char *bstyle = new char[n];
+ strcpy(bstyle,ad->bstyle);
+
+ if (lmp->suffix_enable) {
+ int len = 2 + strlen(bstyle) + strlen(lmp->suffix);
+ char *bsuffix = new char[len];
+ strcpy(bsuffix,bstyle);
+ strcat(bsuffix,"/");
+ strcat(bsuffix,lmp->suffix);
+ ad->bond = force->bond_match(bsuffix);
+ delete [] bsuffix;
+ }
+ // If not set grab regular one instead:
+ if (ad->bond == NULL) ad->bond = force->bond_match(bstyle);
+ if (ad->bond == NULL )
+ error->all(FLERR,"Fix adapt bond style does not exist");
+
+ void *ptr = ad->bond->extract(ad->bparam,ad->bdim);
+
+ if (ptr == NULL)
+ error->all(FLERR,"Fix adapt bond style param not supported");
+
+ // For bond styles you should use a vector
+
+ if (ad->bdim == 1) ad->vector = (double *) ptr;
+
+ if (strcmp(force->bond_style,"hybrid") == 0 ||
+ strcmp(force->bond_style,"hybrid_overlay") == 0)
+ error->all(FLERR,"Fix adapt does not support bond_style hybrid");
+ delete [] bstyle;
+
} else if (ad->which == KSPACE) {
if (force->kspace == NULL)
error->all(FLERR,"Fix adapt kspace style does not exist");
@@ -368,7 +441,7 @@ void FixAdapt::init()
}
}
- // make copy of original pair array values
+ // make copy of original pair/bond array values
for (int m = 0; m < nadapt; m++) {
Adapt *ad = &adapt[m];
@@ -378,7 +451,12 @@ void FixAdapt::init()
ad->array_orig[i][j] = ad->array[i][j];
}else if (ad->which == PAIR && ad->pdim == 0){
ad->scalar_orig = *ad->scalar;
+
+ }else if (ad->which == BOND && ad->bdim == 1){
+ for (i = ad->ilo; i <= ad->ihi; ++i )
+ ad->vector_orig[i] = ad->vector[i];
}
+
}
// fixes that store initial per-atom values
@@ -470,6 +548,18 @@ void FixAdapt::change_settings()
ad->array[i][j] = value;
}
+ // set bond type array values:
+
+ } else if (ad->which == BOND) {
+ if (ad->bdim == 1){
+ if (scaleflag)
+ for (i = ad->ilo; i <= ad->ihi; ++i )
+ ad->vector[i] = value*ad->vector_orig[i];
+ else
+ for (i = ad->ilo; i <= ad->ihi; ++i )
+ ad->vector[i] = value;
+ }
+
// set kspace scale factor
} else if (ad->which == KSPACE) {
@@ -532,6 +622,14 @@ void FixAdapt::change_settings()
}
}
}
+ if (anybond) {
+ for (int m = 0; m < nadapt; ++m ) {
+ Adapt *ad = &adapt[m];
+ if (ad->which == BOND) {
+ ad->bond->reinit();
+ }
+ }
+ }
// reset KSpace charges if charges have changed
@@ -554,6 +652,12 @@ void FixAdapt::restore_settings()
ad->array[i][j] = ad->array_orig[i][j];
}
+ } else if (ad->which == BOND) {
+ if (ad->pdim == 1) {
+ for (int i = ad->ilo; i <= ad->ihi; i++)
+ ad->vector[i] = ad->vector_orig[i];
+ }
+
} else if (ad->which == KSPACE) {
*kspace_scale = 1.0;
@@ -588,6 +692,7 @@ void FixAdapt::restore_settings()
}
if (anypair) force->pair->reinit();
+ if (anybond) force->bond->reinit();
if (chgflag && force->kspace) force->kspace->qsum_qsq();
}
diff --git a/src/fix_adapt.h b/src/fix_adapt.h
index a6d45c78cc4dc5b1f2cca9dc24ccf9fad6a50cb7..6e49f4a284225cffaab614f9a7562eada287bd1c 100644
--- a/src/fix_adapt.h
+++ b/src/fix_adapt.h
@@ -43,7 +43,7 @@ class FixAdapt : public Fix {
private:
int nadapt,resetflag,scaleflag;
- int anypair;
+ int anypair, anybond;
int nlevels_respa;
char *id_fix_diam,*id_fix_chg;
class FixStore *fix_diam,*fix_chg;
@@ -52,12 +52,15 @@ class FixAdapt : public Fix {
int which,ivar;
char *var;
char *pstyle,*pparam;
+ char *bstyle,*bparam;
int ilo,ihi,jlo,jhi;
- int pdim;
+ int pdim,bdim;
double *scalar,scalar_orig;
+ double *vector,*vector_orig;
double **array,**array_orig;
int aparam;
class Pair *pair;
+ class Bond *bond;
};
Adapt *adapt;