From 491d5f341083e176a64013f19b98b8b0d516d634 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Thu, 21 Sep 2017 11:38:59 +0200 Subject: [PATCH 01/44] Added USER-PINNING package --- src/USER-PINNING/README | 13 ++ src/USER-PINNING/fix_rhoKUmbrella.cpp | 246 ++++++++++++++++++++++++++ src/USER-PINNING/fix_rhoKUmbrella.h | 81 +++++++++ 3 files changed, 340 insertions(+) create mode 100644 src/USER-PINNING/README create mode 100644 src/USER-PINNING/fix_rhoKUmbrella.cpp create mode 100644 src/USER-PINNING/fix_rhoKUmbrella.h diff --git a/src/USER-PINNING/README b/src/USER-PINNING/README new file mode 100644 index 0000000000..9f158b1884 --- /dev/null +++ b/src/USER-PINNING/README @@ -0,0 +1,13 @@ +This bias potential is used to study solid-liquid transitions with the interface pinning method. + +Reference: + +[Ulf R. Pedersen, J. Chem. Phys. 139, 104102 (2013)] + +Usage snip: + + fix [fix-name] [groupID] rhoKUmbrella [nx] [ny] [nz] [kappa] [anchor-point] + thermo_style custom step temp pzz pe lz f_umbrella f_umbrella[1] f_umbrella[2] f_umbrella[3] + +where the parameters set the harmonic bias potential U=0.5*kappa*(|rho_k|-anchor-point)^2 +with the wave-vector elements of rho_k to k_x = (2 pi / L_x) * n_x, k_y = (2 pi / L_y) * n_y and k_z = (2 pi / L_z) * n_z. diff --git a/src/USER-PINNING/fix_rhoKUmbrella.cpp b/src/USER-PINNING/fix_rhoKUmbrella.cpp new file mode 100644 index 0000000000..e62eecb31c --- /dev/null +++ b/src/USER-PINNING/fix_rhoKUmbrella.cpp @@ -0,0 +1,246 @@ +/* + fix_rhoK_umbrella.cpp + + A fix to do umbrella sampling on rho(k). + + The usage is as follows: + + fix [name] [groupID] rhoKUmbrella [nx] [ny] [nz] [kappa = spring constant] [rhoK0] + + where k_i = (2 pi / L_i) * n_i + + Written by Ulf Pedersen and Patrick Varilly, 4 Feb 2010 + Tweaked for LAMMPS 15 Jan 2010 version by Ulf Pedersen, 19 Aug 2010 + Tweaked again March 4th 2012 by Ulf Pedersen. + */ + +#include "fix_rhoKUmbrella.h" +#include "error.h" +#include "update.h" +#include "respa.h" +#include "atom.h" +#include "domain.h" + +#include +#include +#include +#include +#include + +using namespace LAMMPS_NS; +using namespace FixConst; + +// Constructor: all the parameters to this fix specified in +// the LAMMPS input get passed in +FixRhoKUmbrella::FixRhoKUmbrella( LAMMPS* inLMP, int inArgc, char** inArgv ) + : Fix( inLMP, inArgc, inArgv ) +{ + // Check arguments + if( inArgc != 8 ) + error->all(FLERR,"Illegal fix rhoKUmbrella command" ); + + // Set up fix flags + scalar_flag = 1; // have compute_scalar + vector_flag = 1; // have compute_vector... + size_vector = 3; // ...with this many components + global_freq = 1; // whose value can be computed at every timestep + //scalar_vector_freq = 1; // OLD lammps: whose value can be computed at every timestep + thermo_energy = 1; // this fix changes system's potential energy + extscalar = 0; // but the deltaPE might not scale with # of atoms + extvector = 0; // neither do the components of the vector + + // Parse fix options + int n[3]; + + n[0] = atoi( inArgv[3] ); + n[1] = atoi( inArgv[4] ); + n[2] = atoi( inArgv[5] ); + + mK[0] = n[0]*(2*M_PI / (domain->boxhi[0] - domain->boxlo[0])); + mK[1] = n[1]*(2*M_PI / (domain->boxhi[1] - domain->boxlo[1])); + mK[2] = n[2]*(2*M_PI / (domain->boxhi[2] - domain->boxlo[2])); + + mKappa = atof( inArgv[6] ); + mRhoK0 = atof( inArgv[7] ); +} + +FixRhoKUmbrella::~FixRhoKUmbrella() +{ +} + +// Methods that this fix implements +// -------------------------------- + +// Tells LAMMPS where this fix should act +int +FixRhoKUmbrella::setmask() +{ + int mask = 0; + + // This fix modifies forces... + mask |= POST_FORCE; + mask |= POST_FORCE_RESPA; + mask |= MIN_POST_FORCE; + + // ...and potential energies + mask |= THERMO_ENERGY; + + return mask; +} + +/*int FixRhoKUmbrella::setmask() +{ + int mask = 0; + mask |= POST_FORCE; + mask |= POST_FORCE_RESPA; + mask |= MIN_POST_FORCE; + return mask; +}*/ + + +// Initializes the fix at the beginning of a run +void +FixRhoKUmbrella::init() +{ + // RESPA boilerplate + if( strcmp( update->integrate_style, "respa" ) == 0 ) + mNLevelsRESPA = ((Respa *) update->integrate)->nlevels; + + // Count the number of affected particles + int nThisLocal = 0; + int *mask = atom->mask; + int nlocal = atom->nlocal; + for( int i = 0; i < nlocal; i++ ) { // Iterate through all atoms on this CPU + if( mask[i] & groupbit ) { // ...only those affected by this fix + nThisLocal++; + } + } + MPI_Allreduce( &nThisLocal, &mNThis, + 1, MPI_INT, MPI_SUM, world ); + mSqrtNThis = sqrt( mNThis ); +} + +// Initial application of the fix to a system (when doing MD) +void +FixRhoKUmbrella::setup( int inVFlag ) +{ + if( strcmp( update->integrate_style, "verlet" ) == 0 ) + post_force( inVFlag ); + else + { + ((Respa *) update->integrate)->copy_flevel_f( mNLevelsRESPA - 1 ); + post_force_respa( inVFlag, mNLevelsRESPA - 1,0 ); + ((Respa *) update->integrate)->copy_f_flevel( mNLevelsRESPA - 1 ); + } +} + +// Initial application of the fix to a system (when doing minimization) +void +FixRhoKUmbrella::min_setup( int inVFlag ) +{ + post_force( inVFlag ); +} + +// Modify the forces calculated in the main force loop of ordinary MD +void +FixRhoKUmbrella::post_force( int inVFlag ) +{ + double **x = atom->x; + double **f = atom->f; + int *mask = atom->mask; + int nlocal = atom->nlocal; + + // Loop over locally-owned atoms affected by this fix and calculate the + // partial rhoK's + mRhoKLocal[0] = 0.0; + mRhoKLocal[1] = 0.0; + + for( int i = 0; i < nlocal; i++ ) { // Iterate through all atoms on this CPU + if( mask[i] & groupbit ) { // ...only those affected by this fix + + // rho_k = sum_i exp( - i k.r_i ) + mRhoKLocal[0] += cos( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); + mRhoKLocal[1] -= sin( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); + } + } + + // Now calculate mRhoKGlobal + MPI_Allreduce( mRhoKLocal, mRhoKGlobal, + 2, MPI_DOUBLE, MPI_SUM, world ); + + // WARNING!!!!! < \sum_{i,j} e^{-ik.(r_i - r_j)} > ~ N, so + // we define rho_k as (1 / sqrt(N)) \sum_i e^{-i k.r_i}, so that + // is intensive. + // + // Don't forget this two years from now when you change the system size!!! + mRhoKGlobal[0] /= mSqrtNThis; + mRhoKGlobal[1] /= mSqrtNThis; + + // We'll need magnitude of rho_k + double rhoK = sqrt( mRhoKGlobal[0]*mRhoKGlobal[0] + + mRhoKGlobal[1]*mRhoKGlobal[1] ); + + for( int i = 0; i < nlocal; i++ ) { // Iterate through all atoms on this CPU + if( mask[i] & groupbit ) { // ...only those affected by this fix + + // Calculate forces + // U = kappa/2 ( |rho_k| - rho_k^0 )^2 + // f_i = -grad_i U = -kappa ( |rho_k| - rho_k^0 ) grad_i |rho_k| + // grad_i |rho_k| = Re( rho_k* (-i k e^{-i k . r_i} / sqrt(N)) ) / |rho_k| + // + // In terms of real and imag parts of rho_k, + // + // Re( rho_k* (-i k e^{-i k . r_i}) ) = + // (- Re[rho_k] * sin( k . r_i ) - Im[rho_k] * cos( k . r_i )) * k + + double sinKRi = sin( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); + double cosKRi = cos( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); + + double prefactor = mKappa * ( rhoK - mRhoK0 ) / rhoK + * (-mRhoKGlobal[0]*sinKRi - mRhoKGlobal[1]*cosKRi) / mSqrtNThis; + f[i][0] -= prefactor * mK[0]; + f[i][1] -= prefactor * mK[1]; + f[i][2] -= prefactor * mK[2]; + } + } +} + +// Forces in RESPA loop +void +FixRhoKUmbrella::post_force_respa( int inVFlag, int inILevel, int inILoop ) +{ + if( inILevel == mNLevelsRESPA - 1 ) + post_force( inVFlag ); +} + +// Forces in minimization loop +void +FixRhoKUmbrella::min_post_force( int inVFlag ) +{ + post_force( inVFlag ); +} + +// Compute the change in the potential energy induced by this fix +double +FixRhoKUmbrella::compute_scalar() +{ + double rhoK = sqrt( mRhoKGlobal[0]*mRhoKGlobal[0] + + mRhoKGlobal[1]*mRhoKGlobal[1] ); + + return 0.5 * mKappa * (rhoK - mRhoK0) * (rhoK - mRhoK0); +} + +// Compute the ith component of the vector +double +FixRhoKUmbrella::compute_vector( int inI ) +{ + if( inI == 0 ) + return mRhoKGlobal[0]; // Real part + else if( inI == 1 ) + return mRhoKGlobal[1]; // Imagniary part + else if( inI == 2 ) + return sqrt( mRhoKGlobal[0]*mRhoKGlobal[0] + + mRhoKGlobal[1]*mRhoKGlobal[1] ); + else + return 12345.0; +} diff --git a/src/USER-PINNING/fix_rhoKUmbrella.h b/src/USER-PINNING/fix_rhoKUmbrella.h new file mode 100644 index 0000000000..d549d31d94 --- /dev/null +++ b/src/USER-PINNING/fix_rhoKUmbrella.h @@ -0,0 +1,81 @@ +/* + fix_rhoK_umbrella.h + + A fix to do umbrella sampling on rho(k). + + The usage is as follows: + + fix [name] [groupID] rhoKUmbrella [kx] [ky] [kz] [kappa = spring constant] [rhoK0] + + where k_i = (2 pi / L_i) * n_i + + Written by Ulf Pedersen and Patrick Varilly, 4 Feb 2010 + Tweaked for LAMMPS 15 Jan 2010 version by Ulf Pedersen, 19 Aug 2010 +*/ + +#ifdef FIX_CLASS + +FixStyle(rhoKUmbrella,FixRhoKUmbrella) + +#else + +#ifndef __FIX_RHOKUMBRELLA__ +#define __FIX_RHOKUMBRELLA__ + +#include "fix.h" + +namespace LAMMPS_NS { + +class FixRhoKUmbrella : public Fix +{ +public: + // Constructor: all the parameters to this fix specified in + // the LAMMPS input get passed in + FixRhoKUmbrella( LAMMPS* inLMP, int inArgc, char** inArgv ); + virtual ~FixRhoKUmbrella(); + + // Methods that this fix implements + // -------------------------------- + + // Tells LAMMPS where this fix should act + int setmask(); + + // Initializes the fix at the beginning of a run + void init(); + + // Initial application of the fix to a system (when doing MD / minimization) + void setup( int inVFlag ); + void min_setup( int inVFlag ); + + // Modify the forces calculated in the main force loop, either when + // doing usual MD, RESPA MD or minimization + void post_force( int inVFlag ); + void post_force_respa( int inVFlag, int inILevel, int inILoop ); + void min_post_force( int inVFlag ); + + // Compute the change in the potential energy induced by this fix + double compute_scalar(); + + // Compute the ith component of the vector associated with this fix + double compute_vector( int inI ); + +private: + // RESPA boilerplate + int mNLevelsRESPA; + + // Defining parameters for this umbrella + double mK[3], mKappa, mRhoK0; + + // Number of particles affected by the fix + int mNThis; + double mSqrtNThis; + + // Real and imaginary parts of rho_k := sum_i exp( - i k . r_i ) + double mRhoKLocal[2], mRhoKGlobal[2]; +}; + +} // namespace LAMMPS_NS + +#endif // __FIX_RHOKUMBRELLA__ +#endif // FIX_CLASS + From 9a9af2ca5e673356a1a4cf66036f36e32663cdca Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Thu, 21 Sep 2017 13:58:51 +0200 Subject: [PATCH 02/44] Delete README --- src/USER-PINNING/README | 13 ------------- 1 file changed, 13 deletions(-) delete mode 100644 src/USER-PINNING/README diff --git a/src/USER-PINNING/README b/src/USER-PINNING/README deleted file mode 100644 index 9f158b1884..0000000000 --- a/src/USER-PINNING/README +++ /dev/null @@ -1,13 +0,0 @@ -This bias potential is used to study solid-liquid transitions with the interface pinning method. - -Reference: - -[Ulf R. Pedersen, J. Chem. Phys. 139, 104102 (2013)] - -Usage snip: - - fix [fix-name] [groupID] rhoKUmbrella [nx] [ny] [nz] [kappa] [anchor-point] - thermo_style custom step temp pzz pe lz f_umbrella f_umbrella[1] f_umbrella[2] f_umbrella[3] - -where the parameters set the harmonic bias potential U=0.5*kappa*(|rho_k|-anchor-point)^2 -with the wave-vector elements of rho_k to k_x = (2 pi / L_x) * n_x, k_y = (2 pi / L_y) * n_y and k_z = (2 pi / L_z) * n_z. From 73708b091c429d1f1a758a42a55a2dd890941a04 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Thu, 21 Sep 2017 16:17:26 +0200 Subject: [PATCH 03/44] Added readme file, and removed old files --- examples/USER/pinning/readme.md | 18 ++++++++++++++++++ 1 file changed, 18 insertions(+) create mode 100644 examples/USER/pinning/readme.md diff --git a/examples/USER/pinning/readme.md b/examples/USER/pinning/readme.md new file mode 100644 index 0000000000..e5563a3056 --- /dev/null +++ b/examples/USER/pinning/readme.md @@ -0,0 +1,18 @@ +This bias potential is used to study solid-liquid transitions with the interface pinning method. + +Reference: + +[Ulf R. Pedersen, J. Chem. Phys. 139, 104102 (2013)] + +Usage snip: + + fix [fix-name] [groupID] rhoKUmbrella [nx] [ny] [nz] [kappa] [anchor-point] + thermo_style custom step temp pzz pe lz f_umbrella f_umbrella[1] f_umbrella[2] f_umbrella[3] + +where the parameters set the harmonic bias potential U=0.5*kappa*(|rho_k|-anchor-point)^2 +with the wave-vector elements of rho_k to k_x = (2 pi / L_x) * n_x, k_y = (2 pi / L_y) * n_y and k_z = (2 pi / L_z) * n_z. + +This package was created by + Ulf R. Pedersen + http://www.urp.dk + ulf AT urp.dk From 3381a4337887c825c98ff0e51776bf5f560add11 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Thu, 21 Sep 2017 16:20:06 +0200 Subject: [PATCH 04/44] Added readme.m --- examples/USER/pinning/readme.md | 52 ++++++++++++++++++++++++++++----- 1 file changed, 45 insertions(+), 7 deletions(-) diff --git a/examples/USER/pinning/readme.md b/examples/USER/pinning/readme.md index e5563a3056..163536ca94 100644 --- a/examples/USER/pinning/readme.md +++ b/examples/USER/pinning/readme.md @@ -1,18 +1,56 @@ -This bias potential is used to study solid-liquid transitions with the interface pinning method. +This package contains a bias potential that is used to study solid-liquid transitions with the interface pinning method. +An interface between a solid and a liquid is simulated by applying a field that bias the system towards two-phase configurations. +This is done by adding a harmonic potential to the Hamiltonian. The bias field couple to an order-parameter of crystallinity Q: -Reference: + U_bias = 0.5*k*(Q-a)^2 -[Ulf R. Pedersen, J. Chem. Phys. 139, 104102 (2013)] +Here, We user long-range order for "crystallinity". Q=rho_k wher rho_k is the collective density field. -Usage snip: +# References +The main reference for the method is + [Ulf R. Pedersen, J. Chem. Phys. 139, 104102 (2013)] - fix [fix-name] [groupID] rhoKUmbrella [nx] [ny] [nz] [kappa] [anchor-point] - thermo_style custom step temp pzz pe lz f_umbrella f_umbrella[1] f_umbrella[2] f_umbrella[3] +Please visit + urp.dk/interface_pinning.htm +for a detailed bibliography. + +# Build +Remember to include the following command when building LAMMPS + make yes-user-pinning + +# Use + + fix [name] [groupID] rhok [nx] [ny] [nz] [kappa] [anchor-point] where the parameters set the harmonic bias potential U=0.5*kappa*(|rho_k|-anchor-point)^2 with the wave-vector elements of rho_k to k_x = (2 pi / L_x) * n_x, k_y = (2 pi / L_y) * n_y and k_z = (2 pi / L_z) * n_z. -This package was created by +# Usage example +In the following we will apply use the interface pinning method for the Lennard-Jones system (trunctaed at 2.5) +at temperature 0.8 and pressure 2.185. This happens to be a coexistence state-point, but we will later show how interface pinning +can be used to determine this. The present directory contains input files, that we will use. + +## Density of crystal +First we will determine the density of the crystal with the following LAMMPS input file +{crystal.lmp} +from the output we get that the average density is 0.9731. We need this density to ensure hydrostatic pressure +when in the crystal slap of a two-phase simulation. + +## Setup two-phase configuration +Next, setup a two-phase configuration using the density determined in the previous step. +{setup.lmp} + + +## Setup two-phase configuration +Finally, we run simulation with the bias field applied. +{pinning.lmp} + +# Contact Ulf R. Pedersen http://www.urp.dk ulf AT urp.dk + +# Cite +Please cite + [Ulf R. Pedersen, J. Chem. Phys. 139, 104102 (2013)] +when using the package for a publication. From fe8244c1c2320b0af83ce739b6a4b55ebbe39237 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Thu, 21 Sep 2017 16:49:46 +0200 Subject: [PATCH 05/44] Rename fix --- .../{fix_rhoKUmbrella.cpp => fix_rhok.cpp} | 35 +++++++++---------- .../{fix_rhoKUmbrella.h => fix_rhok.h} | 16 ++++----- 2 files changed, 25 insertions(+), 26 deletions(-) rename src/USER-PINNING/{fix_rhoKUmbrella.cpp => fix_rhok.cpp} (87%) rename src/USER-PINNING/{fix_rhoKUmbrella.h => fix_rhok.h} (85%) diff --git a/src/USER-PINNING/fix_rhoKUmbrella.cpp b/src/USER-PINNING/fix_rhok.cpp similarity index 87% rename from src/USER-PINNING/fix_rhoKUmbrella.cpp rename to src/USER-PINNING/fix_rhok.cpp index e62eecb31c..5cbc85a0d6 100644 --- a/src/USER-PINNING/fix_rhoKUmbrella.cpp +++ b/src/USER-PINNING/fix_rhok.cpp @@ -5,16 +5,17 @@ The usage is as follows: - fix [name] [groupID] rhoKUmbrella [nx] [ny] [nz] [kappa = spring constant] [rhoK0] + fix [name] [groupID] rhoK [nx] [ny] [nz] [kappa = spring constant] [rhoK0] where k_i = (2 pi / L_i) * n_i Written by Ulf Pedersen and Patrick Varilly, 4 Feb 2010 Tweaked for LAMMPS 15 Jan 2010 version by Ulf Pedersen, 19 Aug 2010 - Tweaked again March 4th 2012 by Ulf Pedersen. + Tweaked again March 4th 2012 by Ulf R. Pedersen, + September 2016 by Ulf R. Pedersen */ -#include "fix_rhoKUmbrella.h" +#include "fix_rhok.h" #include "error.h" #include "update.h" #include "respa.h" @@ -32,7 +33,7 @@ using namespace FixConst; // Constructor: all the parameters to this fix specified in // the LAMMPS input get passed in -FixRhoKUmbrella::FixRhoKUmbrella( LAMMPS* inLMP, int inArgc, char** inArgv ) +FixRhok::FixRhok( LAMMPS* inLMP, int inArgc, char** inArgv ) : Fix( inLMP, inArgc, inArgv ) { // Check arguments @@ -64,7 +65,7 @@ FixRhoKUmbrella::FixRhoKUmbrella( LAMMPS* inLMP, int inArgc, char** inArgv ) mRhoK0 = atof( inArgv[7] ); } -FixRhoKUmbrella::~FixRhoKUmbrella() +FixRhok::~FixRhok() { } @@ -73,7 +74,7 @@ FixRhoKUmbrella::~FixRhoKUmbrella() // Tells LAMMPS where this fix should act int -FixRhoKUmbrella::setmask() +FixRhok::setmask() { int mask = 0; @@ -88,7 +89,7 @@ FixRhoKUmbrella::setmask() return mask; } -/*int FixRhoKUmbrella::setmask() +/*int FixRhok::setmask() { int mask = 0; mask |= POST_FORCE; @@ -100,7 +101,7 @@ FixRhoKUmbrella::setmask() // Initializes the fix at the beginning of a run void -FixRhoKUmbrella::init() +FixRhok::init() { // RESPA boilerplate if( strcmp( update->integrate_style, "respa" ) == 0 ) @@ -122,7 +123,7 @@ FixRhoKUmbrella::init() // Initial application of the fix to a system (when doing MD) void -FixRhoKUmbrella::setup( int inVFlag ) +FixRhok::setup( int inVFlag ) { if( strcmp( update->integrate_style, "verlet" ) == 0 ) post_force( inVFlag ); @@ -136,14 +137,14 @@ FixRhoKUmbrella::setup( int inVFlag ) // Initial application of the fix to a system (when doing minimization) void -FixRhoKUmbrella::min_setup( int inVFlag ) +FixRhok::min_setup( int inVFlag ) { post_force( inVFlag ); } // Modify the forces calculated in the main force loop of ordinary MD void -FixRhoKUmbrella::post_force( int inVFlag ) +FixRhok::post_force( int inVFlag ) { double **x = atom->x; double **f = atom->f; @@ -168,11 +169,9 @@ FixRhoKUmbrella::post_force( int inVFlag ) MPI_Allreduce( mRhoKLocal, mRhoKGlobal, 2, MPI_DOUBLE, MPI_SUM, world ); - // WARNING!!!!! < \sum_{i,j} e^{-ik.(r_i - r_j)} > ~ N, so + // Info: < \sum_{i,j} e^{-ik.(r_i - r_j)} > ~ N, so // we define rho_k as (1 / sqrt(N)) \sum_i e^{-i k.r_i}, so that // is intensive. - // - // Don't forget this two years from now when you change the system size!!! mRhoKGlobal[0] /= mSqrtNThis; mRhoKGlobal[1] /= mSqrtNThis; @@ -207,7 +206,7 @@ FixRhoKUmbrella::post_force( int inVFlag ) // Forces in RESPA loop void -FixRhoKUmbrella::post_force_respa( int inVFlag, int inILevel, int inILoop ) +FixRhok::post_force_respa( int inVFlag, int inILevel, int inILoop ) { if( inILevel == mNLevelsRESPA - 1 ) post_force( inVFlag ); @@ -215,14 +214,14 @@ FixRhoKUmbrella::post_force_respa( int inVFlag, int inILevel, int inILoop ) // Forces in minimization loop void -FixRhoKUmbrella::min_post_force( int inVFlag ) +FixRhok::min_post_force( int inVFlag ) { post_force( inVFlag ); } // Compute the change in the potential energy induced by this fix double -FixRhoKUmbrella::compute_scalar() +FixRhok::compute_scalar() { double rhoK = sqrt( mRhoKGlobal[0]*mRhoKGlobal[0] + mRhoKGlobal[1]*mRhoKGlobal[1] ); @@ -232,7 +231,7 @@ FixRhoKUmbrella::compute_scalar() // Compute the ith component of the vector double -FixRhoKUmbrella::compute_vector( int inI ) +FixRhok::compute_vector( int inI ) { if( inI == 0 ) return mRhoKGlobal[0]; // Real part diff --git a/src/USER-PINNING/fix_rhoKUmbrella.h b/src/USER-PINNING/fix_rhok.h similarity index 85% rename from src/USER-PINNING/fix_rhoKUmbrella.h rename to src/USER-PINNING/fix_rhok.h index d549d31d94..3e0625430f 100644 --- a/src/USER-PINNING/fix_rhoKUmbrella.h +++ b/src/USER-PINNING/fix_rhok.h @@ -1,5 +1,5 @@ /* - fix_rhoK_umbrella.h + fix_rhok.h A fix to do umbrella sampling on rho(k). @@ -15,24 +15,24 @@ #ifdef FIX_CLASS -FixStyle(rhoKUmbrella,FixRhoKUmbrella) +FixStyle(rhok,FixRhok) #else -#ifndef __FIX_RHOKUMBRELLA__ -#define __FIX_RHOKUMBRELLA__ +#ifndef __FIX_RHOK__ +#define __FIX_RHOK__ #include "fix.h" namespace LAMMPS_NS { -class FixRhoKUmbrella : public Fix +class FixRhok : public Fix { public: // Constructor: all the parameters to this fix specified in // the LAMMPS input get passed in - FixRhoKUmbrella( LAMMPS* inLMP, int inArgc, char** inArgv ); - virtual ~FixRhoKUmbrella(); + FixRhok( LAMMPS* inLMP, int inArgc, char** inArgv ); + virtual ~FixRhok(); // Methods that this fix implements // -------------------------------- @@ -76,6 +76,6 @@ private: } // namespace LAMMPS_NS -#endif // __FIX_RHOKUMBRELLA__ +#endif // __FIX_RHOK__ #endif // FIX_CLASS From 4beccf508f9cd072e8ec3467186cc4e513baa30f Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Mon, 25 Sep 2017 18:35:53 +0200 Subject: [PATCH 06/44] Move fix to USER-MISH --- src/{USER-PINNING => USER-MISC}/fix_rhok.cpp | 0 src/{USER-PINNING => USER-MISC}/fix_rhok.h | 0 2 files changed, 0 insertions(+), 0 deletions(-) rename src/{USER-PINNING => USER-MISC}/fix_rhok.cpp (100%) rename src/{USER-PINNING => USER-MISC}/fix_rhok.h (100%) diff --git a/src/USER-PINNING/fix_rhok.cpp b/src/USER-MISC/fix_rhok.cpp similarity index 100% rename from src/USER-PINNING/fix_rhok.cpp rename to src/USER-MISC/fix_rhok.cpp diff --git a/src/USER-PINNING/fix_rhok.h b/src/USER-MISC/fix_rhok.h similarity index 100% rename from src/USER-PINNING/fix_rhok.h rename to src/USER-MISC/fix_rhok.h From b35f2ff8b44bdedca2cb6fa44392ae01a5876674 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Mon, 25 Sep 2017 18:44:24 +0200 Subject: [PATCH 07/44] Example of Interface Pinning Computation --- examples/USER/pinning/crystal.lmp | 35 +++++++++++++++++++ examples/USER/pinning/pinning.lmp | 32 ++++++++++++++++++ examples/USER/pinning/readme.md | 56 ++++++++++++++++--------------- examples/USER/pinning/setup.lmp | 40 ++++++++++++++++++++++ 4 files changed, 136 insertions(+), 27 deletions(-) create mode 100644 examples/USER/pinning/crystal.lmp create mode 100644 examples/USER/pinning/pinning.lmp create mode 100644 examples/USER/pinning/setup.lmp diff --git a/examples/USER/pinning/crystal.lmp b/examples/USER/pinning/crystal.lmp new file mode 100644 index 0000000000..c40e206d46 --- /dev/null +++ b/examples/USER/pinning/crystal.lmp @@ -0,0 +1,35 @@ +units lj +dimension 3 +boundary p p p +atom_style atomic + +# truncated and shifted LJ potential +pair_style lj/cut 2.5 +pair_modify shift yes +lattice fcc 0.9731 +region my_box block 0 8.0 0 8.0 0 20.0 +create_box 1 my_box +region particles block 0 8.0 0 8.0 0 20.0 +create_atoms 1 region particles +pair_coeff 1 1 1.0 1.0 2.5 +pair_modify tail no +pair_modify shift yes +mass 1 1.0 +velocity all create 1.6 1 mom yes rot yes + +# simulation parameters +neighbor 0.6 bin +timestep 0.004 +run_style verlet +fix ensemble all npt temp 0.8 0.8 4.0 aniso 2.185 2.185 8.0 pchain 32 + +# computing long-range order (no bias is added since k=0) +fix bias all rhok 16 0 0 0.0 0.0 + +# output +thermo 50 +thermo_style custom step temp press density f_bias[3] +dump dumpXYZ all xyz 2000 traj.xyz + +# run +run 100000 diff --git a/examples/USER/pinning/pinning.lmp b/examples/USER/pinning/pinning.lmp new file mode 100644 index 0000000000..62e707a1fd --- /dev/null +++ b/examples/USER/pinning/pinning.lmp @@ -0,0 +1,32 @@ +units lj +dimension 3 +boundary p p p +atom_style atomic + +# truncated and shifted LJ potential +pair_style lj/cut 2.5 +pair_modify shift yes +read_data data.halfhalf +pair_coeff 1 1 1.0 1.0 2.5 +mass 1 1.0 + +# simulation parameters +neighbor 0.6 bin +timestep 0.004 +run_style verlet + +velocity all create 0.8 1 mom yes rot yes +fix ensemble all npt temp 0.8 0.8 4.0 z 2.185 2.185 8.0 +fix 100 all momentum 100 linear 1 1 1 + +# harmonic rho_k bias-field +# nx ny nz k a +fix bias all rhok 16 0 0 4.0 26.00 + +# output U_bias rho_k_RE rho_k_IM |rho_k| +thermo_style custom step temp pzz pe lz f_bias f_bias[1] f_bias[2] f_bias[3] +thermo 50 +dump dumpXYZ all xyz 500 traj.xyz + +# run +run 50000 diff --git a/examples/USER/pinning/readme.md b/examples/USER/pinning/readme.md index 163536ca94..3d5c847cb8 100644 --- a/examples/USER/pinning/readme.md +++ b/examples/USER/pinning/readme.md @@ -1,14 +1,13 @@ -This package contains a bias potential that is used to study solid-liquid transitions with the interface pinning method. -An interface between a solid and a liquid is simulated by applying a field that bias the system towards two-phase configurations. -This is done by adding a harmonic potential to the Hamiltonian. The bias field couple to an order-parameter of crystallinity Q: +This package contains a bias potential that can be used to study solid-liquid transitions with the interface pinning method. +This is done by adding a harmonic potential to the Hamiltonian that bias the system towards two-phase configurations. U_bias = 0.5*k*(Q-a)^2 -Here, We user long-range order for "crystallinity". Q=rho_k wher rho_k is the collective density field. +The bias field couple to an order-parameter of crystallinity Q +This implimentation use long-range order: Q=|rho_k|, where rho_k is the collective density field of the wave-vector k. -# References -The main reference for the method is - [Ulf R. Pedersen, J. Chem. Phys. 139, 104102 (2013)] +# Reference +[Ulf R. Pedersen, J. Chem. Phys. 139, 104102 (2013)] Please visit urp.dk/interface_pinning.htm @@ -20,30 +19,34 @@ Remember to include the following command when building LAMMPS # Use - fix [name] [groupID] rhok [nx] [ny] [nz] [kappa] [anchor-point] + fix [name] [groupID] rhok [nx] [ny] [nz] [spring-constant] [anchor-point] -where the parameters set the harmonic bias potential U=0.5*kappa*(|rho_k|-anchor-point)^2 -with the wave-vector elements of rho_k to k_x = (2 pi / L_x) * n_x, k_y = (2 pi / L_y) * n_y and k_z = (2 pi / L_z) * n_z. +include a harmonic bias potential U_bias=0.5*k*(|rho_k|-a)^2 to the force calculation. +The elements of the wave-vector rho_k is k_x = (2 pi / L_x) * n_x, k_y = (2 pi / L_y) * n_y and k_z = (2 pi / L_z) * n_z. -# Usage example -In the following we will apply use the interface pinning method for the Lennard-Jones system (trunctaed at 2.5) -at temperature 0.8 and pressure 2.185. This happens to be a coexistence state-point, but we will later show how interface pinning -can be used to determine this. The present directory contains input files, that we will use. +# The Interface Pinning method for studying melting transitions +We will use the interface pinning method to study melting of the Lennard-Jones system +at temperature 0.8 and pressure 2.185. This is a coexistence state-point, and the method +can be used to show this. The present directory contains the input files: -## Density of crystal -First we will determine the density of the crystal with the following LAMMPS input file -{crystal.lmp} -from the output we get that the average density is 0.9731. We need this density to ensure hydrostatic pressure -when in the crystal slap of a two-phase simulation. + crystal.lmp + setup.lmp + pinning.lmp -## Setup two-phase configuration -Next, setup a two-phase configuration using the density determined in the previous step. -{setup.lmp} +1. First we will determine the density of the crystal with the LAMMPS input file crystal.lmp. + From the output we get that the average density after equbriliation is 0.9731. + We need this density to ensure hydrostatic pressure when in a two-phase simulation. +2. Next, we setup a two-phase configuration using setup.lmp. -## Setup two-phase configuration -Finally, we run simulation with the bias field applied. -{pinning.lmp} +3. Finally, we run a two-phase simulation with the bias-field applied using pinning.lmp. + The last coulmn in the output show |rho_k|. We note that after a equbriliation period + the value fluctuates aroung the anchor point (a) -- showing that this is indeed a coexistence + state point. + +The reference [J. Chem. Phys. 139, 104102 (2013)] gives details on using the method to find coexitence state points, +and the referecee [J. Chem. Phys. 142, 044104 (2015)] show how the crystal growth rate can be computed. +That method have been experienced to be most effective in the slightly super-heated regime above the melting temperature. # Contact Ulf R. Pedersen @@ -52,5 +55,4 @@ Finally, we run simulation with the bias field applied. # Cite Please cite - [Ulf R. Pedersen, J. Chem. Phys. 139, 104102 (2013)] -when using the package for a publication. + [Ulf R. Pedersen, J. Chem. Phys. 139, 104102 (2013)] diff --git a/examples/USER/pinning/setup.lmp b/examples/USER/pinning/setup.lmp new file mode 100644 index 0000000000..4ab9e4498b --- /dev/null +++ b/examples/USER/pinning/setup.lmp @@ -0,0 +1,40 @@ +units lj +dimension 3 +boundary p p p +atom_style atomic + +# truncated and shifted LJ potential +pair_style lj/cut 2.5 +pair_modify shift yes + +# fcc lattice +lattice fcc 0.9731 +region my_box block 0 8.0 0 8.0 0 20.0 +create_box 1 my_box +region particles block 0 8.0 0 8.0 0 20.0 +create_atoms 1 region particles +pair_coeff 1 1 1.0 1.0 2.5 +mass 1 1.0 +change_box all z final 0.0 34 remap units box + +# select particles in one side of the elongated box +region left plane 0 0 10 0 0 1 +group left region left + +velocity left create 6.0 1 mom yes rot yes + +# simulation parameters +neighbor 0.6 bin +timestep 0.004 +run_style verlet +fix ensemble left nve # Note: only move particle in left-hand side +fix langevin left langevin 3.0 0.8 100.0 2017 + +# outout +thermo_style custom step temp pzz pe lz +thermo 100 +dump dumpXYZ all xyz 100 traj.xyz + +# run +run 10000 +write_data data.halfhalf From f1aea57e300e4ed252cd804c0e824843bb7c810f Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Mon, 25 Sep 2017 18:48:21 +0200 Subject: [PATCH 08/44] Update readme.md --- examples/USER/pinning/readme.md | 11 +---------- 1 file changed, 1 insertion(+), 10 deletions(-) diff --git a/examples/USER/pinning/readme.md b/examples/USER/pinning/readme.md index 3d5c847cb8..ee758b6f01 100644 --- a/examples/USER/pinning/readme.md +++ b/examples/USER/pinning/readme.md @@ -1,5 +1,4 @@ -This package contains a bias potential that can be used to study solid-liquid transitions with the interface pinning method. -This is done by adding a harmonic potential to the Hamiltonian that bias the system towards two-phase configurations. +This example demonstrate using a bias potential that can be used to study solid-liquid transitions with the interface pinning method. This is done by adding a harmonic potential to the Hamiltonian that bias the system towards two-phase configurations. U_bias = 0.5*k*(Q-a)^2 @@ -13,10 +12,6 @@ Please visit urp.dk/interface_pinning.htm for a detailed bibliography. -# Build -Remember to include the following command when building LAMMPS - make yes-user-pinning - # Use fix [name] [groupID] rhok [nx] [ny] [nz] [spring-constant] [anchor-point] @@ -52,7 +47,3 @@ That method have been experienced to be most effective in the slightly super-hea Ulf R. Pedersen http://www.urp.dk ulf AT urp.dk - -# Cite -Please cite - [Ulf R. Pedersen, J. Chem. Phys. 139, 104102 (2013)] From 88a882b4575d5a3eb05ea7377e8ba4dbba607645 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Mon, 25 Sep 2017 18:59:23 +0200 Subject: [PATCH 09/44] Added reference to fix rhok implimentation --- src/USER-MISC/README | 1 + 1 file changed, 1 insertion(+) diff --git a/src/USER-MISC/README b/src/USER-MISC/README index 65146abd54..8f37825867 100644 --- a/src/USER-MISC/README +++ b/src/USER-MISC/README @@ -47,6 +47,7 @@ fix imd, Axel Kohlmeyer, akohlmey at gmail.com, 9 Nov 2009 fix ipi, Michele Ceriotti (EPFL Lausanne), michele.ceriotti at gmail.com, 24 Nov 2014 fix nvk, Efrem Braun (UC Berkeley), efrem.braun at gmail.com, https://github.com/lammps/lammps/pull/310 fix pimd, Yuxing Peng (U Chicago), yuxing at uchicago.edu, 24 Nov 2014 +fix rhok, Ulf Pedersen (Roskilde U), ulf AT urp.dk, 25 Sep 2017 fix smd, Axel Kohlmeyer, akohlmey at gmail.com, 19 May 2008 fix ti/spring, Rodrigo Freitas (Unicamp/Brazil), rodrigohb at gmail.com, 7 Nov 2013 fix ttm/mod, Sergey Starikov and Vasily Pisarev (JIHT), pisarevvv at gmail.com, 2 Feb 2015 From 2fda041972013c7519ad04179e64d14bc9d6d588 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Thu, 28 Sep 2017 16:00:16 +0200 Subject: [PATCH 10/44] Update fix_rhok.cpp --- src/USER-MISC/fix_rhok.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/USER-MISC/fix_rhok.cpp b/src/USER-MISC/fix_rhok.cpp index 5cbc85a0d6..7283323f2a 100644 --- a/src/USER-MISC/fix_rhok.cpp +++ b/src/USER-MISC/fix_rhok.cpp @@ -1,7 +1,7 @@ /* fix_rhoK_umbrella.cpp - A fix to do umbrella sampling on rho(k). + A fix to do bias potential on rho(k). The usage is as follows: From 4e1eeca869de469c9564e61b5a2b91e69f095a1e Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Thu, 28 Sep 2017 16:02:00 +0200 Subject: [PATCH 11/44] Update fix_rhok.cpp --- src/USER-MISC/fix_rhok.cpp | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/USER-MISC/fix_rhok.cpp b/src/USER-MISC/fix_rhok.cpp index 7283323f2a..6c75cf6e5d 100644 --- a/src/USER-MISC/fix_rhok.cpp +++ b/src/USER-MISC/fix_rhok.cpp @@ -1,7 +1,7 @@ /* - fix_rhoK_umbrella.cpp + fix_rhok.cpp - A fix to do bias potential on rho(k). + A fix to add harmonic potential that bias |rho(k)|. The usage is as follows: @@ -11,8 +11,8 @@ Written by Ulf Pedersen and Patrick Varilly, 4 Feb 2010 Tweaked for LAMMPS 15 Jan 2010 version by Ulf Pedersen, 19 Aug 2010 - Tweaked again March 4th 2012 by Ulf R. Pedersen, - September 2016 by Ulf R. Pedersen + Tweaked again March 4th 2012 by Ulf R. Pedersen, and + September 2016 by Ulf R. Pedersen. */ #include "fix_rhok.h" From 67e48264d9805d698935903ba9c90e3c4d8d2252 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Thu, 28 Sep 2017 16:02:20 +0200 Subject: [PATCH 12/44] Update fix_rhok.cpp --- src/USER-MISC/fix_rhok.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/USER-MISC/fix_rhok.cpp b/src/USER-MISC/fix_rhok.cpp index 6c75cf6e5d..6f5bd45a1d 100644 --- a/src/USER-MISC/fix_rhok.cpp +++ b/src/USER-MISC/fix_rhok.cpp @@ -12,7 +12,7 @@ Written by Ulf Pedersen and Patrick Varilly, 4 Feb 2010 Tweaked for LAMMPS 15 Jan 2010 version by Ulf Pedersen, 19 Aug 2010 Tweaked again March 4th 2012 by Ulf R. Pedersen, and - September 2016 by Ulf R. Pedersen. + September 2017 by Ulf R. Pedersen. */ #include "fix_rhok.h" From 37e55a825b317d38a2559e47935c7bbcaf6b5f59 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Thu, 28 Sep 2017 16:20:03 +0200 Subject: [PATCH 13/44] Create fix_rhok.txt --- doc/src/fix_rhok.txt | 40 ++++++++++++++++++++++++++++++++++++++++ 1 file changed, 40 insertions(+) create mode 100644 doc/src/fix_rhok.txt diff --git a/doc/src/fix_rhok.txt b/doc/src/fix_rhok.txt new file mode 100644 index 0000000000..58838ef464 --- /dev/null +++ b/doc/src/fix_rhok.txt @@ -0,0 +1,40 @@ +"LAMMPS WWW Site"_lws - "LAMMPS Documentation"_ld - "LAMMPS Commands"_lc :c + +:link(lws,http://lammps.sandia.gov) +:link(ld,Manual.html) +:link(lc,Section_commands.html#comm) + +:line + +fix rhok command :h3 + +fix ID group-ID rhok nx ny nz k a + +ID, group-ID are documented in "fix"_fix.html command +nx,ny,nz = k-vektor of collective density field +k = spring constant of bias potential +a = anchor point of bias potential + +[Examples:] + +fix bias all rhok 16 0 0 4.0 16.0 +fix bias Bs rhok 12 12 0 10.0 32.0 + +[Description:] + +The fix applies an force to atoms +:c,image(Eqs/fix_rhok.jpg) +as described in "(Pedersen)"_#Pedersen. + +[Restrictions:] + +This fix is part of the MISC package. It is only enabled if LAMMPS +was built with that package. See the "Making +LAMMPS"_Section_start.html#start_3 section for more info. + +[Default:] none + +:line + +:link(Pedersen) +[(Pedersen)] Pedersen, J. Chem. Phys., 139, 104102 (2013). From e49f0e396b605ddbd367257e4ba991fdb9b59ba8 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Thu, 28 Sep 2017 16:24:26 +0200 Subject: [PATCH 14/44] Create fix_rhok.txt --- doc/src/Eqs/fix_rhok.txt | 9 +++++++++ 1 file changed, 9 insertions(+) create mode 100644 doc/src/Eqs/fix_rhok.txt diff --git a/doc/src/Eqs/fix_rhok.txt b/doc/src/Eqs/fix_rhok.txt new file mode 100644 index 0000000000..9989a434f8 --- /dev/null +++ b/doc/src/Eqs/fix_rhok.txt @@ -0,0 +1,9 @@ +\documentclass[12pt]{article} + +\begin{document} + +$$ + U_\text{bias} = \frac{k}{2}(|\rho_{k}| - a)^2 +$$ + +\end{document} From fe80c57bdef91e3ff53bc6deb37a0d5d1110d33f Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Thu, 28 Sep 2017 17:01:12 +0200 Subject: [PATCH 15/44] more documentation --- doc/src/Eqs/{fix_rhok.txt => fix_rhok.tex} | 0 doc/src/fix_rhok.txt | 8 +++++--- 2 files changed, 5 insertions(+), 3 deletions(-) rename doc/src/Eqs/{fix_rhok.txt => fix_rhok.tex} (100%) diff --git a/doc/src/Eqs/fix_rhok.txt b/doc/src/Eqs/fix_rhok.tex similarity index 100% rename from doc/src/Eqs/fix_rhok.txt rename to doc/src/Eqs/fix_rhok.tex diff --git a/doc/src/fix_rhok.txt b/doc/src/fix_rhok.txt index 58838ef464..e051b83a40 100644 --- a/doc/src/fix_rhok.txt +++ b/doc/src/fix_rhok.txt @@ -8,22 +8,24 @@ fix rhok command :h3 -fix ID group-ID rhok nx ny nz k a +fix ID group-ID rhok nx ny nz k a :pre ID, group-ID are documented in "fix"_fix.html command nx,ny,nz = k-vektor of collective density field k = spring constant of bias potential -a = anchor point of bias potential +a = anchor point of bias potential :ul [Examples:] fix bias all rhok 16 0 0 4.0 16.0 -fix bias Bs rhok 12 12 0 10.0 32.0 +fix bias Bs rhok 12 12 0 10.0 32.0 :pre [Description:] The fix applies an force to atoms + :c,image(Eqs/fix_rhok.jpg) + as described in "(Pedersen)"_#Pedersen. [Restrictions:] From 75b3f34a582cf5f4aa1d39899fb9ebc7ddc057f8 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Thu, 28 Sep 2017 18:16:06 +0200 Subject: [PATCH 16/44] Update documentation --- doc/src/Eqs/fix_rhok.jpg | Bin 0 -> 12041 bytes doc/src/Eqs/fix_rhok.tex | 8 +++++--- doc/src/fix_rhok.txt | 29 ++++++++++++++++++++++++----- 3 files changed, 29 insertions(+), 8 deletions(-) create mode 100644 doc/src/Eqs/fix_rhok.jpg diff --git a/doc/src/Eqs/fix_rhok.jpg b/doc/src/Eqs/fix_rhok.jpg new file mode 100644 index 0000000000000000000000000000000000000000..db8b4ba47295dde92425b642ce098d41d620b78f GIT binary patch literal 12041 zcmbt)1yo$kvhE(-GXpcYy9I|p0>NQ$cL+{!hd_V?2?P%i+}#Q87Mu_uxVvj`0))UL z=bZcJzxTc?Z@r#gYpTBfs&?&OyL;8{+V?Z}D*#wlN=6C*fj~g?!wnVWHe-CBqU@ER8$l+EDS6xObkp+Y#dxLHVy;_6BCRNhTuVAFc{V&0z!N!AubdK z{p|#V_+W#CjE;buZE0* zh=hs&LVFm6{vH1(;eHOlLVU0RBZ2_{^h@fAW~5CuBip2tv@!^B2(KseIT zf$Dpbt_ze3s#6|x@NCE2|5=K$xi~r>1VVLwdnP(__3P25*K}k|& zr$p>H=ruzIMczc(f%nR2$m&?{73M&p8F5q%-jAv4dfrB`k6PP zwb(yT3RI4M2y)p?oTMI*nFN?#W_Q~n8o``vTulYk9{WC+ZAnCXA$t_#`}bq0n6uDP zHo$U;H(E|w-AKn5e&-(qPCV&DPl#=CH<-BK>EX!Ui>#aK3pCi-6AzYu z8fpNj2K32`%l8lT{~`Tem7HJz0RaI(L_vPoG5+BT5eWeqgaUwJxOh-XT4huSK8FZ3 zFCPIF9VeHOLMvH0WdG*V8d4SMT_P80*R^P_3 z_RKh%nQqx2TjOKLT;~Aav|-8#CSH1BN1uM)>uVJ@EtU?y^UqH5H3WSO?`}6?-<$Qi zu0KsW)JCQiO=)9HwKbm~(%-45k1V!LZcM7U}Lv@FZLJ39a4yFVw~bm-7A? zkl^iVakVpaLh5dM_0q^X!o7RpP48}rxK;5eU(!-RFO*cz_;XSBCp`a{85^UF3D_Nz zmhX3_f(fEc>(7*RM>o7iT*c0t7M-Wb7m07W%nsT0xg>J*rd|Xs^0BgirA879M#kC- zv#_pwRsI$axhH-nWRoz>A~?l5r-+>|7zP3KKzjp+-v1} zd8$@Z-@!Y>6EBONnN0gj8)fbh*X7=Rw{FKx{pGvE%xbbPF9%VuqTNg$+RUh}C!-CuUNsVjOMk?VNv1VJ8PU}a+^->W%_-}Id6 zT;R6%<||CsaKu)ArkfWl!HDEYDoHWy@}i_tcns>zUTQOhPK#gmsQHCJ@=a4?yxY=u z6_ItqHAXL+a6(92&wG_0KBF;AW#x(f!|#q-BMG%{Go;hJJ@lYL&2K2zt!^DPqJ@`kv(5{kv!BD@n6Y5uN}{;AKJn@k_!J908JX=t`yEK<(a((4rKty^+I@6M~@ zM4W6Lq8mb0-(({%@~Pu9F2S38Z_DRpRCE5~N5)I0{0aE)1^BTgbz^HMF6eAr)lIV+ z2m`I9r_gPekKN*m)*6xTfiK=-C(&P5Dm030?}4L3bHMz9c}gPP%WZ?N%&mAR$n|E@ z4XIiryu?rh-=v2?F3E z<56*nLaDh_oT9-rV#a=1m6Ws`%Bub`*)VQi5pfe!XV;EV6dpbaNtdeXvGLVoH4P(2 zw~w(oxt(1n|JprKgdcWKB6zw1$1h z8>0kj%Ui~#HOZFTcc*1=*A$w`{e|PO+Y>nTT)%hVN&Iu;J_NWv{C)5-#Tc{jr?-w~ zN!?acj-!3dN$LhU3;T0V!nG`yJF?3p5Q~Uc%ClmBgsbty@?HihdF95&&|1GDA`0GX zZSQ-zTagi2@^Fz9yX^z}lBd|4O%nx@nV>lUh9v4j0NsaC0 z&?D)(c1fgk0i3faJF@t;LA5C{2WUU3seC|7k78GKx=`lY$8&-90(_oD?<%(Z;8XY+ zS+p-V-0Rz9k9$<1i!3uke2qnTWh9rw=)wGMw;Q`|)o6_#4S0h#U3j}3^*r>$eg^uN zy~fA&R~6GoeTxpaZoN9&n~zRUH+yKg-u%eRv5qBk#{c?O3%lvYPl>EeUZ&H)l*9+DsMA_K(x6WDxAw^|&GNK=wQJ%_L7Dc6O1|odWq#T1V zy|KtPGuUKz^g*3aWZHsMQ;*5kuX>NpypCkL_TjP9?a5E!!2V|Oc(08(p7+~hvk>_r zFN;CXWi39e+x(L_mM`H(*0cI0cISGd{#8^`)5ErlQkpaiOv%vLPX}MUWw^wG4{MVM zta>!e_G7H$NxL>EcgkySk3w~22&owx8G@+BCVzxWuk94|&o$N&2#5QKFZK?HZgD4} zli}Iuze*K|i zHT}&LS07=jZiE|4T#1ykMLuv&oG$c`^zruo{?%F&k#STf)r7}Pl191rC9U15S0%HS zp$>A9%C}V-xy`qV4qCk@mI>J9yG9@v@fdhNq>x+!6AiJ*e1+Q6CtGCQ>LO)~wm_z0*^Qc)Z&(`UZmjw$OgWS58 z>nrN4h)*OPYOc&n6&TcHWWRd~I?`$~@>y(CE!SK^J~TIp@2K4aZhr~B=&aIseMj=A z4U&s4L7dNR3s|swlvGHcqiPmt={T;TRU)Ty)f7DMdFz0v&b%3yKm*yt zA-MRlid1glVy04?Qej-In0@;^G^OAZc*68T<)cnXWh2GL z>ID~SzjMV}L183e=@_cz6tAFgXt1r<4270fU{-ul5pND-93R)R`BB%KTP3vSA$mA9 z?cs1YqfoO$=xY$;6!K(tIv+{1U43U%ThwcY`{jb!+L_LXXvWT%AczgsNZG ztygGk-}DlU_v{hX$ZJiY5 zKvEPFmRwuw8Fa4nxVhELm?lCj8>=|Sw};HSnk{lnRftA9B{A|*7|mI!*9l2ej8Paj zRyTso*1`IyL+cCb=`gABpfA@t?DZo#`=&qS^3(EknMeu4`36<#DxG{%Y@J{>8k77S zGoqX%?DcjH-!HCAALX?5&z#%0Eay7i65foywZYvO%`h#uXZd!+u%N1CIV%r?t}lPj zN&ecQv`%6EK`V&Mz*Ks4(}=B}sQT>leA}#3?v}Oh5vO#L^QJ`A{0s`yJ>sI5d6^nTbqCki31@S6;SvS8RfA=7bT+HAJn3cV zGs%H`pWEUZ_BQq-$VA*Rl6WVvURfJ&>+M38BT&VWl7j{(O9vCaM7h{M(}RAaW3|$N zew4_?ni`lj@O$Qd*P;NUQ=%y}9?YQkZVP9$82)T$qkA=4b8TOFbXX57 z&JRg*cLK7y@AaQ0C9#fs3{*y~5|KkZ8h@oKB0;P$OgFh}_H<^BOt(9@tYmrXH>Y9j zlJ1!SJ^R6Rl~?nF$j!(2mN6PAEBvjI0)jPBa>!@aOB|ZN+`Ez}b}lC6h8mTkBq?(kEkvMrSJncWH*+W=+XP`I(z9rrJ6f-0@p(h*I-) zma^d>P`^`Yz@Z>u3W_TZS%SscEr&8YWp7BvMOMyE1j;$gygoeX(U0msnb2mH;;S51 z)VN3$mp?;Iq}+WX&~gueq;u9yB2?wz>W2_^*~Q$#^EEOeJQ-PD6azD-<%m*re8RD( z)J4~#1?urkhbH|C#Pe(%8U+#Bj%*y})6yjdM~T|T$^(AO*?Dyt*w|>dJjt`dDGKdl z{X+c3$x>$U94>ru(`tLh;{D?M!s2w+!%{o{pHjnLJuFg~d~^v0nujWQHed?$>N?5s z+WtB0kYN!Sz5Tbem!`BQSpR!JfN&pZcPx;aE;Ef3!Pq<~oRG z8dF`)GmlKZcjGm?gA`8$*VIHAjFu^#_ibm}^2>&DZR`9jp3_8iG|3QBMdSL>z@juw zfzdf8XIB@d+P)NuHsc#pj3H`dyHJ_s2+c~~=2MIL#ntlCQ*LL=pLlr}VULnIGk12a z^9y2%O-m>1S7Hs?H4!Xc*2ivbcPO$Sz247jlh4$rfMYK_eFsA*5ZxJ(MD=g>5$0x8Ya8d;to*PLNi$HZAfG;|x?B+aoJK zJ3MpJQBhKYadm#6?!!RFx^|9OY#6yA6ZTP8P$Kd<9IWPkjbOqx-kYr$WKm~7Xc43y z_sW3y>Qm_{_*A>iLH-wwLsFz9QB9MW$1cG(SzD%)qd*~oR0iXBjY}-Ec#j+}Ll%Bi z1c_(^oGrDAD*^7*CVIXh<>lS}Q`EJK4Cxm%G1~EY#awFRtR;?tx0#Pf?5-Qeba?ED zTUmF);*qzmYHhynn=XbB*z}0YxwjM1dC$MTc|K@WJhp05SZqwVqDhCQ*J4a&5fL}f zTw!b>B;K$V>HC)OBTvRq0sqZu`S7vUuJ~=14sG#UFlos__gG5^6w+^sTRp>DfAd|4 zaB|zgyxP2N-fi8u4Vw(|RejyIz|!xHmVWTGdzse4C4s=rk7Z&yx=~ka>&)bG zs)83aRo_=ACL@(Ff|tB&IH)n1x6R2rA`rz1I3#JN;6LUdV1Jpj1j0}zhDZVPpZfNn zx!tkk)W~uLyp_rI$b-*gIFWhJtS1nb8Bh3&KYG9Bi#EAFC0*)>lx^99>{B>MJ3MJm zTMW{+HjXVZqar1FqGit`*<@%&UH|E>v6U9JLqVqQ^%~zk ztLzd^RIB0Zz+NS$iw%q24Q;8rQ;thyvaYz-RPBRwleY@a)PbWqoE}wxlRWV;epx7AgMM@d_F#* ztFQg(i~7ZehbGB7xQFsR+4GeOWe5C4MG;&>YOw>;yV}~XE7b(8tgN+p%6u;I$+JRi z!$<*c(!(JY@1b_bHXb3`%BYEG8NFu~jFQ`@D-w-H#q#IqCJH z0JP*^;Q3#P1y&yg3O?$QFAkZo?{`kfkZq!wwhm9W9=z@y)H^U_2{1X2vbf#9CQI>N zbcMD_ytqjyvFIk-x(B$C6%HipNKCA#VVo^}JmD<#%@o4VVC>R}H8H72$*M;*Z(K;9 zDdfbb#c}QX>DR4P`qeGDILw0TTX*cQn>MLC7k%qJu9BiBxGE}9ow}@Zw9iOVbKW*! zGNU5eFEC4GLEMp5jD#0|JiJn~f*(!=5s`lPe;p(~9R5GNDytYf`GKReC^nCcz~SRV}*L{N{Ve)m}fAgUmG>)gw*Q)wI1fT4puOX+87LcTqkokqok z`c8>LrVd)<%Lsu+@yXV4Y!`AE`IA&rqza_s<7Yr4c<}mMEWwe?R~dt{D*jT!Mrz%7 zOT=u5PZLr?S>MfOneSq-%BRGmdyv#-DPuBm;TUIuRP`#mJ1!q@7vywMF&1u8vNaXi zD?iUc@1&*=xsXMwM8D8Rbh9NM=~J_XA-j{t)qW@HE>GIyKgTy6EM7NvPK&@PEC81v9y}6GOkOkWap{Y*kmPER?LQ1#*nx1@U$LIt*Ss4mLzWSsw+1uWuf`qoKu7Onlr! zfZ(oF$@a)3k-C#f@W~T}SgV^@4)aVHIG5ZrYG8tL2%w}BD+DV!~%QJyJaii}pqdAybe?bQZfGik6YzUy8 zPk9!&aJPWlAYj1Oh)d@EniS2;@zOngDWoqVaCy%$5NlbZEvPigpP;P665oyqvlBe) zNg*-af){kM`?UEZDu>M;h=s}@-HXFvkC?|P7DZsArya6!Co?dEJ#EcI;ONZ7mPPpL z`PUe(8?vN1imp5X^Uo4|)vMdpV~^jP`K7}X-3!n{^B{}_c3$f9d~(_XQq`eJu6JeT zpNQM59OWFtR!z?PzlROSxpk>Yc`hLyYAMfMNUEPXgZ%YDx~hJq8l~qC2g6vnQ7O%p zCWAe{#43|?Rq4>1h{D#QM=G?<_Y@T@Q@RZ_w&9VRS;}l(*8FA1Z<&APEFz_Q;_Fw& zkk#g9$~18I;74!-k5D221*eZ?gNLY}2h-S#23sFpW$y?zCBD{Idaa5kikClylV!)K z8Ns9?PXOjW+{eMCC7*iE@?P%ALWrBaY1&M52m$VhZ_6ow*E&9JgOcBvjhQ^a?6ISE|;?i0sY(@s%Ol4uDDPb(t#QAkRBp)>f zFoK&;Y6qgcAzm(9lfc!|Sr4TnG>X+@7vqTknPvVLR{hQ<@=>Mr7z*mGZLQbx?jgy0U@Y3zt@dPpTeNhWO~vA9awvhyXL=a9Vh^L5 ze3p*(bn?p`;z3el$SZ}a^XSu_qAj5hCaqb+GR>S%Le2L=@D~lAXLn**#F4#UD2Y_4 z!>jbEWo%Ua31Kzf2pU4RqYDkm{IVLLa*i>G7KpBGSlH`(nHQj{-Vgvq^Y%aiI3-2% zX-Mt_eFD(@?XV1lgk+5l>Im;(?xF1otQ7SK^gWC(Q08!UTl}>pE%@(Nc#4NfV~N>N zft6n_9Kr8Nmelvgq-t%s>h*esodVr`j(lq|2-2{E1e`J|ah2gB{l-zdJdLv29czHL zWR)EZ;Y&qG>G7bk9O=6~qH)S2RH}LFc)oI7ihKW1=Hpl( zDKh??KW9%GT^oWZbN*B~3mEK4YN&rw_QE3+nw)KO`}JlF+mZ1Ju+Pg?m#BGEGAX}+ z-I%0~(~~aNY}w^iD7T$4zjL~}Sf`$1*f?cak0=)u@8;jvSv2Wsx6a2%1$5+08*Nf- zDejl~%opwdx>}#)SQUj6kay|XH#{t|EAxSwN~4te*Wh;nSYvns@u$1CY0ry9y#nT~ zqcv<7Rc9L_f-#2cyQ<@zy<(#xqQejA^rsQ52v99=G%mSQd7`c5X5!HYIu@aNo$pw? z#t+7B8eon1UgxwS9=H|?Aw%F<;s&7!acw$*P(7W6I?k;Ad|@hg_o>gB>sj(4W_=>p z>>V)e=#1(vb_b z%j#M*+6!X4G+J3b{8KcXhAYuI|4>D%%EI}d`?=qF0iEoKAvm0QuVHj_XL zm?dQ(jqztP+pEv>=l(-B-;!HT6uJ-C+*PY|da<;GQ$=x{f>ym1WLV)uCKZ;ovR<4l zq%*#iO1S)IQxf`Y=?J-{&KqbcqHM9MQ)V8Dv)Nt+%xhx&$w5;CmwLpY#T7oeMsD-g zGZ8l&7Os9Z+pavgGJGBfE5YGe&QMwy`iNm@C5<}IxIxfig~*{5YTL1 zCf7D6DMGN|toMLsX^~@P0(~xK2)$ztPZFhY6ZOaoWV(cRjkDSYWIedmCDhm>lyT1w zv@q!HW8FiQ>llC7aIx;%x4D!G=SyIbi}uH=6cDy9-Ua52V0A$l0tH9uc8vHC@7f+x zoqv480Q|14f@5T8|ErGse;tA&Mg4o+UxwBn0AkQ-AA~=y;FK{k04@mW@1XAk0YvEw zr}ahu8}k7a{YU@1`GX8c`s3|C4E^s^gonUnkfP8ifq(HSf8UWP06fHo27v#~iSL7o zjs$=xqtKB2A0qiDj{c7T&}nf+A8zE|K+psDH}x-E=6?wSe_Z{s{q6i=N`2u-fDF<@ zBo8zBnR>;ucUK*A$#FK{K?nIK$P7Y{X)s2BFVZU1CZ@B&*l-RbO6KgiIPhK@l}|Prx&5 z-;E%+Z$e{>cKqX(b>D#4_h^7#E=+_D4^oc;BWBM5AQTDmyr9%~P)xK@1X&gM#1fPi zgNjUzZ$mZ@u03FmE*!(+kz~#$Ht( zxIp3&dpnt7JH;{vseG(BJd|emS(Ro`oDzCpR{82O)8z;=7QcH5emgoN1m1* zz2;~7Jnb=sp}KdR%R_&ufZ?4fD(WG)ePwZ5{e`W(BYF!0EX+V&60C>wq^iKwY&t&I zZxC_~z);Y9nKpeLM1P57enA#l$O1H!k-ySAY-hpx}>K*d)5oKgcGDl`Yj4oQ0157xv6xwBy^EKNH?d%9??@Mu!r6Rr| ziR>1ue)Z}qIJLS7Cs7MS(c?V5g+{9S4K;cCc+sdFV3Si{hFN5#aDTSb-Y6z?TXSUG>Wfqt~UF6MLVi=?I_I7PnvIVkW{ZB{tbsn zS=VWJ*@YV#E>?StPQV-y6f}^S$>=$xpJ}#K_JT0VOah`ocje*g+9*Z-GvS$v$-yIf zZJyDMDl(4hB?l>2o03xEPYa|);~$&I)U@bmbUGvWvuex7IbW-2DJrl9YWq3RyXL@n^1%yLh2Tm=0c{&Xjm(JQsr*bdjPjS! z&`eYtmw?oHQFVE*`0|BC3e_yFB;T!6<}VUuo}ioar^A-5kX-8#rVx!gK>y0ZsDEXy zUoaZ?J>?4xZp`Zn|;_OLQf|dZgmX`B+7j&>V~0T30{q}mhSOi z_9+y`el#EnLyI*CPL1gLKGunN4*W8~I2)Dqr6={^5IX*tDFVtwGewnv6z5{wXLl7*M6DpL~&||D#(n7CBeJdbA0V*|8B%M4)tYz5P z))gC@2OY~0g3V-}lhghQWWG0A%Fk|DHSG_1`+}6d8%;A{xPZlA<%g%3~3o0jDI zuY%T*nV6ds1H_j;^aF4dL8Bi5i5$SJx5?|Jg(6u-Q5*igY`C+cEZuXV@i%VwfLOQV zCU4=df>oTVHVPEXcghtQ!r~)Pu$bN|TV5%}iCOSKA%HOump~Q21y9<|XY>1}jZsTO z!KT1x?_DKX?H+oi)sp@g<>$)SF}}pO6j5TcYl)VFa5kAYS=7M^tshRHpWZFpk|D) zS|+IHO_YupGV#w>k2x5Pqk?N4WM+ic*TLnh*%Ci?I-HK0QEW+SCV`mmk^7XE_I{5? zzl0OhF(Jd_X@t^U&dg0-b@OL>-WCPUUCrKZ97Hdo|ENE}QhRvpK**qj)h6e;+xQg} zR{G}Yp6aM3fjyk}EIbvz+AwVxC()v*`4Ha8zhVh|B6SfCF2L$mM(Lyj>M+#(S$xp; zCrUcPz71wd7}db|Yqk}!pUNr@?HNG;jB;&XHPAg^QAj9su~!j{G;Kc}PT`Ml8I~<# zrtX%4^O#yCFF~A)OdiQD&PZ1_8;vg=UZEllq4fT`6Z6-|myo6ya=>Gz7Q2G`+a{Ng z%|*h}_&GIrv`Z5xhumf3WJ(H%Ne-22eZe%kkc2`^-NT`{q8WZX#=`M_Bi&!NG#=)% z=$2g~RAM#JR<~I2++O7ZS+p@$SN~^_p54|qT$$X9)Gq!U7fWQE?nniO$^cz;kVr+{ zOHdNQ8#qAeJka1Ntt|xAB?2R6I#|CdHO-9{#iq>0rz%joJ{BGnc6(nKnuX~ql-xH( zr>a5Oqw7>pBML!X<`lq~HPI}Yp;jeAKbsAnlx*~$BK2FvHIqr&%Y6Ski3E)5T|xD6*!Nw3C~*VGY@d0J69QIa~?jBK7KoXbU*ul0PgtL Ang9R* literal 0 HcmV?d00001 diff --git a/doc/src/Eqs/fix_rhok.tex b/doc/src/Eqs/fix_rhok.tex index 9989a434f8..1496a0c737 100644 --- a/doc/src/Eqs/fix_rhok.tex +++ b/doc/src/Eqs/fix_rhok.tex @@ -2,8 +2,10 @@ \begin{document} -$$ - U_\text{bias} = \frac{k}{2}(|\rho_{k}| - a)^2 -$$ +\begin{eqnarray*} + U &=& \frac{1}{2} K (|\rho_{\vec{k}}| - a)^2 \\ + \rho_{\vec{k}} &=& \sum_i^N \exp(-i\vec{k} \cdot \vec{r}_i )/\sqrt{N} \\ + \vec{k} &=& (2\pi n_x /L_x , 2\pi n_y /L_y , 2\pi n_z/L_z ) +\end{eqnarray*} \end{document} diff --git a/doc/src/fix_rhok.txt b/doc/src/fix_rhok.txt index e051b83a40..162d51e67f 100644 --- a/doc/src/fix_rhok.txt +++ b/doc/src/fix_rhok.txt @@ -11,22 +11,36 @@ fix rhok command :h3 fix ID group-ID rhok nx ny nz k a :pre ID, group-ID are documented in "fix"_fix.html command -nx,ny,nz = k-vektor of collective density field -k = spring constant of bias potential +nx, ny, nz = k-vektor of collective density field +K = spring constant of bias potential a = anchor point of bias potential :ul [Examples:] fix bias all rhok 16 0 0 4.0 16.0 -fix bias Bs rhok 12 12 0 10.0 32.0 :pre +fix bias Bs rhok 12 12 0 10.0 24.0 :pre [Description:] -The fix applies an force to atoms +The fix applies an force to atoms given by the potential :c,image(Eqs/fix_rhok.jpg) -as described in "(Pedersen)"_#Pedersen. +as described in "(Pedersen)"_#Pedersen. +The energy and the value of the +collective density field can be written by including + +# output U_bias rho_k_RE rho_k_IM |rho_k| +thermo_style custom step temp pzz pe lz f_bias f_bias\[1\] f_bias\[2\] f_bias\[3\] :pre + +to the input script. + +This field that bias long-range order +can be used to study crystal-liquid interfaces, +and determine melting temperatures "(Pedersen)"_#Pedersen. +The folder {examples/pinning} of the source code +contains an example of using the interface pinning method +on the Lennard-Jones system. [Restrictions:] @@ -34,9 +48,14 @@ This fix is part of the MISC package. It is only enabled if LAMMPS was built with that package. See the "Making LAMMPS"_Section_start.html#start_3 section for more info. +[Related commands:] + +"thermo_style"_thermo_style.html + [Default:] none :line :link(Pedersen) [(Pedersen)] Pedersen, J. Chem. Phys., 139, 104102 (2013). + From 348c4eb7f30c441b1443aaeaf5317d74bd49e80b Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Thu, 28 Sep 2017 18:18:28 +0200 Subject: [PATCH 17/44] add .cpp and .h to root src --- src/fix_rhok.cpp | 245 +++++++++++++++++++++++++++++++++++++++++++++++ src/fix_rhok.h | 81 ++++++++++++++++ 2 files changed, 326 insertions(+) create mode 100644 src/fix_rhok.cpp create mode 100644 src/fix_rhok.h diff --git a/src/fix_rhok.cpp b/src/fix_rhok.cpp new file mode 100644 index 0000000000..6f5bd45a1d --- /dev/null +++ b/src/fix_rhok.cpp @@ -0,0 +1,245 @@ +/* + fix_rhok.cpp + + A fix to add harmonic potential that bias |rho(k)|. + + The usage is as follows: + + fix [name] [groupID] rhoK [nx] [ny] [nz] [kappa = spring constant] [rhoK0] + + where k_i = (2 pi / L_i) * n_i + + Written by Ulf Pedersen and Patrick Varilly, 4 Feb 2010 + Tweaked for LAMMPS 15 Jan 2010 version by Ulf Pedersen, 19 Aug 2010 + Tweaked again March 4th 2012 by Ulf R. Pedersen, and + September 2017 by Ulf R. Pedersen. + */ + +#include "fix_rhok.h" +#include "error.h" +#include "update.h" +#include "respa.h" +#include "atom.h" +#include "domain.h" + +#include +#include +#include +#include +#include + +using namespace LAMMPS_NS; +using namespace FixConst; + +// Constructor: all the parameters to this fix specified in +// the LAMMPS input get passed in +FixRhok::FixRhok( LAMMPS* inLMP, int inArgc, char** inArgv ) + : Fix( inLMP, inArgc, inArgv ) +{ + // Check arguments + if( inArgc != 8 ) + error->all(FLERR,"Illegal fix rhoKUmbrella command" ); + + // Set up fix flags + scalar_flag = 1; // have compute_scalar + vector_flag = 1; // have compute_vector... + size_vector = 3; // ...with this many components + global_freq = 1; // whose value can be computed at every timestep + //scalar_vector_freq = 1; // OLD lammps: whose value can be computed at every timestep + thermo_energy = 1; // this fix changes system's potential energy + extscalar = 0; // but the deltaPE might not scale with # of atoms + extvector = 0; // neither do the components of the vector + + // Parse fix options + int n[3]; + + n[0] = atoi( inArgv[3] ); + n[1] = atoi( inArgv[4] ); + n[2] = atoi( inArgv[5] ); + + mK[0] = n[0]*(2*M_PI / (domain->boxhi[0] - domain->boxlo[0])); + mK[1] = n[1]*(2*M_PI / (domain->boxhi[1] - domain->boxlo[1])); + mK[2] = n[2]*(2*M_PI / (domain->boxhi[2] - domain->boxlo[2])); + + mKappa = atof( inArgv[6] ); + mRhoK0 = atof( inArgv[7] ); +} + +FixRhok::~FixRhok() +{ +} + +// Methods that this fix implements +// -------------------------------- + +// Tells LAMMPS where this fix should act +int +FixRhok::setmask() +{ + int mask = 0; + + // This fix modifies forces... + mask |= POST_FORCE; + mask |= POST_FORCE_RESPA; + mask |= MIN_POST_FORCE; + + // ...and potential energies + mask |= THERMO_ENERGY; + + return mask; +} + +/*int FixRhok::setmask() +{ + int mask = 0; + mask |= POST_FORCE; + mask |= POST_FORCE_RESPA; + mask |= MIN_POST_FORCE; + return mask; +}*/ + + +// Initializes the fix at the beginning of a run +void +FixRhok::init() +{ + // RESPA boilerplate + if( strcmp( update->integrate_style, "respa" ) == 0 ) + mNLevelsRESPA = ((Respa *) update->integrate)->nlevels; + + // Count the number of affected particles + int nThisLocal = 0; + int *mask = atom->mask; + int nlocal = atom->nlocal; + for( int i = 0; i < nlocal; i++ ) { // Iterate through all atoms on this CPU + if( mask[i] & groupbit ) { // ...only those affected by this fix + nThisLocal++; + } + } + MPI_Allreduce( &nThisLocal, &mNThis, + 1, MPI_INT, MPI_SUM, world ); + mSqrtNThis = sqrt( mNThis ); +} + +// Initial application of the fix to a system (when doing MD) +void +FixRhok::setup( int inVFlag ) +{ + if( strcmp( update->integrate_style, "verlet" ) == 0 ) + post_force( inVFlag ); + else + { + ((Respa *) update->integrate)->copy_flevel_f( mNLevelsRESPA - 1 ); + post_force_respa( inVFlag, mNLevelsRESPA - 1,0 ); + ((Respa *) update->integrate)->copy_f_flevel( mNLevelsRESPA - 1 ); + } +} + +// Initial application of the fix to a system (when doing minimization) +void +FixRhok::min_setup( int inVFlag ) +{ + post_force( inVFlag ); +} + +// Modify the forces calculated in the main force loop of ordinary MD +void +FixRhok::post_force( int inVFlag ) +{ + double **x = atom->x; + double **f = atom->f; + int *mask = atom->mask; + int nlocal = atom->nlocal; + + // Loop over locally-owned atoms affected by this fix and calculate the + // partial rhoK's + mRhoKLocal[0] = 0.0; + mRhoKLocal[1] = 0.0; + + for( int i = 0; i < nlocal; i++ ) { // Iterate through all atoms on this CPU + if( mask[i] & groupbit ) { // ...only those affected by this fix + + // rho_k = sum_i exp( - i k.r_i ) + mRhoKLocal[0] += cos( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); + mRhoKLocal[1] -= sin( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); + } + } + + // Now calculate mRhoKGlobal + MPI_Allreduce( mRhoKLocal, mRhoKGlobal, + 2, MPI_DOUBLE, MPI_SUM, world ); + + // Info: < \sum_{i,j} e^{-ik.(r_i - r_j)} > ~ N, so + // we define rho_k as (1 / sqrt(N)) \sum_i e^{-i k.r_i}, so that + // is intensive. + mRhoKGlobal[0] /= mSqrtNThis; + mRhoKGlobal[1] /= mSqrtNThis; + + // We'll need magnitude of rho_k + double rhoK = sqrt( mRhoKGlobal[0]*mRhoKGlobal[0] + + mRhoKGlobal[1]*mRhoKGlobal[1] ); + + for( int i = 0; i < nlocal; i++ ) { // Iterate through all atoms on this CPU + if( mask[i] & groupbit ) { // ...only those affected by this fix + + // Calculate forces + // U = kappa/2 ( |rho_k| - rho_k^0 )^2 + // f_i = -grad_i U = -kappa ( |rho_k| - rho_k^0 ) grad_i |rho_k| + // grad_i |rho_k| = Re( rho_k* (-i k e^{-i k . r_i} / sqrt(N)) ) / |rho_k| + // + // In terms of real and imag parts of rho_k, + // + // Re( rho_k* (-i k e^{-i k . r_i}) ) = + // (- Re[rho_k] * sin( k . r_i ) - Im[rho_k] * cos( k . r_i )) * k + + double sinKRi = sin( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); + double cosKRi = cos( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); + + double prefactor = mKappa * ( rhoK - mRhoK0 ) / rhoK + * (-mRhoKGlobal[0]*sinKRi - mRhoKGlobal[1]*cosKRi) / mSqrtNThis; + f[i][0] -= prefactor * mK[0]; + f[i][1] -= prefactor * mK[1]; + f[i][2] -= prefactor * mK[2]; + } + } +} + +// Forces in RESPA loop +void +FixRhok::post_force_respa( int inVFlag, int inILevel, int inILoop ) +{ + if( inILevel == mNLevelsRESPA - 1 ) + post_force( inVFlag ); +} + +// Forces in minimization loop +void +FixRhok::min_post_force( int inVFlag ) +{ + post_force( inVFlag ); +} + +// Compute the change in the potential energy induced by this fix +double +FixRhok::compute_scalar() +{ + double rhoK = sqrt( mRhoKGlobal[0]*mRhoKGlobal[0] + + mRhoKGlobal[1]*mRhoKGlobal[1] ); + + return 0.5 * mKappa * (rhoK - mRhoK0) * (rhoK - mRhoK0); +} + +// Compute the ith component of the vector +double +FixRhok::compute_vector( int inI ) +{ + if( inI == 0 ) + return mRhoKGlobal[0]; // Real part + else if( inI == 1 ) + return mRhoKGlobal[1]; // Imagniary part + else if( inI == 2 ) + return sqrt( mRhoKGlobal[0]*mRhoKGlobal[0] + + mRhoKGlobal[1]*mRhoKGlobal[1] ); + else + return 12345.0; +} diff --git a/src/fix_rhok.h b/src/fix_rhok.h new file mode 100644 index 0000000000..3e0625430f --- /dev/null +++ b/src/fix_rhok.h @@ -0,0 +1,81 @@ +/* + fix_rhok.h + + A fix to do umbrella sampling on rho(k). + + The usage is as follows: + + fix [name] [groupID] rhoKUmbrella [kx] [ky] [kz] [kappa = spring constant] [rhoK0] + + where k_i = (2 pi / L_i) * n_i + + Written by Ulf Pedersen and Patrick Varilly, 4 Feb 2010 + Tweaked for LAMMPS 15 Jan 2010 version by Ulf Pedersen, 19 Aug 2010 +*/ + +#ifdef FIX_CLASS + +FixStyle(rhok,FixRhok) + +#else + +#ifndef __FIX_RHOK__ +#define __FIX_RHOK__ + +#include "fix.h" + +namespace LAMMPS_NS { + +class FixRhok : public Fix +{ +public: + // Constructor: all the parameters to this fix specified in + // the LAMMPS input get passed in + FixRhok( LAMMPS* inLMP, int inArgc, char** inArgv ); + virtual ~FixRhok(); + + // Methods that this fix implements + // -------------------------------- + + // Tells LAMMPS where this fix should act + int setmask(); + + // Initializes the fix at the beginning of a run + void init(); + + // Initial application of the fix to a system (when doing MD / minimization) + void setup( int inVFlag ); + void min_setup( int inVFlag ); + + // Modify the forces calculated in the main force loop, either when + // doing usual MD, RESPA MD or minimization + void post_force( int inVFlag ); + void post_force_respa( int inVFlag, int inILevel, int inILoop ); + void min_post_force( int inVFlag ); + + // Compute the change in the potential energy induced by this fix + double compute_scalar(); + + // Compute the ith component of the vector associated with this fix + double compute_vector( int inI ); + +private: + // RESPA boilerplate + int mNLevelsRESPA; + + // Defining parameters for this umbrella + double mK[3], mKappa, mRhoK0; + + // Number of particles affected by the fix + int mNThis; + double mSqrtNThis; + + // Real and imaginary parts of rho_k := sum_i exp( - i k . r_i ) + double mRhoKLocal[2], mRhoKGlobal[2]; +}; + +} // namespace LAMMPS_NS + +#endif // __FIX_RHOK__ +#endif // FIX_CLASS + From d11733d3a0ee2b8b160dc634cf86fdd94bc3df2f Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Thu, 28 Sep 2017 18:28:46 +0200 Subject: [PATCH 18/44] typo in equation --- doc/src/Eqs/fix_rhok.jpg | Bin 12041 -> 18330 bytes doc/src/Eqs/fix_rhok.tex | 2 +- 2 files changed, 1 insertion(+), 1 deletion(-) diff --git a/doc/src/Eqs/fix_rhok.jpg b/doc/src/Eqs/fix_rhok.jpg index db8b4ba47295dde92425b642ce098d41d620b78f..829a866be449ecb0f23dc77f547c5b7d053c0469 100644 GIT binary patch literal 18330 zcmcJ$1z26bnm4-k#@(e@k&U~%ySuwP#i6*nYtiCPacgmh0!53qI7N#TFO*CFGIP#1 z=gd5F?{~BFtd;y;N%F3>la=Jn)6&xhfG#H`D+Pc+AYkVPch=hcKiiwVbih+uRgpP}jfrX8OgM*BQhmVVmkBNRARF1_1^J0UHGg1^a&;p85d{cn~UB914U1fH6Q&7@(&?fCvBqP|vk}uJ7Lg1{Mw) z3JijOregmC{1^FY1we)Z0bq0}bN~STs{C8(|MnefbaXQ8-z$aXa}fEL)M_kC@3tTS zfSMRU>Aa+asR-x)#HMfs#6Zwg03eD29!ut*TH}Z#$NV63V@;s$*`3H$rt-f;c;L{1SNb(f+xXw7A~gW)ace!K*sh?6?jQ1 z-Ie2nqxT;L0rp0l766SsjM%yzZcL4_OwCqG!C}u;5rdl7hGw^2^IBARrsB9aqd28qIvJG zFoIdWjS09yDXXEBxZsQij<@{9W?`D^cW)wSTOJduQC+m@zZB5G(UN*+p1~Ij!b~>= zA(BaQtB7ZfnbAm%TT%1!*@lxoP{&PV2*xIFc)cvsv4bmAV3rhV#jo^Y-gOQmQ~X`{ z?*NL+C=P1OF!}TG0-xlAA!s2XO6ViBI*cR9#Q^gBDNIcgG{&V*dT4G(vgMm$`7uWt zn1>>~`acx_Fe)us$)v_2&$_acB!d#@sU83UIo3dmPbt4&SSWsG3DFCk-oLp2zCo*w z`@IcICjC4AU+-lZsU`mv{TpF-`>p)n1hoEldLsXZ_HSWQMBZnc|64M0&BAXQSt1+h z8M6rh>J?ABxK;uXqj7CS{&fzVVi*Jf$R_QXCaC1~2p1YqoE)z0Ry&nvo|u$1X)r{c z_|a_uNRK_l;2YDBNhZ@qab_q6*9Jf@wXz@cMdx$(eM(Nv8t7QSGJoxZte%=!|N-m7mz{KAs1C%XpDj+x(4Rg8OR)lNEG*G zZ%WS1PVKzpwXyRVjIJg;gkfJ75s{flwJH~`{$ljls&(OK&Rag7Np4&Iuv8erT7;`q z?b;TQj&>1F2?vU>l@D^rasQlr04L;ml^J+FucUwkvajy^&vj@ zBmWPVwS}2-CB+5{oUajOR6rW0K<6d_S z4q^;k1bk{cdM8d9vSh0~5{=k_Rr4$9?GqrW*>`J0Br$$;%l>#i_X;cvAH}{^^9$K5 zNar+_s?o`3nOo2yC)PORt8(G>o7n-QaY2Q8 zhJ8mKIYL+uG0S+x$aRkQXy`CIn;-V_S+~gzNrS~T>xy*>EHS(u+cy;oPy-N0bU$ul z?K2?BFI_zjKG;pk@*#k$!8g9_PUPNNEjXRRty$9)V?R((aly270UsH&PvRTKHsbPl z`RhlrrTQe`(K0#a_4!i0ISW;@p(nHin#R8A_!78hd?3&Yf)O3Awj+OmaCW0CU@pmm z6UxKMi!ki^g<4^RQ8*^_(G^Ti#Whi^(%RZ8`4N>~MJ&y*(sLD-mF>a&c;v)svP>&~ z3>c-cqvS#?gjo{!trmAsUKWRCtNre3v_VrKqu`-el#;fY zaaYu#*Y06DM>aSAIq$vVO5l$bp}Nr4sX^ezFZqZmk3-3Qh8G(nG10sgtAW$?76;}_ z_nTCSaiyApH*b3RcJ}0|in_K;4!_#VGswW3ngo`oH2fs5!@=^hDYMy3G2?0km{YMj z_6PB^$k(WqhEXUJVz*ATJMF(|vGQ@+ADMgy2YALoPfP29RTeM2NybmpvpYTdg&hcp zP)v8jwLaBo7JXV_PIqe?&=00;tVv?N!|fRRsa7jq6`Oza5v|eC^DCsqR-JWoV{n4= z9t(e+@iHzef#x2TgTD1kamy@6>_J!+{FtCnU6x3-flahRc1*iUVCWLja3@)F0@)J~ z`X%rM+?n6@si~y{4u7C%eRTfC5=^69JtK9T^2X=lH}I#s?bbz>RnTdwW{f#D^juDD z)Lp?1wXm%WMN}sU9izC6*nL;muLrV5GbU-7x{gueZ8#e&F}Q6Mr$icT%FybXrii8_ zT5V-LPE8!_MQdYeqCbc^xx&OsviSR!zNIC8NI{%Vm?=;GT;n@IVJdxvK=d4qoW4n~ zJFbpccekfXREky8mkl&eqR>q?7x`3FIx0^nbTq_CO4UpiHN79ZUA}X?ZCPK*?T5dsA;r8CBdjUoKATyVE{*I+Yv*L;O6`$)-gL2gZr{@)=*^)<9csO}7iT<%9Z;2Y zP|B=Gmw>;bHorEpqjW^QZN1yJwHWYg+9l8Ky<_~<9ek$E7e~j5Z`7yk&=!FmX5ZEt zTmn!B{jYPqncIh_)5q)_M0wf3_qF5Mud+WMox^@9w06~!@KfpIOwmC@j>{l_8zr=B z?ZA5i?lTK_%2BNM*0epmwYbHBb~!0@*rRKpr)j0UJ2GSUleavG(9(DUWYrFJh5|Z~ z9%a9XBXjALZ=rSQuA=&i)aw-L7jyY6W-8s&c;%f@aBr`QEl+5gX{kEM|M)y;b(`^V z4!NfCS4!z@V~ zs2>q}DF@kDSP)g(()w_!Pn49wk4n+9cjZ3&9aHnjVw#77T4W(Mr6@+qcgkrj+|X zia%oT>KFW&f}8BoCfxd&m?L#Z5!&M5M@gb?F=979|KxM#SkM5O8^#Z*=uTbq1)pzG zq>RfW{qH#$#(z%L{`>-{|6_c9p2Z1+2$}yIj^{D|FAlBePXG`Y3=M|>kAevM`x^iR zhJwZbFk!IBSh303L{&|LVbMv&)LfDa8z@*fIL%yV&LHgKDkg3rDXD3VeX}@RBIeHS z=|%n9lppGYC3KGI~-amNWCG+%{4x6Er)(!0>*O}o-OiHr#qM^g73Ct`1Zg)t^^tI@@V z;pkdIifRE^r>>Ki@1B6XsLhNr9y;qCdk`61n#g^wMqg1GXNjI(b}Kad=ofBc>B*ba z_lR?2#qqRQH^0U%w{&;75;&PgD$KGhakDaYY3>|{U+D@#hg22J_BC--v&O9Y6@GX* z^`d9^yX3m4@6Of43Tdmi@?&w*?zI17jx*^bkM4_{2Zvl*#|gmBz51gc)QT-ubv;^_FWQ>3WGC!n8EuZ4 z7M+&lViWxnP(ByW`Z#60#8X`XJK5LBJmo}ly_G2~mUcJGj~Gpu>R)!`*edAjOj}4f zW)4Z|p_Arc+LqCt*ARxQs7@Qm%T{(Z%;C;Ab15y6LXMZDMdlqAj}kJu_U6=JBdQj# zWJh)CB*s2LpW~MWrEV4NS5~)rR(-RW6xkM+0S?zu@VXA((5WqA!^Z?%`X3!#-X$~XeRLJy z3eE0!R`B)C6UU~2)dghAT~k`|<590v1h;&~b2u%o;h6TPokXv`2HbNc z#>e?E>Mf`<-O2n+`SoLYE?$vGGvTG%Ximo3de9H`S51qPfbLIMRH~S2t@=k{8InAVgCKx{5V;(SCjds;yr?J| za(bB=VlmNPOh%B%27E6178|HV>4VwCpj_=%MQ^pg!g8dPV_|k7v#PsV@)jz9)-1g- zJ943LZ_<8s;cG>y2djo;8|`^&R-}(?39os)e8WBz_7x*W&?LUUF2B~}t!>-v{zqkI zgOMbDnAK*}(u_?!|FpIWGs~h%-4)^#8fcu$4S8LrGG>}eiiLO&NcK`^O$nOurKc>N zs$A~7*GZSGYyAG#dUj!$;&i)~0^4ZChs!zQ3pq8BZ``B{m z=+-Nwl~=4)Jkz#Evg7Q5z%p&v;gDCxXM~GGDe|0dm+jLck=gYLuXM+`TF|-XC+MR>}F8y5dA0=TBWR#Ts7y zfL}iskEgECzxc4cnQSSaoP{!CL`P4#(!h8TfBZVqws(`zp4wbcJWyx#k?;enpt!4< z*tdFfvCple9ARQ9wH+f=HR}Z2~p$HTNzW!fSR8CwBOAU z1`mYO>aphJMYVh~oZ`^8j6Q4^czvr-3P)Ww;zX0d$0w*S|Jso5MUayJz5oJ@rf6xX zm93UJPV1exr*JNNV4U#z_T2#Y?Z6X2@-kSf#no!%M8mpazk|BZ$EtCjnE=5r8l@l|bhHoe{S<#yD(PqEf{ypK3E${xZyP_et4+au z!LqE|>!t;*aoDH-ID+%t()(jz85vAluX_U8qIn#NG(2|>zkObRC?>#?kQ&;a8Mi;Q zxOnS&og#ZO*nIZb-6~NU=bj-fc}|j09evOo#)6`N<`BZ#6eabmJD>(pkH&)qtj{>M z;&)q-DK%HCX_x30eM8({p7M?iMxB&xxNK8D*9_NlJgK`m!Q3S9W;y1Omd>!bRAJWf z0->c#Cnh0*QBlAzt?OJEL1yj?YPHtldpE4bc>2IQ<4!XwCsh?Cdb_qEiI5cGe23Y< z0#73R(vG-OVC7L~)_-!W5l(<9;g{@Jp)NS=c?w5Pl&UJe_3!p4PJLg?%ASB@rRXJ| zlSB2G_wIOku&poid-OfG#ix|4BXCokDqlm*jr>Y1w6R?lR~^2ZAaq%^>{)z^E^XA% zSG?sypBp7S})RK(Vgct zs&T9mk@b)-EWF>nN;~Ez_Rf%L!TLPDIK~jyAjJp!{37MoXpfXplvGt$G3IDZCVyVj zT;bAx0{G8yaFgxi-cT&1O0!#Q?%~q6%_STf0V(B(pdz=-9B-ajd0Ps6j&ZGw42n0D zHw^;c2+(AmV%cZL1H88a-}pQbvol$+9*}3g1 z(zE$~LCt>r?Qn=Fjs#EW&XY6V73&pvS37=9S)aL#NVjcCREs8{GH$8*>EwCSY{-!OvDYL) zE*hdAoRd=YSZ)QwY2joQy{nZz$BfR6l|8&&T>N!A+%lfN%2T>#XIP;1kfhU2UB0A| z^%KW;LXLM{4af!N))`C@&0@n{&TgaqBz9s zOSViXk)OWoJLGEpMA0&bwU$}OfvZt`Duf)+vbaSozgi7Y2PpEQ<_D7Pop7Waw8iRi))NXr}>5VmcRL(9#P z+fIfoDJ6*TWM+h9w_vxw(eB^G#8i4pBfcPPLU&FK8+%>WlQF_d^1+mD->y;af`Jwk zuftDd^sarO!?0MZ3UqvC^-I~un>2dDXL))oFXAycJp9_MTHF?14TS=H1d95ZJ2oqn z=v`5jY2jR${DBn?{#GF~+y2hM4dbX4^Ib$Wq>j8l34bB9Vd(HEUO#<>)2@Kz6YQ&< zoAMC&+>zQmo=N#UfkEKK$Ii5(Mp6#wh0p7w|3LzUoH^sJSy?>lxC z*L}@eckp3LrK?y`Xmmpx8%|^$j8~-KM!>@jfBk;xuGJfbu8&e4twY>lWL9KkkWcwn z;^Ty$c*zs6{{*ZJu+Pk^qB~|3Gk)4j-SFa<-1?|r>(r%}alHzow8h|Tdk~$Fg^-kG zbbTsJ~kx)H~k2ZDC90xZpA3MyQie(!wzKA#Lx>MH8yV`oxJ z4*k66myBOvx2Js2@H0Q)m44;h(S_GVZmvu7!!^+2-j8mi6ct6|{Q*1Y^_5wZ;>a1I zPCG;`aPk#?LSY+hRcwgrES^R<9Y(e;JC}-bs_&MMX$gnpmW1mXai_y7_rqHwTH0(& z-DRXnmFOXmqc7BF={oEySG{P(J`eroFhT=k=li=O&O7y)w}igq?xrpy5SRK2k!I{$9iD6o9J|^|$tu&#OX^&J&sD)) ztQWNKIak)EtVmnK)BGSK+A?&5R+ryooYuu?koJBM7Zc|E!~b9;#Hz5Lj_i~Y)&-?mJ+zUQM`^9Nm;~6-8cnc06pq8a?$AEfJk0rmH2gUOh54!p(S* z+?*#5Ul2x^J7Q(V2m2r?tN5ZKnc`>h62Rps=v|OO@Zy%y$X3F>j+gWpp zETi^93|_$eilU1uHdRhs8quN3-!k_J(CX=Yl;nz^!R%v$SZ0}iAwWX2#h+Q>Hy$FKMpDQapYm zuRGN@K#3z|61QL)qjuYhD0lVV$3;=Ep;`WUk3r}Rl_M|UGGbokm(rqeuHzl?_kE5i z-r`ct;<%sw_DB7*18iKWm8RI}p>A)EC!|Hzw9P!%KblynypKBaB)E`b4Dj4LsTUh**aF2xmc z-|}S?HBJL$&HcOIUgWZk6;*{ZIU1;cc|~yv-SyMCG(xL&#&?xKsmDj$jwX=W{_Pak zL&A+afdtXaY}L%Le|O8(rmblfq3nEp47_N^gr^ln!K*=8MOm(-G(C=Z9fRYq6JAaO z{;`Zgbr#tl9ReZx0UR+28ZiwO2MoQYSYs-JyQJmOWZ}H51@-t&D9S`H+73T|y@9aS z%Y~Ec5ZyM?zTdQlPX%R z4tUZF!zpa_1;VwK6J`4~nq}-T8eTX@aggjm)cbNB`~` zTe+Je$ovSA`=Up)>8p>cb8VLAsv5uc&!t_OGS&u~w=%L4Y%%<2I?MTW?q&NT%HR*l z`M1->nG3eJe9)fj(m$zOOtX>onoSy9jk;=9mYtHWpr*b;7kB-U?BEcJKqzjGz~-O; zfyx!Zs*N)-cOCnjlW&qX$4WYdjT*1CjBQ$vTjLZ(p73tAf5N`qM32plnK zD@cV-vUL)00U9>pHj@hv$9)ZGCo#VW`qBoXBbQBLad^kX#BG{9&&SfBwI+qwhbfYK zC)rR!i^=FHl~L<1m*zoBn%;`h`ta6@G9C9w5v5ouwL`kz?0?rLlUXT^&c>1Nk|9|l zCc2V$$EZjgDA;^#8NIxCpJn&%D6pkNs})~NS<|*?xLF->Eq^;@?~W?8Hg_G)V)EBv zQ&Uk_>%e71rVC4$kX30Q|3nJffr{&O8cat*WVit8NR8M*ADqA%5~-vn%LWDpP6HGb z?xtq7?QuM%WtxDoHaZj+(+U@_9RyeJ9Xhhcb}?h{y^jx_-RP8=eAS)e1iuok?VG)1 zyf1{~ZY_x=x%hS_2}Km$yc%rxs6x~QkD5Yp4D!-Yy zlF3d46U@zoiGen2R#LbWC1$Xi)*DAD=zCja7t*Keag9ZzZ75cFD<~tOP-89(DA~Mt zJ~{fxNmK@td2e!>xMc-?-DGI!;%?UT{}gRM8;u8pfrFO~}$v)RLh$v&xav#S8Z=j?ckjvO9Bd zDP<2D{h=WzcF67D4e9Kea6lPIGqP`Dg(MWuCGB2jj9CNprsa#Y$NS<@Z_XuTK?$=y zC=EmwHym$@AVDLxLLsw`V;}j$_(Ip1cVoVph^-o<#UxnDeSPbWd!GK7GQHvEwc|rk zWLn*3OgQpcrvRC$S)-kt!wD7H;Frq9orJ?QNw3*EZ4>If&?msWSb6fe6z@uC>0LXm z{uS{jvShAT3*~#)KE;K@#NC-LcFPfl^#k~(FF8`NC;BD45H5-_6s}(QiHmF2x+7zD zQs81a-~uljDX+y2xFfTPylcREX=RSYm%J15$t`#h`|LFyL_D83wPE{D_MH@ad|mzA z#f4(zAL}%X`|0mpOJt=b9`;gG%}9}Ik9r9wKDOk?CV!*0mYIpwK{~aRSB6z_5Op%Q z6JLLUd&}yZ+#=&7o7V{X=*eBSID2M^nrfW0EY&vL!CcZtY?W;?Q?)9&TP#P7s65E^ z+N@6hb;+=II6irdNrca|ulT9T%e_yZjz6c`zm+8fwU=JaMrtp2r{jaV4#Z3kEI+BD z#XQ{H*L@1C@|hZE=@x`tnk@mTREFnoVUKXBP;nGVy!*~88G-i`24Zp?#u0W~xaWEv z!oj7)>6IV72PYx4#gvgzBX$J3K#XgufD2laLP08ZI^Rrv^91aIWPCz~QDP^{GO$+3 zS<>&i@-uUYFxyCwIcPEzeRE}N8Zil}x>Lh*i)oa84q3Aelzj{f>?j`@6<1;*x1*K9 zO3gvuH=)1;$hPVCMWffO<&ei+-m z2K0V?#_=f1{TG>-jAYXA0?h9nGc?o$Q|r0y6>$uj;*R|M3Es_Kp0)C3@ zVU(UluS_xOec$a`&?J3fPJAJhk#Zx3$TpPOZ7x%7b~D7>VlAOV+B;sFbcB#3$fGu6Kk#%#0Xs<&RWo1hv&E0> z*J1lB9d4Py-vc?`gfr|wrDcOTec##l1+15k_9_NFU z#^T9_j+i?Wn5K*`r3|jnZ+~|#&wWgJ&2cf(_Uu{WBFbCKQaE@2J8ey;_g%O&)HT>q zI|WQOTaq@E=1c?;8_RRxxIqV5LqO}zdVwX6j=8JNp56fFt<09G<5v3WF~b(r3JH_s zwCjwK5#|;co;};T6N{N&%=_LnCz>(~Rb{HFgf7In9syf2t7N#I$_!zYBePs88kv!j z!EGUra1d-g$#IVfwZnE-V%=ae5KCJ_AT#!}M#1ciN3EeZ3u%Sqq81U^4k7eJuY23W z7v2*!zMkTX>Yc3wrC^9)_>92<)yYR2^c!Mgn8&x?w;1>C3{(A-V{c9NWh!Hq#E-mO za%deY`?GCT z2tT4$;EUOssY@M#0-Z zMhoc|EsU$@2vy!Bglp5xT=IPa;3;)vY-f{d%EYz|I_!GRq$YB{>k{I}N@t7XJK<9% zsG}mWF&R>3;E9q6w3Zj>z~{lRpWP|v)_G?#DQa4(&0m%g}q{AS% z>M=(%hLBfLar~P&!9A&0L*>?&ZO&i*!r$hIoja!^NuqJhBYNRt^J%0r;ZQ7O<&jFA z^$I5urXP<;+4JdyKgJrpNh)cQb<_1pE@$LfU9;(jnc}U{VauA3WDJ94VZ(`VHcaU} z1FsGYUM-NlG&2%hpPk!irrU7M z^wtol^G;EK31@f0!O>$%MU^Zjui1djuJQ*3@LaR%ymO?tq9E0xT(}=Fc6p&FnluG=>{tFZiw##Gj_h9CR*HN*;O z&%ROSFbEH8NsIRHQBI zkX`A6BpYUA!HQbO?HQRBu#3j|n(hgxQ1~GyAI^zCU2&u@hP}rCQD^vxoT8rsWoq2` zU6MA0TF*kY#4;-91S1(I7hTaPGB{GgwU+bBr#Bp?MunNJ40t~{jav}`9L7M`NGZ}0 zOXUY4r&%bLoTqwSO4*+sXQSa!TZrg} zi%HM*g6Aqk3;Tej7ke#YV)C*wGi%feRiZ1jaKm(ls9MqAsGdyJ9vf+*GKi^x1Ez{B zh_5J^9wSvQQnkDQit+`PT4o3$G%ho2ki0_ZS=J3^h7JiiBJy-M%r~e>Be|YKEP>gS zO+_gVI=;`-;k>a3*cQRD^gEE)YUl{vLb>xaDz*6)=*^5XYWOShwG63+ItT$8T@*S<1`;Amiuo4_0?0y#JV($o_$)%2Bue@_@kjiF2SESJP^3wJ z<9|}0ujmMWasCJL|26S{UEv?E4W{hx);-7PCOpUIY|8I~pfN@Nz!CnypG}AG2hTGJ z0Ro7Sh9Ceb@T4iC&k-674+(n4|9Bip2+yV<{gXN5Sqk*0w%?{mde;1R4Ei(i%6-J4pVuK7&qb>cSG7{Jc2hT-b2_ zA8Rw{Pk?7Lnp%r(8jS^dDVoQG8kZ%AL!UBmIyl`{3;S&toBpk0i) z#xE_dZL-75`i%ox9!BMir_CJ^8k0CQeK*_HaL%fk4ks!ZnUwy$)|`4s!SiGyY&fSL z`3dO$IBp=$h$b~`L?mXexC9%x?rYUCXkmgLhG~cn*eDRGxP@p0t!XmNOUbG8L@98* zRvd(SOv-_oTF~_lke<=nGezRw@*23h6QK)oG^(G846>_*#tPi&+Ycue9#4jY=60o$ ze+()Q4kAI8NW)VfI)58FamAN>iQ#EBZh=dW%5H1u%sSczLH*S4K~B6^?aJ1Ez~NO+ zf&xDq>n4M`;jLiydCs1_zb@#^^e2nKwA=zLYk@U>rCq$!c!&&WvD4z?64Yp>J`U3%8r!(rfc%En0OUz$o2w-mPeevgWkt{p{iY&0?YBdA<7E6Hq<& z8hWJyp%JS_W33L=icx$DhJJHMdAaGeA%j8@jPtAxbeVC8_%I*0a+2wT@N&pd_*C$G z>?{TW*GR630;vSLoYAO-Y%)iMw8^h7Plp9hb%_d;{q=@@kzN>0Z>Y{PDaUB+A<}Xs z*q`V!MO4ZBLw)qJG3>KYIvB|x)3tXv2ubYwRdkX0-$*qDbQ;^vzocM+v;9k-kMc(C5<-Wpi{ z>&K_X@|xzj_i7n;<{ z*qf7?MhRCq2f_zvxzz`WL}f&v;mhz-NVR&~f`1{0yRnz2zxQ<(2qYqqTj|9SIZRuv}} zbT|Btqj+=o_07u~nm5N5LRI2s+-K1Gl~%e1AkFi5$uC>(>3PIZ$G|dxPYEu-5xH}I zulM~+(*a-{7p-0QF+25lTq#|p)kN|WHn*$ zIQMz7668MYGtnlP(Pl!47}}QCotBg-?^SN{@Jiy^V3|@w@oo&8-w}cOnzrSLgvKqNUQ90K9u=qbNwgv@7F*Zho^acnK0UV=Bw_v`sND)a{H!=&=~hM$@}UBmoI)}x z)CMNpWcx~f6r0qgo{}{}9jdU~vpmh$zotN%m~e}_y+xk*hmFM(5X4YPO#q!9FKQGC zJzA|A9vl^vNL@)OsuUAMR0(LRz$ny`ay3IL^l(E_U(Q}{Mk%Y=j1-1Q=-?*waR>Pm zgO4512C=R$Op~BNK(9y9=5tho~_&iUmr`)ozB~F<~)~h zz>f`}QCWQ<8z%i706ELU>idvuH}L-C*V-5nAqo)$7 ziad`ZgP#d`yb2PbFV7KvG%{NI=;nEg?fH&s5W51mZ~r@~o9<`KS?;r*Q)WXu8%iG> zHA!ToeiQ`ooep5bb~_>`+HgP`1!YiPUSnA>U7UHcJFO3GBiRYNCPsk_`ys4~fz`g) z6j@T^JvLpfs|{hZG$~FqRBETBY2r-JhAA_b*p$u}3H%C?PM%-2>Z5Jzk7D3O$P*M;_ zAYyC7fuh3Z#Hge`d|}#X2u|LJqgxCW4+<_AN9Jx+3D9{!Rp&#Hrgtn0#SXFkbU zV-=f09$W=GRAI{CB~jrC@Z0Ff+v}%F<_PG6`c!5`Zn-b|K|0)VcHvQ z?H*Ai-+=;0?ci8zhsy#va&i;B>96)dvyJD*!W%aoRq1PD_AWuDp(vkW`$lo#`)%e1MQMju|b-Vz(&ocv*2-3ytWxpa33B^FU}(P zD!k?ih=jXkSzhiIg|3)GlD*@j<@7&ydx!(nP~S16XXO-XZgSVI>E1?!$kCHZkQ^rU zat_f}OBnl`brRQ^UW^h~hBTefxf*4ohQ{7;XX!D+vnx8fW#T|y+c1JSf1D_kVw%}4 z(7nLd!1^-Idj}COsywaH4sAIW(uHba(ef`og=Ghhbk<_yc@l$i=UVuQ|rcZWgLU1pGUGfnasuRuu3 z1`-&b7kwC!!~3YMQ!_wbZdpY~wF{S0N+Pb-EC~}8QA7sE(b(^0o$*|fC-+i-vi$d? zFX3otV?RTZu0u4?tdqZf#vEIM*p`1s)?+X=MP-4g95S}Hd%b%z168OTN8a8|;4Bqc zDjh2BOB7u`JF7q@q9^5gQZQ`^!S@Ns4Q=QixbQJ{k2t3B4p3ghj%>y(@D)h7uFp{F zqfiufHk&a=f_Difk5qMR8R&8!jv_|%GsPJ?n+xF(q!S<@T-s55-p(*Z1kbW(I>Lt1 zk`ByW?hwaczf7I|S;B29{(P==$N^O+7j9oDgX~H-BIWJ`7W5Bh=OdY{H;x& zXLFV>O$bKj!p~8V*&d<`&NmwO1g{l1?AQb0G6t41mDu=K@l1JsYosX-NjSD57@_JI z=2kIv@8KnCU?vrIe+IFm?|Tjt#6^*&so+;(pU!k7yi2v?EDNFhw&~b4l|~8#k+7s} zQCo25&H1`!G4ikc0&Nez{=~bT#=hh#1NByHys4`C!d1P0uNR#3esg{q1>zP5*oFib zI58w_x|oD<8xn)25Mc>u<}K8*{@$av(5zD$pkL9~T`AhhkdQK9~Rn zu5JmLB0*Z)V^Sbo^7qsFk$8k}4iC7B^=!zE0&rK?4#b|%f^Q=*k$Z~2{=(fxn=;o` zw?QHxo%62yDSP|a7>Lybq~W*RDW{yyFSxU7NT>(0P;y7ZeP+Js7F;PIb8 z#YkShox{afr;!|P8_rbfc0!lsl0j@PfbXkJhk2DU!$3a|s90kUH&DJMm?y%9&HjmR zZMPE{C8Sb!cMA#y2|i6Q9SIeMR3s8E4Q-}x(P?_+ytxw|RW zTok@wSImN2$4ia=_8HQBY~+9h8Yh+mh>>g(qrl->BY7m4t^_+_(jJxU(NyNFvOh%= zhU2UZV(?w?-Gw9=&f9}H(H`m>l*$|`DeizFKLpWfeuymSxcqfU`#qUDN3tc_ z3p4f|Lt^aEq4712ASs$Mq#C%!V~{8*#(*)ooFqp`)B({<@gWT&A7BKaBe!UDIjteh zvIXhhlHnl#=y*#(s`q-K!$J}jT2qFWIQuh^_;^vZ0R#>lva7Z*$kF``?$P6Q$p z*y=E!p_V;?f+YCOx@1c7Y*7jK1o>M~fL%!#00X6wLQRd%f;|joz(~4vlnHIHPX2FA z!$v$hj}{aP*yfQb^Ans)!nd$si@euII(w{fZoB-5i2t#EN&@^H)Ej_BQV%^PC~z5B z_i3W=3s?w4L4C-zOPqLFvM7AuKWKKf>9Fa+%?^tY(P4464u{XZMZNGq0L6bWcweqc Q$`BQ5uJjzddRqQ}0B&iL?f?J) literal 12041 zcmbt)1yo$kvhE(-GXpcYy9I|p0>NQ$cL+{!hd_V?2?P%i+}#Q87Mu_uxVvj`0))UL z=bZcJzxTc?Z@r#gYpTBfs&?&OyL;8{+V?Z}D*#wlN=6C*fj~g?!wnVWHe-CBqU@ER8$l+EDS6xObkp+Y#dxLHVy;_6BCRNhTuVAFc{V&0z!N!AubdK z{p|#V_+W#CjE;buZE0* zh=hs&LVFm6{vH1(;eHOlLVU0RBZ2_{^h@fAW~5CuBip2tv@!^B2(KseIT zf$Dpbt_ze3s#6|x@NCE2|5=K$xi~r>1VVLwdnP(__3P25*K}k|& zr$p>H=ruzIMczc(f%nR2$m&?{73M&p8F5q%-jAv4dfrB`k6PP zwb(yT3RI4M2y)p?oTMI*nFN?#W_Q~n8o``vTulYk9{WC+ZAnCXA$t_#`}bq0n6uDP zHo$U;H(E|w-AKn5e&-(qPCV&DPl#=CH<-BK>EX!Ui>#aK3pCi-6AzYu z8fpNj2K32`%l8lT{~`Tem7HJz0RaI(L_vPoG5+BT5eWeqgaUwJxOh-XT4huSK8FZ3 zFCPIF9VeHOLMvH0WdG*V8d4SMT_P80*R^P_3 z_RKh%nQqx2TjOKLT;~Aav|-8#CSH1BN1uM)>uVJ@EtU?y^UqH5H3WSO?`}6?-<$Qi zu0KsW)JCQiO=)9HwKbm~(%-45k1V!LZcM7U}Lv@FZLJ39a4yFVw~bm-7A? zkl^iVakVpaLh5dM_0q^X!o7RpP48}rxK;5eU(!-RFO*cz_;XSBCp`a{85^UF3D_Nz zmhX3_f(fEc>(7*RM>o7iT*c0t7M-Wb7m07W%nsT0xg>J*rd|Xs^0BgirA879M#kC- zv#_pwRsI$axhH-nWRoz>A~?l5r-+>|7zP3KKzjp+-v1} zd8$@Z-@!Y>6EBONnN0gj8)fbh*X7=Rw{FKx{pGvE%xbbPF9%VuqTNg$+RUh}C!-CuUNsVjOMk?VNv1VJ8PU}a+^->W%_-}Id6 zT;R6%<||CsaKu)ArkfWl!HDEYDoHWy@}i_tcns>zUTQOhPK#gmsQHCJ@=a4?yxY=u z6_ItqHAXL+a6(92&wG_0KBF;AW#x(f!|#q-BMG%{Go;hJJ@lYL&2K2zt!^DPqJ@`kv(5{kv!BD@n6Y5uN}{;AKJn@k_!J908JX=t`yEK<(a((4rKty^+I@6M~@ zM4W6Lq8mb0-(({%@~Pu9F2S38Z_DRpRCE5~N5)I0{0aE)1^BTgbz^HMF6eAr)lIV+ z2m`I9r_gPekKN*m)*6xTfiK=-C(&P5Dm030?}4L3bHMz9c}gPP%WZ?N%&mAR$n|E@ z4XIiryu?rh-=v2?F3E z<56*nLaDh_oT9-rV#a=1m6Ws`%Bub`*)VQi5pfe!XV;EV6dpbaNtdeXvGLVoH4P(2 zw~w(oxt(1n|JprKgdcWKB6zw1$1h z8>0kj%Ui~#HOZFTcc*1=*A$w`{e|PO+Y>nTT)%hVN&Iu;J_NWv{C)5-#Tc{jr?-w~ zN!?acj-!3dN$LhU3;T0V!nG`yJF?3p5Q~Uc%ClmBgsbty@?HihdF95&&|1GDA`0GX zZSQ-zTagi2@^Fz9yX^z}lBd|4O%nx@nV>lUh9v4j0NsaC0 z&?D)(c1fgk0i3faJF@t;LA5C{2WUU3seC|7k78GKx=`lY$8&-90(_oD?<%(Z;8XY+ zS+p-V-0Rz9k9$<1i!3uke2qnTWh9rw=)wGMw;Q`|)o6_#4S0h#U3j}3^*r>$eg^uN zy~fA&R~6GoeTxpaZoN9&n~zRUH+yKg-u%eRv5qBk#{c?O3%lvYPl>EeUZ&H)l*9+DsMA_K(x6WDxAw^|&GNK=wQJ%_L7Dc6O1|odWq#T1V zy|KtPGuUKz^g*3aWZHsMQ;*5kuX>NpypCkL_TjP9?a5E!!2V|Oc(08(p7+~hvk>_r zFN;CXWi39e+x(L_mM`H(*0cI0cISGd{#8^`)5ErlQkpaiOv%vLPX}MUWw^wG4{MVM zta>!e_G7H$NxL>EcgkySk3w~22&owx8G@+BCVzxWuk94|&o$N&2#5QKFZK?HZgD4} zli}Iuze*K|i zHT}&LS07=jZiE|4T#1ykMLuv&oG$c`^zruo{?%F&k#STf)r7}Pl191rC9U15S0%HS zp$>A9%C}V-xy`qV4qCk@mI>J9yG9@v@fdhNq>x+!6AiJ*e1+Q6CtGCQ>LO)~wm_z0*^Qc)Z&(`UZmjw$OgWS58 z>nrN4h)*OPYOc&n6&TcHWWRd~I?`$~@>y(CE!SK^J~TIp@2K4aZhr~B=&aIseMj=A z4U&s4L7dNR3s|swlvGHcqiPmt={T;TRU)Ty)f7DMdFz0v&b%3yKm*yt zA-MRlid1glVy04?Qej-In0@;^G^OAZc*68T<)cnXWh2GL z>ID~SzjMV}L183e=@_cz6tAFgXt1r<4270fU{-ul5pND-93R)R`BB%KTP3vSA$mA9 z?cs1YqfoO$=xY$;6!K(tIv+{1U43U%ThwcY`{jb!+L_LXXvWT%AczgsNZG ztygGk-}DlU_v{hX$ZJiY5 zKvEPFmRwuw8Fa4nxVhELm?lCj8>=|Sw};HSnk{lnRftA9B{A|*7|mI!*9l2ej8Paj zRyTso*1`IyL+cCb=`gABpfA@t?DZo#`=&qS^3(EknMeu4`36<#DxG{%Y@J{>8k77S zGoqX%?DcjH-!HCAALX?5&z#%0Eay7i65foywZYvO%`h#uXZd!+u%N1CIV%r?t}lPj zN&ecQv`%6EK`V&Mz*Ks4(}=B}sQT>leA}#3?v}Oh5vO#L^QJ`A{0s`yJ>sI5d6^nTbqCki31@S6;SvS8RfA=7bT+HAJn3cV zGs%H`pWEUZ_BQq-$VA*Rl6WVvURfJ&>+M38BT&VWl7j{(O9vCaM7h{M(}RAaW3|$N zew4_?ni`lj@O$Qd*P;NUQ=%y}9?YQkZVP9$82)T$qkA=4b8TOFbXX57 z&JRg*cLK7y@AaQ0C9#fs3{*y~5|KkZ8h@oKB0;P$OgFh}_H<^BOt(9@tYmrXH>Y9j zlJ1!SJ^R6Rl~?nF$j!(2mN6PAEBvjI0)jPBa>!@aOB|ZN+`Ez}b}lC6h8mTkBq?(kEkvMrSJncWH*+W=+XP`I(z9rrJ6f-0@p(h*I-) zma^d>P`^`Yz@Z>u3W_TZS%SscEr&8YWp7BvMOMyE1j;$gygoeX(U0msnb2mH;;S51 z)VN3$mp?;Iq}+WX&~gueq;u9yB2?wz>W2_^*~Q$#^EEOeJQ-PD6azD-<%m*re8RD( z)J4~#1?urkhbH|C#Pe(%8U+#Bj%*y})6yjdM~T|T$^(AO*?Dyt*w|>dJjt`dDGKdl z{X+c3$x>$U94>ru(`tLh;{D?M!s2w+!%{o{pHjnLJuFg~d~^v0nujWQHed?$>N?5s z+WtB0kYN!Sz5Tbem!`BQSpR!JfN&pZcPx;aE;Ef3!Pq<~oRG z8dF`)GmlKZcjGm?gA`8$*VIHAjFu^#_ibm}^2>&DZR`9jp3_8iG|3QBMdSL>z@juw zfzdf8XIB@d+P)NuHsc#pj3H`dyHJ_s2+c~~=2MIL#ntlCQ*LL=pLlr}VULnIGk12a z^9y2%O-m>1S7Hs?H4!Xc*2ivbcPO$Sz247jlh4$rfMYK_eFsA*5ZxJ(MD=g>5$0x8Ya8d;to*PLNi$HZAfG;|x?B+aoJK zJ3MpJQBhKYadm#6?!!RFx^|9OY#6yA6ZTP8P$Kd<9IWPkjbOqx-kYr$WKm~7Xc43y z_sW3y>Qm_{_*A>iLH-wwLsFz9QB9MW$1cG(SzD%)qd*~oR0iXBjY}-Ec#j+}Ll%Bi z1c_(^oGrDAD*^7*CVIXh<>lS}Q`EJK4Cxm%G1~EY#awFRtR;?tx0#Pf?5-Qeba?ED zTUmF);*qzmYHhynn=XbB*z}0YxwjM1dC$MTc|K@WJhp05SZqwVqDhCQ*J4a&5fL}f zTw!b>B;K$V>HC)OBTvRq0sqZu`S7vUuJ~=14sG#UFlos__gG5^6w+^sTRp>DfAd|4 zaB|zgyxP2N-fi8u4Vw(|RejyIz|!xHmVWTGdzse4C4s=rk7Z&yx=~ka>&)bG zs)83aRo_=ACL@(Ff|tB&IH)n1x6R2rA`rz1I3#JN;6LUdV1Jpj1j0}zhDZVPpZfNn zx!tkk)W~uLyp_rI$b-*gIFWhJtS1nb8Bh3&KYG9Bi#EAFC0*)>lx^99>{B>MJ3MJm zTMW{+HjXVZqar1FqGit`*<@%&UH|E>v6U9JLqVqQ^%~zk ztLzd^RIB0Zz+NS$iw%q24Q;8rQ;thyvaYz-RPBRwleY@a)PbWqoE}wxlRWV;epx7AgMM@d_F#* ztFQg(i~7ZehbGB7xQFsR+4GeOWe5C4MG;&>YOw>;yV}~XE7b(8tgN+p%6u;I$+JRi z!$<*c(!(JY@1b_bHXb3`%BYEG8NFu~jFQ`@D-w-H#q#IqCJH z0JP*^;Q3#P1y&yg3O?$QFAkZo?{`kfkZq!wwhm9W9=z@y)H^U_2{1X2vbf#9CQI>N zbcMD_ytqjyvFIk-x(B$C6%HipNKCA#VVo^}JmD<#%@o4VVC>R}H8H72$*M;*Z(K;9 zDdfbb#c}QX>DR4P`qeGDILw0TTX*cQn>MLC7k%qJu9BiBxGE}9ow}@Zw9iOVbKW*! zGNU5eFEC4GLEMp5jD#0|JiJn~f*(!=5s`lPe;p(~9R5GNDytYf`GKReC^nCcz~SRV}*L{N{Ve)m}fAgUmG>)gw*Q)wI1fT4puOX+87LcTqkokqok z`c8>LrVd)<%Lsu+@yXV4Y!`AE`IA&rqza_s<7Yr4c<}mMEWwe?R~dt{D*jT!Mrz%7 zOT=u5PZLr?S>MfOneSq-%BRGmdyv#-DPuBm;TUIuRP`#mJ1!q@7vywMF&1u8vNaXi zD?iUc@1&*=xsXMwM8D8Rbh9NM=~J_XA-j{t)qW@HE>GIyKgTy6EM7NvPK&@PEC81v9y}6GOkOkWap{Y*kmPER?LQ1#*nx1@U$LIt*Ss4mLzWSsw+1uWuf`qoKu7Onlr! zfZ(oF$@a)3k-C#f@W~T}SgV^@4)aVHIG5ZrYG8tL2%w}BD+DV!~%QJyJaii}pqdAybe?bQZfGik6YzUy8 zPk9!&aJPWlAYj1Oh)d@EniS2;@zOngDWoqVaCy%$5NlbZEvPigpP;P665oyqvlBe) zNg*-af){kM`?UEZDu>M;h=s}@-HXFvkC?|P7DZsArya6!Co?dEJ#EcI;ONZ7mPPpL z`PUe(8?vN1imp5X^Uo4|)vMdpV~^jP`K7}X-3!n{^B{}_c3$f9d~(_XQq`eJu6JeT zpNQM59OWFtR!z?PzlROSxpk>Yc`hLyYAMfMNUEPXgZ%YDx~hJq8l~qC2g6vnQ7O%p zCWAe{#43|?Rq4>1h{D#QM=G?<_Y@T@Q@RZ_w&9VRS;}l(*8FA1Z<&APEFz_Q;_Fw& zkk#g9$~18I;74!-k5D221*eZ?gNLY}2h-S#23sFpW$y?zCBD{Idaa5kikClylV!)K z8Ns9?PXOjW+{eMCC7*iE@?P%ALWrBaY1&M52m$VhZ_6ow*E&9JgOcBvjhQ^a?6ISE|;?i0sY(@s%Ol4uDDPb(t#QAkRBp)>f zFoK&;Y6qgcAzm(9lfc!|Sr4TnG>X+@7vqTknPvVLR{hQ<@=>Mr7z*mGZLQbx?jgy0U@Y3zt@dPpTeNhWO~vA9awvhyXL=a9Vh^L5 ze3p*(bn?p`;z3el$SZ}a^XSu_qAj5hCaqb+GR>S%Le2L=@D~lAXLn**#F4#UD2Y_4 z!>jbEWo%Ua31Kzf2pU4RqYDkm{IVLLa*i>G7KpBGSlH`(nHQj{-Vgvq^Y%aiI3-2% zX-Mt_eFD(@?XV1lgk+5l>Im;(?xF1otQ7SK^gWC(Q08!UTl}>pE%@(Nc#4NfV~N>N zft6n_9Kr8Nmelvgq-t%s>h*esodVr`j(lq|2-2{E1e`J|ah2gB{l-zdJdLv29czHL zWR)EZ;Y&qG>G7bk9O=6~qH)S2RH}LFc)oI7ihKW1=Hpl( zDKh??KW9%GT^oWZbN*B~3mEK4YN&rw_QE3+nw)KO`}JlF+mZ1Ju+Pg?m#BGEGAX}+ z-I%0~(~~aNY}w^iD7T$4zjL~}Sf`$1*f?cak0=)u@8;jvSv2Wsx6a2%1$5+08*Nf- zDejl~%opwdx>}#)SQUj6kay|XH#{t|EAxSwN~4te*Wh;nSYvns@u$1CY0ry9y#nT~ zqcv<7Rc9L_f-#2cyQ<@zy<(#xqQejA^rsQ52v99=G%mSQd7`c5X5!HYIu@aNo$pw? z#t+7B8eon1UgxwS9=H|?Aw%F<;s&7!acw$*P(7W6I?k;Ad|@hg_o>gB>sj(4W_=>p z>>V)e=#1(vb_b z%j#M*+6!X4G+J3b{8KcXhAYuI|4>D%%EI}d`?=qF0iEoKAvm0QuVHj_XL zm?dQ(jqztP+pEv>=l(-B-;!HT6uJ-C+*PY|da<;GQ$=x{f>ym1WLV)uCKZ;ovR<4l zq%*#iO1S)IQxf`Y=?J-{&KqbcqHM9MQ)V8Dv)Nt+%xhx&$w5;CmwLpY#T7oeMsD-g zGZ8l&7Os9Z+pavgGJGBfE5YGe&QMwy`iNm@C5<}IxIxfig~*{5YTL1 zCf7D6DMGN|toMLsX^~@P0(~xK2)$ztPZFhY6ZOaoWV(cRjkDSYWIedmCDhm>lyT1w zv@q!HW8FiQ>llC7aIx;%x4D!G=SyIbi}uH=6cDy9-Ua52V0A$l0tH9uc8vHC@7f+x zoqv480Q|14f@5T8|ErGse;tA&Mg4o+UxwBn0AkQ-AA~=y;FK{k04@mW@1XAk0YvEw zr}ahu8}k7a{YU@1`GX8c`s3|C4E^s^gonUnkfP8ifq(HSf8UWP06fHo27v#~iSL7o zjs$=xqtKB2A0qiDj{c7T&}nf+A8zE|K+psDH}x-E=6?wSe_Z{s{q6i=N`2u-fDF<@ zBo8zBnR>;ucUK*A$#FK{K?nIK$P7Y{X)s2BFVZU1CZ@B&*l-RbO6KgiIPhK@l}|Prx&5 z-;E%+Z$e{>cKqX(b>D#4_h^7#E=+_D4^oc;BWBM5AQTDmyr9%~P)xK@1X&gM#1fPi zgNjUzZ$mZ@u03FmE*!(+kz~#$Ht( zxIp3&dpnt7JH;{vseG(BJd|emS(Ro`oDzCpR{82O)8z;=7QcH5emgoN1m1* zz2;~7Jnb=sp}KdR%R_&ufZ?4fD(WG)ePwZ5{e`W(BYF!0EX+V&60C>wq^iKwY&t&I zZxC_~z);Y9nKpeLM1P57enA#l$O1H!k-ySAY-hpx}>K*d)5oKgcGDl`Yj4oQ0157xv6xwBy^EKNH?d%9??@Mu!r6Rr| ziR>1ue)Z}qIJLS7Cs7MS(c?V5g+{9S4K;cCc+sdFV3Si{hFN5#aDTSb-Y6z?TXSUG>Wfqt~UF6MLVi=?I_I7PnvIVkW{ZB{tbsn zS=VWJ*@YV#E>?StPQV-y6f}^S$>=$xpJ}#K_JT0VOah`ocje*g+9*Z-GvS$v$-yIf zZJyDMDl(4hB?l>2o03xEPYa|);~$&I)U@bmbUGvWvuex7IbW-2DJrl9YWq3RyXL@n^1%yLh2Tm=0c{&Xjm(JQsr*bdjPjS! z&`eYtmw?oHQFVE*`0|BC3e_yFB;T!6<}VUuo}ioar^A-5kX-8#rVx!gK>y0ZsDEXy zUoaZ?J>?4xZp`Zn|;_OLQf|dZgmX`B+7j&>V~0T30{q}mhSOi z_9+y`el#EnLyI*CPL1gLKGunN4*W8~I2)Dqr6={^5IX*tDFVtwGewnv6z5{wXLl7*M6DpL~&||D#(n7CBeJdbA0V*|8B%M4)tYz5P z))gC@2OY~0g3V-}lhghQWWG0A%Fk|DHSG_1`+}6d8%;A{xPZlA<%g%3~3o0jDI zuY%T*nV6ds1H_j;^aF4dL8Bi5i5$SJx5?|Jg(6u-Q5*igY`C+cEZuXV@i%VwfLOQV zCU4=df>oTVHVPEXcghtQ!r~)Pu$bN|TV5%}iCOSKA%HOump~Q21y9<|XY>1}jZsTO z!KT1x?_DKX?H+oi)sp@g<>$)SF}}pO6j5TcYl)VFa5kAYS=7M^tshRHpWZFpk|D) zS|+IHO_YupGV#w>k2x5Pqk?N4WM+ic*TLnh*%Ci?I-HK0QEW+SCV`mmk^7XE_I{5? zzl0OhF(Jd_X@t^U&dg0-b@OL>-WCPUUCrKZ97Hdo|ENE}QhRvpK**qj)h6e;+xQg} zR{G}Yp6aM3fjyk}EIbvz+AwVxC()v*`4Ha8zhVh|B6SfCF2L$mM(Lyj>M+#(S$xp; zCrUcPz71wd7}db|Yqk}!pUNr@?HNG;jB;&XHPAg^QAj9su~!j{G;Kc}PT`Ml8I~<# zrtX%4^O#yCFF~A)OdiQD&PZ1_8;vg=UZEllq4fT`6Z6-|myo6ya=>Gz7Q2G`+a{Ng z%|*h}_&GIrv`Z5xhumf3WJ(H%Ne-22eZe%kkc2`^-NT`{q8WZX#=`M_Bi&!NG#=)% z=$2g~RAM#JR<~I2++O7ZS+p@$SN~^_p54|qT$$X9)Gq!U7fWQE?nniO$^cz;kVr+{ zOHdNQ8#qAeJka1Ntt|xAB?2R6I#|CdHO-9{#iq>0rz%joJ{BGnc6(nKnuX~ql-xH( zr>a5Oqw7>pBML!X<`lq~HPI}Yp;jeAKbsAnlx*~$BK2FvHIqr&%Y6Ski3E)5T|xD6*!Nw3C~*VGY@d0J69QIa~?jBK7KoXbU*ul0PgtL Ang9R* diff --git a/doc/src/Eqs/fix_rhok.tex b/doc/src/Eqs/fix_rhok.tex index 1496a0c737..a468dfedc9 100644 --- a/doc/src/Eqs/fix_rhok.tex +++ b/doc/src/Eqs/fix_rhok.tex @@ -4,7 +4,7 @@ \begin{eqnarray*} U &=& \frac{1}{2} K (|\rho_{\vec{k}}| - a)^2 \\ - \rho_{\vec{k}} &=& \sum_i^N \exp(-i\vec{k} \cdot \vec{r}_i )/\sqrt{N} \\ + \rho_{\vec{k}} &=& \sum_j^N \exp(-i\vec{k} \cdot \vec{r}_j )/\sqrt{N} \\ \vec{k} &=& (2\pi n_x /L_x , 2\pi n_y /L_y , 2\pi n_z/L_z ) \end{eqnarray*} From 10d1741e7f0894393a29cf179b9d0247882e1ecf Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Thu, 28 Sep 2017 18:38:25 +0200 Subject: [PATCH 19/44] Update fix_rhok.txt --- doc/src/fix_rhok.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/doc/src/fix_rhok.txt b/doc/src/fix_rhok.txt index 162d51e67f..8d3e5dcb17 100644 --- a/doc/src/fix_rhok.txt +++ b/doc/src/fix_rhok.txt @@ -8,7 +8,7 @@ fix rhok command :h3 -fix ID group-ID rhok nx ny nz k a :pre +fix ID group-ID rhok nx ny nz K a :pre ID, group-ID are documented in "fix"_fix.html command nx, ny, nz = k-vektor of collective density field From 0f52dd7c5f29e7a72548cd91d8c6ab1366b01782 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Thu, 28 Sep 2017 18:41:06 +0200 Subject: [PATCH 20/44] Update fix_rhok.h --- src/USER-MISC/fix_rhok.h | 11 ----------- 1 file changed, 11 deletions(-) diff --git a/src/USER-MISC/fix_rhok.h b/src/USER-MISC/fix_rhok.h index 3e0625430f..784bf1d14a 100644 --- a/src/USER-MISC/fix_rhok.h +++ b/src/USER-MISC/fix_rhok.h @@ -1,16 +1,5 @@ /* fix_rhok.h - - A fix to do umbrella sampling on rho(k). - - The usage is as follows: - - fix [name] [groupID] rhoKUmbrella [kx] [ky] [kz] [kappa = spring constant] [rhoK0] - - where k_i = (2 pi / L_i) * n_i - - Written by Ulf Pedersen and Patrick Varilly, 4 Feb 2010 - Tweaked for LAMMPS 15 Jan 2010 version by Ulf Pedersen, 19 Aug 2010 */ #ifdef FIX_CLASS From 285a123c904d45a27bdc49e2510f4780fce73e3c Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Thu, 28 Sep 2017 18:42:15 +0200 Subject: [PATCH 21/44] Update fix_rhok.cpp --- src/USER-MISC/fix_rhok.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/USER-MISC/fix_rhok.cpp b/src/USER-MISC/fix_rhok.cpp index 6f5bd45a1d..344b4726a8 100644 --- a/src/USER-MISC/fix_rhok.cpp +++ b/src/USER-MISC/fix_rhok.cpp @@ -5,7 +5,7 @@ The usage is as follows: - fix [name] [groupID] rhoK [nx] [ny] [nz] [kappa = spring constant] [rhoK0] + fix [name] [groupID] rhok [nx] [ny] [nz] [K] [a] where k_i = (2 pi / L_i) * n_i From c66ddf9ac0e8c0c256a14220a80395ab581db31e Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Fri, 29 Sep 2017 08:03:02 +0200 Subject: [PATCH 22/44] Update fix_rhok.txt --- doc/src/fix_rhok.txt | 22 +++++++++------------- 1 file changed, 9 insertions(+), 13 deletions(-) diff --git a/doc/src/fix_rhok.txt b/doc/src/fix_rhok.txt index 8d3e5dcb17..5046303329 100644 --- a/doc/src/fix_rhok.txt +++ b/doc/src/fix_rhok.txt @@ -18,7 +18,9 @@ a = anchor point of bias potential :ul [Examples:] fix bias all rhok 16 0 0 4.0 16.0 -fix bias Bs rhok 12 12 0 10.0 24.0 :pre +fix 1 all npt temp 0.8 0.8 4.0 z 2.2 2.2 8.0 +# output U_bias rho_k_RE rho_k_IM |rho_k| +thermo_style custom step temp pzz lz f_bias f_bias\[1\] f_bias\[2\] f_bias\[3\] :pre [Description:] @@ -27,20 +29,14 @@ The fix applies an force to atoms given by the potential :c,image(Eqs/fix_rhok.jpg) as described in "(Pedersen)"_#Pedersen. -The energy and the value of the -collective density field can be written by including -# output U_bias rho_k_RE rho_k_IM |rho_k| -thermo_style custom step temp pzz pe lz f_bias f_bias\[1\] f_bias\[2\] f_bias\[3\] :pre - -to the input script. - -This field that bias long-range order -can be used to study crystal-liquid interfaces, +This field, that bias configurations with long-range order, +can be used to study crystal-liquid interfaces and determine melting temperatures "(Pedersen)"_#Pedersen. -The folder {examples/pinning} of the source code -contains an example of using the interface pinning method -on the Lennard-Jones system. + +An example of using the interface pinning method +is located in the folder +{examples/USER/pinning} of the source code to LAMMPS. [Restrictions:] From fe14eeccac974c48058c6c63490e4a84eca2d7a8 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 10:54:22 +0200 Subject: [PATCH 23/44] Delete fix_rhok.h --- src/fix_rhok.h | 81 -------------------------------------------------- 1 file changed, 81 deletions(-) delete mode 100644 src/fix_rhok.h diff --git a/src/fix_rhok.h b/src/fix_rhok.h deleted file mode 100644 index 3e0625430f..0000000000 --- a/src/fix_rhok.h +++ /dev/null @@ -1,81 +0,0 @@ -/* - fix_rhok.h - - A fix to do umbrella sampling on rho(k). - - The usage is as follows: - - fix [name] [groupID] rhoKUmbrella [kx] [ky] [kz] [kappa = spring constant] [rhoK0] - - where k_i = (2 pi / L_i) * n_i - - Written by Ulf Pedersen and Patrick Varilly, 4 Feb 2010 - Tweaked for LAMMPS 15 Jan 2010 version by Ulf Pedersen, 19 Aug 2010 -*/ - -#ifdef FIX_CLASS - -FixStyle(rhok,FixRhok) - -#else - -#ifndef __FIX_RHOK__ -#define __FIX_RHOK__ - -#include "fix.h" - -namespace LAMMPS_NS { - -class FixRhok : public Fix -{ -public: - // Constructor: all the parameters to this fix specified in - // the LAMMPS input get passed in - FixRhok( LAMMPS* inLMP, int inArgc, char** inArgv ); - virtual ~FixRhok(); - - // Methods that this fix implements - // -------------------------------- - - // Tells LAMMPS where this fix should act - int setmask(); - - // Initializes the fix at the beginning of a run - void init(); - - // Initial application of the fix to a system (when doing MD / minimization) - void setup( int inVFlag ); - void min_setup( int inVFlag ); - - // Modify the forces calculated in the main force loop, either when - // doing usual MD, RESPA MD or minimization - void post_force( int inVFlag ); - void post_force_respa( int inVFlag, int inILevel, int inILoop ); - void min_post_force( int inVFlag ); - - // Compute the change in the potential energy induced by this fix - double compute_scalar(); - - // Compute the ith component of the vector associated with this fix - double compute_vector( int inI ); - -private: - // RESPA boilerplate - int mNLevelsRESPA; - - // Defining parameters for this umbrella - double mK[3], mKappa, mRhoK0; - - // Number of particles affected by the fix - int mNThis; - double mSqrtNThis; - - // Real and imaginary parts of rho_k := sum_i exp( - i k . r_i ) - double mRhoKLocal[2], mRhoKGlobal[2]; -}; - -} // namespace LAMMPS_NS - -#endif // __FIX_RHOK__ -#endif // FIX_CLASS - From 4dcc49ebe2c397911bd7392baafbdbc685f080c6 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 10:56:07 +0200 Subject: [PATCH 24/44] Delete fix_rhok.cpp --- src/fix_rhok.cpp | 245 ----------------------------------------------- 1 file changed, 245 deletions(-) delete mode 100644 src/fix_rhok.cpp diff --git a/src/fix_rhok.cpp b/src/fix_rhok.cpp deleted file mode 100644 index 6f5bd45a1d..0000000000 --- a/src/fix_rhok.cpp +++ /dev/null @@ -1,245 +0,0 @@ -/* - fix_rhok.cpp - - A fix to add harmonic potential that bias |rho(k)|. - - The usage is as follows: - - fix [name] [groupID] rhoK [nx] [ny] [nz] [kappa = spring constant] [rhoK0] - - where k_i = (2 pi / L_i) * n_i - - Written by Ulf Pedersen and Patrick Varilly, 4 Feb 2010 - Tweaked for LAMMPS 15 Jan 2010 version by Ulf Pedersen, 19 Aug 2010 - Tweaked again March 4th 2012 by Ulf R. Pedersen, and - September 2017 by Ulf R. Pedersen. - */ - -#include "fix_rhok.h" -#include "error.h" -#include "update.h" -#include "respa.h" -#include "atom.h" -#include "domain.h" - -#include -#include -#include -#include -#include - -using namespace LAMMPS_NS; -using namespace FixConst; - -// Constructor: all the parameters to this fix specified in -// the LAMMPS input get passed in -FixRhok::FixRhok( LAMMPS* inLMP, int inArgc, char** inArgv ) - : Fix( inLMP, inArgc, inArgv ) -{ - // Check arguments - if( inArgc != 8 ) - error->all(FLERR,"Illegal fix rhoKUmbrella command" ); - - // Set up fix flags - scalar_flag = 1; // have compute_scalar - vector_flag = 1; // have compute_vector... - size_vector = 3; // ...with this many components - global_freq = 1; // whose value can be computed at every timestep - //scalar_vector_freq = 1; // OLD lammps: whose value can be computed at every timestep - thermo_energy = 1; // this fix changes system's potential energy - extscalar = 0; // but the deltaPE might not scale with # of atoms - extvector = 0; // neither do the components of the vector - - // Parse fix options - int n[3]; - - n[0] = atoi( inArgv[3] ); - n[1] = atoi( inArgv[4] ); - n[2] = atoi( inArgv[5] ); - - mK[0] = n[0]*(2*M_PI / (domain->boxhi[0] - domain->boxlo[0])); - mK[1] = n[1]*(2*M_PI / (domain->boxhi[1] - domain->boxlo[1])); - mK[2] = n[2]*(2*M_PI / (domain->boxhi[2] - domain->boxlo[2])); - - mKappa = atof( inArgv[6] ); - mRhoK0 = atof( inArgv[7] ); -} - -FixRhok::~FixRhok() -{ -} - -// Methods that this fix implements -// -------------------------------- - -// Tells LAMMPS where this fix should act -int -FixRhok::setmask() -{ - int mask = 0; - - // This fix modifies forces... - mask |= POST_FORCE; - mask |= POST_FORCE_RESPA; - mask |= MIN_POST_FORCE; - - // ...and potential energies - mask |= THERMO_ENERGY; - - return mask; -} - -/*int FixRhok::setmask() -{ - int mask = 0; - mask |= POST_FORCE; - mask |= POST_FORCE_RESPA; - mask |= MIN_POST_FORCE; - return mask; -}*/ - - -// Initializes the fix at the beginning of a run -void -FixRhok::init() -{ - // RESPA boilerplate - if( strcmp( update->integrate_style, "respa" ) == 0 ) - mNLevelsRESPA = ((Respa *) update->integrate)->nlevels; - - // Count the number of affected particles - int nThisLocal = 0; - int *mask = atom->mask; - int nlocal = atom->nlocal; - for( int i = 0; i < nlocal; i++ ) { // Iterate through all atoms on this CPU - if( mask[i] & groupbit ) { // ...only those affected by this fix - nThisLocal++; - } - } - MPI_Allreduce( &nThisLocal, &mNThis, - 1, MPI_INT, MPI_SUM, world ); - mSqrtNThis = sqrt( mNThis ); -} - -// Initial application of the fix to a system (when doing MD) -void -FixRhok::setup( int inVFlag ) -{ - if( strcmp( update->integrate_style, "verlet" ) == 0 ) - post_force( inVFlag ); - else - { - ((Respa *) update->integrate)->copy_flevel_f( mNLevelsRESPA - 1 ); - post_force_respa( inVFlag, mNLevelsRESPA - 1,0 ); - ((Respa *) update->integrate)->copy_f_flevel( mNLevelsRESPA - 1 ); - } -} - -// Initial application of the fix to a system (when doing minimization) -void -FixRhok::min_setup( int inVFlag ) -{ - post_force( inVFlag ); -} - -// Modify the forces calculated in the main force loop of ordinary MD -void -FixRhok::post_force( int inVFlag ) -{ - double **x = atom->x; - double **f = atom->f; - int *mask = atom->mask; - int nlocal = atom->nlocal; - - // Loop over locally-owned atoms affected by this fix and calculate the - // partial rhoK's - mRhoKLocal[0] = 0.0; - mRhoKLocal[1] = 0.0; - - for( int i = 0; i < nlocal; i++ ) { // Iterate through all atoms on this CPU - if( mask[i] & groupbit ) { // ...only those affected by this fix - - // rho_k = sum_i exp( - i k.r_i ) - mRhoKLocal[0] += cos( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); - mRhoKLocal[1] -= sin( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); - } - } - - // Now calculate mRhoKGlobal - MPI_Allreduce( mRhoKLocal, mRhoKGlobal, - 2, MPI_DOUBLE, MPI_SUM, world ); - - // Info: < \sum_{i,j} e^{-ik.(r_i - r_j)} > ~ N, so - // we define rho_k as (1 / sqrt(N)) \sum_i e^{-i k.r_i}, so that - // is intensive. - mRhoKGlobal[0] /= mSqrtNThis; - mRhoKGlobal[1] /= mSqrtNThis; - - // We'll need magnitude of rho_k - double rhoK = sqrt( mRhoKGlobal[0]*mRhoKGlobal[0] - + mRhoKGlobal[1]*mRhoKGlobal[1] ); - - for( int i = 0; i < nlocal; i++ ) { // Iterate through all atoms on this CPU - if( mask[i] & groupbit ) { // ...only those affected by this fix - - // Calculate forces - // U = kappa/2 ( |rho_k| - rho_k^0 )^2 - // f_i = -grad_i U = -kappa ( |rho_k| - rho_k^0 ) grad_i |rho_k| - // grad_i |rho_k| = Re( rho_k* (-i k e^{-i k . r_i} / sqrt(N)) ) / |rho_k| - // - // In terms of real and imag parts of rho_k, - // - // Re( rho_k* (-i k e^{-i k . r_i}) ) = - // (- Re[rho_k] * sin( k . r_i ) - Im[rho_k] * cos( k . r_i )) * k - - double sinKRi = sin( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); - double cosKRi = cos( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); - - double prefactor = mKappa * ( rhoK - mRhoK0 ) / rhoK - * (-mRhoKGlobal[0]*sinKRi - mRhoKGlobal[1]*cosKRi) / mSqrtNThis; - f[i][0] -= prefactor * mK[0]; - f[i][1] -= prefactor * mK[1]; - f[i][2] -= prefactor * mK[2]; - } - } -} - -// Forces in RESPA loop -void -FixRhok::post_force_respa( int inVFlag, int inILevel, int inILoop ) -{ - if( inILevel == mNLevelsRESPA - 1 ) - post_force( inVFlag ); -} - -// Forces in minimization loop -void -FixRhok::min_post_force( int inVFlag ) -{ - post_force( inVFlag ); -} - -// Compute the change in the potential energy induced by this fix -double -FixRhok::compute_scalar() -{ - double rhoK = sqrt( mRhoKGlobal[0]*mRhoKGlobal[0] - + mRhoKGlobal[1]*mRhoKGlobal[1] ); - - return 0.5 * mKappa * (rhoK - mRhoK0) * (rhoK - mRhoK0); -} - -// Compute the ith component of the vector -double -FixRhok::compute_vector( int inI ) -{ - if( inI == 0 ) - return mRhoKGlobal[0]; // Real part - else if( inI == 1 ) - return mRhoKGlobal[1]; // Imagniary part - else if( inI == 2 ) - return sqrt( mRhoKGlobal[0]*mRhoKGlobal[0] - + mRhoKGlobal[1]*mRhoKGlobal[1] ); - else - return 12345.0; -} From 40ae6f215bb7cbd65eb29f7694c068f88c0626dd Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 10:58:21 +0200 Subject: [PATCH 25/44] add the usual LAMMPS copyright --- src/USER-MISC/fix_rhok.cpp | 30 ++++++++++++++---------------- 1 file changed, 14 insertions(+), 16 deletions(-) diff --git a/src/USER-MISC/fix_rhok.cpp b/src/USER-MISC/fix_rhok.cpp index 344b4726a8..d6e516ac43 100644 --- a/src/USER-MISC/fix_rhok.cpp +++ b/src/USER-MISC/fix_rhok.cpp @@ -1,19 +1,17 @@ -/* - fix_rhok.cpp - - A fix to add harmonic potential that bias |rho(k)|. - - The usage is as follows: - - fix [name] [groupID] rhok [nx] [ny] [nz] [K] [a] - - where k_i = (2 pi / L_i) * n_i - - Written by Ulf Pedersen and Patrick Varilly, 4 Feb 2010 - Tweaked for LAMMPS 15 Jan 2010 version by Ulf Pedersen, 19 Aug 2010 - Tweaked again March 4th 2012 by Ulf R. Pedersen, and - September 2017 by Ulf R. Pedersen. - */ +/* ---------------------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + http://lammps.sandia.gov, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + Copyright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ + +/* ---------------------------------------------------------------------- + Contributing author: Ulf R. Pedersen, ulf@urp.dk +------------------------------------------------------------------------- */ #include "fix_rhok.h" #include "error.h" From 8f79f5ddb9daaee3c6cbc908785a8a978c4f689f Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 10:59:46 +0200 Subject: [PATCH 26/44] add the LAMMPS copyright --- src/USER-MISC/fix_rhok.h | 13 ++++++++++--- 1 file changed, 10 insertions(+), 3 deletions(-) diff --git a/src/USER-MISC/fix_rhok.h b/src/USER-MISC/fix_rhok.h index 784bf1d14a..145b3e7c1e 100644 --- a/src/USER-MISC/fix_rhok.h +++ b/src/USER-MISC/fix_rhok.h @@ -1,6 +1,13 @@ -/* - fix_rhok.h -*/ +/* -*- c++ -*- ---------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + http://lammps.sandia.gov, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + Copyright (2003) Sandia Corporation. Under the terms of Contract + DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains + certain rights in this software. This software is distributed under + the GNU General Public License. + See the README file in the top-level LAMMPS directory. +------------------------------------------------------------------------- */ #ifdef FIX_CLASS From 8a2cf5ce8e778db9c29bc4a1205fd215c296914f Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 11:13:10 +0200 Subject: [PATCH 27/44] reformatting and use citeme class --- src/USER-MISC/fix_rhok.cpp | 35 +++++++++++++++++++++-------------- 1 file changed, 21 insertions(+), 14 deletions(-) diff --git a/src/USER-MISC/fix_rhok.cpp b/src/USER-MISC/fix_rhok.cpp index d6e516ac43..df8b8bd667 100644 --- a/src/USER-MISC/fix_rhok.cpp +++ b/src/USER-MISC/fix_rhok.cpp @@ -15,6 +15,7 @@ #include "fix_rhok.h" #include "error.h" +#include "citeme.h" #include "update.h" #include "respa.h" #include "atom.h" @@ -29,14 +30,30 @@ using namespace LAMMPS_NS; using namespace FixConst; -// Constructor: all the parameters to this fix specified in -// the LAMMPS input get passed in +static const char cite_flow_gauss[] = + "Bias on the collective density field (fix rhok):\n\n" + "@Article{pedersen_jcp139_104102_2013,\n" + "title = {Direct calculation of the solid-liquid Gibbs free energy difference in a single equilibrium simulation},\n" + "volume = {139},\n" + "number = {10},\n" + "url = {http://aip.scitation.org/doi/10.1063/1.4818747},\n" + "doi = {10.1063/1.4818747},\n" + "urldate = {2017-10-03},\n" + "journal = {J. Chem. Phys.},\n" + "author = {Pedersen, Ulf R.},\n" + "year = {2013},\n" + "pages = {104102}\n" + "}\n\n"; + FixRhok::FixRhok( LAMMPS* inLMP, int inArgc, char** inArgv ) : Fix( inLMP, inArgc, inArgv ) { + + if (lmp->citeme) lmp->citeme->add(cite_flow_gauss); + // Check arguments if( inArgc != 8 ) - error->all(FLERR,"Illegal fix rhoKUmbrella command" ); + error->all(FLERR,"Illegal fix rhoKUmbrella command" ); // Set up fix flags scalar_flag = 1; // have compute_scalar @@ -87,16 +104,6 @@ FixRhok::setmask() return mask; } -/*int FixRhok::setmask() -{ - int mask = 0; - mask |= POST_FORCE; - mask |= POST_FORCE_RESPA; - mask |= MIN_POST_FORCE; - return mask; -}*/ - - // Initializes the fix at the beginning of a run void FixRhok::init() @@ -151,7 +158,7 @@ FixRhok::post_force( int inVFlag ) // Loop over locally-owned atoms affected by this fix and calculate the // partial rhoK's - mRhoKLocal[0] = 0.0; + mRhoKLocal[0] = 0.0; mRhoKLocal[1] = 0.0; for( int i = 0; i < nlocal; i++ ) { // Iterate through all atoms on this CPU From e6d687faaca53a6dc03d17be1178ba3ccb189e41 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 11:15:49 +0200 Subject: [PATCH 28/44] Rename crystal.lmp to in.crystal --- examples/USER/pinning/{crystal.lmp => in.crystal} | 0 1 file changed, 0 insertions(+), 0 deletions(-) rename examples/USER/pinning/{crystal.lmp => in.crystal} (100%) diff --git a/examples/USER/pinning/crystal.lmp b/examples/USER/pinning/in.crystal similarity index 100% rename from examples/USER/pinning/crystal.lmp rename to examples/USER/pinning/in.crystal From 7e8bbe84814ba2545c3cb5f3943c1737fb09d97c Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 11:16:29 +0200 Subject: [PATCH 29/44] Rename pinning.lmp to in.pinning --- examples/USER/pinning/{pinning.lmp => in.pinning} | 0 1 file changed, 0 insertions(+), 0 deletions(-) rename examples/USER/pinning/{pinning.lmp => in.pinning} (100%) diff --git a/examples/USER/pinning/pinning.lmp b/examples/USER/pinning/in.pinning similarity index 100% rename from examples/USER/pinning/pinning.lmp rename to examples/USER/pinning/in.pinning From 245bf745521d77564b3d45ba0d20eb2cf179488c Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 11:16:54 +0200 Subject: [PATCH 30/44] Rename setup.lmp to in.setup --- examples/USER/pinning/{setup.lmp => in.setup} | 0 1 file changed, 0 insertions(+), 0 deletions(-) rename examples/USER/pinning/{setup.lmp => in.setup} (100%) diff --git a/examples/USER/pinning/setup.lmp b/examples/USER/pinning/in.setup similarity index 100% rename from examples/USER/pinning/setup.lmp rename to examples/USER/pinning/in.setup From 35cc79597261b1a96a1c25fc987705135182bcd1 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 11:17:43 +0200 Subject: [PATCH 31/44] Update in.setup --- examples/USER/pinning/in.setup | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/examples/USER/pinning/in.setup b/examples/USER/pinning/in.setup index 4ab9e4498b..0febd3c4f6 100644 --- a/examples/USER/pinning/in.setup +++ b/examples/USER/pinning/in.setup @@ -1,11 +1,11 @@ -units lj -dimension 3 -boundary p p p -atom_style atomic +units lj +dimension 3 +boundary p p p +atom_style atomic # truncated and shifted LJ potential -pair_style lj/cut 2.5 -pair_modify shift yes +pair_style lj/cut 2.5 +pair_modify shift yes # fcc lattice lattice fcc 0.9731 From 1e790fbafef82af620934d09bfece96e1adafed7 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 11:44:23 +0200 Subject: [PATCH 32/44] Update readme.md Now use the standard namings. Corrected spelling errors. --- examples/USER/pinning/readme.md | 73 +++++++++++++++++++++------------ 1 file changed, 46 insertions(+), 27 deletions(-) diff --git a/examples/USER/pinning/readme.md b/examples/USER/pinning/readme.md index ee758b6f01..2071b63b4a 100644 --- a/examples/USER/pinning/readme.md +++ b/examples/USER/pinning/readme.md @@ -1,49 +1,68 @@ -This example demonstrate using a bias potential that can be used to study solid-liquid transitions with the interface pinning method. This is done by adding a harmonic potential to the Hamiltonian that bias the system towards two-phase configurations. +This example demonstrate using a bias potential that can be used to study solid-liquid transitions +with the interface pinning method. This is done by adding a harmonic potential to the Hamiltonian +that bias the system towards two-phase configurations. - U_bias = 0.5*k*(Q-a)^2 + U_bias = 0.5*K*(Q-a)^2 -The bias field couple to an order-parameter of crystallinity Q -This implimentation use long-range order: Q=|rho_k|, where rho_k is the collective density field of the wave-vector k. +The bias field couple to an order-parameter of crystallinity Q. The implementation use long-range order: + + Q=|rho_k|, + +where rho_k is the collective density field of the wave-vector k. +For future reference we note that the structure factor S(k) is given by the variance of the collective density field: + + S(k)=|rho_k|^2. # Reference -[Ulf R. Pedersen, J. Chem. Phys. 139, 104102 (2013)] -Please visit +It is recommended to get familiar with the interface pinning method by reading: + + [Ulf R. Pedersen, J. Chem. Phys. 139, 104102 (2013)] + +A detailed bibliography is provided at + urp.dk/interface_pinning.htm -for a detailed bibliography. -# Use +# Use of rhok fix - fix [name] [groupID] rhok [nx] [ny] [nz] [spring-constant] [anchor-point] +For this example we will be using the rhok fix. -include a harmonic bias potential U_bias=0.5*k*(|rho_k|-a)^2 to the force calculation. -The elements of the wave-vector rho_k is k_x = (2 pi / L_x) * n_x, k_y = (2 pi / L_y) * n_y and k_z = (2 pi / L_z) * n_z. + fix [name] [groupID] rhok [nx] [ny] [nz] [K] [a] -# The Interface Pinning method for studying melting transitions -We will use the interface pinning method to study melting of the Lennard-Jones system +This fix include a harmonic bias potential U_bias=0.5*K*(|rho_k|-a)^2 to the force calculation. +The elements of the wave-vector k is given by the nx, ny and nz input: + + k_x = (2 pi / L_x) * n_x, k_y = (2 pi / L_y) * n_y and k_z = (2 pi / L_z) * n_z. + +We will use a k vector that correspond to a Bragg peak. + +# The Interface Pinning method for studying melting transitions of the Lennard-Jones (LJ) system + +We will use the interface pinning method to study melting of the LJ system at temperature 0.8 and pressure 2.185. This is a coexistence state-point, and the method -can be used to show this. The present directory contains the input files: +can be used to show this. The present directory contains the input files that we will use: - crystal.lmp - setup.lmp - pinning.lmp + in.crystal + in.setup + in.pinning -1. First we will determine the density of the crystal with the LAMMPS input file crystal.lmp. - From the output we get that the average density after equbriliation is 0.9731. - We need this density to ensure hydrostatic pressure when in a two-phase simulation. +1. First we will determine the density of the crystal with the LAMMPS input file in.crystal. + From the output we get that the average density after equilibration is 0.9731. + We need this density to ensure hydrostatic pressure when in a two-phase simulation. -2. Next, we setup a two-phase configuration using setup.lmp. +2. Next, we setup a two-phase configuration using in.setup. -3. Finally, we run a two-phase simulation with the bias-field applied using pinning.lmp. - The last coulmn in the output show |rho_k|. We note that after a equbriliation period - the value fluctuates aroung the anchor point (a) -- showing that this is indeed a coexistence - state point. +3. Finally, we run a two-phase simulation with the bias-field applied using in.pinning. + The last column in the output show |rho_k|. We note that after a equilibration period + the value fluctuates around the anchor point (a) -- showing that this is indeed a coexistence + state point. -The reference [J. Chem. Phys. 139, 104102 (2013)] gives details on using the method to find coexitence state points, -and the referecee [J. Chem. Phys. 142, 044104 (2015)] show how the crystal growth rate can be computed. +The reference [J. Chem. Phys. 139, 104102 (2013)] gives details on using the method to find coexistence state points, +and the reference [J. Chem. Phys. 142, 044104 (2015)] show how the crystal growth rate can be computed from fluctuations. That method have been experienced to be most effective in the slightly super-heated regime above the melting temperature. # Contact + Ulf R. Pedersen http://www.urp.dk ulf AT urp.dk From e44f370d4949f196c6abf33d35a3138410540929 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 11:44:52 +0200 Subject: [PATCH 33/44] Update readme.md --- examples/USER/pinning/readme.md | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/examples/USER/pinning/readme.md b/examples/USER/pinning/readme.md index 2071b63b4a..8a804a5127 100644 --- a/examples/USER/pinning/readme.md +++ b/examples/USER/pinning/readme.md @@ -13,7 +13,7 @@ For future reference we note that the structure factor S(k) is given by the vari S(k)=|rho_k|^2. -# Reference +## Reference It is recommended to get familiar with the interface pinning method by reading: @@ -23,7 +23,7 @@ A detailed bibliography is provided at urp.dk/interface_pinning.htm -# Use of rhok fix +## Use of rhok fix For this example we will be using the rhok fix. From 250ef9f837bece6c60c28bf2a3ccd57ef5815665 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 11:46:08 +0200 Subject: [PATCH 34/44] Update readme.md --- examples/USER/pinning/readme.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/examples/USER/pinning/readme.md b/examples/USER/pinning/readme.md index 8a804a5127..c7d5d792c3 100644 --- a/examples/USER/pinning/readme.md +++ b/examples/USER/pinning/readme.md @@ -21,7 +21,7 @@ It is recommended to get familiar with the interface pinning method by reading: A detailed bibliography is provided at - urp.dk/interface_pinning.htm + ## Use of rhok fix From 5c59eb637bca8e9bacf8753840a01af51283e2a2 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 11:49:57 +0200 Subject: [PATCH 35/44] Update readme.md --- examples/USER/pinning/readme.md | 12 +++++++----- 1 file changed, 7 insertions(+), 5 deletions(-) diff --git a/examples/USER/pinning/readme.md b/examples/USER/pinning/readme.md index c7d5d792c3..6e53e62d2f 100644 --- a/examples/USER/pinning/readme.md +++ b/examples/USER/pinning/readme.md @@ -1,3 +1,5 @@ +# The Interface Pinning method for studying solid-liquid transitions + This example demonstrate using a bias potential that can be used to study solid-liquid transitions with the interface pinning method. This is done by adding a harmonic potential to the Hamiltonian that bias the system towards two-phase configurations. @@ -13,7 +15,7 @@ For future reference we note that the structure factor S(k) is given by the vari S(k)=|rho_k|^2. -## Reference +### Reference It is recommended to get familiar with the interface pinning method by reading: @@ -23,7 +25,7 @@ A detailed bibliography is provided at -## Use of rhok fix +### Use of rhok fix For this example we will be using the rhok fix. @@ -36,9 +38,9 @@ The elements of the wave-vector k is given by the nx, ny and nz input: We will use a k vector that correspond to a Bragg peak. -# The Interface Pinning method for studying melting transitions of the Lennard-Jones (LJ) system +## Example: the Lennard-Jones (LJ) model -We will use the interface pinning method to study melting of the LJ system +We will use the interface pinning method to study melting of the LJ model at temperature 0.8 and pressure 2.185. This is a coexistence state-point, and the method can be used to show this. The present directory contains the input files that we will use: @@ -61,7 +63,7 @@ The reference [J. Chem. Phys. 139, 104102 (2013)] gives details on using the met and the reference [J. Chem. Phys. 142, 044104 (2015)] show how the crystal growth rate can be computed from fluctuations. That method have been experienced to be most effective in the slightly super-heated regime above the melting temperature. -# Contact +## Contact Ulf R. Pedersen http://www.urp.dk From fd8f5f8f9e486b064c0f552e9466e41f0c039ec7 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 11:52:08 +0200 Subject: [PATCH 36/44] Update readme.md --- examples/USER/pinning/readme.md | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/examples/USER/pinning/readme.md b/examples/USER/pinning/readme.md index 6e53e62d2f..9df85432ad 100644 --- a/examples/USER/pinning/readme.md +++ b/examples/USER/pinning/readme.md @@ -19,7 +19,7 @@ For future reference we note that the structure factor S(k) is given by the vari It is recommended to get familiar with the interface pinning method by reading: - [Ulf R. Pedersen, J. Chem. Phys. 139, 104102 (2013)] + [Ulf R. Pedersen, JCP 139, 104102 (2013)](http://dx.doi.org/10.1063/1.4818747) A detailed bibliography is provided at @@ -59,8 +59,8 @@ can be used to show this. The present directory contains the input files that we the value fluctuates around the anchor point (a) -- showing that this is indeed a coexistence state point. -The reference [J. Chem. Phys. 139, 104102 (2013)] gives details on using the method to find coexistence state points, -and the reference [J. Chem. Phys. 142, 044104 (2015)] show how the crystal growth rate can be computed from fluctuations. +The reference [JCP 139, 104102 (2013)](http://dx.doi.org/10.1063/1.4818747) gives details on using the method to find coexistence state points, +and the reference [JCP 142, 044104 (2015)](http://dx.doi.org/10.1063/1.4818747) show how the crystal growth rate can be computed from fluctuations. That method have been experienced to be most effective in the slightly super-heated regime above the melting temperature. ## Contact From e9b07a7a109260575a81cdc2f1344c5c15ccfe20 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 11:52:48 +0200 Subject: [PATCH 37/44] Update readme.md --- examples/USER/pinning/readme.md | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/examples/USER/pinning/readme.md b/examples/USER/pinning/readme.md index 9df85432ad..cd9a33167d 100644 --- a/examples/USER/pinning/readme.md +++ b/examples/USER/pinning/readme.md @@ -65,6 +65,6 @@ That method have been experienced to be most effective in the slightly super-hea ## Contact - Ulf R. Pedersen - http://www.urp.dk + Ulf R. Pedersen; + ; ulf AT urp.dk From 74dcf0bf9bfccce311466cb454658b8eeaf4ebe0 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 11:54:46 +0200 Subject: [PATCH 38/44] Update readme.md --- examples/USER/pinning/readme.md | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/examples/USER/pinning/readme.md b/examples/USER/pinning/readme.md index cd9a33167d..4a652cb1c3 100644 --- a/examples/USER/pinning/readme.md +++ b/examples/USER/pinning/readme.md @@ -1,8 +1,8 @@ # The Interface Pinning method for studying solid-liquid transitions -This example demonstrate using a bias potential that can be used to study solid-liquid transitions -with the interface pinning method. This is done by adding a harmonic potential to the Hamiltonian -that bias the system towards two-phase configurations. +In this example we will use the interface pinnig method to study a solid-liquid transition. +This is done by adding a harmonic potential to the Hamiltonian +that bias the system towards two-phase configurations: U_bias = 0.5*K*(Q-a)^2 From 30aaa7e47bb58ef6ed3defd599b779c448510bcc Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 12:00:43 +0200 Subject: [PATCH 39/44] Update readme.md --- examples/USER/pinning/readme.md | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/examples/USER/pinning/readme.md b/examples/USER/pinning/readme.md index 4a652cb1c3..5e9392ba47 100644 --- a/examples/USER/pinning/readme.md +++ b/examples/USER/pinning/readme.md @@ -25,6 +25,10 @@ A detailed bibliography is provided at +and a brief introduction can be found at YouTube: + + [![Interface Pinning](http://img.youtube.com/vi/F_79JZNdyoQ/0.jpg)](http://www.youtube.com/watch?v=F_79JZNdyoQ) + ### Use of rhok fix For this example we will be using the rhok fix. From 5f527091b8a6b157f04cd16aba75ba52ff045eb6 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Tue, 3 Oct 2017 12:02:01 +0200 Subject: [PATCH 40/44] Update readme.md --- examples/USER/pinning/readme.md | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/examples/USER/pinning/readme.md b/examples/USER/pinning/readme.md index 5e9392ba47..4e011255fc 100644 --- a/examples/USER/pinning/readme.md +++ b/examples/USER/pinning/readme.md @@ -15,7 +15,7 @@ For future reference we note that the structure factor S(k) is given by the vari S(k)=|rho_k|^2. -### Reference +### About the method It is recommended to get familiar with the interface pinning method by reading: @@ -29,7 +29,7 @@ and a brief introduction can be found at YouTube: [![Interface Pinning](http://img.youtube.com/vi/F_79JZNdyoQ/0.jpg)](http://www.youtube.com/watch?v=F_79JZNdyoQ) -### Use of rhok fix +### Implimentation in LAMMPS For this example we will be using the rhok fix. From f07719e924b22da61dc7bffae1b5ec286a3cf980 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Tue, 3 Oct 2017 10:08:38 -0400 Subject: [PATCH 41/44] make fix rhok examples more consistent with LAMMPS conventions: - move example folder to examples/USER/misc/ - comment out writing of trajectory files - reduce run length (for easier testing for regressions) - record example outputs for 1 and 4 MPI processes - rename readme.md to README.md for visibility --- .../readme.md => misc/pinning/README.md} | 0 examples/USER/{ => misc}/pinning/in.crystal | 7 +- examples/USER/{ => misc}/pinning/in.pinning | 7 +- examples/USER/{ => misc}/pinning/in.setup | 7 +- .../misc/pinning/log.22Sep2017.crystal.g++.1 | 187 ++++++++++++++++++ .../misc/pinning/log.22Sep2017.crystal.g++.4 | 187 ++++++++++++++++++ .../misc/pinning/log.22Sep2017.pinning.g++.1 | 186 +++++++++++++++++ .../misc/pinning/log.22Sep2017.pinning.g++.4 | 186 +++++++++++++++++ .../misc/pinning/log.22Sep2017.setup.g++.1 | 141 +++++++++++++ .../misc/pinning/log.22Sep2017.setup.g++.4 | 141 +++++++++++++ 10 files changed, 1040 insertions(+), 9 deletions(-) rename examples/USER/{pinning/readme.md => misc/pinning/README.md} (100%) rename examples/USER/{ => misc}/pinning/in.crystal (85%) rename examples/USER/{ => misc}/pinning/in.pinning (87%) rename examples/USER/{ => misc}/pinning/in.setup (89%) create mode 100644 examples/USER/misc/pinning/log.22Sep2017.crystal.g++.1 create mode 100644 examples/USER/misc/pinning/log.22Sep2017.crystal.g++.4 create mode 100644 examples/USER/misc/pinning/log.22Sep2017.pinning.g++.1 create mode 100644 examples/USER/misc/pinning/log.22Sep2017.pinning.g++.4 create mode 100644 examples/USER/misc/pinning/log.22Sep2017.setup.g++.1 create mode 100644 examples/USER/misc/pinning/log.22Sep2017.setup.g++.4 diff --git a/examples/USER/pinning/readme.md b/examples/USER/misc/pinning/README.md similarity index 100% rename from examples/USER/pinning/readme.md rename to examples/USER/misc/pinning/README.md diff --git a/examples/USER/pinning/in.crystal b/examples/USER/misc/pinning/in.crystal similarity index 85% rename from examples/USER/pinning/in.crystal rename to examples/USER/misc/pinning/in.crystal index c40e206d46..55e9e59a06 100644 --- a/examples/USER/pinning/in.crystal +++ b/examples/USER/misc/pinning/in.crystal @@ -29,7 +29,8 @@ fix bias all rhok 16 0 0 0.0 0.0 # output thermo 50 thermo_style custom step temp press density f_bias[3] -dump dumpXYZ all xyz 2000 traj.xyz +# dump dumpXYZ all xyz 2000 traj.xyz -# run -run 100000 +# NOTE: this is cut short to 5000 steps for demonstration purposes +# run 100000 +run 5000 diff --git a/examples/USER/pinning/in.pinning b/examples/USER/misc/pinning/in.pinning similarity index 87% rename from examples/USER/pinning/in.pinning rename to examples/USER/misc/pinning/in.pinning index 62e707a1fd..607335cf1a 100644 --- a/examples/USER/pinning/in.pinning +++ b/examples/USER/misc/pinning/in.pinning @@ -26,7 +26,8 @@ fix bias all rhok 16 0 0 4.0 26.00 # output U_bias rho_k_RE rho_k_IM |rho_k| thermo_style custom step temp pzz pe lz f_bias f_bias[1] f_bias[2] f_bias[3] thermo 50 -dump dumpXYZ all xyz 500 traj.xyz +# dump dumpXYZ all xyz 500 traj.xyz -# run -run 50000 +# NOTE: run reduced for demonstration purposes +# run 50000 +run 5000 diff --git a/examples/USER/pinning/in.setup b/examples/USER/misc/pinning/in.setup similarity index 89% rename from examples/USER/pinning/in.setup rename to examples/USER/misc/pinning/in.setup index 0febd3c4f6..649b0f534c 100644 --- a/examples/USER/pinning/in.setup +++ b/examples/USER/misc/pinning/in.setup @@ -33,8 +33,9 @@ fix langevin left langevin 3.0 0.8 100.0 2017 # outout thermo_style custom step temp pzz pe lz thermo 100 -dump dumpXYZ all xyz 100 traj.xyz +# dump dumpXYZ all xyz 100 traj.xyz -# run -run 10000 +# run reduced for demonstration purposes +# run 10000 +run 5000 write_data data.halfhalf diff --git a/examples/USER/misc/pinning/log.22Sep2017.crystal.g++.1 b/examples/USER/misc/pinning/log.22Sep2017.crystal.g++.1 new file mode 100644 index 0000000000..05fadb5c03 --- /dev/null +++ b/examples/USER/misc/pinning/log.22Sep2017.crystal.g++.1 @@ -0,0 +1,187 @@ +LAMMPS (22 Sep 2017) + using 1 OpenMP thread(s) per MPI task +units lj +dimension 3 +boundary p p p +atom_style atomic + +# truncated and shifted LJ potential +pair_style lj/cut 2.5 +pair_modify shift yes +lattice fcc 0.9731 +Lattice spacing in x,y,z = 1.6019 1.6019 1.6019 +region my_box block 0 8.0 0 8.0 0 20.0 +create_box 1 my_box +Created orthogonal box = (0 0 0) to (12.8152 12.8152 32.0379) + 1 by 1 by 1 MPI processor grid +region particles block 0 8.0 0 8.0 0 20.0 +create_atoms 1 region particles +Created 5120 atoms +pair_coeff 1 1 1.0 1.0 2.5 +pair_modify tail no +pair_modify shift yes +mass 1 1.0 +velocity all create 1.6 1 mom yes rot yes + +# simulation parameters +neighbor 0.6 bin +timestep 0.004 +run_style verlet +fix ensemble all npt temp 0.8 0.8 4.0 aniso 2.185 2.185 8.0 pchain 32 + +# computing long-range order (no bias is added since k=0) +fix bias all rhok 16 0 0 0.0 0.0 + +# output +thermo 50 +thermo_style custom step temp press density f_bias[3] +# dump dumpXYZ all xyz 2000 traj.xyz + +# NOTE: this is cut short to 5000 steps for demonstration purposes +# run 100000 +run 5000 +Neighbor list info ... + update every 1 steps, delay 10 steps, check yes + max neighbors/atom: 2000, page size: 100000 + master list distance cutoff = 3.1 + ghost atom cutoff = 3.1 + binsize = 1.55, bins = 9 9 21 + 1 neighbor lists, perpetual/occasional/extra = 1 0 0 + (1) pair lj/cut, perpetual + attributes: half, newton on + pair build: half/bin/atomonly/newton + stencil: half/bin/3d/newton + bin: standard +Per MPI rank memory allocation (min/avg/max) = 4.523 | 4.523 | 4.523 Mbytes +Step Temp Press Density f_bias[3] + 0 1.6 -2.7568106 0.9731 71.554175 + 50 0.78457786 3.1029192 0.97362639 54.327705 + 100 0.85528971 2.4670259 0.97213457 55.189308 + 150 0.85241818 2.3210306 0.9698027 56.138125 + 200 0.82301385 2.3448692 0.96708227 55.735326 + 250 0.83076383 2.0890816 0.96425763 55.320625 + 300 0.81602823 2.0118796 0.96173925 54.095736 + 350 0.81084006 1.9122192 0.95979392 54.526429 + 400 0.80776593 1.8502174 0.95869117 54.434901 + 450 0.80694697 1.8435873 0.95851085 53.20809 + 500 0.81384248 1.8111331 0.95917305 53.419395 + 550 0.81027072 1.9222272 0.96056019 54.36723 + 600 0.81199582 2.0291945 0.96248486 54.888582 + 650 0.82507964 2.0706462 0.96467227 55.807137 + 700 0.832562 2.1471442 0.9668913 56.721267 + 750 0.83358138 2.2674672 0.968984 56.723838 + 800 0.83477542 2.3658275 0.97072603 56.234689 + 850 0.84722921 2.3506233 0.97189674 56.262424 + 900 0.83526965 2.4532068 0.97248856 56.219103 + 950 0.83174583 2.4763958 0.97249527 56.409813 + 1000 0.83022557 2.4334341 0.97194093 55.890858 + 1050 0.83208978 2.3478416 0.97092452 54.934691 + 1100 0.82789545 2.272404 0.9696152 54.90894 + 1150 0.82678617 2.1798046 0.96819776 54.927782 + 1200 0.8088841 2.1960256 0.96687735 54.914327 + 1250 0.81512784 2.0736261 0.96579008 53.927291 + 1300 0.81271067 2.0297138 0.96504188 54.289698 + 1350 0.8201767 1.9493976 0.96464115 55.342131 + 1400 0.80880489 2.0016987 0.96468463 55.757758 + 1450 0.8114196 2.0282699 0.96514115 55.865676 + 1500 0.81085664 2.0838361 0.96591869 56.553425 + 1550 0.81257075 2.1283157 0.96694549 56.921544 + 1600 0.82617645 2.1017986 0.96817075 56.858808 + 1650 0.82616141 2.1885582 0.96941073 56.717917 + 1700 0.81634174 2.2996967 0.97047447 56.453745 + 1750 0.82447573 2.2924266 0.97128663 56.916813 + 1800 0.83610432 2.236456 0.97178453 56.400752 + 1850 0.82479203 2.3103493 0.97197318 55.891368 + 1900 0.82298992 2.3059289 0.97181084 55.680563 + 1950 0.82098556 2.2801003 0.97138609 55.754406 + 2000 0.8181203 2.2480175 0.97078591 55.801363 + 2050 0.82822293 2.1208884 0.97004107 55.687 + 2100 0.7976818 2.2711199 0.96930169 55.459844 + 2150 0.81817848 2.0680351 0.96860201 56.514731 + 2200 0.80707457 2.1112141 0.96810519 55.504308 + 2250 0.81651111 2.0077603 0.96781161 55.635702 + 2300 0.80634534 2.0662241 0.96777177 56.051086 + 2350 0.80892831 2.0619333 0.96799037 56.548711 + 2400 0.82454203 1.9585394 0.9684672 56.695235 + 2450 0.81517178 2.075283 0.96921622 56.613082 + 2500 0.80969595 2.1624581 0.97010528 56.57516 + 2550 0.80862964 2.2088622 0.97100774 57.072594 + 2600 0.81468816 2.2293973 0.97192868 56.879212 + 2650 0.82063107 2.2244887 0.97269715 55.454502 + 2700 0.81691618 2.2789954 0.97319841 54.421943 + 2750 0.8141787 2.2981247 0.97340453 54.469921 + 2800 0.81973871 2.2422136 0.9733278 55.959235 + 2850 0.82037399 2.201016 0.97302727 56.685826 + 2900 0.80650164 2.2672955 0.9726128 56.574395 + 2950 0.81752783 2.1317541 0.97207545 56.809412 + 3000 0.80836945 2.1461483 0.97151192 57.205206 + 3050 0.80785109 2.1189056 0.97103049 57.418763 + 3100 0.79835058 2.146416 0.97069705 57.329383 + 3150 0.79792089 2.1388267 0.97051679 57.279852 + 3200 0.79934603 2.1049562 0.97046851 56.351494 + 3250 0.79523232 2.1549779 0.97063956 56.00356 + 3300 0.8004458 2.1145975 0.97096375 55.725509 + 3350 0.79772742 2.166292 0.97143785 55.558075 + 3400 0.80621087 2.1309217 0.97198456 55.816704 + 3450 0.80540626 2.1727557 0.97263267 55.671283 + 3500 0.80867606 2.1905129 0.97321538 55.390086 + 3550 0.80917896 2.2144872 0.97370472 55.742085 + 3600 0.80930722 2.2288938 0.974093 56.23064 + 3650 0.80390523 2.2777327 0.97431886 56.084731 + 3700 0.79620093 2.3143541 0.97435103 55.942797 + 3750 0.80252393 2.2564638 0.9741875 56.042055 + 3800 0.78981264 2.3156481 0.9739121 55.971352 + 3850 0.80391951 2.1804938 0.97351088 55.855858 + 3900 0.81268129 2.0855818 0.97308521 56.288315 + 3950 0.7958182 2.175259 0.97273088 56.140141 + 4000 0.80054484 2.1163279 0.97243129 56.366818 + 4050 0.79760187 2.105362 0.97225308 56.684619 + 4100 0.79283424 2.1357603 0.972206 56.203341 + 4150 0.79543088 2.1036951 0.97227608 56.606315 + 4200 0.79410999 2.1402049 0.97253758 56.277478 + 4250 0.7985469 2.1285154 0.97293622 56.356076 + 4300 0.79700387 2.1470614 0.97337091 56.722298 + 4350 0.80479321 2.1403244 0.97384674 57.212574 + 4400 0.79505512 2.224463 0.97434415 56.561877 + 4450 0.78346648 2.3347865 0.97478611 56.681362 + 4500 0.79811284 2.259123 0.97510069 57.365929 + 4550 0.80015561 2.2345254 0.97523653 57.34799 + 4600 0.79648318 2.2651869 0.97525975 57.502318 + 4650 0.80524865 2.1943025 0.97507638 57.702488 + 4700 0.80397778 2.1758629 0.97478268 57.162107 + 4750 0.78914913 2.2470191 0.9744625 56.849565 + 4800 0.79324889 2.2028993 0.97408817 57.572344 + 4850 0.78993209 2.181763 0.97373372 57.683552 + 4900 0.79041263 2.1604768 0.97348692 56.922312 + 4950 0.79741332 2.1105901 0.97332545 57.488932 + 5000 0.7891178 2.163416 0.97328963 57.365252 +Loop time of 33.6467 on 1 procs for 5000 steps with 5120 atoms + +Performance: 51357.258 tau/day, 148.603 timesteps/s +99.7% CPU use with 1 MPI tasks x 1 OpenMP threads + +MPI task timing breakdown: +Section | min time | avg time | max time |%varavg| %total +--------------------------------------------------------------- +Pair | 24.699 | 24.699 | 24.699 | 0.0 | 73.41 +Neigh | 2.8894 | 2.8894 | 2.8894 | 0.0 | 8.59 +Comm | 0.34907 | 0.34907 | 0.34907 | 0.0 | 1.04 +Output | 0.0056 | 0.0056 | 0.0056 | 0.0 | 0.02 +Modify | 5.5718 | 5.5718 | 5.5718 | 0.0 | 16.56 +Other | | 0.1319 | | | 0.39 + +Nlocal: 5120 ave 5120 max 5120 min +Histogram: 1 0 0 0 0 0 0 0 0 0 +Nghost: 7594 ave 7594 max 7594 min +Histogram: 1 0 0 0 0 0 0 0 0 0 +Neighs: 323081 ave 323081 max 323081 min +Histogram: 1 0 0 0 0 0 0 0 0 0 + +Total # of neighbors = 323081 +Ave neighs/atom = 63.1018 +Neighbor list builds = 248 +Dangerous builds = 0 + +Please see the log.cite file for references relevant to this simulation + +Total wall time: 0:00:33 diff --git a/examples/USER/misc/pinning/log.22Sep2017.crystal.g++.4 b/examples/USER/misc/pinning/log.22Sep2017.crystal.g++.4 new file mode 100644 index 0000000000..cec9c69aff --- /dev/null +++ b/examples/USER/misc/pinning/log.22Sep2017.crystal.g++.4 @@ -0,0 +1,187 @@ +LAMMPS (22 Sep 2017) + using 1 OpenMP thread(s) per MPI task +units lj +dimension 3 +boundary p p p +atom_style atomic + +# truncated and shifted LJ potential +pair_style lj/cut 2.5 +pair_modify shift yes +lattice fcc 0.9731 +Lattice spacing in x,y,z = 1.6019 1.6019 1.6019 +region my_box block 0 8.0 0 8.0 0 20.0 +create_box 1 my_box +Created orthogonal box = (0 0 0) to (12.8152 12.8152 32.0379) + 1 by 1 by 4 MPI processor grid +region particles block 0 8.0 0 8.0 0 20.0 +create_atoms 1 region particles +Created 5120 atoms +pair_coeff 1 1 1.0 1.0 2.5 +pair_modify tail no +pair_modify shift yes +mass 1 1.0 +velocity all create 1.6 1 mom yes rot yes + +# simulation parameters +neighbor 0.6 bin +timestep 0.004 +run_style verlet +fix ensemble all npt temp 0.8 0.8 4.0 aniso 2.185 2.185 8.0 pchain 32 + +# computing long-range order (no bias is added since k=0) +fix bias all rhok 16 0 0 0.0 0.0 + +# output +thermo 50 +thermo_style custom step temp press density f_bias[3] +# dump dumpXYZ all xyz 2000 traj.xyz + +# NOTE: this is cut short to 5000 steps for demonstration purposes +# run 100000 +run 5000 +Neighbor list info ... + update every 1 steps, delay 10 steps, check yes + max neighbors/atom: 2000, page size: 100000 + master list distance cutoff = 3.1 + ghost atom cutoff = 3.1 + binsize = 1.55, bins = 9 9 21 + 1 neighbor lists, perpetual/occasional/extra = 1 0 0 + (1) pair lj/cut, perpetual + attributes: half, newton on + pair build: half/bin/atomonly/newton + stencil: half/bin/3d/newton + bin: standard +Per MPI rank memory allocation (min/avg/max) = 3.23 | 3.23 | 3.23 Mbytes +Step Temp Press Density f_bias[3] + 0 1.6 -2.7568106 0.9731 71.554175 + 50 0.78457786 3.1029192 0.97362639 54.327705 + 100 0.85528971 2.4670259 0.97213457 55.189308 + 150 0.85241818 2.3210306 0.9698027 56.138125 + 200 0.82301385 2.3448692 0.96708227 55.735326 + 250 0.83076383 2.0890816 0.96425763 55.320625 + 300 0.81602823 2.0118796 0.96173925 54.095736 + 350 0.81084006 1.9122192 0.95979392 54.526429 + 400 0.80776593 1.8502174 0.95869117 54.434901 + 450 0.80694697 1.8435873 0.95851085 53.20809 + 500 0.81384248 1.8111331 0.95917305 53.419395 + 550 0.81027072 1.9222272 0.96056019 54.36723 + 600 0.81199582 2.0291945 0.96248486 54.888582 + 650 0.82507964 2.0706462 0.96467227 55.807137 + 700 0.832562 2.1471442 0.9668913 56.721267 + 750 0.83358138 2.2674672 0.968984 56.723838 + 800 0.83477542 2.3658275 0.97072603 56.234689 + 850 0.84722921 2.3506233 0.97189674 56.262424 + 900 0.83526965 2.4532068 0.97248856 56.219103 + 950 0.83174583 2.4763958 0.97249527 56.409813 + 1000 0.83022557 2.4334341 0.97194093 55.890858 + 1050 0.83208978 2.3478416 0.97092452 54.934691 + 1100 0.82789545 2.272404 0.9696152 54.90894 + 1150 0.82678617 2.1798046 0.96819776 54.927782 + 1200 0.8088841 2.1960256 0.96687735 54.914327 + 1250 0.81512784 2.0736261 0.96579008 53.927291 + 1300 0.81271067 2.0297138 0.96504188 54.289698 + 1350 0.8201767 1.9493976 0.96464115 55.342131 + 1400 0.80880489 2.0016987 0.96468463 55.757758 + 1450 0.8114196 2.0282699 0.96514115 55.865676 + 1500 0.81085664 2.0838361 0.96591869 56.553425 + 1550 0.81257075 2.1283157 0.96694549 56.921544 + 1600 0.82617645 2.1017986 0.96817075 56.858808 + 1650 0.82616141 2.1885582 0.96941073 56.717917 + 1700 0.81634174 2.2996967 0.97047447 56.453745 + 1750 0.82447573 2.2924266 0.97128663 56.916813 + 1800 0.83610432 2.236456 0.97178453 56.400752 + 1850 0.824792 2.3103491 0.97197318 55.891368 + 1900 0.82298989 2.3059287 0.97181084 55.680562 + 1950 0.82098545 2.2801009 0.97138609 55.754404 + 2000 0.81812031 2.2480166 0.97078591 55.801371 + 2050 0.82822262 2.1208887 0.97004108 55.687001 + 2100 0.79768162 2.2711186 0.9693017 55.459852 + 2150 0.81817874 2.0680317 0.96860202 56.514744 + 2200 0.80707412 2.1112032 0.96810521 55.504308 + 2250 0.81650921 2.0077757 0.96781164 55.635717 + 2300 0.80634656 2.066186 0.96777181 56.051088 + 2350 0.80893174 2.0619084 0.96799042 56.548711 + 2400 0.82453783 1.9585503 0.96846727 56.695111 + 2450 0.81517275 2.0752617 0.96921631 56.614046 + 2500 0.80969622 2.1624476 0.9701054 56.574846 + 2550 0.80861922 2.2089505 0.97100787 57.072334 + 2600 0.81468888 2.2293754 0.97192875 56.879416 + 2650 0.82061239 2.2245462 0.97269723 55.442015 + 2700 0.81687473 2.2792015 0.97319852 54.420301 + 2750 0.81416567 2.2982988 0.97340467 54.469427 + 2800 0.81978563 2.2418723 0.97332803 55.965451 + 2850 0.82069759 2.1988948 0.97302752 56.686807 + 2900 0.80631184 2.2684466 0.97261407 56.585682 + 2950 0.81759744 2.1312328 0.97207888 56.812431 + 3000 0.80748056 2.152676 0.97151807 57.178849 + 3050 0.80789237 2.118162 0.97103728 57.433724 + 3100 0.79882523 2.1414744 0.97070338 57.34686 + 3150 0.79803949 2.1359043 0.97052875 57.382544 + 3200 0.79170386 2.1548392 0.97049349 56.465806 + 3250 0.78848813 2.1990144 0.97067557 55.929088 + 3300 0.79820555 2.1304609 0.97101444 55.624487 + 3350 0.79250565 2.1971235 0.97149233 55.933615 + 3400 0.80584844 2.1417239 0.97206083 55.85922 + 3450 0.80685744 2.1640501 0.97266047 55.135963 + 3500 0.80751888 2.1858277 0.97318703 55.407581 + 3550 0.79882754 2.2796452 0.97363149 55.392366 + 3600 0.80219171 2.2715765 0.97392571 55.867887 + 3650 0.79061794 2.3492866 0.97410985 56.0192 + 3700 0.8058483 2.2327904 0.97411924 56.491303 + 3750 0.79460746 2.2941868 0.97397764 55.929912 + 3800 0.80447478 2.2018009 0.97367627 55.663208 + 3850 0.80355335 2.17638 0.97333164 55.637261 + 3900 0.80388417 2.1531434 0.9729647 56.03794 + 3950 0.79557409 2.1853318 0.9726503 56.132348 + 4000 0.79547396 2.1457051 0.97235244 55.552675 + 4050 0.8058384 2.0637678 0.97213346 56.185416 + 4100 0.7976931 2.1028246 0.97208255 56.050347 + 4150 0.79555522 2.115473 0.97216375 56.868136 + 4200 0.79324134 2.1510383 0.97246129 56.462635 + 4250 0.80788167 2.0534887 0.97287821 55.650788 + 4300 0.79389865 2.2019815 0.97337765 55.596846 + 4350 0.79786309 2.1851119 0.97389825 57.000921 + 4400 0.79986518 2.1997541 0.97443778 57.551564 + 4450 0.8063901 2.1893874 0.97493151 57.236138 + 4500 0.80005802 2.250364 0.97533075 57.341358 + 4550 0.79707443 2.2995576 0.97557554 57.338713 + 4600 0.79869949 2.2807889 0.97563277 57.084504 + 4650 0.79694427 2.2673215 0.97544638 57.025663 + 4700 0.79023986 2.2884131 0.97511483 57.131188 + 4750 0.79566823 2.2215519 0.97464304 57.045676 + 4800 0.78936986 2.2268037 0.97410626 57.384178 + 4850 0.79025913 2.1836718 0.973616 57.78438 + 4900 0.80138424 2.0657609 0.9732124 57.888266 + 4950 0.77853735 2.207944 0.97296347 57.312213 + 5000 0.79115984 2.1035893 0.97285578 57.109472 +Loop time of 9.53489 on 4 procs for 5000 steps with 5120 atoms + +Performance: 181229.223 tau/day, 524.390 timesteps/s +99.1% CPU use with 4 MPI tasks x 1 OpenMP threads + +MPI task timing breakdown: +Section | min time | avg time | max time |%varavg| %total +--------------------------------------------------------------- +Pair | 6.312 | 6.4238 | 6.5139 | 3.1 | 67.37 +Neigh | 0.72062 | 0.73538 | 0.74531 | 1.1 | 7.71 +Comm | 0.52697 | 0.64152 | 0.78688 | 14.1 | 6.73 +Output | 0.0028393 | 0.0029888 | 0.0033851 | 0.4 | 0.03 +Modify | 1.6249 | 1.669 | 1.7253 | 2.9 | 17.50 +Other | | 0.06221 | | | 0.65 + +Nlocal: 1280 ave 1289 max 1266 min +Histogram: 1 0 0 0 0 1 0 0 0 2 +Nghost: 3346.25 ave 3379 max 3331 min +Histogram: 1 2 0 0 0 0 0 0 0 1 +Neighs: 80701.8 ave 81534 max 79755 min +Histogram: 1 0 1 0 0 0 0 0 1 1 + +Total # of neighbors = 322807 +Ave neighs/atom = 63.0482 +Neighbor list builds = 248 +Dangerous builds = 0 + +Please see the log.cite file for references relevant to this simulation + +Total wall time: 0:00:09 diff --git a/examples/USER/misc/pinning/log.22Sep2017.pinning.g++.1 b/examples/USER/misc/pinning/log.22Sep2017.pinning.g++.1 new file mode 100644 index 0000000000..c2aaa9a581 --- /dev/null +++ b/examples/USER/misc/pinning/log.22Sep2017.pinning.g++.1 @@ -0,0 +1,186 @@ +LAMMPS (22 Sep 2017) + using 1 OpenMP thread(s) per MPI task +units lj +dimension 3 +boundary p p p +atom_style atomic + +# truncated and shifted LJ potential +pair_style lj/cut 2.5 +pair_modify shift yes +read_data data.halfhalf + orthogonal box = (0 0 0) to (12.8152 12.8152 34) + 1 by 1 by 1 MPI processor grid + reading atoms ... + 5120 atoms + reading velocities ... + 5120 velocities +pair_coeff 1 1 1.0 1.0 2.5 +mass 1 1.0 + +# simulation parameters +neighbor 0.6 bin +timestep 0.004 +run_style verlet + +velocity all create 0.8 1 mom yes rot yes +fix ensemble all npt temp 0.8 0.8 4.0 z 2.185 2.185 8.0 +fix 100 all momentum 100 linear 1 1 1 + +# harmonic rho_k bias-field +# nx ny nz k a +fix bias all rhok 16 0 0 4.0 26.00 + +# output U_bias rho_k_RE rho_k_IM |rho_k| +thermo_style custom step temp pzz pe lz f_bias f_bias[1] f_bias[2] f_bias[3] +thermo 50 +dump dumpXYZ all xyz 500 traj.xyz + +# NOTE: run reduced for demonstration purposes +# run 50000 +run 5000 +Neighbor list info ... + update every 1 steps, delay 10 steps, check yes + max neighbors/atom: 2000, page size: 100000 + master list distance cutoff = 3.1 + ghost atom cutoff = 3.1 + binsize = 1.55, bins = 9 9 22 + 1 neighbor lists, perpetual/occasional/extra = 1 0 0 + (1) pair lj/cut, perpetual + attributes: half, newton on + pair build: half/bin/atomonly/newton + stencil: half/bin/3d/newton + bin: standard +Per MPI rank memory allocation (min/avg/max) = 5.723 | 5.723 | 5.723 Mbytes +Step Temp Pzz PotEng Lz f_bias f_bias[1] f_bias[2] f_bias[3] + 0 0.8 5.1566801 -4.8556711 34 179.52419 35.473155 -0.2832763 35.474286 + 50 1.072533 3.8158392 -5.2704532 34.024206 0.010596224 25.927135 -0.063106738 25.927212 + 100 1.1093231 3.6703116 -5.3380965 34.094814 1.8552612 26.958236 -0.51403326 26.963136 + 150 1.1080721 3.8202601 -5.3568368 34.207473 0.39188605 26.420755 -1.0759575 26.442655 + 200 1.1435287 3.3445987 -5.4365298 34.354119 3.0758718 27.239465 -0.19115251 27.240135 + 250 1.1203046 3.4669456 -5.4293867 34.511473 0.8543814 26.652785 -0.20818214 26.653598 + 300 1.1012709 3.4583154 -5.4281803 34.664509 2.4837156 27.097674 0.9518628 27.114387 + 350 1.0439632 3.8953869 -5.368619 34.810399 0.55385719 26.518391 0.64525272 26.52624 + 400 1.0083878 4.0523864 -5.3418278 34.957669 0.34806057 26.404011 0.83368604 26.417169 + 450 0.9675244 4.310087 -5.3089468 35.114208 0.7921285 26.607512 1.077889 26.629336 + 500 0.94605476 4.1050097 -5.3062273 35.284018 0.87757754 26.639125 1.1140858 26.662411 + 550 0.92662323 3.9299063 -5.3062927 35.458565 1.3746441 26.773494 1.7256603 26.829049 + 600 0.89723165 3.7683555 -5.289725 35.629881 0.46692943 26.372973 2.4135502 26.483182 + 650 0.90612566 3.1098837 -5.3267851 35.788537 0.032662126 25.918784 3.2982102 26.127793 + 700 0.9152508 2.6527976 -5.3597799 35.923343 0.014621588 25.834591 3.6093492 26.085503 + 750 0.90156356 2.3474851 -5.3545938 36.031813 0.75225637 26.307204 4.0247215 26.613293 + 800 0.89748513 1.9825103 -5.3610202 36.113888 0.33402511 26.261326 2.7858039 26.408672 + 850 0.89496343 1.8906342 -5.3673514 36.169424 0.85431557 26.534648 2.5150347 26.653573 + 900 0.89463983 1.5654217 -5.3753283 36.20181 1.5689239 26.764737 2.5474794 26.885699 + 950 0.88663832 1.4399476 -5.3703322 36.209971 0.044436903 25.818418 1.2963356 25.850941 + 1000 0.87407208 1.485718 -5.3572665 36.195386 1.4405611 26.828072 1.0520795 26.848693 + 1050 0.87580489 1.163155 -5.3647269 36.160279 0.15319559 26.234791 1.4845964 26.276763 + 1100 0.86978111 1.3743181 -5.3594907 36.104958 1.1313537 25.19895 1.5711793 25.247885 + 1150 0.86987861 1.3212927 -5.3628503 36.035486 0.039865678 25.841762 0.93898962 25.858816 + 1200 0.87142486 1.3293818 -5.3676854 35.954411 0.16827389 25.70952 -0.14639427 25.709936 + 1250 0.87582265 1.3203803 -5.3764058 35.86575 0.25946652 25.639682 0.082696867 25.639815 + 1300 0.87371627 1.4680294 -5.375151 35.772824 0.17697069 25.701417 0.2397926 25.702535 + 1350 0.88617453 1.5923057 -5.3954912 35.681046 0.00049155526 25.973634 -0.74521794 25.984323 + 1400 0.87809636 1.5821707 -5.3850722 35.594706 0.024050814 26.107395 -0.34393685 26.10966 + 1450 0.87912192 1.7820174 -5.3885842 35.514722 0.20999914 25.667238 -0.66933655 25.675964 + 1500 0.88293618 2.0295275 -5.3963602 35.443445 0.60232374 25.376395 -1.9501461 25.451218 + 1550 0.90012542 1.9476472 -5.4249456 35.382791 0.4488038 26.448928 -1.1452474 26.473711 + 1600 0.89155063 2.2462603 -5.4153432 35.332095 0.039621687 26.138157 -0.36825239 26.140751 + 1650 0.8942624 2.343747 -5.4233433 35.294954 0.0089980332 26.064277 0.38189192 26.067075 + 1700 0.90047841 2.451289 -5.4376312 35.27234 0.86985171 26.646438 0.83408084 26.659489 + 1750 0.87586052 2.6381221 -5.4067182 35.264564 6.346204 27.652722 2.6699692 27.78132 + 1800 0.87392582 2.6338176 -5.4109056 35.270073 0.046414129 26.016188 2.6651053 26.152339 + 1850 0.86540415 2.5434301 -5.4058587 35.285902 0.054615472 26.074279 2.1799787 26.165251 + 1900 0.87043082 2.5776772 -5.4216997 35.309062 0.68978148 26.38648 3.2614091 26.587274 + 1950 0.86281992 2.3107762 -5.4188978 35.338501 0.0072672577 25.736893 3.2375012 25.93972 + 2000 0.85905576 2.2894047 -5.4215995 35.36787 0.095633435 26.072085 2.7685848 26.21867 + 2050 0.85793751 2.2382039 -5.4279351 35.395213 0.13602344 25.598457 2.6881027 25.739209 + 2100 0.85585253 2.0765811 -5.4324511 35.418877 0.0059888115 25.754128 3.1436222 25.945279 + 2150 0.86701057 1.8449875 -5.4562208 35.436124 0.097328618 25.413697 4.3268293 25.7794 + 2200 0.85168154 1.9024923 -5.4395776 35.44246 0.20764576 25.094788 5.4406104 25.677784 + 2250 0.8429719 1.870335 -5.4320586 35.438363 0.34419961 24.998478 5.4475709 25.585151 + 2300 0.84176891 1.7100228 -5.4351472 35.422863 0.76036958 24.697018 5.8629967 25.383409 + 2350 0.84601588 1.8539039 -5.4456629 35.395979 0.38437531 25.647986 6.4163366 26.438392 + 2400 0.84637647 1.6299091 -5.4498948 35.36125 0.074236719 24.995872 7.8269968 26.192661 + 2450 0.85650449 1.6828907 -5.4683101 35.316669 0.3671827 25.280669 7.7040329 26.428476 + 2500 0.84963707 1.7305222 -5.4605394 35.265508 0.1406965 25.236741 7.2780025 26.265232 + 2550 0.84084365 1.8758368 -5.4497083 35.208725 0.33937687 24.544376 7.2334512 25.588067 + 2600 0.85317342 1.7781674 -5.4702734 35.149747 0.60378248 24.046307 8.3370138 25.450554 + 2650 0.85487644 2.0065374 -5.4747643 35.090431 0.22483651 24.937101 8.4669004 26.335288 + 2700 0.84550083 1.9363031 -5.4628401 35.034349 0.43442577 24.250196 7.9943738 25.533939 + 2750 0.85843419 2.0473138 -5.484528 34.980671 0.45959294 24.17438 8.179356 25.520629 + 2800 0.86047607 2.0754522 -5.4899966 34.932466 0.00038123477 24.619856 8.3153434 25.986194 + 2850 0.86375793 2.2751324 -5.4977459 34.892337 0.0016455263 24.927259 7.289789 25.971316 + 2900 0.84438986 2.3790377 -5.4721407 34.863512 1.2372354 25.819445 7.132603 26.786523 + 2950 0.8551438 2.2721926 -5.4925958 34.84473 1.5405388 25.956466 6.976385 26.87765 + 3000 0.83737707 2.4009609 -5.4707188 34.834171 0.28507766 25.643879 6.1778846 26.377543 + 3050 0.84923235 2.4187994 -5.4938573 34.830836 0.036512025 25.139252 7.1457857 26.135115 + 3100 0.83872396 2.3811576 -5.4838787 34.833673 0.246984 24.21358 8.4588719 25.648586 + 3150 0.83957817 2.3901421 -5.4913118 34.84163 0.20477984 24.309852 10.088243 26.319984 + 3200 0.84283033 2.17292 -5.5025459 34.853975 1.3367154 24.581685 10.72011 26.817531 + 3250 0.84002379 2.1247709 -5.5044955 34.866106 0.11434509 24.463842 9.4874246 26.239108 + 3300 0.83311101 2.1492058 -5.5000847 34.875625 0.0053284993 23.815298 10.560222 26.051616 + 3350 0.83216701 1.9587594 -5.5043446 34.881623 0.58985562 23.934253 11.475462 26.543073 + 3400 0.82396039 2.1914951 -5.4971506 34.881199 0.098206955 23.393402 10.82936 25.778407 + 3450 0.83483253 1.9783612 -5.5182844 34.877327 7.6571212e-05 23.675355 10.761012 26.006188 + 3500 0.82712062 1.9718522 -5.5111818 34.869214 0.014836125 23.314122 11.312845 25.913872 + 3550 0.8342927 1.9114357 -5.5259968 34.855179 1.4050442 22.442758 11.377192 25.161834 + 3600 0.82631637 1.9836457 -5.5176244 34.835738 0.084637609 23.413286 10.824194 25.794285 + 3650 0.82425697 1.9218541 -5.5178548 34.811901 0.11000071 22.788707 12.022258 25.765478 + 3700 0.82491437 1.9624493 -5.521738 34.782417 0.034984027 23.011433 12.384217 26.132257 + 3750 0.82758167 2.0856442 -5.5283493 34.748872 0.001362163 23.030662 12.122144 26.026098 + 3800 0.81891108 1.9858824 -5.5177774 34.714618 0.17075993 23.21344 12.345683 26.292199 + 3850 0.83392227 2.1631514 -5.5426333 34.681146 0.82106473 22.510204 11.678329 25.359272 + 3900 0.82230654 2.0017132 -5.5276756 34.650221 0.48735732 23.444809 12.339117 26.493638 + 3950 0.81929288 2.1749936 -5.5256673 34.61976 0.089219805 23.540062 11.527925 26.211211 + 4000 0.83415169 2.0446791 -5.5506187 34.591266 0.15593937 23.742282 11.26508 26.279231 + 4050 0.82362522 2.1998083 -5.5375157 34.563164 0.25405351 23.913834 11.081011 26.356408 + 4100 0.82589505 2.3074345 -5.543718 34.537763 0.080213125 24.03253 10.435108 26.200266 + 4150 0.83855297 2.2424199 -5.5658171 34.517758 0.62913338 23.974257 8.5079223 25.439138 + 4200 0.82522111 2.2622619 -5.5493275 34.502472 1.8756517 25.754617 7.9996898 26.968414 + 4250 0.82083124 2.4135193 -5.5464932 34.4919 1.1217436 25.8944 6.7070444 26.748914 + 4300 0.83059704 2.1375109 -5.5653245 34.487366 0.53623038 26.05979 4.9072346 26.517798 + 4350 0.82755047 2.1159821 -5.5650889 34.484506 0.10017723 25.405936 4.3532342 25.776195 + 4400 0.83192877 2.180851 -5.5759565 34.480909 0.053664012 25.993034 2.9844338 26.163805 + 4450 0.81860572 2.2333381 -5.5602138 34.477183 0.037864077 25.792233 1.9038859 25.862406 + 4500 0.82821762 2.1142023 -5.5788682 34.474784 0.088221344 26.20329 0.59417897 26.210025 + 4550 0.8205154 2.0896984 -5.5715531 34.472405 0.016076192 26.083166 -0.58187024 26.089655 + 4600 0.81294948 2.2274108 -5.5642678 34.469014 0.033774986 25.869616 0.14951307 25.870048 + 4650 0.80890532 2.1556346 -5.5622407 34.465277 0.67402048 25.413229 -0.56341819 25.419474 + 4700 0.82070227 1.9852605 -5.583747 34.460206 0.052623237 26.158394 -0.44673492 26.162209 + 4750 0.81451857 2.1097726 -5.5779782 34.451438 0.12221733 25.733718 -0.9911436 25.752798 + 4800 0.81300453 2.0211325 -5.5790076 34.439504 0.34536082 26.358606 -1.7335167 26.415548 + 4850 0.82035497 1.9489595 -5.5929886 34.424097 0.70899626 26.575865 -1.0191012 26.595397 + 4900 0.8127066 2.1312269 -5.584271 34.405998 0.087959314 26.185217 -1.1329105 26.209713 + 4950 0.81252621 2.1094866 -5.5866296 34.387869 0.79067667 26.564722 -1.8456354 26.628759 + 5000 0.80575936 2.1875995 -5.579054 34.370679 0.031787364 26.027557 -2.2666774 26.12607 +Loop time of 32.2397 on 1 procs for 5000 steps with 5120 atoms + +Performance: 53598.557 tau/day, 155.088 timesteps/s +99.4% CPU use with 1 MPI tasks x 1 OpenMP threads + +MPI task timing breakdown: +Section | min time | avg time | max time |%varavg| %total +--------------------------------------------------------------- +Pair | 22.967 | 22.967 | 22.967 | 0.0 | 71.24 +Neigh | 2.9914 | 2.9914 | 2.9914 | 0.0 | 9.28 +Comm | 0.37485 | 0.37485 | 0.37485 | 0.0 | 1.16 +Output | 0.064337 | 0.064337 | 0.064337 | 0.0 | 0.20 +Modify | 5.7143 | 5.7143 | 5.7143 | 0.0 | 17.72 +Other | | 0.1281 | | | 0.40 + +Nlocal: 5120 ave 5120 max 5120 min +Histogram: 1 0 0 0 0 0 0 0 0 0 +Nghost: 7962 ave 7962 max 7962 min +Histogram: 1 0 0 0 0 0 0 0 0 0 +Neighs: 296101 ave 296101 max 296101 min +Histogram: 1 0 0 0 0 0 0 0 0 0 + +Total # of neighbors = 296101 +Ave neighs/atom = 57.8322 +Neighbor list builds = 283 +Dangerous builds = 0 + +Please see the log.cite file for references relevant to this simulation + +Total wall time: 0:00:32 diff --git a/examples/USER/misc/pinning/log.22Sep2017.pinning.g++.4 b/examples/USER/misc/pinning/log.22Sep2017.pinning.g++.4 new file mode 100644 index 0000000000..c7a77fb9cc --- /dev/null +++ b/examples/USER/misc/pinning/log.22Sep2017.pinning.g++.4 @@ -0,0 +1,186 @@ +LAMMPS (22 Sep 2017) + using 1 OpenMP thread(s) per MPI task +units lj +dimension 3 +boundary p p p +atom_style atomic + +# truncated and shifted LJ potential +pair_style lj/cut 2.5 +pair_modify shift yes +read_data data.halfhalf + orthogonal box = (0 0 0) to (12.8152 12.8152 34) + 1 by 1 by 4 MPI processor grid + reading atoms ... + 5120 atoms + reading velocities ... + 5120 velocities +pair_coeff 1 1 1.0 1.0 2.5 +mass 1 1.0 + +# simulation parameters +neighbor 0.6 bin +timestep 0.004 +run_style verlet + +velocity all create 0.8 1 mom yes rot yes +fix ensemble all npt temp 0.8 0.8 4.0 z 2.185 2.185 8.0 +fix 100 all momentum 100 linear 1 1 1 + +# harmonic rho_k bias-field +# nx ny nz k a +fix bias all rhok 16 0 0 4.0 26.00 + +# output U_bias rho_k_RE rho_k_IM |rho_k| +thermo_style custom step temp pzz pe lz f_bias f_bias[1] f_bias[2] f_bias[3] +thermo 50 +dump dumpXYZ all xyz 500 traj.xyz + +# NOTE: run reduced for demonstration purposes +# run 50000 +run 5000 +Neighbor list info ... + update every 1 steps, delay 10 steps, check yes + max neighbors/atom: 2000, page size: 100000 + master list distance cutoff = 3.1 + ghost atom cutoff = 3.1 + binsize = 1.55, bins = 9 9 22 + 1 neighbor lists, perpetual/occasional/extra = 1 0 0 + (1) pair lj/cut, perpetual + attributes: half, newton on + pair build: half/bin/atomonly/newton + stencil: half/bin/3d/newton + bin: standard +Per MPI rank memory allocation (min/avg/max) = 4.023 | 4.027 | 4.03 Mbytes +Step Temp Pzz PotEng Lz f_bias f_bias[1] f_bias[2] f_bias[3] + 0 0.8 5.1872229 -4.8361269 34 152.02869 34.711006 -0.72709593 34.71862 + 50 1.0819371 3.9250728 -5.2655842 34.02563 0.51385908 26.505979 0.2187864 26.506882 + 100 1.1182271 3.5397251 -5.3331401 34.100753 2.1059904 27.025883 -0.12127124 27.026156 + 150 1.1121434 3.7845571 -5.3440494 34.213993 1.7206575 26.894862 -1.3261751 26.927539 + 200 1.1446439 3.4114364 -5.4199119 34.358914 2.383615 27.054401 -1.4211151 27.091699 + 250 1.1115073 3.5603047 -5.3988013 34.517397 0.60819391 26.51404 -1.4089688 26.55145 + 300 1.0828478 3.7411116 -5.3842818 34.673948 0.73987465 26.528178 -2.062382 26.608225 + 350 1.0342597 3.9206217 -5.3384367 34.825105 0.99965014 26.657737 -1.6211152 26.706983 + 400 1.0064356 3.9929044 -5.324003 34.97579 0.41927007 26.402623 -1.7087432 26.457859 + 450 0.96799277 4.2764255 -5.2947175 35.133839 0.77739461 26.503753 -2.5217998 26.623456 + 500 0.94691076 4.1007962 -5.2922124 35.301893 0.17015805 26.212252 -2.0421698 26.291683 + 550 0.93675297 3.7404088 -5.3056917 35.474301 0.56247039 26.335926 -3.205722 26.530316 + 600 0.92707577 3.5439822 -5.3176094 35.641282 0.04054693 25.679461 -3.0301039 25.857615 + 650 0.91828226 3.1833648 -5.3278237 35.794766 0.8427054 26.4003 -3.6331211 26.649117 + 700 0.9084826 2.8374306 -5.3327944 35.928138 1.5635222 26.605971 -3.8575939 26.884173 + 750 0.91219559 2.46172 -5.3548299 36.039156 6.3772911 27.350725 -4.8971146 27.785678 + 800 0.90000337 2.2187716 -5.3499181 36.126451 4.9080838 27.085156 -5.1291678 27.566538 + 850 0.9003432 1.8634244 -5.3614648 36.189019 0.0081092188 25.497333 -5.4038153 26.063676 + 900 0.89741573 1.5697398 -5.3660799 36.226074 0.011155479 25.312724 -6.2574069 26.074684 + 950 0.88871454 1.4427205 -5.3604669 36.237407 0.3287078 25.659237 -6.232896 26.405406 + 1000 0.88606353 1.3088636 -5.3626576 36.226015 0.30068168 24.554896 -7.2832017 25.612262 + 1050 0.88527541 1.3194263 -5.3666773 36.19311 0.10646314 24.514921 -7.9419424 25.76928 + 1100 0.87522001 1.2852124 -5.3556811 36.143056 0.13675329 24.865981 -8.446822 26.261489 + 1150 0.8805978 1.246973 -5.3671716 36.0781 0.00043275463 24.187039 -9.4985495 25.98529 + 1200 0.85711495 1.376588 -5.3346243 36.002427 0.47623639 23.691349 -9.4648541 25.512026 + 1250 0.88116805 1.3562001 -5.3731036 35.919289 0.32797055 23.322103 -10.54448 25.595049 + 1300 0.87178482 1.5046564 -5.3610798 35.831278 0.17704849 24.190231 -10.314689 26.29753 + 1350 0.87022621 1.6830825 -5.3603618 35.743318 0.0052854997 23.731157 -10.747465 26.051408 + 1400 0.89019669 1.6144812 -5.3921986 35.659687 1.4152796 22.8393 -10.551347 25.158787 + 1450 0.88852819 1.7587964 -5.3918592 35.580319 0.63560961 23.599033 -12.195 26.563742 + 1500 0.89029085 1.8772498 -5.3966098 35.509232 0.20895386 23.055083 -12.703366 26.323229 + 1550 0.88639722 2.2284824 -5.3933288 35.449043 0.44413965 22.448774 -12.156068 25.528757 + 1600 0.88816451 2.2167704 -5.3994757 35.401661 0.12210235 23.108351 -12.44643 26.247085 + 1650 0.89154791 2.3397824 -5.4086923 35.365815 0.4820208 23.090699 -12.984179 26.490928 + 1700 0.88518032 2.5351236 -5.4041601 35.343757 0.080806002 22.749825 -12.99762 26.201005 + 1750 0.86848721 2.5527491 -5.3851928 35.336433 0.045102165 22.357111 -13.564328 26.15017 + 1800 0.88501061 2.5215825 -5.4169341 35.340849 0.27488483 22.086584 -14.408273 26.370732 + 1850 0.8716061 2.5809558 -5.4045854 35.355038 0.042909785 21.270956 -14.695278 25.853525 + 1900 0.85672517 2.4836326 -5.3902797 35.375469 0.72877764 21.639909 -15.474764 26.603646 + 1950 0.85133731 2.3141629 -5.3902573 35.398523 0.0016908803 21.106617 -15.132733 25.970924 + 2000 0.86152109 2.1562002 -5.4132601 35.419851 0.371016 21.325237 -15.614625 26.430706 + 2050 0.86243551 2.019931 -5.4220349 35.436069 0.017935421 20.4131 -16.255418 26.094698 + 2100 0.87417672 1.8083823 -5.4464117 35.445091 0.18429432 19.75625 -17.365705 26.303558 + 2150 0.85872128 1.7608768 -5.4293103 35.44341 0.91209166 20.149648 -17.480387 26.675312 + 2200 0.86615373 1.8372778 -5.4458315 35.430616 0.10151993 18.559234 -17.885469 25.7747 + 2250 0.85053605 1.7198437 -5.4272104 35.408688 0.96154548 17.200861 -18.562206 25.306622 + 2300 0.85400281 1.7939644 -5.4364682 35.377708 0.12283263 18.759325 -18.358539 26.247823 + 2350 0.85495278 1.5856029 -5.4417321 35.337987 0.20564329 18.967923 -18.248149 26.320658 + 2400 0.84606771 1.7782708 -5.4315646 35.287411 0.10063977 19.185527 -17.878215 26.224321 + 2450 0.85210051 1.8190391 -5.4432116 35.232321 0.69988647 19.268861 -18.325448 26.59156 + 2500 0.85304715 1.7466204 -5.4470889 35.175245 0.0048314937 18.09176 -18.74157 26.04915 + 2550 0.85401123 1.8601945 -5.4509309 35.115748 0.99467901 17.170045 -18.574587 25.294777 + 2600 0.85778606 1.974012 -5.4586742 35.058013 0.0026599702 17.438966 -19.333395 26.036469 + 2650 0.8521239 2.0606329 -5.4526006 35.003616 0.091056354 17.16363 -19.244738 25.786627 + 2700 0.85918482 2.0766792 -5.4658947 34.954171 0.89590606 15.77108 -19.822153 25.330707 + 2750 0.85786577 2.225549 -5.4667773 34.911468 0.26577575 15.769018 -21.128817 26.364538 + 2800 0.86764664 2.2325018 -5.4849414 34.877604 0.47167555 14.950515 -20.675229 25.514369 + 2850 0.85209564 2.3434319 -5.465734 34.852715 2.7350296 13.51553 -20.829996 24.830592 + 2900 0.85757283 2.3512971 -5.4786051 34.836138 0.14816492 14.06033 -21.545946 25.727819 + 2950 0.86098926 2.3480431 -5.4890615 34.826408 0.26401534 13.381395 -22.714827 26.363329 + 3000 0.85413421 2.3243973 -5.4844129 34.823242 0.024244334 12.739486 -22.538687 25.889899 + 3050 0.85015323 2.5479266 -5.4844303 34.825228 0.4463147 12.990582 -21.975063 25.527605 + 3100 0.8530523 2.3643505 -5.495343 34.834883 0.12144265 12.844293 -22.321989 25.753583 + 3150 0.85098478 2.2521299 -5.4990526 34.848419 0.33194916 12.747856 -23.126671 26.4074 + 3200 0.84391449 2.2650057 -5.495222 34.862626 0.031888328 12.788845 -22.782174 26.12627 + 3250 0.84807155 2.1715388 -5.5080873 34.877548 0.082426694 13.316219 -22.09441 25.796989 + 3300 0.83028242 2.242889 -5.4878846 34.89175 1.1334975 14.326678 -22.593363 26.752827 + 3350 0.82924001 2.0324002 -5.4924558 34.903232 0.35473989 14.354166 -22.181868 26.421153 + 3400 0.83032841 2.0003371 -5.4997142 34.908733 0.041677437 14.528378 -21.735998 26.144356 + 3450 0.82908891 1.8683902 -5.5029185 34.907936 0.02365857 15.069507 -21.053887 25.891237 + 3500 0.82842914 1.9165344 -5.5064218 34.898681 0.17663531 15.27043 -20.674834 25.702817 + 3550 0.82735822 1.98221 -5.5088197 34.88272 1.5607134 14.915228 -20.208431 25.116622 + 3600 0.82642915 1.8422766 -5.5110752 34.8611 1.1861112 15.312314 -20.051953 25.229899 + 3650 0.82556781 1.9351408 -5.5130349 34.833406 1.018872 16.152478 -19.454871 25.286252 + 3700 0.82360651 1.9791184 -5.5128431 34.802021 0.14080727 16.907104 -19.401616 25.734663 + 3750 0.83017793 1.9855734 -5.5253254 34.768644 0.15311334 16.969506 -19.331958 25.723311 + 3800 0.82362926 2.1029656 -5.5179624 34.734178 0.10807487 17.892584 -18.542426 25.76754 + 3850 0.82313508 2.0781681 -5.5196175 34.70093 0.13343085 19.072706 -17.28778 25.741707 + 3900 0.83643385 2.0570262 -5.5421224 34.669761 0.00022792038 19.551677 -17.1548 26.010675 + 3950 0.82346174 2.0842322 -5.5252757 34.640849 0.0093759386 20.892792 -15.590263 26.068469 + 4000 0.83485868 2.1196669 -5.5451736 34.612396 0.31198053 21.630258 -15.126984 26.394956 + 4050 0.82729429 2.2033274 -5.5365945 34.585721 0.53752252 21.283533 -14.011497 25.481578 + 4100 0.82040242 2.1757309 -5.5292269 34.562271 0.36031984 22.047609 -12.961927 25.575548 + 4150 0.81932521 2.285666 -5.5307807 34.542102 0.84343149 22.486289 -11.70555 25.350604 + 4200 0.83819319 2.231174 -5.5625532 34.526447 0.47190752 23.311855 -12.57189 26.485751 + 4250 0.82542274 2.1874789 -5.5472057 34.513795 0.70518398 23.411553 -12.614639 26.593795 + 4300 0.81971158 2.241167 -5.5424504 34.503969 0.26707612 23.089805 -12.727793 26.365429 + 4350 0.83255377 2.1295532 -5.5657895 34.496326 0.072548591 23.003138 -12.52181 26.190458 + 4400 0.8128474 2.3327845 -5.5402264 34.490126 0.0013023434 23.020811 -12.029795 25.974482 + 4450 0.82013491 2.3069915 -5.5554953 34.488039 0.041123896 23.632908 -11.178674 26.143394 + 4500 0.81411544 2.2247193 -5.5509183 34.488014 0.54440601 23.010678 -10.938506 25.478269 + 4550 0.82814624 2.1142779 -5.5763482 34.487885 0.1518945 23.696817 -11.351972 26.275585 + 4600 0.82929492 2.090881 -5.5823492 34.486698 0.0045520899 23.538527 -10.929741 25.952292 + 4650 0.81061417 1.9818043 -5.5584018 34.484038 0.012526806 23.993543 -10.219174 26.079142 + 4700 0.81816105 1.9605811 -5.5735005 34.476764 1.2079835 25.151166 -9.1888856 26.777169 + 4750 0.81657042 2.0064313 -5.5744795 34.465784 1.2045017 25.487486 -8.2063886 26.776048 + 4800 0.81789335 2.0838696 -5.5796632 34.451996 0.27642542 24.647157 -7.023095 25.62823 + 4850 0.80649339 1.9892413 -5.5654796 34.436067 0.024697945 25.09823 -6.3492244 25.888874 + 4900 0.81673441 2.0125635 -5.5835037 34.416236 0.0011188576 25.446818 -5.2182483 25.976348 + 4950 0.82250033 1.9770391 -5.5946082 34.394723 0.72696707 26.37002 -3.5122842 26.602896 + 5000 0.80762758 2.075517 -5.5746076 34.371696 0.12796344 26.102184 -2.8094827 26.252946 +Loop time of 10.3394 on 4 procs for 5000 steps with 5120 atoms + +Performance: 167127.370 tau/day, 483.586 timesteps/s +99.0% CPU use with 4 MPI tasks x 1 OpenMP threads + +MPI task timing breakdown: +Section | min time | avg time | max time |%varavg| %total +--------------------------------------------------------------- +Pair | 5.029 | 6.0128 | 7.1918 | 35.0 | 58.15 +Neigh | 0.65673 | 0.75825 | 0.88597 | 10.3 | 7.33 +Comm | 0.43982 | 1.5284 | 2.4112 | 60.5 | 14.78 +Output | 0.022835 | 0.023039 | 0.023453 | 0.2 | 0.22 +Modify | 1.7294 | 1.9472 | 2.5687 | 25.7 | 18.83 +Other | | 0.06978 | | | 0.67 + +Nlocal: 1280 ave 1404 max 1214 min +Histogram: 2 0 0 1 0 0 0 0 0 1 +Nghost: 3521.25 ave 3581 max 3426 min +Histogram: 1 0 0 0 0 0 1 1 0 1 +Neighs: 73872.2 ave 87973 max 64161 min +Histogram: 1 1 0 0 1 0 0 0 0 1 + +Total # of neighbors = 295489 +Ave neighs/atom = 57.7127 +Neighbor list builds = 278 +Dangerous builds = 0 + +Please see the log.cite file for references relevant to this simulation + +Total wall time: 0:00:10 diff --git a/examples/USER/misc/pinning/log.22Sep2017.setup.g++.1 b/examples/USER/misc/pinning/log.22Sep2017.setup.g++.1 new file mode 100644 index 0000000000..8606d4ed94 --- /dev/null +++ b/examples/USER/misc/pinning/log.22Sep2017.setup.g++.1 @@ -0,0 +1,141 @@ +LAMMPS (22 Sep 2017) + using 1 OpenMP thread(s) per MPI task +units lj +dimension 3 +boundary p p p +atom_style atomic + +# truncated and shifted LJ potential +pair_style lj/cut 2.5 +pair_modify shift yes + +# fcc lattice +lattice fcc 0.9731 +Lattice spacing in x,y,z = 1.6019 1.6019 1.6019 +region my_box block 0 8.0 0 8.0 0 20.0 +create_box 1 my_box +Created orthogonal box = (0 0 0) to (12.8152 12.8152 32.0379) + 1 by 1 by 1 MPI processor grid +region particles block 0 8.0 0 8.0 0 20.0 +create_atoms 1 region particles +Created 5120 atoms +pair_coeff 1 1 1.0 1.0 2.5 +mass 1 1.0 +change_box all z final 0.0 34 remap units box + orthogonal box = (0 0 0) to (12.8152 12.8152 34) + +# select particles in one side of the elongated box +region left plane 0 0 10 0 0 1 +group left region left +2688 atoms in group left + +velocity left create 6.0 1 mom yes rot yes + +# simulation parameters +neighbor 0.6 bin +timestep 0.004 +run_style verlet +fix ensemble left nve # Note: only move particle in left-hand side +fix langevin left langevin 3.0 0.8 100.0 2017 + +# outout +thermo_style custom step temp pzz pe lz +thermo 100 +# dump dumpXYZ all xyz 100 traj.xyz + +# run reduced for demonstration purposes +# run 10000 +run 5000 +Neighbor list info ... + update every 1 steps, delay 10 steps, check yes + max neighbors/atom: 2000, page size: 100000 + master list distance cutoff = 3.1 + ghost atom cutoff = 3.1 + binsize = 1.55, bins = 9 9 22 + 1 neighbor lists, perpetual/occasional/extra = 1 0 0 + (1) pair lj/cut, perpetual + attributes: half, newton on + pair build: half/bin/atomonly/newton + stencil: half/bin/3d/newton + bin: standard +Per MPI rank memory allocation (min/avg/max) = 4.524 | 4.524 | 4.524 Mbytes +Step Temp Pzz PotEng Lz + 0 3.1494433 -3.4735106 -6.8707307 34 + 100 1.7727555 6.5330255 -4.8035477 34 + 200 1.7462368 7.0070325 -4.7646426 34 + 300 1.7564888 6.6190123 -4.7894637 34 + 400 1.7641186 6.609684 -4.8064772 34 + 500 1.7383511 6.7304936 -4.7708095 34 + 600 1.731708 6.8574656 -4.7612918 34 + 700 1.7332167 6.6530919 -4.7670014 34 + 800 1.7487537 6.5644963 -4.7907458 34 + 900 1.7353648 6.7115188 -4.7772149 34 + 1000 1.728878 6.4175719 -4.7797216 34 + 1100 1.7471488 6.5346083 -4.813376 34 + 1200 1.7188149 6.2502104 -4.7822235 34 + 1300 1.7151194 6.792534 -4.7781701 34 + 1400 1.7406603 6.6639592 -4.8170174 34 + 1500 1.7090537 6.4677579 -4.770701 34 + 1600 1.7014954 6.2853535 -4.7679742 34 + 1700 1.7064354 6.4352857 -4.7812978 34 + 1800 1.7169971 6.5808758 -4.799426 34 + 1900 1.6822712 6.3746758 -4.7522464 34 + 2000 1.7126546 6.534969 -4.8091595 34 + 2100 1.7086108 6.4679932 -4.8146664 34 + 2200 1.6974952 6.3802129 -4.8052505 34 + 2300 1.6868035 6.4009243 -4.7935769 34 + 2400 1.7107125 6.2318869 -4.8358765 34 + 2500 1.660241 6.4891487 -4.7661183 34 + 2600 1.6801816 6.1988356 -4.8024291 34 + 2700 1.6940298 6.1328187 -4.8290053 34 + 2800 1.6755061 6.4150693 -4.8145473 34 + 2900 1.6749928 6.4248792 -4.8213509 34 + 3000 1.6310737 6.6491291 -4.7673027 34 + 3100 1.6559915 6.2726719 -4.8109181 34 + 3200 1.6574579 5.7132029 -4.8189484 34 + 3300 1.6816136 5.7697439 -4.8652811 34 + 3400 1.6489483 6.4463349 -4.8247812 34 + 3500 1.6557974 5.9763333 -4.8383712 34 + 3600 1.6215459 6.2806534 -4.7954657 34 + 3700 1.6484987 6.0671609 -4.8470777 34 + 3800 1.6473922 5.8688108 -4.8555351 34 + 3900 1.6435957 5.930425 -4.8562076 34 + 4000 1.6514434 6.1962122 -4.872998 34 + 4100 1.6138337 6.4808124 -4.8219373 34 + 4200 1.6215239 5.9467966 -4.8412146 34 + 4300 1.6129295 5.9377323 -4.8414596 34 + 4400 1.6020549 6.1104301 -4.8395939 34 + 4500 1.6047738 6.0816222 -4.8538151 34 + 4600 1.6053565 6.183466 -4.8686817 34 + 4700 1.6088152 5.7416542 -4.894114 34 + 4800 1.5954309 5.694319 -4.8840198 34 + 4900 1.5582564 6.1199614 -4.8429998 34 + 5000 1.5786672 5.8813574 -4.8907344 34 +Loop time of 28.3867 on 1 procs for 5000 steps with 5120 atoms + +Performance: 60873.483 tau/day, 176.139 timesteps/s +99.4% CPU use with 1 MPI tasks x 1 OpenMP threads + +MPI task timing breakdown: +Section | min time | avg time | max time |%varavg| %total +--------------------------------------------------------------- +Pair | 22.269 | 22.269 | 22.269 | 0.0 | 78.45 +Neigh | 4.7222 | 4.7222 | 4.7222 | 0.0 | 16.64 +Comm | 0.40821 | 0.40821 | 0.40821 | 0.0 | 1.44 +Output | 0.0042329 | 0.0042329 | 0.0042329 | 0.0 | 0.01 +Modify | 0.88231 | 0.88231 | 0.88231 | 0.0 | 3.11 +Other | | 0.1005 | | | 0.35 + +Nlocal: 5120 ave 5120 max 5120 min +Histogram: 1 0 0 0 0 0 0 0 0 0 +Nghost: 7768 ave 7768 max 7768 min +Histogram: 1 0 0 0 0 0 0 0 0 0 +Neighs: 297167 ave 297167 max 297167 min +Histogram: 1 0 0 0 0 0 0 0 0 0 + +Total # of neighbors = 297167 +Ave neighs/atom = 58.0404 +Neighbor list builds = 474 +Dangerous builds = 246 +write_data data.halfhalf +Total wall time: 0:00:28 diff --git a/examples/USER/misc/pinning/log.22Sep2017.setup.g++.4 b/examples/USER/misc/pinning/log.22Sep2017.setup.g++.4 new file mode 100644 index 0000000000..14088f7c95 --- /dev/null +++ b/examples/USER/misc/pinning/log.22Sep2017.setup.g++.4 @@ -0,0 +1,141 @@ +LAMMPS (22 Sep 2017) + using 1 OpenMP thread(s) per MPI task +units lj +dimension 3 +boundary p p p +atom_style atomic + +# truncated and shifted LJ potential +pair_style lj/cut 2.5 +pair_modify shift yes + +# fcc lattice +lattice fcc 0.9731 +Lattice spacing in x,y,z = 1.6019 1.6019 1.6019 +region my_box block 0 8.0 0 8.0 0 20.0 +create_box 1 my_box +Created orthogonal box = (0 0 0) to (12.8152 12.8152 32.0379) + 1 by 1 by 4 MPI processor grid +region particles block 0 8.0 0 8.0 0 20.0 +create_atoms 1 region particles +Created 5120 atoms +pair_coeff 1 1 1.0 1.0 2.5 +mass 1 1.0 +change_box all z final 0.0 34 remap units box + orthogonal box = (0 0 0) to (12.8152 12.8152 34) + +# select particles in one side of the elongated box +region left plane 0 0 10 0 0 1 +group left region left +2688 atoms in group left + +velocity left create 6.0 1 mom yes rot yes + +# simulation parameters +neighbor 0.6 bin +timestep 0.004 +run_style verlet +fix ensemble left nve # Note: only move particle in left-hand side +fix langevin left langevin 3.0 0.8 100.0 2017 + +# outout +thermo_style custom step temp pzz pe lz +thermo 100 +# dump dumpXYZ all xyz 100 traj.xyz + +# run reduced for demonstration purposes +# run 10000 +run 5000 +Neighbor list info ... + update every 1 steps, delay 10 steps, check yes + max neighbors/atom: 2000, page size: 100000 + master list distance cutoff = 3.1 + ghost atom cutoff = 3.1 + binsize = 1.55, bins = 9 9 22 + 1 neighbor lists, perpetual/occasional/extra = 1 0 0 + (1) pair lj/cut, perpetual + attributes: half, newton on + pair build: half/bin/atomonly/newton + stencil: half/bin/3d/newton + bin: standard +Per MPI rank memory allocation (min/avg/max) = 3.23 | 3.23 | 3.23 Mbytes +Step Temp Pzz PotEng Lz + 0 3.1494433 -3.4735106 -6.8707307 34 + 100 1.7914373 6.4805818 -4.8420353 34 + 200 1.740256 6.6108149 -4.7672571 34 + 300 1.7663827 6.5188941 -4.8103672 34 + 400 1.7440644 6.5156543 -4.7769467 34 + 500 1.7471724 6.5208992 -4.7843928 34 + 600 1.7320106 6.6557835 -4.7654637 34 + 700 1.6839043 6.7689759 -4.7045352 34 + 800 1.7216746 6.66436 -4.7601673 34 + 900 1.7342542 6.3242367 -4.7790803 34 + 1000 1.7338566 6.5803438 -4.7854529 34 + 1100 1.7328856 6.3846366 -4.7902625 34 + 1200 1.7546906 6.5048137 -4.8213443 34 + 1300 1.7163891 6.3903221 -4.7665145 34 + 1400 1.7011627 6.5782672 -4.7517875 34 + 1500 1.7105234 6.5811813 -4.7677748 34 + 1600 1.7334403 6.5032837 -4.8067749 34 + 1700 1.7252102 6.5443871 -4.8058994 34 + 1800 1.721958 6.3378188 -4.8150073 34 + 1900 1.6797892 6.6780506 -4.7538618 34 + 2000 1.7001774 6.3578192 -4.7894018 34 + 2100 1.7127784 6.3219105 -4.8161059 34 + 2200 1.696825 6.536793 -4.7946902 34 + 2300 1.6704578 6.7186933 -4.7609628 34 + 2400 1.6772498 6.3432817 -4.7778471 34 + 2500 1.7073862 6.2153226 -4.8299181 34 + 2600 1.6951557 6.4397257 -4.8156787 34 + 2700 1.6845984 6.0123544 -4.8136864 34 + 2800 1.6550565 6.2489392 -4.7829639 34 + 2900 1.6892315 6.158499 -4.8423004 34 + 3000 1.6814436 6.07976 -4.8400696 34 + 3100 1.6387025 6.330166 -4.7878978 34 + 3200 1.6747855 6.0767043 -4.8481995 34 + 3300 1.6508768 6.2749233 -4.8181888 34 + 3400 1.6426364 6.3934935 -4.8223824 34 + 3500 1.6576512 6.0638185 -4.8559078 34 + 3600 1.6444173 6.1376573 -4.8463113 34 + 3700 1.6480039 5.9943705 -4.8601776 34 + 3800 1.6467212 6.0556591 -4.8722719 34 + 3900 1.6271804 6.116738 -4.8547278 34 + 4000 1.6158134 5.9089534 -4.8477829 34 + 4100 1.6388157 5.9890465 -4.8920284 34 + 4200 1.6182368 6.0639887 -4.8724963 34 + 4300 1.647633 5.6333906 -4.9267536 34 + 4400 1.5856411 6.2675475 -4.8471239 34 + 4500 1.5773417 6.1789163 -4.8469057 34 + 4600 1.6181445 5.7988068 -4.922419 34 + 4700 1.5876712 5.7398111 -4.8853849 34 + 4800 1.5708353 6.2204997 -4.8718872 34 + 4900 1.5514708 5.9782256 -4.8523812 34 + 5000 1.553347 5.9286523 -4.86582 34 +Loop time of 8.10259 on 4 procs for 5000 steps with 5120 atoms + +Performance: 213265.164 tau/day, 617.087 timesteps/s +99.2% CPU use with 4 MPI tasks x 1 OpenMP threads + +MPI task timing breakdown: +Section | min time | avg time | max time |%varavg| %total +--------------------------------------------------------------- +Pair | 5.2964 | 5.6236 | 5.8982 | 9.0 | 69.40 +Neigh | 1.0562 | 1.1907 | 1.3257 | 8.7 | 14.70 +Comm | 0.43963 | 0.98786 | 1.5968 | 42.5 | 12.19 +Output | 0.0023124 | 0.004741 | 0.0090873 | 4.0 | 0.06 +Modify | 0.018652 | 0.22213 | 0.39884 | 36.4 | 2.74 +Other | | 0.07357 | | | 0.91 + +Nlocal: 1280 ave 1337 max 1211 min +Histogram: 1 0 0 0 0 1 1 0 0 1 +Nghost: 3416.25 ave 3549 max 3297 min +Histogram: 2 0 0 0 0 0 0 0 1 1 +Neighs: 74269.8 ave 77932 max 69612 min +Histogram: 1 0 0 0 0 1 0 1 0 1 + +Total # of neighbors = 297079 +Ave neighs/atom = 58.0232 +Neighbor list builds = 474 +Dangerous builds = 247 +write_data data.halfhalf +Total wall time: 0:00:08 From 251f28049a3d43c3516df9d500bd9357c89a9818 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Tue, 3 Oct 2017 10:10:38 -0400 Subject: [PATCH 42/44] make fix rhok more similar to other LAMMPS sources - re-indent to 2 blanks - white space cleanup - use force->numeric() and force->inumeric() instead of atof() and atoi() - include system headers before local LAMMPS headers --- src/USER-MISC/fix_rhok.cpp | 209 ++++++++++++++++++------------------- src/USER-MISC/fix_rhok.h | 30 +++--- 2 files changed, 117 insertions(+), 122 deletions(-) diff --git a/src/USER-MISC/fix_rhok.cpp b/src/USER-MISC/fix_rhok.cpp index df8b8bd667..58b0e95a97 100644 --- a/src/USER-MISC/fix_rhok.cpp +++ b/src/USER-MISC/fix_rhok.cpp @@ -13,24 +13,24 @@ Contributing author: Ulf R. Pedersen, ulf@urp.dk ------------------------------------------------------------------------- */ -#include "fix_rhok.h" -#include "error.h" -#include "citeme.h" -#include "update.h" -#include "respa.h" -#include "atom.h" -#include "domain.h" - #include #include #include -#include #include +#include "fix_rhok.h" +#include "atom.h" +#include "domain.h" +#include "error.h" +#include "force.h" +#include "respa.h" +#include "update.h" +#include "citeme.h" + using namespace LAMMPS_NS; using namespace FixConst; -static const char cite_flow_gauss[] = +static const char cite_fix_rhok[] = "Bias on the collective density field (fix rhok):\n\n" "@Article{pedersen_jcp139_104102_2013,\n" "title = {Direct calculation of the solid-liquid Gibbs free energy difference in a single equilibrium simulation},\n" @@ -48,40 +48,35 @@ static const char cite_flow_gauss[] = FixRhok::FixRhok( LAMMPS* inLMP, int inArgc, char** inArgv ) : Fix( inLMP, inArgc, inArgv ) { - - if (lmp->citeme) lmp->citeme->add(cite_flow_gauss); - + + if (lmp->citeme) lmp->citeme->add(cite_fix_rhok); + // Check arguments - if( inArgc != 8 ) - error->all(FLERR,"Illegal fix rhoKUmbrella command" ); + if( inArgc != 8 ) + error->all(FLERR,"Illegal fix rhoKUmbrella command" ); // Set up fix flags scalar_flag = 1; // have compute_scalar vector_flag = 1; // have compute_vector... size_vector = 3; // ...with this many components global_freq = 1; // whose value can be computed at every timestep - //scalar_vector_freq = 1; // OLD lammps: whose value can be computed at every timestep thermo_energy = 1; // this fix changes system's potential energy extscalar = 0; // but the deltaPE might not scale with # of atoms extvector = 0; // neither do the components of the vector // Parse fix options - int n[3]; - - n[0] = atoi( inArgv[3] ); - n[1] = atoi( inArgv[4] ); - n[2] = atoi( inArgv[5] ); - - mK[0] = n[0]*(2*M_PI / (domain->boxhi[0] - domain->boxlo[0])); - mK[1] = n[1]*(2*M_PI / (domain->boxhi[1] - domain->boxlo[1])); - mK[2] = n[2]*(2*M_PI / (domain->boxhi[2] - domain->boxlo[2])); - - mKappa = atof( inArgv[6] ); - mRhoK0 = atof( inArgv[7] ); -} + int n[3]; -FixRhok::~FixRhok() -{ + n[0] = force->inumeric(FLERR,inArgv[3]); + n[1] = force->inumeric(FLERR,inArgv[4]); + n[2] = force->inumeric(FLERR,inArgv[5]); + + mK[0] = n[0]*(2*M_PI / (domain->boxhi[0] - domain->boxlo[0])); + mK[1] = n[1]*(2*M_PI / (domain->boxhi[1] - domain->boxlo[1])); + mK[2] = n[2]*(2*M_PI / (domain->boxhi[2] - domain->boxlo[2])); + + mKappa = force->numeric(FLERR,inArgv[6]); + mRhoK0 = force->numeric(FLERR,inArgv[7]); } // Methods that this fix implements @@ -110,20 +105,20 @@ FixRhok::init() { // RESPA boilerplate if( strcmp( update->integrate_style, "respa" ) == 0 ) - mNLevelsRESPA = ((Respa *) update->integrate)->nlevels; - - // Count the number of affected particles - int nThisLocal = 0; - int *mask = atom->mask; - int nlocal = atom->nlocal; - for( int i = 0; i < nlocal; i++ ) { // Iterate through all atoms on this CPU - if( mask[i] & groupbit ) { // ...only those affected by this fix - nThisLocal++; - } - } - MPI_Allreduce( &nThisLocal, &mNThis, - 1, MPI_INT, MPI_SUM, world ); - mSqrtNThis = sqrt( mNThis ); + mNLevelsRESPA = ((Respa *) update->integrate)->nlevels; + + // Count the number of affected particles + int nThisLocal = 0; + int *mask = atom->mask; + int nlocal = atom->nlocal; + for( int i = 0; i < nlocal; i++ ) { // Iterate through all atoms on this CPU + if( mask[i] & groupbit ) { // ...only those affected by this fix + nThisLocal++; + } + } + MPI_Allreduce( &nThisLocal, &mNThis, + 1, MPI_INT, MPI_SUM, world ); + mSqrtNThis = sqrt( mNThis ); } // Initial application of the fix to a system (when doing MD) @@ -131,13 +126,13 @@ void FixRhok::setup( int inVFlag ) { if( strcmp( update->integrate_style, "verlet" ) == 0 ) - post_force( inVFlag ); + post_force( inVFlag ); else - { - ((Respa *) update->integrate)->copy_flevel_f( mNLevelsRESPA - 1 ); - post_force_respa( inVFlag, mNLevelsRESPA - 1,0 ); - ((Respa *) update->integrate)->copy_f_flevel( mNLevelsRESPA - 1 ); - } + { + ((Respa *) update->integrate)->copy_flevel_f( mNLevelsRESPA - 1 ); + post_force_respa( inVFlag, mNLevelsRESPA - 1,0 ); + ((Respa *) update->integrate)->copy_f_flevel( mNLevelsRESPA - 1 ); + } } // Initial application of the fix to a system (when doing minimization) @@ -159,54 +154,54 @@ FixRhok::post_force( int inVFlag ) // Loop over locally-owned atoms affected by this fix and calculate the // partial rhoK's mRhoKLocal[0] = 0.0; - mRhoKLocal[1] = 0.0; + mRhoKLocal[1] = 0.0; for( int i = 0; i < nlocal; i++ ) { // Iterate through all atoms on this CPU - if( mask[i] & groupbit ) { // ...only those affected by this fix + if( mask[i] & groupbit ) { // ...only those affected by this fix - // rho_k = sum_i exp( - i k.r_i ) - mRhoKLocal[0] += cos( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); - mRhoKLocal[1] -= sin( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); - } + // rho_k = sum_i exp( - i k.r_i ) + mRhoKLocal[0] += cos( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); + mRhoKLocal[1] -= sin( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); + } } - - // Now calculate mRhoKGlobal - MPI_Allreduce( mRhoKLocal, mRhoKGlobal, - 2, MPI_DOUBLE, MPI_SUM, world ); - - // Info: < \sum_{i,j} e^{-ik.(r_i - r_j)} > ~ N, so - // we define rho_k as (1 / sqrt(N)) \sum_i e^{-i k.r_i}, so that - // is intensive. - mRhoKGlobal[0] /= mSqrtNThis; - mRhoKGlobal[1] /= mSqrtNThis; - - // We'll need magnitude of rho_k - double rhoK = sqrt( mRhoKGlobal[0]*mRhoKGlobal[0] - + mRhoKGlobal[1]*mRhoKGlobal[1] ); - - for( int i = 0; i < nlocal; i++ ) { // Iterate through all atoms on this CPU - if( mask[i] & groupbit ) { // ...only those affected by this fix - - // Calculate forces - // U = kappa/2 ( |rho_k| - rho_k^0 )^2 - // f_i = -grad_i U = -kappa ( |rho_k| - rho_k^0 ) grad_i |rho_k| - // grad_i |rho_k| = Re( rho_k* (-i k e^{-i k . r_i} / sqrt(N)) ) / |rho_k| - // - // In terms of real and imag parts of rho_k, - // - // Re( rho_k* (-i k e^{-i k . r_i}) ) = - // (- Re[rho_k] * sin( k . r_i ) - Im[rho_k] * cos( k . r_i )) * k - double sinKRi = sin( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); - double cosKRi = cos( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); - - double prefactor = mKappa * ( rhoK - mRhoK0 ) / rhoK - * (-mRhoKGlobal[0]*sinKRi - mRhoKGlobal[1]*cosKRi) / mSqrtNThis; - f[i][0] -= prefactor * mK[0]; - f[i][1] -= prefactor * mK[1]; - f[i][2] -= prefactor * mK[2]; - } - } + // Now calculate mRhoKGlobal + MPI_Allreduce( mRhoKLocal, mRhoKGlobal, + 2, MPI_DOUBLE, MPI_SUM, world ); + + // Info: < \sum_{i,j} e^{-ik.(r_i - r_j)} > ~ N, so + // we define rho_k as (1 / sqrt(N)) \sum_i e^{-i k.r_i}, so that + // is intensive. + mRhoKGlobal[0] /= mSqrtNThis; + mRhoKGlobal[1] /= mSqrtNThis; + + // We'll need magnitude of rho_k + double rhoK = sqrt( mRhoKGlobal[0]*mRhoKGlobal[0] + + mRhoKGlobal[1]*mRhoKGlobal[1] ); + + for( int i = 0; i < nlocal; i++ ) { // Iterate through all atoms on this CPU + if( mask[i] & groupbit ) { // ...only those affected by this fix + + // Calculate forces + // U = kappa/2 ( |rho_k| - rho_k^0 )^2 + // f_i = -grad_i U = -kappa ( |rho_k| - rho_k^0 ) grad_i |rho_k| + // grad_i |rho_k| = Re( rho_k* (-i k e^{-i k . r_i} / sqrt(N)) ) / |rho_k| + // + // In terms of real and imag parts of rho_k, + // + // Re( rho_k* (-i k e^{-i k . r_i}) ) = + // (- Re[rho_k] * sin( k . r_i ) - Im[rho_k] * cos( k . r_i )) * k + + double sinKRi = sin( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); + double cosKRi = cos( mK[0]*x[i][0] + mK[1]*x[i][1] + mK[2]*x[i][2] ); + + double prefactor = mKappa * ( rhoK - mRhoK0 ) / rhoK + * (-mRhoKGlobal[0]*sinKRi - mRhoKGlobal[1]*cosKRi) / mSqrtNThis; + f[i][0] -= prefactor * mK[0]; + f[i][1] -= prefactor * mK[1]; + f[i][2] -= prefactor * mK[2]; + } + } } // Forces in RESPA loop @@ -214,7 +209,7 @@ void FixRhok::post_force_respa( int inVFlag, int inILevel, int inILoop ) { if( inILevel == mNLevelsRESPA - 1 ) - post_force( inVFlag ); + post_force( inVFlag ); } // Forces in minimization loop @@ -228,9 +223,9 @@ FixRhok::min_post_force( int inVFlag ) double FixRhok::compute_scalar() { - double rhoK = sqrt( mRhoKGlobal[0]*mRhoKGlobal[0] - + mRhoKGlobal[1]*mRhoKGlobal[1] ); - + double rhoK = sqrt( mRhoKGlobal[0]*mRhoKGlobal[0] + + mRhoKGlobal[1]*mRhoKGlobal[1] ); + return 0.5 * mKappa * (rhoK - mRhoK0) * (rhoK - mRhoK0); } @@ -238,13 +233,13 @@ FixRhok::compute_scalar() double FixRhok::compute_vector( int inI ) { - if( inI == 0 ) - return mRhoKGlobal[0]; // Real part - else if( inI == 1 ) - return mRhoKGlobal[1]; // Imagniary part - else if( inI == 2 ) - return sqrt( mRhoKGlobal[0]*mRhoKGlobal[0] - + mRhoKGlobal[1]*mRhoKGlobal[1] ); - else - return 12345.0; + if( inI == 0 ) + return mRhoKGlobal[0]; // Real part + else if( inI == 1 ) + return mRhoKGlobal[1]; // Imagniary part + else if( inI == 2 ) + return sqrt( mRhoKGlobal[0]*mRhoKGlobal[0] + + mRhoKGlobal[1]*mRhoKGlobal[1] ); + else + return 12345.0; } diff --git a/src/USER-MISC/fix_rhok.h b/src/USER-MISC/fix_rhok.h index 145b3e7c1e..c950c08b1d 100644 --- a/src/USER-MISC/fix_rhok.h +++ b/src/USER-MISC/fix_rhok.h @@ -15,8 +15,8 @@ FixStyle(rhok,FixRhok) #else -#ifndef __FIX_RHOK__ -#define __FIX_RHOK__ +#ifndef LMP_FIX_RHOK_H +#define LMP_FIX_RHOK_H #include "fix.h" @@ -28,8 +28,8 @@ public: // Constructor: all the parameters to this fix specified in // the LAMMPS input get passed in FixRhok( LAMMPS* inLMP, int inArgc, char** inArgv ); - virtual ~FixRhok(); - + virtual ~FixRhok() {}; + // Methods that this fix implements // -------------------------------- @@ -51,23 +51,23 @@ public: // Compute the change in the potential energy induced by this fix double compute_scalar(); - - // Compute the ith component of the vector associated with this fix + + // Compute the ith component of the vector associated with this fix double compute_vector( int inI ); - + private: // RESPA boilerplate int mNLevelsRESPA; // Defining parameters for this umbrella - double mK[3], mKappa, mRhoK0; - - // Number of particles affected by the fix - int mNThis; - double mSqrtNThis; - - // Real and imaginary parts of rho_k := sum_i exp( - i k . r_i ) - double mRhoKLocal[2], mRhoKGlobal[2]; + double mK[3], mKappa, mRhoK0; + + // Number of particles affected by the fix + int mNThis; + double mSqrtNThis; + + // Real and imaginary parts of rho_k := sum_i exp( - i k . r_i ) + double mRhoKLocal[2], mRhoKGlobal[2]; }; } // namespace LAMMPS_NS From 1bb7af9ef91a2f74198f8251ac6a960cc9a86163 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Tue, 3 Oct 2017 10:45:08 -0400 Subject: [PATCH 43/44] integrate fix rhok into the LAMMPS source code management and documentation system --- doc/src/Section_commands.txt | 1 + doc/src/fixes.txt | 1 + doc/src/lammps.book | 1 + src/.gitignore | 2 ++ 4 files changed, 5 insertions(+) diff --git a/doc/src/Section_commands.txt b/doc/src/Section_commands.txt index 0d46a01424..4ee0d8b51b 100644 --- a/doc/src/Section_commands.txt +++ b/doc/src/Section_commands.txt @@ -728,6 +728,7 @@ package"_Section_start.html#start_3. "qtb"_fix_qtb.html, "reax/c/bonds"_fix_reax_bonds.html, "reax/c/species"_fix_reaxc_species.html, +"rhok"_fix_rhok.html, "rx"_fix_rx.html, "saed/vtk"_fix_saed_vtk.html, "shardlow"_fix_shardlow.html, diff --git a/doc/src/fixes.txt b/doc/src/fixes.txt index 7000a66c51..e6aad1c530 100644 --- a/doc/src/fixes.txt +++ b/doc/src/fixes.txt @@ -124,6 +124,7 @@ Fixes :h1 fix_reaxc_species fix_recenter fix_restrain + fix_rhok fix_rigid fix_rx fix_saed_vtk diff --git a/doc/src/lammps.book b/doc/src/lammps.book index 86dfe78af3..9f2e270452 100644 --- a/doc/src/lammps.book +++ b/doc/src/lammps.book @@ -253,6 +253,7 @@ fix_reaxc_species.html fix_recenter.html fix_restrain.html fix_rigid.html +fix_rhok.html fix_rx.html fix_saed_vtk.html fix_setforce.html diff --git a/src/.gitignore b/src/.gitignore index 1571065b72..504db2cdde 100644 --- a/src/.gitignore +++ b/src/.gitignore @@ -503,6 +503,8 @@ /fix_reaxc_bonds.h /fix_reaxc_species.cpp /fix_reaxc_species.h +/fix_rhok.cpp +/fix_rhok.h /fix_rigid.cpp /fix_rigid.h /fix_rigid_nh.cpp From 957263431aafcad7f5a8734358d99971d1036358 Mon Sep 17 00:00:00 2001 From: "Ulf R. Pedersen" Date: Wed, 4 Oct 2017 09:38:43 +0200 Subject: [PATCH 44/44] Ensure consistency with documentation --- examples/USER/misc/pinning/in.pinning | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/examples/USER/misc/pinning/in.pinning b/examples/USER/misc/pinning/in.pinning index 607335cf1a..0c220f480b 100644 --- a/examples/USER/misc/pinning/in.pinning +++ b/examples/USER/misc/pinning/in.pinning @@ -20,7 +20,7 @@ fix ensemble all npt temp 0.8 0.8 4.0 z 2.185 2.185 8.0 fix 100 all momentum 100 linear 1 1 1 # harmonic rho_k bias-field -# nx ny nz k a +# nx ny nz K a fix bias all rhok 16 0 0 4.0 26.00 # output U_bias rho_k_RE rho_k_IM |rho_k|