From 39e8123a02d8c211a09e421e77d22faa4d0aa1fa Mon Sep 17 00:00:00 2001 From: athomps Date: Tue, 13 Oct 2015 15:53:41 +0000 Subject: [PATCH 1/4] added tfac_insert keyword git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@14109 f3b2605a-c512-4ea7-a41b-209d697bcdaa --- doc/fix_gcmc.txt | 22 +++++++++++++++++++--- 1 file changed, 19 insertions(+), 3 deletions(-) diff --git a/doc/fix_gcmc.txt b/doc/fix_gcmc.txt index 7fb0ff6966..e4a0df29cc 100644 --- a/doc/fix_gcmc.txt +++ b/doc/fix_gcmc.txt @@ -23,7 +23,7 @@ T = temperature of the ideal gas reservoir (temperature units) :l mu = chemical potential of the ideal gas reservoir (energy units) :l translate = maximum Monte Carlo translation distance (length units) :l zero or more keyword/value pairs may be appended to args :l -keyword = {mol}, {region}, {maxangle}, {pressure}, {fugacity_coeff}, {full_energy}, {charge}, {group}, {grouptype}, or {intra_energy} +keyword = {mol}, {region}, {maxangle}, {pressure}, {fugacity_coeff}, {full_energy}, {charge}, {group}, {grouptype}, {intra_energy}, or {tfac_insert} {mol} value = template-ID template-ID = ID of molecule template specified in a separate "molecule"_molecule.html command {shake} value = fix-ID @@ -40,7 +40,8 @@ keyword = {mol}, {region}, {maxangle}, {pressure}, {fugacity_coeff}, {full_energ {grouptype} values = type group-ID type = atom type (int) group-ID = group-ID for inserted atoms (string) - {intra_energy} value = intramolecular energy (energy units) :pre + {intra_energy} value = intramolecular energy (energy units) + {tfac_insert} value = scale up/down temperature of inserted atoms (unitless) :pre :ule [Examples:] @@ -131,7 +132,9 @@ distance. See the "neighbor"_neighbor.html command for details. When an atom or molecule is to be inserted, its coordinates are chosen at a random position within the current simulation cell or region, and new atom velocities are randomly chosen from -the specified temperature distribution given by T. Relative +the specified temperature distribution given by T. The effective +temperature for new atom velocities can be increased or decreased +using the optional keyword {tfac_insert} (see below). Relative coordinates for atoms in a molecule are taken from the template molecule provided by the user. The center of mass of the molecule is placed at the insertion point. The orientation of the molecule @@ -217,6 +220,19 @@ deleted. For molecules that have a non-zero intramolecular energy, this will ensure roughly the same behavior whether or not the {full_energy} option is used. +Inserted atoms and molecules are assigned random velocities based on the +specified temperature T. Because the relative velocity of +all atoms in the molecule is zero, this may result in inserted molecules +that are systematically too cold. In addition, the intramolecular potential +energy of the inserted molecule may casue the kinetic energy +of the molecule to quickly increase or decrease after insertion. +The {tfac_insert} keyword allows the user to counteract these effects +by changing the temperature used to assign velocities to +inserted atoms and molecules by a constant factor. For a +particular application, some experimentation may be required +to find a value of {tfac_insert} that results in inserted molecules that +equilibrate quickly to the correct temperature. + Some fixes have an associated potential energy. Examples of such fixes include: "efield"_fix_efield.html, "gravity"_fix_gravity.html, "addforce"_fix_addforce.html, "langevin"_fix_langevin.html, From 537e951c6d8d4f339832a60a291cefa3b4f5102e Mon Sep 17 00:00:00 2001 From: athomps Date: Tue, 13 Oct 2015 15:54:51 +0000 Subject: [PATCH 2/4] added tfac_insert keyword and changed molecule insertion velocities PRNG git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@14110 f3b2605a-c512-4ea7-a41b-209d697bcdaa --- src/MC/fix_gcmc.cpp | 20 ++++++++++++-------- src/MC/fix_gcmc.h | 1 + 2 files changed, 13 insertions(+), 8 deletions(-) diff --git a/src/MC/fix_gcmc.cpp b/src/MC/fix_gcmc.cpp index a9a070c771..8482d453c9 100644 --- a/src/MC/fix_gcmc.cpp +++ b/src/MC/fix_gcmc.cpp @@ -240,6 +240,7 @@ void FixGCMC::options(int narg, char **arg) grouptypes = NULL; grouptypebits = NULL; energy_intra = 0.0; + tfac_insert = 1.0; int iarg = 0; while (iarg < narg) { @@ -330,6 +331,10 @@ void FixGCMC::options(int narg, char **arg) if (iarg+2 > narg) error->all(FLERR,"Illegal fix gcmc command"); energy_intra = force->numeric(FLERR,arg[iarg+1]); iarg += 2; + } else if (strcmp(arg[iarg],"tfac_insert") == 0) { + if (iarg+2 > narg) error->all(FLERR,"Illegal fix gcmc command"); + tfac_insert = force->numeric(FLERR,arg[iarg+1]); + iarg += 2; } else error->all(FLERR,"Illegal fix gcmc command"); } } @@ -530,7 +535,6 @@ void FixGCMC::init() onemols[imol]->compute_mass(); onemols[imol]->compute_com(); gas_mass = onemols[imol]->masstotal; - printf("gas_mass = %g\n",gas_mass); for (int i = 0; i < onemols[imol]->natoms; i++) { onemols[imol]->x[i][0] -= onemols[imol]->com[0]; onemols[imol]->x[i][1] -= onemols[imol]->com[1]; @@ -566,7 +570,7 @@ void FixGCMC::init() double lambda = sqrt(force->hplanck*force->hplanck/ (2.0*MY_PI*gas_mass*force->mvv2e* force->boltz*reservoir_temperature)); - sigma = sqrt(force->boltz*reservoir_temperature/gas_mass/force->mvv2e); + sigma = sqrt(force->boltz*reservoir_temperature*tfac_insert/gas_mass/force->mvv2e); zz = exp(beta*chemical_potential)/(pow(lambda,3.0)); if (pressure_flag) zz = pressure*fugacity_coeff*beta/force->nktv2p; @@ -1272,9 +1276,9 @@ void FixGCMC::attempt_molecule_insertion() int nlocalprev = atom->nlocal; double vnew[3]; - vnew[0] = random_unequal->gaussian()*sigma; - vnew[1] = random_unequal->gaussian()*sigma; - vnew[2] = random_unequal->gaussian()*sigma; + vnew[0] = random_equal->gaussian()*sigma; + vnew[1] = random_equal->gaussian()*sigma; + vnew[2] = random_equal->gaussian()*sigma; for (int i = 0; i < natoms_per_molecule; i++) { if (procflag[i]) { @@ -1897,9 +1901,9 @@ void FixGCMC::attempt_molecule_insertion_full() MathExtra::quat_to_mat(quat,rotmat); double vnew[3]; - vnew[0] = random_unequal->gaussian()*sigma; - vnew[1] = random_unequal->gaussian()*sigma; - vnew[2] = random_unequal->gaussian()*sigma; + vnew[0] = random_equal->gaussian()*sigma; + vnew[1] = random_equal->gaussian()*sigma; + vnew[2] = random_equal->gaussian()*sigma; for (int i = 0; i < natoms_per_molecule; i++) { double xtmp[3]; diff --git a/src/MC/fix_gcmc.h b/src/MC/fix_gcmc.h index 0c9a7a4878..437a6c7610 100644 --- a/src/MC/fix_gcmc.h +++ b/src/MC/fix_gcmc.h @@ -97,6 +97,7 @@ class FixGCMC : public Fix { int max_region_attempts; double gas_mass; double reservoir_temperature; + double tfac_insert; double chemical_potential; double displace; double max_rotation_angle; From d357a2f44a002fc45d6cecf4085ddacbbd5862a9 Mon Sep 17 00:00:00 2001 From: athomps Date: Tue, 13 Oct 2015 16:11:13 +0000 Subject: [PATCH 3/4] added tfac_insert keyword git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@14111 f3b2605a-c512-4ea7-a41b-209d697bcdaa --- doc/fix_gcmc.txt | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/doc/fix_gcmc.txt b/doc/fix_gcmc.txt index e4a0df29cc..e062be7bd1 100644 --- a/doc/fix_gcmc.txt +++ b/doc/fix_gcmc.txt @@ -120,7 +120,14 @@ until the atom or molecule center-of-mass is inside the specified region. If used with "fix_nvt"_fix_nvt.html, the temperature of the imaginary reservoir, T, should be set to be equivalent to the target temperature used in "fix_nvt"_fix_nvt.html. Otherwise, the imaginary reservoir -will not be in thermal equilibrium with the simulation cell. +will not be in thermal equilibrium with the simulation cell. Also, +it is important that the temperature used by fix nvt be dynamic, +which can be achieved as follows: + +compute mdtemp mdatoms temp +compute_modify mdtemp dynamic yes +fix mdnvt mdatoms nvt temp 300.0 300.0 10.0 +fix_modify mdnvt temp mdtemp :pre Note that neighbor lists are re-built every timestep that this fix is invoked, so you should not set N to be too small. However, periodic From b437d74c75c042831acf151251e4a56819aa9cf8 Mon Sep 17 00:00:00 2001 From: sjplimp Date: Tue, 13 Oct 2015 20:07:00 +0000 Subject: [PATCH 4/4] git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@14112 f3b2605a-c512-4ea7-a41b-209d697bcdaa --- src/Make.py | 124 ++++++++++++++++++++++++++++++++++++---------------- 1 file changed, 87 insertions(+), 37 deletions(-) diff --git a/src/Make.py b/src/Make.py index 774afe2a82..8fd4e7601a 100755 --- a/src/Make.py +++ b/src/Make.py @@ -5,7 +5,7 @@ # Syntax: Make.py -h (for help) # Notes: needs python 2.7 (not Python 3) -import sys,os,commands,re,copy +import sys,os,commands,re,copy,subprocess # switch abbrevs # switch classes = created class for each switch @@ -332,11 +332,11 @@ class Actions: make.addvar("CCFLAGS","-restrict") if final["user-omp"]: + if compile_check(compiler,"-fopenmp",1): + make.addvar("CCFLAGS","-fopenmp") + make.addvar("LINKFLAGS","-fopenmp") if compile_check(compiler,"-restrict",0): make.addvar("CCFLAGS","-restrict") - if compile_check(compiler,"-fopenmp",1): - make.addvar("CCFLAGS","-fopenmp") - make.addvar("LINKFLAGS","-fopenmp") if final["user-intel"]: if intel.mode == "cpu": @@ -496,7 +496,7 @@ class Actions: def clean(self): str = "cd %s; make clean-auto" % dir.src commands.getoutput(str) - if verbose: print "Performed make clean-auto" + print "Performed make clean-auto" # build LAMMPS using Makefile.auto and -j setting # invoke self.file() first, to test makefile compile/link @@ -516,10 +516,15 @@ class Actions: print "Created src/STUBS/libmpi_stubs.a" if jmake: str = "cd %s; make -j %d auto" % (dir.src,jmake.n) else: str = "cd %s; make auto" % dir.src - txt = commands.getoutput(str) - if verbose: print txt + + # if verbose, print output as build proceeds, else only print if fails + + if verbose: subprocess.call(str,shell=True) + else: + try: subprocess.check_output(str,stderr=subprocess.STDOUT,shell=True) + except Exception as e: print e.output + if not os.path.isfile("%s/lmp_auto" % dir.src): - if not verbose: print txt error('Unsuccessful "make auto"') elif not output: print "Created src/lmp_auto" @@ -968,11 +973,16 @@ class ATC: commands.getoutput("cd %s; make -f Makefile.auto clean" % libdir) if jmake: str = "cd %s; make -j %d -f Makefile.auto" % (libdir,jmake.n) else: str = "cd %s; make -f Makefile.auto" % libdir - txt = commands.getoutput(str) - if verbose: print txt + + # if verbose, print output as build proceeds, else only print if fails + + if verbose: subprocess.call(str,shell=True) + else: + try: subprocess.check_output(str,stderr=subprocess.STDOUT,shell=True) + except Exception as e: print e.output + if not os.path.isfile("%s/libatc.a" % libdir) or \ not os.path.isfile("%s/Makefile.lammps" % libdir): - if not verbose: print txt error("Unsuccessful build of lib/atc library") else: print "Created lib/atc library" @@ -1015,11 +1025,16 @@ class AWPMD: commands.getoutput("cd %s; make -f Makefile.auto clean" % libdir) if jmake: str = "cd %s; make -j %d -f Makefile.auto" % (libdir,jmake.n) else: str = "cd %s; make -f Makefile.auto" % libdir - txt = commands.getoutput(str) - if verbose: print txt + + # if verbose, print output as build proceeds, else only print if fails + + if verbose: subprocess.call(str,shell=True) + else: + try: subprocess.check_output(str,stderr=subprocess.STDOUT,shell=True) + except Exception as e: print e.output + if not os.path.isfile("%s/libawpmd.a" % libdir) or \ not os.path.isfile("%s/Makefile.lammps" % libdir): - if not verbose: print txt error("Unsuccessful build of lib/awpmd library") else: print "Created lib/awpmd library" @@ -1062,11 +1077,16 @@ class COLVARS: commands.getoutput("cd %s; make -f Makefile.auto clean" % libdir) if jmake: str = "cd %s; make -j %d -f Makefile.auto" % (libdir,jmake.n) else: str = "cd %s; make -f Makefile.auto" % libdir - txt = commands.getoutput(str) - if verbose: print txt + + # if verbose, print output as build proceeds, else only print if fails + + if verbose: subprocess.call(str,shell=True) + else: + try: subprocess.check_output(str,stderr=subprocess.STDOUT,shell=True) + except Exception as e: print e.output + if not os.path.isfile("%s/libcolvars.a" % libdir) or \ not os.path.isfile("%s/Makefile.lammps" % libdir): - if not verbose: print txt error("Unsuccessful build of lib/colvars library") else: print "Created lib/colvars library" @@ -1115,11 +1135,16 @@ class CUDA: (libdir,jmake.n,n,self.arch) else: str = str = "cd %s; make precision=%d arch=%s" % \ (libdir,n,self.arch) - txt = commands.getoutput(str) - if verbose: print txt + + # if verbose, print output as build proceeds, else only print if fails + + if verbose: subprocess.call(str,shell=True) + else: + try: subprocess.check_output(str,stderr=subprocess.STDOUT,shell=True) + except Exception as e: print e.output + if not os.path.isfile("%s/liblammpscuda.a" % libdir) or \ not os.path.isfile("%s/Makefile.lammps" % libdir): - if not verbose: print txt error("Unsuccessful build of lib/cuda library") else: print "Created lib/cuda library" @@ -1185,11 +1210,16 @@ class GPU: commands.getoutput("cd %s; make -f Makefile.auto clean" % libdir) if jmake: str = "cd %s; make -j %d -f Makefile.auto" % (libdir,jmake.n) else: str = "cd %s; make -f Makefile.auto" % libdir - txt = commands.getoutput(str) - if verbose: print txt + + # if verbose, print output as build proceeds, else only print if fails + + if verbose: subprocess.call(str,shell=True) + else: + try: subprocess.check_output(str,stderr=subprocess.STDOUT,shell=True) + except Exception as e: print e.output + if not os.path.isfile("%s/libgpu.a" % libdir) or \ not os.path.isfile("%s/Makefile.lammps" % libdir): - if not verbose: print txt error("Unsuccessful build of lib/gpu library") else: print "Created lib/gpu library" @@ -1232,11 +1262,16 @@ class MEAM: commands.getoutput("cd %s; make -f Makefile.auto clean" % libdir) # do not use -j for MEAM build, parallel build does not work str = "cd %s; make -f Makefile.auto" % libdir - txt = commands.getoutput(str) - if verbose: print txt + + # if verbose, print output as build proceeds, else only print if fails + + if verbose: subprocess.call(str,shell=True) + else: + try: subprocess.check_output(str,stderr=subprocess.STDOUT,shell=True) + except Exception as e: print e.output + if not os.path.isfile("%s/libmeam.a" % libdir) or \ not os.path.isfile("%s/Makefile.lammps" % libdir): - if not verbose: print txt error("Unsuccessful build of lib/meam library") else: print "Created lib/meam library" @@ -1279,11 +1314,16 @@ class POEMS: commands.getoutput("cd %s; make -f Makefile.auto clean" % libdir) if jmake: str = "cd %s; make -j %d -f Makefile.auto" % (libdir,jmake.n) else: str = "cd %s; make -f Makefile.auto" % libdir - txt = commands.getoutput(str) - if verbose: print txt + + # if verbose, print output as build proceeds, else only print if fails + + if verbose: subprocess.call(str,shell=True) + else: + try: subprocess.check_output(str,stderr=subprocess.STDOUT,shell=True) + except Exception as e: print e.output + if not os.path.isfile("%s/libpoems.a" % libdir) or \ not os.path.isfile("%s/Makefile.lammps" % libdir): - if not verbose: print txt error("Unsuccessful build of lib/poems library") else: print "Created lib/poems library" @@ -1326,11 +1366,16 @@ class QMMM: commands.getoutput("cd %s; make -f Makefile.auto clean" % libdir) if jmake: str = "cd %s; make -j %d -f Makefile.auto" % (libdir,jmake.n) else: str = "cd %s; make -f Makefile.auto" % libdir - txt = commands.getoutput(str) - if verbose: print txt + + # if verbose, print output as build proceeds, else only print if fails + + if verbose: subprocess.call(str,shell=True) + else: + try: subprocess.check_output(str,stderr=subprocess.STDOUT,shell=True) + except Exception as e: print e.output + if not os.path.isfile("%s/libqmmm.a" % libdir) or \ not os.path.isfile("%s/Makefile.lammps" % libdir): - if not verbose: print txt error("Unsuccessful build of lib/qmmm library") else: print "Created lib/qmmm library" @@ -1373,11 +1418,16 @@ class REAX: commands.getoutput("cd %s; make -f Makefile.auto clean" % libdir) if jmake: str = "cd %s; make -j %d -f Makefile.auto" % (libdir,jmake.n) else: str = "cd %s; make -f Makefile.auto" % libdir - txt = commands.getoutput(str) - if verbose: print txt + + # if verbose, print output as build proceeds, else only print if fails + + if verbose: subprocess.call(str,shell=True) + else: + try: subprocess.check_output(str,stderr=subprocess.STDOUT,shell=True) + except Exception as e: print e.output + if not os.path.isfile("%s/libreax.a" % libdir) or \ not os.path.isfile("%s/Makefile.lammps" % libdir): - if not verbose: print txt error("Unsuccessful build of lib/reax library") else: print "Created lib/reax library" @@ -1573,10 +1623,10 @@ class Fft: change FFT settings in makefile mode is required, all other args are optional removes all current FFT variable settings - mode = none or fftw or fftw3 of ... + mode = none or fftw or fftw3 or ... adds -DFFT_MODE setting lib = name of FFT library to link with (def is libname = mode) - adds -lliblib setting, e.g. -llibfftw3 + adds -llib{libname} setting, e.g. -llibfftw3 dir = home dir for include and library files (def = none) adds -Idir/include and -Ldir/lib settings if set, overrides idir and ldir args