diff --git a/.github/workflows/build-and-test-with-lammps-release-branch.yml b/.github/workflows/build-and-test-with-lammps-release-branch.yml new file mode 100644 index 0000000..fdb48e9 --- /dev/null +++ b/.github/workflows/build-and-test-with-lammps-release-branch.yml @@ -0,0 +1,56 @@ +name: build-and-test-with-lammps-release-branch # workflow + +# Triggers for the jobs in this workflow +on: + push: + branches: [ "master", "develop", "*/*" ] + pull_request: + branches: [ "master", "develop" ] + +# Steps to build lammps with octp and run tests +jobs: + build-and-test: + + runs-on: ubuntu-22.04 # GitHub hosted runner (machine) on which lammps will be built + + steps: + - name: Checkout OCTP + uses: actions/checkout@v4 + with: + path: octp + + - name: Checkout LAMMPS release branch + uses: actions/checkout@v4 + with: + repository: lammps/lammps + ref: release + path: lammps + + - name: Copy OCTP files to LAMMPS src directory + run: | + cp --verbose octp/src/* lammps/src/ + + - name: Install OpenMPI in runner + run: | + sudo apt-get install openmpi-bin openmpi-common openmpi-doc libopenmpi-dev + echo "mpi version" + mpicc --showme:version + + - name: Build LAMMPS with OCTP plugin + run: | + cd lammps/src + make yes-asphere + make yes-body + make yes-class2 + make yes-dipole + make yes-granular + make yes-kspace + make yes-manybody + make yes-molecule + make yes-rigid + make yes-shock + make -j 8 mpi + if ! [ -f lmp_mpi ]; then + echo "::error:: lammps binary does not exist" + exit 1 + fi \ No newline at end of file diff --git a/.github/workflows/build-and-test-with-lammps-stable-release.yml b/.github/workflows/build-and-test-with-lammps-stable-release.yml new file mode 100644 index 0000000..792fd61 --- /dev/null +++ b/.github/workflows/build-and-test-with-lammps-stable-release.yml @@ -0,0 +1,61 @@ +name: build-and-test-with-lammps-stable-release # workflow + +# Triggers for the jobs in this workflow +on: + push: + branches: [ "master", "develop", "*/*" ] + pull_request: + branches: [ "master", "develop" ] + +# Variables +env: + # tag or commit ID of the lammps release with which we want to build and test + LAMMPS_RELEASE: 'stable_23Jun2022_update4' + +# Steps to build lammps with octp and run tests +jobs: + build-and-test: + + runs-on: ubuntu-22.04 # GitHub hosted runner (machine) on which lammps will be built + + steps: + - name: Checkout OCTP + uses: actions/checkout@v4 + with: + path: octp + + - name: Checkout LAMMPS latest stable_* tagged release + uses: actions/checkout@v4 + with: + repository: lammps/lammps + ref: ${{ env.LAMMPS_RELEASE }} + path: lammps + + - name: Copy OCTP files to LAMMPS src directory + run: | + cp --verbose octp/src/* lammps/src/ + + - name: Install OpenMPI in runner + run: | + sudo apt-get install openmpi-bin openmpi-common openmpi-doc libopenmpi-dev + echo "mpi version" + mpicc --showme:version + + - name: Build LAMMPS with OCTP plugin + run: | + cd lammps/src + make yes-asphere + make yes-body + make yes-class2 + make yes-dipole + make yes-granular + make yes-kspace + make yes-manybody + make yes-molecule + make yes-rigid + make yes-shock + make -j 8 mpi + if ! [ -f lmp_mpi ]; then + echo "::error:: lammps binary does not exist" + exit 1 + fi \ No newline at end of file diff --git a/.gitignore b/.gitignore index e915029..feb5ffe 100644 --- a/.gitignore +++ b/.gitignore @@ -2,4 +2,25 @@ # This .gitignore file was automatically created by Microsoft(R) Visual Studio. ################################################################################ -/.vs +# visual studio code settings +.vs*/ + +# Precompiled Headers +*.gch +*.pch + +# Compiled Object files +*.o +*.obj + +# Compiled Dynamic libraries +*.so +*.dll + +# Compiled Static libraries +*.a +*.lib + +# Executables +*.exe +*.out diff --git a/.vscode/settings.json b/.vscode/settings.json deleted file mode 100644 index d8cb326..0000000 --- a/.vscode/settings.json +++ /dev/null @@ -1,5 +0,0 @@ -{ - "files.associations": { - "string": "cpp" - } -} \ No newline at end of file diff --git a/Comment on thermal conductivity computation b/Comment on thermal conductivity computation deleted file mode 100644 index 0492587..0000000 --- a/Comment on thermal conductivity computation +++ /dev/null @@ -1,5 +0,0 @@ -The OCTP plugin uses the stress per atom computation of LAMMPS to yield the MSD of heat flux from which the thermal -conductivity can be computed. However, the OCTP user should be aware of the fact that the atomic stress approximation, as -implemented in LAMMPS, may produce erroneous results in the case of heat flux when applied to systems with many-body -interactions (e.g., angles, torsions, improper torsions). This results in incorrect thermal conductivity values. More -details can be found in the work by Surblys et al. (PHYSICAL REVIEW E 99, 051301, 2019 - DOI: 10.1103/PhysRevE.99.051301). diff --git a/compute_rdf_ext.h b/compute_rdf_ext.h deleted file mode 100644 index d37d592..0000000 --- a/compute_rdf_ext.h +++ /dev/null @@ -1,109 +0,0 @@ -/* ---------------------------------------------------------------------- - compute rdf/ext is a child class of "compute", developed based on two - classes of "compute msd" and "compute rdf", provided by LAMMPS. - This command is distributed under the GNU General Public License. -------------------------------------------------------------------------- */ - -/* -*- 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 COMPUTE_CLASS - -ComputeStyle(rdf/ext,ComputeRDFExtend) - -#else - -#ifndef LMP_COMPUTE_RDF__Extend_H -#define LMP_COMPUTE_RDF__Extend_H - -#include "compute.h" - -namespace LAMMPS_NS { - -class ComputeRDFExtend : public Compute { - public: - ComputeRDFExtend(class LAMMPS *, int, char **); - virtual ~ComputeRDFExtend(); - void init(); - virtual double compute_scalar(); // remember what this returns - - protected: - #define MAXGROUPS 32 - #define SQR(x) ((x)*(x)) - #define CUBE(x) ((x)*(x)*(x)) - // The parameters related to MPI and gathering the positions - int me; // This is for the definition of rank of processor in LAMMPS - int nprocs; // This is for the number of processors - int *tmprecvcnts; // Number of atoms per each processor [array] - int *recvcnts; // Number of atoms per each processor [array] - int *displs; // Displacement array for receiving the correct data set from each processor [array] - double *sendbuff; - double *recvbuff; - - // New static arrays or dynamic arrays - double **TmpPos; - int **Groups; - - // RDF - int binsize; - double r, r2, sf, pairs; - double rin , rout, vin, vout; - int whichbin; - double NumberOfParticles; - double ***Gg; - double ***Gglocal; - double boxsize; - double maxdist; - double hBox; - double Delta; - double Delta_1; - double Cube_deltaboxsize; - double Sqr_deltaboxsize; - double Ggt; - double dr[3]; - // RDF Correction - double ***PartSum; - double SphereVolFrac; - double gr_correction; - // The things in common with the other one - double timeinterval; - int count ; - int WriteFileEvery; - FILE *FilePtr; // The file handle for writing the calculated RDF - int samplerate; - int realatom; - char *filename; // filename - - // Finding corresponding groups to each atom - int tmpgroup[MAXGROUPS][3]; // group_id ; groupbits ; number_of_particles - int tmpnumgroup ; - int tmpfoundgroup ; - int groupatom1, groupatom2; -}; - -} - -#endif -#endif - -/* ERROR/WARNING messages: - -E: Illegal ... command - -Self-explanatory. - -E: Not a cubic simulation box for rdf/ext - -Self-explanatory. - -*/ diff --git a/fix_ordern.h b/fix_ordern.h deleted file mode 100644 index efb5766..0000000 --- a/fix_ordern.h +++ /dev/null @@ -1,212 +0,0 @@ -/* ---------------------------------------------------------------------- - fix ordern is a child class of "fix", developed based on two classes - of "fix ave/time" and "fix ave/correlate/long", provided by LAMMPS. - This command is distributed under the GNU General Public License. -------------------------------------------------------------------------- */ - -/* -*- 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 - -FixStyle(ordern,FixOrderN) - -#else - -#ifndef LMP_FIX_ORDERN_H -#define LMP_FIX_ORDERN_H - -#include "fix.h" - -#include - -namespace LAMMPS_NS { - -class FixOrderN : public Fix { - public: - FixOrderN(class LAMMPS *, int, char **); // constructor - ~FixOrderN(); // destructor - int setmask(); // the code is needed at the end of step - void init(); // initialization - void setup(int); - void end_of_step(); - void invoke_scalar(bigint); - - void write_restart(FILE *); - void restart(char *); - - private: - - int me; // the ID of the local processor - int mode; // which type of transport properties (defined as enum) - int nvalues; // the number of input values (it must be 1) - int nfreq; // the frequency of writing files to disk - int icompute; // the ID of compute for the tranpsort property - int nrows; // number of data passed to fix ordern via compute - int startstep; // the first timestep that it starts sampling - int tnb; // total number of blocks (different time scales) - int tnbe; // total number of block elements (elemets of similar timescales) - int vecsize; // number of vector elements constructed from nrows (exp. diff) - int sampsize; // number of vector elements for constructing sample arrays - int flag_Dxyz; // if components of diffusivities in x, y, & z should be written - int flag_TCconv; // if convective components of thermal conductivity should be written - int restart_continue; // checks to not sample twice in a timestep due to restarting - - bigint nvalid; // the next(current) timestep to call end_of_step() - bigint nvalid_last; // the previous timestep that end_of_step() was called - - double deltat; // timeinterval (nevery*dt) - double time; // correlation time (only for writting data) - double boltz; // Boltzmann constant - double nktv2p; // conversion factors - double volume; // volume of the simulation box - double coef; // conversion coefficient - - - long filepos1; // the location of header lines for file 1 - long filepos2; // the location of header lines for file 2 - - double *recdata; // data passed from compute to fix; - - char *filename1; // filename for self-diffusion, shear viscosity, and conductivity - char *filename2; // filename for onsager coefficients and bulk viscosity - - char *format; // the default format of writing data to files (%g) - char *format_user; // user-defined format of writing data to files - char *title; // optional text for the first line of the text - char *idcompute; // the ID of compute - - FILE *fp1; // output file 1 (self-diffusion, shear viscosity, thermal conductivity) - FILE *fp2; // output file 2 (onsager coefficients, bulk viscosity) - - - // General variables - int count; // how many cycles have been passed - int icount; // how many nevery, used for integration - // ORDER-N ALGORITHM variables - double **nsamp; // total number of samples - double ***oldint; // The lowest bound of the integral - double *rint; // running integral - int *nbe; // (BlockLength) number of active elements of blocks - int cnb; // (NumberOfBlock) current active number of blocks - int cnbe; // (CurrentBlockLength) current active number of elemets of the block - - - // DIFFUSIVITY vairables - double ***PosC_ii; - double ***PosC_iix; - double ***PosC_iiy; - double ***PosC_iiz; - double ****PosC_ij; - double ****PosC_ijx; - double ****PosC_ijy; - double ****PosC_ijz; - double ****PosCorrSum; - int **groupinfo; - int **atomingroup; - int atomgroup; // the ID of each group - int tngroup; // total # of defined groups (tngroup >= ngroup) - int ngroup; // total # of groups - int tnatom; // total # of atoms in system - int natom; // total # of atoms in groups - int sortID; // A virtual ID of the atom in a for loop (to tnatom) - int ID; // the real ID of an atom on all cores (to tnatom) - int atomID; // the ID of an atom inside a group (to natom) - int atommask; // the group mask of an atom - - // VISCOSITY/THERMCOND vairables - double ***samp; // samples of MSDs - double *data; // accumulated integrals - double sumP, numP, avgP; // hydrostatic pressure (sum/num = avg) - double *simpf0, *simpf1; // Simpson's rule of integration - double dist; // (integral) The integral we need to sample, i.e., "accint-oldint" - - // PRIVATE METHODS - bigint nextvalid(); - void integrate(); - void write_diffusivity(); - void write_viscosity(); - void write_thermcond(); - -}; - -} - -#endif -#endif - -/* ERROR/WARNING messages: - -E: Illegal ... command - -Self-explanatory. Check the input script syntax and compare to the -documentation for the command. You can use -echo screen as a -command-line option when running LAMMPS to see the offending line. - -E: Incorrect number of inputs for fix ordern command - -There should be only one compute ID. - -E: All components of the vector are required for fix ordern command - -No use of brackets [] - -E: Illegal fix ordern command: illegal numbers for fix ordern command - -Self-explanatory. - -E: Illegal fix ordern command: nevery is not a factor of nfreq - -In case of viscosity and thermal conductivity, 2*nevery should be a factor of nfreq - -E: Illegal fix ordern command: nevery is not a factor of start - -Self-explanatory. - -E: No compute ID for fix ordern command - -Self-explanatory. - -E: No global compute vector is computed for fix ordern command - -Self-explanatory. - -E: Input vector for fix ordern command has variable size - -Self-explanatory. - -E: Cannot open fix ordern files - -Self-explanatory. - -E: Error in writing file header for fix ordern command - -Self-explanatory. - -E: Invalid timestep reset for fix order - -Self-explanatory. - -E: restart is not possible with different parameters - -Self-explanatory. - -E: restart error: change in the number of group - -Self-explanatory. - -E: restart error: change in the number of atoms - -Self-explanatory. - -*/ diff --git a/readme.md b/readme.md index 1918a32..303372f 100644 --- a/readme.md +++ b/readme.md @@ -1,10 +1,57 @@ -On-the-fly Calculation of Transport Properties of Fluids using the order-n algorithm in Equilibrium Molecular Dynamics. More information on how this plugin works and how this plugin can be used with LAMMPS is thoroughly explained in the work by Jamali et al. J. Chem. Inf. Model. 2019, 59, 4, 1290-1294. If this package is used, please cite the following article: -Jamali SH, Wolff L, Becker TM, de Groen M, Ramdin M, Hartkamp R, Bardow A, Vlugt TJH, Moultos OA. OCTP: A Tool for On-the-Fly Calculation of Transport Properties of Fluids with the Order-n Algorithm in LAMMPS. Journal of Chemical Information and Modeling. 2019, 59, 4, 1290-1294 (DOI: 10.1021/acs.jcim.8b00939). +# OCTP Plugin for LAMMPS + +On-the-fly Calculation of Transport Properties of Fluids (OCTP) using the order-n algorithm in Equilibrium Molecular Dynamics. More information on how this plugin works and how this plugin can be used with [LAMMPS](https://www.lammps.org/) is thoroughly explained in the work by [Jamali et al. J. Chem. Inf. Model. 2019, 59, 4, 1290-1294](https://pubs.acs.org/doi/10.1021/acs.jcim.8b00939). This plugin consists of two new compute commands (position and rdf/ext) and one new fix command (ordern) for LAMMPS. For more information on the use of these new commands as well as a LAMMPS input file for the calcualation of transport properties for a water-methanol mixture, the user is referred to the work of Jamali et al. -"compute position" computes a global vector of 5*N components, where N is the total number of atoms in the simulation box. This vector includes the position (x,y,z), the atom ID, and the group mask of each atom. This compute can be used by "fix ordern" to calculate the self-diffusion and Maxwell-Stefan diffusion coefficients of a multicomponent mixture (see J. Chem. Theory Comput., 2018, 14 (5), pp 2667–2677). The components of the mixture should be pre-defined via the command group. The only restriction on the calculation of diffusion coefficients is that no atom should belong to more than 1 group. +## Compatibility + +| OCTP branch | LAMMPS | Status | +|-------------|--------|--------| +| [master](https://github.com/yiquintero/octp/tree/master) | [stable_23Jun2022_update4](https://github.com/lammps/lammps/releases/tag/stable_23Jun2022_update4) | [![Build and Test](https://github.com/yiquintero/octp/actions/workflows/build-and-test-with-lammps-stable-release.yml/badge.svg?branch=master)](https://github.com/yiquintero/octp/actions/workflows/build-and-test-with-lammps-stable-release.yml) +| [master](https://github.com/yiquintero/octp/tree/master) | [release branch](https://github.com/lammps/lammps/tree/release) | [![Build and Test](https://github.com/yiquintero/octp/actions/workflows/build-and-test-with-lammps-release-branch.yml/badge.svg?branch=master)](https://github.com/yiquintero/octp/actions/workflows/build-and-test-with-lammps-release-branch.yml) +| [develop](https://github.com/yiquintero/octp/tree/develop) | [stable_23Jun2022_update4](https://github.com/lammps/lammps/releases/tag/stable_23Jun2022_update4) | [![Build and Test](https://github.com/yiquintero/octp/actions/workflows/build-and-test-with-lammps-stable-release.yml/badge.svg?branch=develop)](https://github.com/yiquintero/octp/actions/workflows/build-and-test-with-lammps-stable-release.yml) +| [develop](https://github.com/yiquintero/octp/tree/develop) | [release branch](https://github.com/lammps/lammps/tree/release) | [![Build and Test](https://github.com/yiquintero/octp/actions/workflows/build-and-test-with-lammps-release-branch.yml/badge.svg?branch=develop)](https://github.com/yiquintero/octp/actions/workflows/build-and-test-with-lammps-release-branch.yml) + +## Installation + +It is not possible to use the OCTP plugin with pre-built LAMMPS binaries or installations. To use any of the OCTP commands, the source code from LAMMPS and the plugin need to be compiled together - i.e. LAMMPS needs to be built from source. On Linux: + +```bash +# Clone the OCTP plugin repo and checkout the master branch +git clone https://github.com/yiquintero/octp.git +cd octp +git checkout master +cd .. + +# Clone the LAMMPS repo and checkout the latest stable release +# A list of releases can be found in https://github.com/lammps/lammps/releases +git clone https://github.com/lammps/lammps.git +cd lammps +git checkout stable_23Jun2022_update4 + +# Copy the plugin's source code to the LAMMPS source directory +cp ../octp/src/* src/ + +# Build LAMMPS +cd src/ +make mpi +``` + +More information on how to build LAMMPS from source can be found in the official [LAMMPS documentation](https://docs.lammps.org/Build.html). + +## OCTP Commands + +- **compute position** computes a global vector of 5*N components, where N is the total number of atoms in the simulation box. This vector includes the position (x,y,z), the atom ID, and the group mask of each atom. This compute can be used by "fix ordern" to calculate the self-diffusion and Maxwell-Stefan diffusion coefficients of a multicomponent mixture (see J. Chem. Theory Comput., 2018, 14 (5), pp 2667–2677). The components of the mixture should be pre-defined via the command group. The only restriction on the calculation of diffusion coefficients is that no atom should belong to more than 1 group. + +- **compute rdf/ext** is an extension of "compute rdf", which is already available in LAMMPS. This compute calculates the RDF up to sqrt(2)/2 of the box length. In addition, RDFs are corrected due to the finite-size of the simulation box are computed. For the finite-size correction, the method of van der Vegt et al. (see J. Chem. Theory Comput., 2013, 9 (3), pp 1347–1355 and J. Phys. Chem. B, 2018, 122 (21), pp 5515–5526) has been implemented. Similar to "compute position", the mixture should be pre-defined via the command group. The restriction on the calculation of RDFs is that no atom should belong to more than 1 group. This compute command directly generates output files; however, in order to invoke it during the simulation, a fix command (e.g., fix ave/time) should be specified. + +- **fix ordern** uses the modified order-n algorithm by Dubbeldam et al. (Mol. Simul., 2009, 35 (12–13), pp 1084–1097) to compute mean-square displacement (MSD) for the self-diffusion coefficient, Maxwell-Stefan diffusion coefficient, shear viscosity, bulk viscosity, and thermal conductivity. Each transport property can be calculated from the slope of the time-MSD plot, where a linear relation between MSD and time is found. As the input to this fix command, the ID of "compute position", "compute pressure", and "compute heat/flux" should be assigned. This fix generates output files containing MSD as a function of time. The definition of reported MSD can be found at pages S14, S17, and S19 of the Supporting Information of the article "OCTP: A Tool for On-the-fly Calculation of Transport Properties of Fluids with the Order-n Algorithm in LAMMPS" (J. Chem. Inf. Model., 2019). It should be noted that the reported MSD must be divided by a factor explicitly mentioned at the top of each output file. + +The thermal conductivity computation may produce erroneous results in the case of heat flux when applied to systems with many-body interactions. Read the [Wiki page](https://github.com/yiquintero/octp/wiki/Thermal-Conductivity-Computation) for more information. + +## Citation -"compute rdf/ext" is an extension of "compute rdf", which is already available in LAMMPS. This compute calculates the RDF up to sqrt(2)/2 of the box length. In addition, RDFs are corrected due to the finite-size of the simulation box are computed. For the finite-size correction, the method of van der Vegt et al. (see J. Chem. Theory Comput., 2013, 9 (3), pp 1347–1355 and J. Phys. Chem. B, 2018, 122 (21), pp 5515–5526) has been implemented. Similar to "compute position", the mixture should be pre-defined via the command group. The restriction on the calculation of RDFs is that no atom should belong to more than 1 group. This compute command directly generates output files; however, in order to invoke it during the simulation, a fix command (e.g., fix ave/time) should be specified. +If this package is used, please cite the following article: -"fix ordern" uses the modified order-n algorithm by Dubbeldam et al. (Mol. Simul., 2009, 35 (12–13), pp 1084–1097) to compute mean-square displacement (MSD) for the self-diffusion coefficient, Maxwell-Stefan diffusion coefficient, shear viscosity, bulk viscosity, and thermal conductivity. Each transport property can be calculated from the slope of the time-MSD plot, where a linear relation between MSD and time is found. As the input to this fix command, the ID of "compute position", "compute pressure", and "compute heat/flux" should be assigned. This fix generates output files containing MSD as a function of time. The definition of reported MSD can be found at pages S14, S17, and S19 of the Supporting Information of the article "OCTP: A Tool for On-the-fly Calculation of Transport Properties of Fluids with the Order-n Algorithm in LAMMPS" (J. Chem. Inf. Model., 2019). It should be noted that the reported MSD must be divided by a factor explicitly mentioned at the top of each output file. \ No newline at end of file +Jamali SH, Wolff L, Becker TM, de Groen M, Ramdin M, Hartkamp R, Bardow A, Vlugt TJH, Moultos OA. OCTP: A Tool for On-the-Fly Calculation of Transport Properties of Fluids with the Order-n Algorithm in LAMMPS. Journal of Chemical Information and Modeling. 2019, 59, 4, 1290-1294 (DOI: 10.1021/acs.jcim.8b00939). \ No newline at end of file diff --git a/compute_position.cpp b/src/compute_position.cpp similarity index 87% rename from compute_position.cpp rename to src/compute_position.cpp index 14770f7..b80d037 100755 --- a/compute_position.cpp +++ b/src/compute_position.cpp @@ -27,10 +27,12 @@ using namespace LAMMPS_NS; /* ---------------------------------------------------------------------- */ -ComputePosition::ComputePosition(LAMMPS *lmp, int narg, char **arg) : - Compute(lmp, narg, arg) +ComputePosition::ComputePosition(LAMMPS *lmp, int narg, char **arg) + : Compute(lmp, narg, arg) { - if (narg != 3) error->all(FLERR,"Illegal compute position command"); + if (narg != 3) { + error->all(FLERR,"Illegal compute position command"); + } vector_flag = 1; size_vector = atom->natoms*5; @@ -66,8 +68,7 @@ ComputePosition::~ComputePosition() void ComputePosition::init() { - for (int ii = 0; ii < size_vector; ii++) - { + for (int ii = 0; ii < size_vector; ii++) { vector[ii]=-1.0; } } @@ -90,11 +91,9 @@ void ComputePosition::compute_vector() int xbox,ybox,zbox; double xtmp, ytmp, ztmp; - if (domain->triclinic == 0) - { - for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit) - { + if (domain->triclinic == 0) { + for (int i = 0; i < nlocal; i++) { + if (mask[i] & groupbit) { xbox = (image[i] & IMGMASK) - IMGMAX; ybox = (image[i] >> IMGBITS & IMGMASK) - IMGMAX; zbox = (image[i] >> IMG2BITS) - IMGMAX; @@ -107,10 +106,11 @@ void ComputePosition::compute_vector() sendbuff[5*i+3] = (double) ((atom->tag[i]) - 1); sendbuff[5*i+4] = (double) mask[i]; } - } else { - for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit) - { + } + } + else { + for (int i = 0; i < nlocal; i++) { + if (mask[i] & groupbit) { xbox = (image[i] & IMGMASK) - IMGMAX; ybox = (image[i] >> IMGBITS & IMGMASK) - IMGMAX; zbox = (image[i] >> IMG2BITS) - IMGMAX; @@ -123,20 +123,24 @@ void ComputePosition::compute_vector() sendbuff[5*i+3] = (double) ((atom->tag[i]) - 1); sendbuff[5*i+4] = (double) mask[i]; } + } } + // Sending the position of all atoms from all cores (0, 1, ..., n) to all cores - for (int ii = 0; ii < nprocs; ii++) + for (int ii = 0; ii < nprocs; ii++) { tmprecvcnts[ii] = (ii == me) ? 5*nlocal : 0 ; + } MPI_Allreduce(tmprecvcnts, recvcnts, nprocs, MPI_INT , MPI_SUM, world); - for (int ii = 0; ii < nprocs; ii++) - { + for (int ii = 0; ii < nprocs; ii++) { displs[ii] = 0; } - for (int ii = 1; ii < nprocs; ii++) - for (int jj = ii; jj < nprocs ; jj++) - { + + for (int ii = 1; ii < nprocs; ii++) { + for (int jj = ii; jj < nprocs ; jj++) { displs[jj] += recvcnts[ii-1]; } + } + MPI_Allgatherv(sendbuff,5*nlocal, MPI_DOUBLE, vector, recvcnts, displs, MPI_DOUBLE, world); } diff --git a/compute_position.h b/src/compute_position.h similarity index 65% rename from compute_position.h rename to src/compute_position.h index f8a6501..55750fa 100755 --- a/compute_position.h +++ b/src/compute_position.h @@ -18,39 +18,38 @@ ------------------------------------------------------------------------- */ #ifdef COMPUTE_CLASS - ComputeStyle(position,ComputePosition) - #else - -#ifndef LMP_COMPUTE_POSITION_H +#ifndef LMP_COMPUTE_POSITION_H #define LMP_COMPUTE_POSITION_H #include "compute.h" -namespace LAMMPS_NS { - -class ComputePosition : public Compute { - public: - ComputePosition(class LAMMPS *, int, char **); - virtual ~ComputePosition(); - void init(); - virtual void compute_vector(); - - protected: - int me; // This is for the definition of rank of processor in LAMMPS - int nprocs; // This is for the number of processors - int *tmprecvcnts; // Number of atoms per each processor - int *recvcnts; // Number of atoms per each processor - int *displs; // Displacement array for receiving the correct data set from each processor - double *sendbuff; // The sending array from all processors - -}; +namespace LAMMPS_NS +{ + class ComputePosition : public Compute + { + public: + ComputePosition(class LAMMPS *, int, char **); + virtual ~ComputePosition(); + + virtual void compute_vector(); + void init(); + + protected: + int me; // this is for the definition of rank of processor in LAMMPS + int nprocs; // number of processors + int *tmprecvcnts; // number of atoms per processor + int *recvcnts; // number of atoms per processor + int *displs; // displacement array for receiving the correct data set from each processor + double *sendbuff; // sending array from all processors + }; } #endif #endif + /* ERROR/WARNING messages: E: Illegal ... command diff --git a/compute_rdf_ext.cpp b/src/compute_rdf_ext.cpp similarity index 75% rename from compute_rdf_ext.cpp rename to src/compute_rdf_ext.cpp index 5319342..b88b2ee 100644 --- a/compute_rdf_ext.cpp +++ b/src/compute_rdf_ext.cpp @@ -29,10 +29,12 @@ using namespace LAMMPS_NS; /* ---------------------------------------------------------------------- */ -ComputeRDFExtend::ComputeRDFExtend(LAMMPS *lmp, int narg, char **arg) : - Compute(lmp, narg, arg) +ComputeRDFExtend::ComputeRDFExtend(LAMMPS *lmp, int narg, char **arg) + : Compute(lmp, narg, arg) { - if (narg < 3) error->all(FLERR,"Illegal compute rdf/ext command"); + if (narg < 3) { + error->all(FLERR,"Illegal compute rdf/ext command"); + } scalar_flag = 1; //size_scalar = 4; @@ -48,23 +50,33 @@ ComputeRDFExtend::ComputeRDFExtend(LAMMPS *lmp, int narg, char **arg) : int iarg = 3; while (iarg < narg) { if (strcmp(arg[iarg],"Nbin") == 0) { - if (iarg+2 > narg) error->all(FLERR,"Illegal compute rdf/ext command"); + if (iarg+2 > narg) { + error->all(FLERR,"Illegal compute rdf/ext command"); + } binsize = utils::inumeric(FLERR,arg[iarg+1],false,lmp); iarg += 2; - } else if (strcmp(arg[iarg],"Nwrite") == 0) { - if (iarg+2 > narg) error->all(FLERR,"Illegal compute rdf/ext command"); + } + else if (strcmp(arg[iarg],"Nwrite") == 0) { + if (iarg+2 > narg) { + error->all(FLERR,"Illegal compute rdf/ext command"); + } WriteFileEvery = utils::inumeric(FLERR,arg[iarg+1],false,lmp); iarg += 2; - } else if (strcmp(arg[iarg],"file") == 0) { - if (iarg+2 > narg) error->all(FLERR,"Illegal fix ordern command"); + } + else if (strcmp(arg[iarg],"file") == 0) { + if (iarg+2 > narg) { + error->all(FLERR,"Illegal fix ordern command"); + } delete [] filename; filename = new char[strlen(arg[iarg+1])+1]; strcpy(filename,arg[iarg+1]); iarg += 2; - } else error->all(FLERR,"Illegal compute rdf/ext command"); + } + else { + error->all(FLERR,"Illegal compute rdf/ext command"); + } } - MPI_Comm_rank(world, &me); // Assigning the rank of a molecule for each core MPI_Comm_size(world, &nprocs); // The number of processors tmprecvcnts = new int[nprocs]; @@ -84,7 +96,6 @@ ComputeRDFExtend::ComputeRDFExtend(LAMMPS *lmp, int narg, char **arg) : ComputeRDFExtend::~ComputeRDFExtend() { - delete [] tmprecvcnts; delete [] recvcnts; delete [] displs; @@ -101,21 +112,18 @@ ComputeRDFExtend::~ComputeRDFExtend() /* ---------------------------------------------------------------------- */ void ComputeRDFExtend::init() -{ - +{ tmpnumgroup = 0; count = 0; - for (int ii=0; iixprd; hBox = 0.5 * boxsize ; @@ -124,9 +132,6 @@ void ComputeRDFExtend::init() Delta_1 = 1.0 / Delta; Cube_deltaboxsize = CUBE(Delta/boxsize); Sqr_deltaboxsize = SQR(Delta/boxsize); - - - } /* ---------------------------------------------------------------------- */ @@ -151,7 +156,7 @@ double ComputeRDFExtend::compute_scalar() int *currentgroupbit = group->bitmask; // The bitmask of a group if (domain->triclinic == 0) { - for (int i = 0; i < nlocal; i++) + for (int i = 0; i < nlocal; i++) { if (mask[i] & groupbit) { sendbuff[5*i] = x[i][0]; sendbuff[5*i+1] = x[i][1]; @@ -159,37 +164,47 @@ double ComputeRDFExtend::compute_scalar() sendbuff[5*i+3] = (double) ((atom->tag[i]) - 1); sendbuff[5*i+4] = (double) mask[i]; } - } else error->all(FLERR,"Not a cubic simulation box for rdf/ext"); - + } + } + else { + error->all(FLERR,"Not a cubic simulation box for rdf/ext"); + } + // Sending data from all cores (0, 1, ..., n) to all cores - for (int ii = 0; ii < nprocs; ii++) + for (int ii = 0; ii < nprocs; ii++) { tmprecvcnts[ii] = (ii == me) ? 5*nlocal : 0 ; + } MPI_Allreduce(tmprecvcnts, recvcnts, nprocs, MPI_INT , MPI_SUM, world); - for (int ii = 0; ii < nprocs; ii++) - for (int jj = ii; jj < nprocs ; jj++) - { - if (ii == 0) displs[jj] = 0; - else displs[jj] += recvcnts[ii-1]; + + for (int ii = 0; ii < nprocs; ii++) { + for (int jj = ii; jj < nprocs ; jj++) { + if (ii == 0) { + displs[jj] = 0; + } + else { + displs[jj] += recvcnts[ii-1]; + } } + } MPI_Allgatherv(sendbuff,5*nlocal, MPI_DOUBLE, recvbuff, recvcnts, displs, MPI_DOUBLE, world); + // Definining the group number of each molecules - for (int ii = 0; ii < natom; ii++) + for (int ii = 0; ii < natom; ii++) { realatom = (int) (recvbuff[5*ii+3]+0.1); TmpPos[realatom][0] = recvbuff[5*ii]; TmpPos[realatom][1] = recvbuff[5*ii+1]; TmpPos[realatom][2] = recvbuff[5*ii+2]; + // Finding corresponding groups to each atom at the first time - if (count == 0) // only run during the first time step + if (count == 0) // only run during the first time step { TmpPos[realatom][3] = recvbuff[5*ii+4]+0.1; tmpfoundgroup = 0; - if (tmpnumgroup > 0 ) // First try to find the group number if the group is available + if (tmpnumgroup > 0 ) // First try to find the group number if the group is available { - for (int jj = 1; jj <= tmpnumgroup; jj++) - { - if (((int) TmpPos[realatom][3]) & tmpgroup[jj][1]) - { + for (int jj = 1; jj <= tmpnumgroup; jj++) { + if (((int) TmpPos[realatom][3]) & tmpgroup[jj][1]) { Groups[realatom][0] = jj; Groups[realatom][1] = tmpgroup[jj][1]; Groups[realatom][2] = tmpgroup[jj][0]; @@ -201,10 +216,8 @@ double ComputeRDFExtend::compute_scalar() } if (tmpfoundgroup == 0) // If couldn't find any, tries to find the new group { - for (int jj = 1; jj < MAXGROUPS; jj++) - { - if ( ((int) TmpPos[realatom][3]) & group->bitmask[jj]) - { + for (int jj = 1; jj < MAXGROUPS; jj++) { + if ( ((int) TmpPos[realatom][3]) & group->bitmask[jj]) { tmpnumgroup++; tmpfoundgroup = 1; tmpgroup[tmpnumgroup][0] = jj; @@ -216,8 +229,8 @@ double ComputeRDFExtend::compute_scalar() break; } } - if (tmpfoundgroup == 0) // Finally, this atom does not belong to the available groups - { + + if (tmpfoundgroup == 0) { // Finally, this atom does not belong to the available groups Groups[realatom][0] = -1; Groups[realatom][1] = -1; Groups[realatom][2] = -1; @@ -225,26 +238,29 @@ double ComputeRDFExtend::compute_scalar() } } } + + // CALCULATING RDF - if (count == 0) - { + if (count == 0) { samplerate = (update->ntimestep); - } else if (count == 1) { + } + else if (count == 1) { samplerate = (update->ntimestep) - samplerate; timeinterval = (double) (samplerate)*(update->dt); } + Ggt+=1.0; - for(int i=0;itag[i]) - 1; groupatom1 = Groups[realatom][0]; - if (groupatom1>0) + if (groupatom1 > 0) { - for(int j=realatom+1;j0) ) - { + if (groupatom2 > 0) { dr[0] = TmpPos[realatom][0] - TmpPos[j][0]; dr[1] = TmpPos[realatom][1] - TmpPos[j][1]; dr[2] = TmpPos[realatom][2] - TmpPos[j][2]; @@ -255,8 +271,7 @@ double ComputeRDFExtend::compute_scalar() r2 = sqrt( dr[0]*dr[0] + dr[1]*dr[1] + dr[2]*dr[2] ); - if (r2 < maxdist) - { + if (r2 < maxdist) { whichbin = (int) (r2*Delta_1); Gglocal[whichbin][groupatom1][groupatom2] += 1; } @@ -264,57 +279,54 @@ double ComputeRDFExtend::compute_scalar() } } } + + // Writting in the file - if ( ((count > 0) && (count%WriteFileEvery == 0)) ) + if (count > 0 && count%WriteFileEvery == 0) { int sizeofGg = binsize * MAXGROUPS * MAXGROUPS ; // Everything should be passed // pass the address of the first element of a multidimensional array MPI_Reduce(&Gglocal[0][0][0] , &Gg[0][0][0] , sizeofGg , MPI_DOUBLE , MPI_SUM , 0 , world); - if ( me == 0) + + if (me == 0) { // Start RDF postprocessing (Only the first binsize data should be processed) - for ( int i = 0; i0) - { + + for (int i = 0; i0) { PartSum[i][j][k] = PartSum[i-1][j][k] + Gg[i][j][k]; - } else { + } + else { PartSum[i][j][k] = Gg[i][j][k]; } } } } + // Open a file and write the header FilePtr = fopen(filename,"w"); fprintf(FilePtr , "# Radius \t"); - for (int j = 1 ; j <= tmpnumgroup ; j++) - { - for (int k = j ; k <= tmpnumgroup ; k++) - { + for (int j = 1; j <= tmpnumgroup; j++) { + for (int k = j; k <= tmpnumgroup; k++) { // PRINT g(r), g(r)*correction) - fprintf(FilePtr , "finiterdf__%s_%s\t",group->names[tmpgroup[j][0]],group->names[tmpgroup[k][0]]); - fprintf(FilePtr , "rdf__%s_%s\t",group->names[tmpgroup[j][0]],group->names[tmpgroup[k][0]]); + fprintf(FilePtr, "finiterdf__%s_%s\t", group->names[tmpgroup[j][0]], group->names[tmpgroup[k][0]]); + fprintf(FilePtr, "rdf__%s_%s\t", group->names[tmpgroup[j][0]], group->names[tmpgroup[k][0]]); } } - fprintf(FilePtr , "\n"); + fprintf(FilePtr, "\n"); // Write all the bins - for ( int i = 0 ; i < binsize ; i++ ) + for (int i = 0; i < binsize; i++) { r = (i+0.5)*Delta; rin = 1.0*i*Delta; @@ -328,25 +340,23 @@ double ComputeRDFExtend::compute_scalar() fprintf(FilePtr,"%g ",r); NumberOfParticles = 0; - for (int j = 1 ; j <= tmpnumgroup ; j++) - { - for (int k = j ; k <= tmpnumgroup ; k++) - { - if ( j != k ) - { + for (int j = 1; j <= tmpnumgroup; j++) { + for (int k = j; k <= tmpnumgroup; k++) { + if ( j != k ) { pairs = 1.0 * sf * tmpgroup[j][2] * tmpgroup[k][2] ; gr_correction = tmpgroup[k][2] * (1.0 - SphereVolFrac) / (tmpgroup[k][2] - PartSum[i][j][k]/(tmpgroup[j][2]*Ggt)); - } else { + } + else { pairs = 0.5 * sf * tmpgroup[j][2] * tmpgroup[k][2] ; gr_correction = tmpgroup[k][2] * (1.0 - SphereVolFrac) / (tmpgroup[k][2] - 2.0*PartSum[i][j][k]/(tmpgroup[j][2]*Ggt)-1.0); } // PRINT g(r), g(r)*correction) - fprintf(FilePtr , "%g\t",Gg[i][j][k]/pairs); - fprintf(FilePtr , "%g\t",gr_correction*Gg[i][j][k]/pairs); + fprintf(FilePtr, "%g\t", Gg[i][j][k]/pairs); + fprintf(FilePtr, "%g\t", gr_correction*Gg[i][j][k]/pairs); } NumberOfParticles += tmpgroup[j][2] ; } - fprintf(FilePtr , "\n"); + fprintf(FilePtr, "\n"); } fclose(FilePtr); } diff --git a/src/compute_rdf_ext.h b/src/compute_rdf_ext.h new file mode 100644 index 0000000..37f4c2d --- /dev/null +++ b/src/compute_rdf_ext.h @@ -0,0 +1,111 @@ +/* ---------------------------------------------------------------------- + compute rdf/ext is a child class of "compute", developed based on two + classes of "compute msd" and "compute rdf", provided by LAMMPS. + This command is distributed under the GNU General Public License. +------------------------------------------------------------------------- */ + +/* -*- 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 COMPUTE_CLASS +ComputeStyle(rdf/ext,ComputeRDFExtend) +#else +#ifndef LMP_COMPUTE_RDF__Extend_H +#define LMP_COMPUTE_RDF__Extend_H + +#include "compute.h" + +namespace LAMMPS_NS +{ + class ComputeRDFExtend : public Compute + { + public: + ComputeRDFExtend(class LAMMPS *, int, char **); + virtual ~ComputeRDFExtend(); + + virtual double compute_scalar(); // remember what this returns + void init(); + + protected: + #define MAXGROUPS 32 + #define SQR(x) ((x)*(x)) + #define CUBE(x) ((x)*(x)*(x)) + + // parameters related to MPI and gathering the positions + int me; // This is for the definition of rank of processor in LAMMPS + int nprocs; // number of processors + int *tmprecvcnts; // number of atoms per each processor [array] + int *recvcnts; // number of atoms per each processor [array] + int *displs; // displacement array for receiving the correct data set from each processor [array] + double *sendbuff; + double *recvbuff; + + // New static arrays or dynamic arrays + double **TmpPos; + int **Groups; + + // RDF + int binsize; + double r, r2, sf, pairs; + double rin , rout, vin, vout; + int whichbin; + double NumberOfParticles; + double ***Gg; + double ***Gglocal; + double boxsize; + double maxdist; + double hBox; + double Delta; + double Delta_1; + double Cube_deltaboxsize; + double Sqr_deltaboxsize; + double Ggt; + double dr[3]; + + // RDF Correction + double ***PartSum; + double SphereVolFrac; + double gr_correction; + + // The things in common with the other one + double timeinterval; + int count ; + int WriteFileEvery; + FILE *FilePtr; // file handle for writing the calculated RDF + int samplerate; + int realatom; + char *filename; // filename + + // Finding corresponding groups to each atom + int tmpgroup[MAXGROUPS][3]; // group_id ; groupbits ; number_of_particles + int tmpnumgroup ; + int tmpfoundgroup ; + int groupatom1, groupatom2; + }; +} + +#endif +#endif + + +/* ERROR/WARNING messages: + +E: Illegal ... command + +Self-explanatory. + +E: Not a cubic simulation box for rdf/ext + +Self-explanatory. + +*/ diff --git a/fix_ordern.cpp b/src/fix_ordern.cpp similarity index 100% rename from fix_ordern.cpp rename to src/fix_ordern.cpp diff --git a/src/fix_ordern.h b/src/fix_ordern.h new file mode 100644 index 0000000..79af39f --- /dev/null +++ b/src/fix_ordern.h @@ -0,0 +1,209 @@ +/* ---------------------------------------------------------------------- + fix ordern is a child class of "fix", developed based on two classes + of "fix ave/time" and "fix ave/correlate/long", provided by LAMMPS. + This command is distributed under the GNU General Public License. +------------------------------------------------------------------------- */ + +/* -*- 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 +FixStyle(ordern,FixOrderN) +#else +#ifndef LMP_FIX_ORDERN_H +#define LMP_FIX_ORDERN_H + +#include "fix.h" +#include + +namespace LAMMPS_NS +{ + class FixOrderN : public Fix + { + public: + FixOrderN(class LAMMPS *, int, char **); // constructor + ~FixOrderN(); // destructor + + int setmask(); // the code is needed at the end of step + void init(); // initialization + void setup(int); + void end_of_step(); + void invoke_scalar(bigint); + + void write_restart(FILE *); + void restart(char *); + + private: + + int me; // the ID of the local processor + int mode; // which type of transport properties (defined as enum) + int nvalues; // the number of input values (it must be 1) + int nfreq; // the frequency of writing files to disk + int icompute; // the ID of compute for the tranpsort property + int nrows; // number of data passed to fix ordern via compute + int startstep; // the first timestep that it starts sampling + int tnb; // total number of blocks (different time scales) + int tnbe; // total number of block elements (elemets of similar timescales) + int vecsize; // number of vector elements constructed from nrows (exp. diff) + int sampsize; // number of vector elements for constructing sample arrays + int flag_Dxyz; // if components of diffusivities in x, y, & z should be written + int flag_TCconv; // if convective components of thermal conductivity should be written + int restart_continue; // checks to not sample twice in a timestep due to restarting + + bigint nvalid; // the next(current) timestep to call end_of_step() + bigint nvalid_last; // the previous timestep that end_of_step() was called + + double deltat; // timeinterval (nevery*dt) + double time; // correlation time (only for writting data) + double boltz; // Boltzmann constant + double nktv2p; // conversion factors + double volume; // volume of the simulation box + double coef; // conversion coefficient + + long filepos1; // the location of header lines for file 1 + long filepos2; // the location of header lines for file 2 + + double *recdata; // data passed from compute to fix; + + char *filename1; // filename for self-diffusion, shear viscosity, and conductivity + char *filename2; // filename for onsager coefficients and bulk viscosity + + char *format; // the default format of writing data to files (%g) + char *format_user; // user-defined format of writing data to files + char *title; // optional text for the first line of the text + char *idcompute; // the ID of compute + + FILE *fp1; // output file 1 (self-diffusion, shear viscosity, thermal conductivity) + FILE *fp2; // output file 2 (onsager coefficients, bulk viscosity) + + // General variables + int count; // how many cycles have been passed + int icount; // how many nevery, used for integration + + // ORDER-N ALGORITHM variables + double **nsamp; // total number of samples + double ***oldint; // The lowest bound of the integral + double *rint; // running integral + int *nbe; // (BlockLength) number of active elements of blocks + int cnb; // (NumberOfBlock) current active number of blocks + int cnbe; // (CurrentBlockLength) current active number of elemets of the block + + // DIFFUSIVITY vairables + double ***PosC_ii; + double ***PosC_iix; + double ***PosC_iiy; + double ***PosC_iiz; + double ****PosC_ij; + double ****PosC_ijx; + double ****PosC_ijy; + double ****PosC_ijz; + double ****PosCorrSum; + int **groupinfo; + int **atomingroup; + int atomgroup; // the ID of each group + int tngroup; // total # of defined groups (tngroup >= ngroup) + int ngroup; // total # of groups + int tnatom; // total # of atoms in system + int natom; // total # of atoms in groups + int sortID; // A virtual ID of the atom in a for loop (to tnatom) + int ID; // the real ID of an atom on all cores (to tnatom) + int atomID; // the ID of an atom inside a group (to natom) + int atommask; // the group mask of an atom + + // VISCOSITY/THERMCOND vairables + double ***samp; // samples of MSDs + double *data; // accumulated integrals + double sumP, numP, avgP; // hydrostatic pressure (sum/num = avg) + double *simpf0, *simpf1; // Simpson's rule of integration + double dist; // (integral) The integral we need to sample, i.e., "accint-oldint" + + // PRIVATE METHODS + bigint nextvalid(); + void integrate(); + void write_diffusivity(); + void write_viscosity(); + void write_thermcond(); + + }; + +} + +#endif +#endif + + +/* ERROR/WARNING messages: + +E: Illegal ... command + +Self-explanatory. Check the input script syntax and compare to the +documentation for the command. You can use -echo screen as a +command-line option when running LAMMPS to see the offending line. + +E: Incorrect number of inputs for fix ordern command + +There should be only one compute ID. + +E: All components of the vector are required for fix ordern command + +No use of brackets [] + +E: Illegal fix ordern command: illegal numbers for fix ordern command + +Self-explanatory. + +E: Illegal fix ordern command: nevery is not a factor of nfreq + +In case of viscosity and thermal conductivity, 2*nevery should be a factor of nfreq + +E: Illegal fix ordern command: nevery is not a factor of start + +Self-explanatory. + +E: No compute ID for fix ordern command + +Self-explanatory. + +E: No global compute vector is computed for fix ordern command + +Self-explanatory. + +E: Input vector for fix ordern command has variable size + +Self-explanatory. + +E: Cannot open fix ordern files + +Self-explanatory. + +E: Error in writing file header for fix ordern command + +Self-explanatory. + +E: Invalid timestep reset for fix order + +Self-explanatory. + +E: restart is not possible with different parameters + +Self-explanatory. + +E: restart error: change in the number of group + +Self-explanatory. + +E: restart error: change in the number of atoms + +Self-explanatory. + +*/